The BDF method ( English B ackward D ifferentiation F ormulas ) are linear multi-step process for the numerical solution of initial value problems of ordinary differential equations :
-
.
The following is calculated for an approximate solution at the intermediate points :
![y (x)](https://wikimedia.org/api/rest_v1/media/math/render/svg/e871993bfd131a8b0c3591c26084cf8171a74dcd)
![x_ {i}](https://wikimedia.org/api/rest_v1/media/math/render/svg/e87000dd6142b81d041896a30fe58f0c3acb2158)
-
.
The methods were introduced in 1952 by Charles Francis Curtiss and Joseph Oakland Hirschfelder and have been widely used as solvers for stiff initial value problems since the publication of C. William Gear in 1971 .
description
In contrast to the Adams-Moulton method , in the BDF method the right side is not approximated by an interpolation polynomial , instead a polynomial with (maximum) degree is constructed , which interpolates the last approximations to the solution as well as the unknown value :
![q_ {k}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f27215e46abcad60f100434d2c8003310580af95)
![k](https://wikimedia.org/api/rest_v1/media/math/render/svg/c3c9a2c7b599b37105512c5d570edc034056dd40)
![k](https://wikimedia.org/api/rest_v1/media/math/render/svg/c3c9a2c7b599b37105512c5d570edc034056dd40)
![{\ displaystyle y_ {n}, \ dotsc, y_ {n + k-1}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/8cdcf9cdcd5aa1035e1497eea667d6336361d406)
![{\ displaystyle y_ {n + k}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/98caca6a96b374a5068ef079467994418ef5181f)
-
.
In addition, one demands that the interpolation polynomial solves the given differential equation at the point , i.e. that it holds
![q_ {k}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f27215e46abcad60f100434d2c8003310580af95)
![{\ displaystyle x_ {n + k}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/2688bfda2bd1e7803ebe6481f0925ed7c7da61e2)
-
,
and thus receives a non-linear system of equations for determining the implicitly given value .
![{\ displaystyle y_ {n + k}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/98caca6a96b374a5068ef079467994418ef5181f)
Lagrange representation
One possibility for the representation of the interpolation polynomial is the Lagrange representation . The Lagrange base polynomials with the support points are defined by
![q_ {k}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f27215e46abcad60f100434d2c8003310580af95)
![k + 1](https://wikimedia.org/api/rest_v1/media/math/render/svg/552a558062ed4c0486297b5b5531c5ee044dbd9b)
![{\ displaystyle x_ {n}, \ dots, x_ {n + k}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/5b3a72d38c672dbd2511a121466703af5ecb286e)
![{\ displaystyle l_ {j} (x_ {n + i}) = \ delta _ {ji} = {\ begin {cases} 1 & {\ text {if}} j = i, \\ 0 & {\ text {if} } j \ neq i. \ end {cases}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/9faa18d8d5671a034a15a87e168a7d4dfe968283)
where is the Kronecker delta . Because of this, the representation follows directly
![{\ displaystyle \ delta _ {ji}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/e75d9441b1f6290aabb853b44d5337fdf884b5ee)
![{\ displaystyle q_ {k} (x_ {n + i}) = \ sum _ {j = 0} ^ {k} l_ {j} (x_ {n + j}) y_ {n + j} = y_ {n + i}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/315e8c54a783de393eafec8155aa2763c0a186e4)
-
.
The requirement now gives the linear recursion formula for the BDF method:
![{\ displaystyle q_ {k} '(x_ {n + k}) = f (x_ {n + k}, y_ {n + k})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f2397e90ff9eb192cc49b09c46eb20492c5475d1)
-
,
where the coefficients are given by
![\ alpha _ {j}](https://wikimedia.org/api/rest_v1/media/math/render/svg/293a364991ab1ee55c25b0f60fd9e52af7b7dbde)
-
.
Alternative Lagrange representation
Alternatively, consider the Lagrange base polynomials defined by
![{\ displaystyle L_ {j} (- s) = \ delta _ {js} = {\ begin {cases} 1 & {\ text {falls}} j = s, \\ 0 & {\ text {if}} j \ neq s. \ end {cases}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/caf285360ae16519028d6ed0268c745c9a0aa53d)
The illustration follows
-
.
This is the distance between the support points and the constant step size of the method. With the requirement of being here
![{\ displaystyle h = x_ {i + 1} -x_ {i}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/2bcbeb15d1bffb9c684d57b6e44ab13758a8e44b)
![{\ displaystyle q_ {k} '(x_ {n + k}) = f (x_ {n + k}, y_ {n + k})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f2397e90ff9eb192cc49b09c46eb20492c5475d1)
![{\ displaystyle q_ {k} '(x) = {\ frac {d} {dx}} q_ {k} (x) = {\ frac {d} {d (sh)}} q_ {k} (x_ { n + k} + sh) = {\ frac {1} {h}} {\ frac {d} {ds}} q_ {k} (x_ {n + k} + sh)}](https://wikimedia.org/api/rest_v1/media/math/render/svg/8befc37a99a0a395198bbe21d8d2ecb77fdb484c)
applies, one now obtains for the calculation of the coefficients
![{\ displaystyle \ alpha _ {j} = L_ {kj} '(0), \ quad j = 0, \ dots, k}](https://wikimedia.org/api/rest_v1/media/math/render/svg/968dc9e9dba9dd040973af2bb1e5d3fa3c1c4d76)
and with it the recursion formula
![{\ displaystyle \ sum _ {j = 0} ^ {k} L_ {kj} '(0) y_ {n + j} = hf (x_ {n + k}, y_ {n + k})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/46d020313c9cb758ecfffb546c4cdcb9ad816605)
Newton representation
The Newtonian representation of the interpolation polynomial uses backward differences, which are recursively defined by
![q_ {k}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f27215e46abcad60f100434d2c8003310580af95)
![{\ displaystyle \ nabla ^ {0} y_ {i} = y_ {i}, \ quad \ nabla ^ {j + 1} y_ {i} = \ nabla ^ {j} y_ {i} - \ nabla ^ {j } y_ {i-1}.}](https://wikimedia.org/api/rest_v1/media/math/render/svg/77a92f413524e5e8ceca2ee785753fb340f1e89e)
This can be written as
![q_ {k}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f27215e46abcad60f100434d2c8003310580af95)
-
.
This formula leads because of for on the representation
![{\ displaystyle \ left. {\ frac {d} {ds}} (- 1) ^ {j} {\ binom {-s} {j}} \ right | _ {s = 0} = {\ frac {1 } {j}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/b28788c166a60488964a859c354db96c6b6bbd69)
![j = 1, \ dots, k](https://wikimedia.org/api/rest_v1/media/math/render/svg/17b9f8d6a741f04e9280f1aa8ad318057beff36f)
![{\ displaystyle \ sum \ limits _ {j = 1} ^ {k} {\ frac {1} {j}} \ nabla ^ {j} y_ {n + k} = hf (x_ {n + k}, y_ {n + k})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/61833937acbac248fbd4e0c21382177daafd9586)
the BDF process.
Calculation formulas
All the representations of the calculation formulas considered above are equivalent, since they only used different types of representation of the unique interpolation polynomial . For the implicit calculation formulas of the BDF (k) method are:
![q_ {k}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f27215e46abcad60f100434d2c8003310580af95)
![{\ displaystyle k \ leq 6}](https://wikimedia.org/api/rest_v1/media/math/render/svg/de8af239d3c448147f46d694042b7e7484941c79)
- BDF (1) - implicit Euler method:
![{\ displaystyle y_ {n + 1} -y_ {n} = hf (x_ {n + 1}, y_ {n + 1})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/b9e2136b1a66f4bddbc4bd8b128522f8a2fb89c7)
![{\ displaystyle 3y_ {n + 2} -4y_ {n + 1} + y_ {n} = 2hf (x_ {n + 2}, y_ {n + 2})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a5ae4648a3322d653b3248848e8d59a174f4e807)
![{\ displaystyle 11y_ {n + 3} -18y_ {n + 2} + 9y_ {n + 1} -2y_ {n} = 6hf (x_ {n + 3}, y_ {n + 3})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/e304adb88bc3ee2d652df196b3169c6de39adc71)
![{\ displaystyle 25y_ {n + 4} -48y_ {n + 3} + 36y_ {n + 2} -16y_ {n + 1} + 3y_ {n} = 12hf (x_ {n + 4}, y_ {n + 4 })}](https://wikimedia.org/api/rest_v1/media/math/render/svg/96c680a168e8669c185be5f4d2b7b00f8a6c1974)
![{\ displaystyle 137y_ {n + 5} -300y_ {n + 4} + 300y_ {n + 3} -200y_ {n + 2} + 75y_ {n + 1} -12y_ {n} = 60hf (x_ {n + 5 }, y_ {n + 5})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/c79a3e7501260705e4248a37a527f879a42a51b1)
![{\ displaystyle 147y_ {n + 6} -360y_ {n + 5} + 450y_ {n + 4} -400y_ {n + 3} + 225y_ {n + 2} -72y_ {n + 1} + 10y_ {n} = 60hf (x_ {n + 6}, y_ {n + 6})}](https://wikimedia.org/api/rest_v1/media/math/render/svg/1f82b3515386c8e3b17b5d65a66f1cd8f229009d)
properties
The BDF methods are all implicit because the unknown value goes into the equation. BDF (k) has exactly the consistency order k . The BDF (1) method is the implicit Euler method . This and BDF (2) are A-stable , the higher-order methods are A ( ) -stable, with the opening angle decreasing with the higher order. In particular, BDF (2) is very popular for the calculation of stiff differential equations due to its optimal properties with regard to the second Dahlquist barrier . For k <6 , the methods are stable and consistent and therefore also convergent. The greatest incentive of the BDF methods is their large stability areas, which is why they are suitable for use in solving rigid initial value problems. For k> 6 the methods are unstable.
![\alpha](https://wikimedia.org/api/rest_v1/media/math/render/svg/b79333175c8b3f0840bfb4ec41b8072c83ea88d3)
![\alpha](https://wikimedia.org/api/rest_v1/media/math/render/svg/b79333175c8b3f0840bfb4ec41b8072c83ea88d3)
- Stability areas of the BDF process
literature
- E. Hairer, Syvert P. Nørsett, Gerhard Wanner: Solving Ordinary Differential Equations I, Nonstiff Problems , Springer Verlag, ISBN 3540566708
- E. Hairer, G. Wanner: Solving Ordinary Differential Equations II, Stiff problems , Springer Verlag, ISBN 3-540-60452-9
- HR Schwarz, N. Köckler: Numerical Mathematics , Teubner (2004)
- Curtiss, Hirschfelder Integration of stiff equations , Proc. Nat. Acad. Sci. USA, Vol. 38, 1952, 235-243.
Web links