Timestep splitting

The numerical jet model [9-11] is based on the numerical solution of the time-dependent, compressible flow conservation equations for total mass, energy, momentum, and chemical species number densities, with appropriate in-flow/outfiow open-boundary conditions and an ideal gas equation of state. In the reactive simulations, multispecies temperature-dependent diffusion and thermal conduction processes [11, 12] are calculated explicitly using central difference approximations and coupled to chemical kinetics and convection using timestep-splitting techniques [13]. Global models for hydrogen [14] and propane chemistry [15] have been used in the 3D, time-dependent reactive jet simulations. Extensive comparisons with laboratory experiments have been reported for non-reactive jets [9, 16] validation of the reactive/diffusive models is discussed in [14]. 

The second approach is a fractional-step method we call asymptotic timestep-splitting. It is developed by consideration of the specific physics of the problem being solved. Stiffness in the governing equations can be handled "asymptotically" as well as implicitly. The individual terms, including those which lead to the stiff behavior, are solved as independently and accurately as possible. Examples of such methods include the Selected Asymptotic Integration Method (4,5) for kinetics problems and the asymptotic slow flow algorithm for hydrodynamic problems where the sound speed is so fast that the pressure is essentially constant (6, 2). ... 

In contrast, the asymptotic approach puts minimal strain on the computer but demands more of the modeller. The convergence of the computed solutions is usually easy to test with respect to spatial and temporal resolution, but situations exist where reducing the timestep can make an asymptotic treatment of a "stiff" phenomenon less accurate rather than more accurate. This follows because the disparity of time scales between fast and slow phenomena is often exploited in the asymptotic approach rather than tolerated. Furthermore, the non-convergence of any particular solution is often easier to spot in timestep splitting with asymptotics because the manner of degradation is usually catastrophic. In kinetics calculations, lack of conservation of mass or atoms signals inaccuracy rather clearly. 

Equation (25.85), the basis of every atmospheric model, is a set of time-dependent, nonlinear, coupled partial differential equations. Several methods have been proposed for their solution including global finite differences, operator splitting, finite element methods, spectral methods, and the method of lines (Oran and Boris 1987). Operator splitting, also called the fractional step method or timestep splitting, allows significant flexibility and is used in most atmospheric chemical transport models. 

The skeletal LN procedure is a dual timestep scheme, At, Atm, of two practical tasks (a) constructing the Hessian H in system (17) every Atm interval, and (b) solving system (17), where R is given by eq. (3), at the timestep At by procedure (23) outlined for LIN above. When a force-splitting procedure is also applied to LN, a value At > Atm is used to update the slow forces less often than the linearized model. A suitable frequency for the linearization is 1-3 fs (the smaller value is used for water systems), and the appropriate inner timestep is 0.5 fs, as in LIN. This inner timestep parallels the update frequency of the fast motions in force splitting approaches, and the linearization frequency Atm) is analogous to the medium timestep used in such three-class schemes (see below). 

A detailed examination of LN behavior is available [88] for the blocked alanine model, the proteins BPTI and lysozyme, and a large water system, compared to reference Langevin trajectories, in terms of energetic, geometric, and dynamic behavior. The middle timestep in LN can be considered an adjustable quantity (when force splitting is used), whose value does not significantly affect performance but does affect accuracy with respect to the reference trajectories. For example, we have used Atm = 3 fs for the proteins in vacuum, but 1 fs for the water system, where librational motions are rapid. 

The many approaches to the challenging timestep problem in biomolecular dynamics have achieved success with similar final schemes. However, the individual routes taken to produce these methods — via implicit integration, harmonic approximation, other separating frameworks, and/or force splitting into frequency classes — have been quite different. Each path has encountered different problems along the way which only increased our understanding of the numerical, computational, and accuracy issues involved. This contribution reported on our experiences in this quest. LN has its roots in LIN, which... 

Avery powerful method for the solution of system 8-68, widely in use, is the so-called alternating-direction implicit method, which is based on the idea of splitting each timestep... 

The qualitative criterion for the validity of (25.97) or (25.100) is that the concentration values must not change too quickly over the applied timestep from any of the individual processes. The error because of the splitting becomes zero as At —> 0, but unfortunately the maximum allowable value of At cannot be easily determined a priori. The allowable At depends strongly on the simulated processes and scenario. For example, for the application of (25.101) in an urban airshed, McRae et al. (1982a) recommended a value of At smaller than 10 min. 

So the total cost to advance over one timestep using splitting would be proportional to 2/V3, which is considerably less than the N1 cost for the fully implicit case. If each onedimensional method is stable, then the overall splitting technique is usually stable. Splitting does not require exclusive use of finite difference methods for example, finite clement techniques could be used to solve the above one-dimensional problems. 

In Section 25.6.1 we discussed finite difference schemes for the solution of the onedimensional diffusion equation. This explicit scheme of (25.93) is stable only if At < (Ax)2/2. If K is not equal to unity, the corresponding stability criterion is KAt < (Ax)2/2. Therefore explicit schemes cannot be used efficiently because stability considerations dictate relatively small timesteps. Thus implicit methods are used for the diffusion aspect of a problem. In higher dimensions this requirement implies that splitting will almost certainly have to be used. 

To understand why this is the case, we consider a step of the primitive splitting algorithm applied to a simple system consisting of a fixed body at the origin and a single body (mass=l) moving under an external potential and in occasional contact with the fixed body, both bodies having radius 0.5. Let us consider three possibilities (i) there is no collision in the timestep, (ii) there is a collision only at... 

We may use a splitting approach to develop discretization algorithms for other forms of dynamics. In Chap.4, we included holonomic constraints to the ODEs defining constant-energy Hamiltonian dynamics, as a way to increase the stability of our algorithms and increase the usable timestep. We can add constraints to the Langevin dynamics SDKs similarly as... 

We illustrate this numerically in Fig. 7.13, by applying the splitting scheme (7.48) using various timesteps and frictions. Each gridpoint marks an independent Simula-... 

Figure 4.2. Example 1 - Splitting dimension and timestep behavior for block-Schur decomposition at each discretization step (41 steps)...

Figure 4.5. Example 2 - Splitting dimension and timestep behavior.

---