# Dormand–Prince 5(4): How to update the stepsize and make accept/reject decision?

https://en.wikipedia.org/wiki/Dormand–Prince_method

I want to implement the Dormand-Prince 4(5) version to solve Initial Value problems. Using regular notation I have $$A$$ matrix and the $$c,b,\hat{b}$$ vectors.

I also have the error tolerence of $$ETOL$$ absolute error tolerence $$ATOL$$ and relative error tolerance $$RTOL$$

At each time step I compute $$x_n$$ by use of $$b$$ and $$\hat{x}_n$$ by us of $$\hat{b}$$ and estimate the error by as the difference $$e_n=x_n-\hat{x}_n$$

The way I understand this method, when the error is approximated then I have to do the following to things:

• Decide whether to accept the step just computed or not
• update the stepsize

Which exact decision rules/algorithm should I use to do the those steps?

the book I am using (Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations by Petzold and Ascher) does not explain this part.

Side Note:

And I assume there is a way to update the stepsize regarding what kind of RUnge Kutta method you use.

• I kindly ask people not to use use Scientific papers when answering since I have very limited experience in reading these papers/articles and I have already looked at some paper to find answer to this problem – k.dkhk May 2 '19 at 7:45
• There are simple methods, but there are also a lot of options as this is still an active area of research. It's hard to know which method you are looking for. What is your goal? – David Ketcheson May 2 '19 at 12:21
• @DavidKetcheson My goal is to replicate the results from Matllabs ode45 and SciPy's rk45. I get that there might be more than one option, in that case I am interested in the most common way for dopri54 (what is the consensus). I could have asked the question. If you needed to implement this dopir5, how would YOU have decided the step sizes – k.dkhk May 3 '19 at 17:30

The basic idea is that

1. You use the estimated error given to you (cheaply) by the embedded methods;
2. You use a metric to define acceptance using a user-defined relative and absolute tolerance;
3. Based on the order properties of the code and this metric, a new step size is computed;
4. You avoid large differences in step sizes from one step to the next by limiting the maximum increase and maximum decrease.

For two approximations to the solution $$y_{1i}$$ (using the lower order) and $$\hat{y}_{1i}$$ (using the higher order), we want the error to satisfy (componentwise) $$|y_{1i} - \hat{y}_{1i}| \leq sc_i|, \quad sc_i = Atol_i + \max\left(|y_{0i}|,|y_{1i}|\right)\cdot Rtol_i.$$

As a measure of the error, one takes $$err = \sqrt{\frac{1}{n}\sum_{i=1}^n\left(\frac{y_{1i}-\hat{y}_{1i}}{sc_i}\right)^{2}}$$ and this value is compared to one. From the expected error behaviour $$err \propto C\cdot h^{q+1}$$ the optimal step size is obtained as $$h_{\text{opt}} = h\cdot\left(\frac{1}{err}\right)^{\frac{1}{q+1}}$$ where $$h$$ is the previous step size.

Typically, one would add a kind of safety factor so that the probability of having selected a good step size is increased. Typically this safety factor $$fac$$ is taken $$fac = 0.8$$ or some power of the order like $$fac = \left(0.25\right)^{\frac{1}{q+1}}$$ so that $$h_{\text{opt}} = fac\cdot h\cdot\left(\frac{1}{err}\right)^{\frac{1}{q+1}}$$ .

As a second measure of safety, one wants to prohibit the stepsize to vary too much from one step to the next. Typically, one limits the new stepsize to be less than a factor five of the previous stepsize. And in some codes, a stepsize increase is forbidden if the previous step was rejected (so if in a previous step, $$h$$ was decreased, it can only increase again if a step with that step size is accepted).

The initial step size for the algorithm can be set by the user or it could be estimated using a procedure based on the order of the method. Hindmarsch and Shampine have developed methods to do so.

• Do note that this is close to but not the same as the algorithm used in dopri5. – Chris Rackauckas May 3 '19 at 7:17
• Thank you very much. I will implement this method and comment if I need more. This was more or less what I was looking for! – k.dkhk May 3 '19 at 17:31