6
$\begingroup$

Suppose I wish to derive a numerical integrator for the first-order ODE $$x'(t)=F(x(t)).$$ By differentiating both sides of the expression in $t$, I can write a second-order relation also satisfied by $x(t)$: $$x''(t)=\frac{d}{dt}\left(x'(t)\right)=\frac{d}{dt}F(x(t))=DF(x(t))\cdot x'(t)=DF(x(t))\cdot F(x(t)),$$ by the chain rule, where $DF$ denotes the Jacobian of $F$.

Using Taylor series, I could write $$x(t+h)=x(t)+F(x(t))\cdot h + \frac12 DF(x(t))\cdot F(x(t))\cdot h^2 + O(h^3).$$ If we omit the $O(h^3)$ term, we get quadratic-error explicit time integrator for our ODE.

Does this time integration strategy where we know derivatives of $F(\cdot)$ have a name? What are its drawbacks?

Improve this question
$\endgroup$
3
  • $\begingroup$ I'm not sure if I understand 100%, but as per my comment below, I use a numerical technique that I call the "Taylor Series Method" (it is also used for manual calculations as you show), although with a name like that it is very hard to search for! Also the old autodiff.org site seems to be defunct (listing software like FADBAD, C-XSC, Tides)! If you want to search, the method seems to be originally by a guy called Ramon E. Moore. $\endgroup$
    – m4r35n357
    Commented 15 hours ago
  • $\begingroup$ Addendum: The method is also related to something called the Differential Transform researchgate.net/publication/…, but is numerical throughout. $\endgroup$
    – m4r35n357
    Commented 15 hours ago
  • $\begingroup$ If you apply this to a PDE problem, then you get Lax-Wendroff Method. There are also what are called two-derivative Runge-Kutta methods which use not only $F(x)$ but its derivative also. Earliest paper is likely: Chan, R.P.K., Tsai, A.Y.J.: On explicit two-derivative Runge–Kutta methods. Numer. Algorithms 53, 171–194 (2010). There is quite an explosion of work on these kind of methods, especially for conservation laws, flow problems. $\endgroup$
    – cfdlab
    Commented 14 hours ago

2 Answers 2

Reset to default
10
$\begingroup$

Not too surprisingly, it is referred to as a "Taylor series method." The drawback is that, in general, it may not be easy to compute the derivatives of $F$, particularly as you add more terms to the Taylor series to improve the accuracy.

But the method can be quite attractive if $F$ has a specific form, say a polynomial in $t$.

Other approaches to implementation of the method include using automatic differentiation methods or symbolic calculus systems to calculate the higher-order derivatives.

Improve this answer
$\endgroup$
6
  • 1
    $\begingroup$ Thanks! This is indeed the Google search term I was missing ;-) . I can't find many documents that take this method seriously --- it seems to be among the zoo of relatively unused numerical integrators. I can't tell if it's just because of the cost of computing DF, or if the integrator itself isn't useful. Maybe now that automatic differentiation is more common (e.g. using libraries like Pytorch) we should revisit... $\endgroup$
    – Justin Solomon
    Commented Jul 8 at 21:00
  • 2
    $\begingroup$ In many real applications, $F$ is very complicated. Computing its derivatives is not always easy. Worse, $F$ is in many cases no differentiable -- say if a part of the device you're modeling switches on or off at specific times or based on conditions on $x$. $\endgroup$
    – Wolfgang Bangerth
    Commented Jul 9 at 3:25
  • $\begingroup$ This makes sense! I work on a lot of graphics problems where $F$ is likely to be fairly clean/differentiable --- maybe I'll play with these! $\endgroup$
    – Justin Solomon
    Commented Jul 9 at 15:11
  • $\begingroup$ One other potential pitfall is that, although using derivative information can lead to a higher order scheme, the error coefficient may grow quite large and the higher-order scheme may require smaller step sizes than the lower order scheme to hit some desired error tolerance. $\endgroup$
    – whpowell96
    Commented Jul 11 at 16:25
  • $\begingroup$ I must interject here that the statement "it may not be easy to compute the derivatives of F" applies only to explicit ("by hand") calculations. If you are doing this numerically (as I do in my ODE project) then there are simple (but apparently not very well known) recurrence relations that can be used to calculate derivatives arbitrarily (up to degrees of hundreds and more!) for compositions of many common (and uncommon) functions. $\endgroup$
    – m4r35n357
    Commented 17 hours ago
0
$\begingroup$

To follow up from my comments, for which I can't easily find any references.

If you have a Linux box, a Mac, or a VM of either, you can try a Lorenz system simulation using GNU bc.

define void output (x, y, z, t) { print x, " ", y, " ", z, " ", t, " _ _ _ 0.0\n" }

define horner (n, j[], h) {
    auto i, s
    for (i = n; i > -1; i--) s = s * h + j[i]
    return s
}

define mul (u[], v[], k) {
    auto j, p
    for (j = 0; j < k + 1; j++) p += u[j] * v[k - j]
    return p
}

scale=15; n=12; h=0.01; steps=10000; x=-15.8; y=-18.48; z=35.64; sigma=10; rho=28; beta=8.0/3.0

output(x, y, z, 0.0)
for (step = 1; step < steps + 1; step++) {
    x[0] = x; y[0] = y; z[0] = z
    for (k = 0; k < n; k++) {
        x[k + 1] = sigma * (y[k] - x[k]) / (k + 1)
        y[k + 1] = (rho * x[k] - mul(x[], z[], k) - y[k]) / (k + 1)
        z[k + 1] = (mul(x[], y[], k) - beta * z[k]) / (k + 1)
    }
    x = horner(n, x[], h); y = horner(n, y[], h); z = horner(n, z[], h)
    output(x, y, z, h * step)
}
quit

Copy this to a file called lorenz.bc, then run:

BC_LINE_LENGTH=0 /usr/bin/bc lorenz.bc

The parameter "n" is the order of the Taylor series, experiment with it! Note that all derivatives are calculated, not approximated - step size is not used at all in the derivative calculations (the "k" loop), only the final "integration" step (like an Euler step, but using Horner's method). The "scale" parameter determines precision (in decimal places), and like "n" can be set to arbitrary values.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.