16.1.2. Modified Damped Newton’s Method for 0-D Reactors

Newton’s method determines a sequence of iterations or approximate solutions that approach the true solution. For the sake of notation simplicity, we call these approximate solution vectors, . When an arbitrary is substituted into the governing equations, the equations generally do not equal zero unless also represents the true solution; in general the equations equal a residual vector . The goal then is to find such that

(16–1)

For 0-D homogeneous reactors, this vector may be, for example:

(16–2)

where:

(16–3)

where and are the gas and electron temperatures, the ’s are the gas species mass fractions, the ’s are the surface site fractions within each surface phase, the ’s are the mole fractions of species in the bulk phases, and the ’s are the surface site densities of those surface phases whose site densities are allowed to change. The exact components of the solution vector may vary depending on problem type and reactor model. The total number of unknowns, , is defined by Equation 8–153 . The corresponding vector is composed of a corresponding set of residuals of the gas energy equation, the electron energy equation, the species equations, and the continuity equation for the surface site densities.

Provided the initial estimate of the solution is sufficiently good, Newton’s method produces a sequence of iterates that converges to the solution of the nonlinear equations . The purest form of the Newton algorithm,

(16–4)

is usually difficult to implement in practice. Evaluating the Jacobian matrices, , is time consuming, and convergence typically requires a very good initial estimate . Twopnt employs the usual remedies for these problems. First, Twopnt retains the Jacobian matrix through several iteration steps, rather than updating at every iteration as indicated in Equation 16–10 . Thus, the Jacobian used at the current iteration, , may be based on a solution that is several iterations old. The user can specify in the input the maximum number of iterations Twopnt performs before calculating a new Jacobian. Second, the advancement of the iterate to is damped by a factor . The modified Newton algorithm is then,

(16–5)

where , and

(16–6)

Internally, Twopnt will generate a new Jacobian whenever convergence with the old Jacobian fails. While Equation 16–7 correctly indicates the relation between the new and the old iterate, Twopnt does not compute the inverse of the Jacobian matrix, but rather solves a system of linear equations, , for the undamped vector, .

The Twopnt algorithm determines the damping parameter and the need for a new Jacobian based on several criteria designed to keep the iteration stable and within solution bounds. To accept a new solution iterate , Twopnt requires that the undamped steps decrease in magnitude, [121] that is,

(16–7)

If the solution fails this criterion, Twopnt rejects it and retries the step with half the damping parameter or a new Jacobian matrix. Twopnt adjusts the damping parameter to ensure that the evolving solution always remains within the solution bounds. Examples of physical bounds that are imposed on the solution are: the temperature(s) must be positive and the species mass fractions must be between zero and one. For the steady-state solver, it is actually possible to modify the lower bounds placed on the species mass fractions to be slightly negative (species "floor" value). Allowing some species mole fractions to become slightly negative during iteration sometimes enhances the convergence rate, especially when the solution composition has species mass fractions that vary over many orders of magnitude.

The Newton iteration procedure continues along these lines until the user-defined convergence criteria are met. Convergence requires that the maximum norm of the undamped correction vector has been reduced to within user-defined absolute and relative tolerances, that is,

(16–8)

where is the absolute tolerance and is the relative tolerance. The relative tolerance roughly indicates the number of significant digits in the converged solution, and the absolute tolerance serves to exclude solution components smaller than from the relative convergence criterion. Typically, the absolute tolerance should be smaller than the maximum mass fraction of any species of interest. The user should be particularly careful in specifying for plasma solutions, since the small electron mass will result in a much smaller mass fraction for electrons than for other important species. The relative tolerance should generally be in the range of 10-3 to 10-4.

If damping does not produce a suitable correction vector, Twopnt computes a new Jacobian. In the case when a new Jacobian has just been computed, and a damped Newton step still cannot produce a suitable correction vector, Twopnt begins to take time steps. The time-stepping procedure is described in the section following the description of the Jacobian Matrix below.