---
title: "Benchmarks"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Benchmarks}
  %\VignetteEncoding{UTF-8}
  %\VignetteEngine{knitr::rmarkdown}
editor_options: 
  chunk_output_type: console
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```


```{r setup}
library(earthtide)
library(bench)


eval_chunks <- TRUE # may not want to run on CRAN because of threads and running time
```

This vignette describes a few ways to speed up the computation of Earth tides and in some cases reduce memory consumption.  The examples below are kept small to minimize computation time for CRAN, but the methods can scale to larger problems.


The following techniques are presented below:
- Irregular time steps
- Change wave catalog
- Change wave amplitude cutoff
- Change how often astronomical parameters are updated
- Use parallel computation
- Interpolations


# Irregular time steps

Some times you may not need to predict at regular time steps. Irregular time steps are allowed, however, the \code{astro_update} parameter should be set to 1L if you are not using a regular time series.

```{r irregular, eval = eval_chunks}
set.seed(123)
tms <- as.POSIXct("1990-01-01", tz = "UTC") + 0:(900)
indices <- sort(sample(0:900, 100, replace = FALSE)) 

wave_groups <- data.frame(start = 0, end = 8)

check_fun <- function(target, current) (all.equal(target, current, check.attributes = FALSE))

bench::mark(
  et <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups
  )[indices, ],
  et_irregular <- calc_earthtide(
    utc = tms[indices],
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups
  ), check = check_fun, iterations = 1
)

```


# Catalog

Using a catalog with fewer waves will be faster. Here we compare ksm04 and hw95s.

```{r catalog, eval = eval_chunks}

tms <- as.POSIXct("1990-01-01", tz = "UTC") + 0:(900)

wave_groups <- data.frame(start = 0, end = 8)

bench::mark(
  et <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups
  ),
  et_catalog <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "hw95s",
    wave_groups = wave_groups
  ), check = FALSE, iterations = 1
)

```


# cutoff parameter

Increasing the cutoff will decrease the number of waves and thus the speed
increases. Results will not be the same.

```{r cutoff, eval = eval_chunks}

tms <- as.POSIXct("1990-01-01", tz = "UTC") + 0:(1800)

wave_groups <- data.frame(start = 0, end = 8)

bench::mark(
  et <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups
  ),
  et_cutoff <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-5,
    catalog = "ksm04",
    wave_groups = wave_groups
  ), check = FALSE, iterations = 1
)

```



# astro_update parameter

Increasing the \code{astro_update} parameter leads to an approximation that may 
speed up computation. Results will not be exactly the same but can be very close 
as in the following example.  The default is that parameters are updated for 
every time-step.

```{r astro_update, eval = eval_chunks}

tms <- as.POSIXct("1990-01-01", tz = "UTC") + 0:(900)

wave_groups <- data.frame(start = 0, end = 8)

bench::mark(
  et <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups
  ),
  et_astro <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups,
    astro_update = 30L
  ), iterations = 1
)


```

# n_thread parameter

Adjust the number of threads used for parallel computation. This should result
in equivalent values.

```{r threads, eval = eval_chunks}

tms <- as.POSIXct("1990-01-01", tz = "UTC") + 0:(900)

wave_groups <- data.frame(start = 0, end = 8)

bench::mark(
  et <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups
  ),
  et_threads <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups,
    n_thread = 10L
  ), iterations = 1
)


```

# Predict and interpolate

For one second output you can predict every minute and interpolate.  
Interpolation is done via \code{stat::spline} which achieves good accuracy with 
larger approximations. The number of samples skipped may need to be adjusted 
depending on the size of your time step. 
Results will not be the exactly the same but can be very close as in 
the following example.

```{r predictinterp, eval = eval_chunks}

tms <- as.POSIXct("1990-01-01", tz = "UTC") + 0:(900)
tms_interp <- as.POSIXct("1990-01-01", tz = "UTC") + seq(0, 900, 180)

wave_groups <- data.frame(start = 0, end = 8)

bench::mark(
  et <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups
  ),
  et_interp <- calc_earthtide(
    utc = tms_interp,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups,
    utc_interp = tms
  ), iterations = 1
)

  
```


# Combination of the above

We will use a larger dataset to compare  approximation methods. In general, 
interpolation will give the best speed-up to accuracy if your time-steps are 
small (seconds).

```{r combination, eval = eval_chunks}

tms <- as.POSIXct("1990-01-01", tz = "UTC") + 0:(86400)
tms_interp <- as.POSIXct("1990-01-01", tz = "UTC") + seq(0, 86400, 180)

wave_groups <- data.frame(start = 0, end = 8)

bench::mark(
  et_astro_threads <- calc_earthtide(
    utc = tms,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups,
    astro_update = 60L,
    n_thread = 10L
  ),
  et_interp_threads <- calc_earthtide(
    utc = tms_interp,
    do_predict = TRUE,
    method = c("tidal_potential", "lod_tide", "pole_tide"),
    latitude = 52.3868,
    longitude = 9.7144,
    elevation = 110,
    gravity = 9.8127,
    cutoff = 1.0e-10,
    catalog = "ksm04",
    wave_groups = wave_groups,
    utc_interp = tms,
    n_thread = 10L
  ), iterations = 1
)

  
```

