# Convergence Simulations

The convergence simulation type is useful for deriving order of convergence estimates
from a group of simulations. This object will automatically assemble error vectors
into a more useful manner and provide plotting functionality. Convergence estimates
are also given by pair-wise estimates.

One can automatically have DifferentialEquations.jl perform the error analysis by
passing a `ConvergenceSimulation` a vector of solutions, or using one of the provided
`test_convergence` functions. These will give order of convergence estimates and
provide plotting functionality. This requires that the true solution was provided
in the problem definition.

`ConvergenceSimulation`s can either be created by passing the constructor the
appropriate solution array or by using one of the provided `test_convergence` functions.

## The ConvergenceSimulation Type

A type which holds the data from a convergence simulation.

### Fields

  - `solutions::Array{<:DESolution}`: Holds the solutions. For the Monte-Carlo method of
    `test_convergence` these hold a single representative trajectory unless
    `retain_solutions = true` is passed — see
    [Memory use of Monte-Carlo studies](@ref).

  - `errors`: Dictionary of the error calculations. Can contain:
    
      + `h1Errors`: Vector of the H1 errors.
      + `l2Errors`: Vector of the L2 errors.
      + `maxErrors`: Vector of the nodal maximum errors.
      + `node2Errors`: Vector of the nodal l2 errors.
  - `N`: The number of simulations.
  - `auxdata`: Auxiliary data of the convergence simulation. Entries can include:
    
      + `dts`: The dt's in the simulations.
      + `dxs`: The dx's in the simulations.
      + `μs`: The CFL μ's in the simulations.
      + `νs`: The CFL ν's in the simulations.
  - `𝒪est`: Dictionary of order estimates. Can contain:
    
      + `ConvEst_h1`: The H1 error order of convergence estimate for the convergence simulation.
        Generated via `log2(error[i+1]/error[i])`. Thus only valid if generated by halving/doubling
        the dt/dx. If alternate scaling, modify by dividing of log(base,ConvEst_h1)
      + `ConvEst_l2`: The L2 error order of convergence estimate for the convergence simulation.
        Generated via `log2(error[i+1]/error[i])`. Thus only valid if generated by halving/doubling
        the dt/dx. If alternate scaling, modify by dividing of log(base,ConvEst_l2)
      + `ConvEst_max`: The nodal maximum error order of convergence estimate for the convergence simulation.
        Generated via `log2(error[i+1]/error[i])`. Thus only valid if generated by halving/doubling
        the dt/dx. If alternate scaling, modify by dividing of log(base,ConvEst_max)
      + `ConvEst_node2`: The nodal l2 error order of convergence estimate for the convergence simulation.
        Generated via `log2(error[i+1]/error[i])`. Thus only valid if generated by halving/doubling
        the dt/dx. If alternate scaling, modify by dividing of log(base,ConvEst_node2)
  - `convergence_axis`: The axis along which convergence is calculated. For example, if
    we calculate the dt convergence, convergence_axis is the dts used in the calculation.

## Plot Functions

The plot functionality is provided by a Plots.jl recipe. What is plotted is a
line series for each calculated error along the convergence axis. To plot a
convergence simulation, simply use:

```julia
plot(sim::ConvergenceSimulation)
```

All of the functionality (keyword arguments) provided by Plots.jl are able to
be used in this command. Please see the Plots.jl documentation for more information.

## ODE

`test_convergence(dts::AbstractArray,prob::AbstractODEProblem)`

Tests the order of the time convergence of the given algorithm on the given problem
solved over the given dts. Keyword arguments are passed to the ODE solver.

## SDE

`test_convergence(dts::AbstractArray,prob::AbstractSDEProblem)`

Tests the strong order time convergence of the given algorithm on the given problem
solved over the given dts. Keyword arguments are passed to the ODE solver. Except:

  - `numMonte`: The number of simulations for each dt. Default is 10000.

### Order Estimation

`calc𝒪estimates(error::Vector{Number})``

Computes the pairwise convergence estimate for a convergence test done by
halving/doubling stepsizes via

log2(error[i+1]/error[i])

Returns the mean of the convergence estimates.

## Memory use of Monte-Carlo studies

For an SDE or RODE problem, `test_convergence` runs an ensemble per step size. A full
solution per trajectory carries the solver cache, the noise process, the problem and its
interpolation, so keeping them makes a study over many trajectories retain far more than
the error estimates it reports — at hundreds of thousands of trajectories, tens of
gigabytes.

They are therefore **not** kept by default. Each trajectory is reduced to a
[`ConvergenceTrajectory`](@ref) by the ensemble's `output_func` as it is solved, so the
full solutions never coexist, and once the errors have been computed the ensemble is
stripped to a single representative trajectory. `errors`, `weak_errors`, `error_means`
and `𝒪est` come from the full ensemble and are unaffected — what changes is that
`solutions[i].u` holds one entry rather than `trajectories` of them.

Pass `retain_solutions = true` when the trajectories themselves are the point, such as
comparing two algorithms path by path. `sim[i, j]` forwards into `solutions[i][j]`, so
it too reaches only the representative trajectory unless the solutions are retained.

The reduction is skipped where it cannot apply, and the solutions retained: when you
supply your own `EnsembleProblem` (set its `output_func` instead), when `expected_value`
is given, and when `weak_timeseries_errors` or `weak_dense_errors` is set.

## API

```@docs
DiffEqDevTools.ConvergenceSimulation
DiffEqDevTools.ConvergenceTrajectory
DiffEqDevTools.test_convergence
DiffEqDevTools.analyticless_test_convergence
```
