Dense Error Estimation calcs in explicit.c Stepper

[gregjewi]

Hi All,

I've scoured the literature and couldn't find anything on how to calculate error for dense outputs. I assume calculating dense error is to maintain error thresholds for the interpolated values? Looking at explicit.c lines 177 - 195, the stepper calculates errors using an array d[ ] that has error coefficients. Where do these come from? What's the process? Where along the interpolation (what theta,) is this error calculated for? Any texts that describe this process that I can read? Relevant code below.

I'm looking at lines 177 - 195.

    //Check the dense error (in inf norm) to determine if the step can be accepted
    double err_d;
    dcopy(temp_k[0], sum, 0, link_i->dim);
    dscal(h * d[0], sum, 0, link_i->dim);

    for (unsigned int i = 1; i < num_stages; i++)
        daxpy(h * d[i], temp_k[i], sum, 0, link_i->dim);

    for (unsigned int i = 0; i < dim; i++)
        temp[i] = max(fabs(new_y[i]), fabs(y_0[i])) * error->reltol_dense[i] + error->abstol_dense[i];

    err_d = nrminf2(sum, temp, 0, link_i->dim);

    double value_d = pow(1.0 / err_d, 1.0 / meth->d_order);

and below are the coefficients from dopri5_dense.c file:

static const double d[] = { .610351562499951, 0.0, -2.105795148247852, 18.310546874999346, -25.185639003536881, 20.749496981890658, -12.378961267605213 };
//d[0] = .610351562499951;
//d[1] = 0.0;
//d[2] = -2.105795148247852;
//d[3] = 18.310546874999346;
//d[4] = -25.185639003536881;
//d[5] = 20.749496981890658;
//d[6] = -12.378961267605213;
method->d = d;

Thanks for your responses
G