Issues integrating Verner 9(8) with EmbeddedRungeKutta

I am trying to write a RK9(8) method using Verner’s coefficients and extending the EmbeddedRungeKutta class, but I am having some issues. From my understanding, all I needed to do was read in the coefficients, set up the error estimation, and have a custom interpolator to go with it, but my results are off when compared with DOP853 already provided in orekit. My diagnosis of the issue so far is that the error estimate always return estimates less than 1, and therefore, the step size is always accepted. Any help would be greatly appreciated. A caviar to this is that I have mostly been adjusting the DormandPrince Integrator to suite my purpose.

This is my error estimate below, which is basically a copy of the DormandPrince853 with a few adjustments (especially to the coefficients).

protected double estimateError(double[][] yDotK, double[] y0, double[] y1, double h
		) {
	
	timerErr.start();
		
    final StepsizeHelper helper = getStepSizeHelper();
    double error = 0;

    for (int j = 0; j < helper.getMainSetDimension(); ++j) 
    {
        final double errSum = B1 * yDotK[0][j] +  B8 * yDotK[7][j] +
                              B9 * yDotK[8][j] +  B10 * yDotK[9][j] +
                              B11 * yDotK[10][j] +  B12 * yDotK[11][j] +
                              B13 * yDotK[12][j] +  B14 * yDotK[13][j] +
                              B15 * yDotK[14][j] +  B16 * yDotK[15][j];

        final double tol = helper.getTolerance(j, FastMath.max(FastMath.abs(y0[j]), FastMath.abs(y1[j])));
        final double ratio  = errSum / tol;
        error += ratio * ratio;
        
    }
    timerErr.stop();
    
    if (error <=0)
    {
    	error = 1.0;
    }

    return FastMath.abs(h) * error / FastMath.sqrt(helper.getMainSetDimension());
}

A loooooong time ago (between 2002 and 2004), I wrote a small application that could compute and check the order of Runge-Kutta type methods (both the internal weights and the error estimation weights). The methods would be defined using xml files. The application could also create a TeX file with the analytical equations for any order (beware there are thousands of such equations, so when you generate equations for a high order, you get a TeX file that translates to a PDF with several hundreds of pages, so it is probably useless).

Here is the source code for this old project:
rkcheck.tar.bz2 (39.0 KB)

Beware it is something I have not touched in more than 20 years, so it is old, and I am not even sure it compiles with recent Java versions, but it could help you. I would suggest to create an xml file with your coefficients and to see if it really gives you the expected order for integration and for error estimation. If it outputs an order that is below what you expect, it probably means there are errors in the coefficients.

1 Like

Thanks a lot for the help Luc. I suspected that the issue was with the way I was implementing the coefficients, and the application should go a long way in helping me resolve that.