Starting BDF

Among computing professionals solving ode, the usual practice has been what might be called the rational start, see Fig. 4.5. This starts with the method BI, which can be regarded as 2-point BDF, to generate the first new point, then uses 3-point BDF to generate the next, then 4-point, and so on, until the desired k has been reached, and continues from there. Inevitably, the first few points will then have errors of a lower order than later points. This does a little better in terms of accuracy than the simple start (without the correction). 

We allow ourselves a short digression here, in order to make a special point. There are two ways of presenting an error in a numerical solution of a differential equation. The usual way is to refer to the error in the quantity computed at each new time interval that is, the difference between the numerical approximation and the underlying exact solution. Another way is to compute, for each calculated value, the time at which that value is exact, and to express the error as a time shift, the difference between the calculated time and the time at that iteration number. It is called a time shift because in many kinds of simulations dealt with in this book. 

Feldberg and Goldstein [14] extended BDF to exponentially expanding time intervals, by using a general method for computing the coefficients for any time sequence. This also yielded good results. 

Finally, Lambert [22] describes a high-order start for general multi-level methods, based on Taylor expansions using higher derivatives. This seems less practical to use as, for example, KW. 

This does a little better in terms of accuracy than the simple start (without the correction). 

---

The 3-point BDF KW start is actually very simple to implement, requiring only a 2 X 2 system whose solution (for yi and j2) is easily expressed, and so it could be feasible for use in pde. However, the table shows that it results in no better errors than simpH- or the rational start, so it does not recommend itself. It is interesting to note, regarding the error orders, that both simpH- and rational show an order close to 2, regardless of the BDF order, meaning that with these starts, BDF using more than three points is no improvement over three-point BDF. The only start that enables the full accuracy of higher BDF orders is the KW start, which follows the BDF order. 

For simultaneous solution of (16), however, the equivalent set of DAEs (and the problem index) changes over the time domain as different constraints are active. Therefore, reformulation strategies cannot be applied since the active sets are unknown a priori. Instead, we need to determine a maximum index for (16) and apply a suitable discretization, if it exists. Moreover, BDF and other linear multistep methods are also not appropriate for (16), since they are not self-starting. Therefore, implicit Runge-Kutta (IRK) methods, including orthogonal collocation, need to be considered. 

In this way, the coefficients for any y((n) can be calculated. Table A.l in Appendix A shows them all, as whole numbers m/3j, where m is the multiplier mentioned above. For each n, the Table shows forward differences (at index 1), backward derivatives (at index n) and derivatives applying at points between the two ends. For n up to 6, all possible forms are included, as they will be needed later, while for n = 7, only the forward and backward formulae are shown, as only these are needed. In case the reader wonders why all this is of interest the forms y[(n) will be used to approximate the current in simulations (see the next section) the backward forms y n(n) will be used in the section on the BDF method in Chaps. 4 and 9, and the intermediate forms shown in the Table will be used for the Kimble White (high-order) start of the BDF method, also described in these chapters. The coefficients have a long history. Collatz [169] derived some of them in 1935 and presents more of them in [170]. Bickley tabulated a number of them in 1941 [88], The three-point current approximation, essentially y((3) in the present notation, was first used in electrochemistry by Randles [460] (preempted by two years by Eyres et al. [225] for heat flow simulations), then by Heinze et al. [301], and schemes of up to seven-point were provided in [133]. 

This works rather well with odes, as also seen in Fig. 4.7. For use with pdes, however, it may be considered too much trouble to program, especially as there are easier options, for example extrapolation, which produce results that are just as good. Also, if BDF is nonetheless chosen, it was found in Sect. 4.8.1 and proved mathematically in Appendix B, that the simple start with a simple time correction produces about equally good results for much less effort. 

There is no problem with varying time intervals with two-level simulation methods, but with a method like BDF, there is the problem that one needs multipoint time derivatives calculated from unequally spaced points in time. Feldberg and Goldstein [236] show how to do this and even show how to apply the Feldbergian correction of half a time interval in this case, that becomes necessary when using the simple start for BDF, described in Chap. 4 (see also the consistency proof for this procedure in Appendix B). 

The BDF method has been described. One starts with the system such as (8.9), and goes on from there as described. This was first suggested by Richt-myer in 1957 [470J, who suggested the three-point variant, and was first used in electrochemistry by Mocak and Feldberg [402] and later refined to variable time intervals by Feldberg and Goldstein [236]. These workers call it FIRM,... 

There have been attempts to improve the performance of BDF, which is normally limited by the second-order (in the spatial interval H) discretisation of the spatial derivative. Higher-order spatial second derivatives have been tried out in connection with BDF [152,154], They can only work as intended if a high-order start is used, such as the KW start as described in Sect. 4.8.1. This start was not found to be efficient in [154], but it may be that a technique other than the one used there, such as Numerov (see Chap. 9), which does not produce banded matrices, will make the use of KW efficient and thus interesting. For this reason, the KW start is described below. 

Now for the KW start for BDF. The description in [154] will be followed here. First of all, (8.47) is rewritten in ode form for the whole system, replacing the left-hand side by the time derivative and the right-hand side by the general matrix form... 

The DuFort-Frankel scheme has apparently been dropped in favour of more interesting schemes such as BDF, which can be driven to higher orders, and for which the start-up problem has been overcome (Chap. 4). 

Some experiments show that the 2(2) forms are sufficient here, the 2(3) forms not leading to further improvement in accuracy. This is no doubt because the three-point BDF algorithm used, started with a BI step, is second order accurate in time, so a third-order form cannot improve the accuracy. A higher-order algorithm, such as ROWDA3 as used by Bieniasz [108] would make the higher 2(3) form more useful. 

Wu and White [577] have described a new method that is reminiscent of the earher work of Kimble and White [338] but makes use of the Hermitian method (that is, using derivatives) to achieve higher-order solutions for several concentration rows at a time. They also suggest, but do not demonstrate, the use of their new scheme as a possible start-up for BDF. The reader is referred to their paper for details. 

This behaviour is typical of CN and known since its inception [185] and is due to the so-called A-stability of the method (see below). As mentioned above, however, there are various means of damping out the oscillations, discussed in Sect. 8.5.1, Chap. 8. The attraction of CN is, of course, that it presents no difficulties with start-up, as does BDF. 

Table 14.1. Error orders (and errors at t 1) for some BDF starts, for the ode...

The two methods that stand out in terms of efficiency and convenience are BDF and extrapolation. Both require minimal programming effort, and can be extended to higher-order spatial derivatives. However, in the case of BDF, a limit is encountered. For the most convenient start-up methods such as the simple or the rational start, the accuracy from BDF is limited to 0(8T2). This means for one thing that one need not go beyond 3-point BDF (which is 0 8T2) in itself), but that no marked improvement can be gained from higher-order spatial derivative approximations, because there will then be a mismatch between the accuracy orders with respect to the time and spatial intervals. 

BDF Chaps. 4 and 8. L-stable and (for small number of levels) non-oscillatory method that can be driven to higher orders with respect to the time interval. Realistic starting strategies reduce the order to 0(5T2), so the three-point variant is recommended, using the simple start with subsequent subtraction of half a time interval (Sect. 4.8.1), which in this case is justified. Reasonably easy to program. 

---