Matrix Newton iteration

Compute a new set of values of the T) tear variables by solving simultaneously the set of N energy-balance equations (13-72), which are nonlinear in the temperatures that determine the enthalpy values. When linearized by a Newton iterative procedure, a tridiagonal-matrix equation that is solved by the Thomas gorithm is obtained. If we set gj equal to Eq. (13-72), i.e., its residual, the hnearized equations to be solved simultaneously are... 

Steady-state solutions are found by iterative solution of the nonlinear residual equations R(a,P) = 0 using Newton s methods, as described elsewhere (28). Contributions to the Jacobian matrix are formed explicitly in terms of the finite element coefficients for the interface shape and the field variables. Special matrix software (31) is used for Gaussian elimination of the linear equation sets which result at each Newton iteration. This software accounts for the special "arrow structure of the Jacobian matrix and computes an LU-decomposition of the matrix so that qu2usi-Newton iteration schemes can be used for additional savings. 

In solving the underlying model problem, the Jacobian matrix is an iteration matrix used in a modified Newton iteration. Thus it usually doesn t need to be computed too accurately or updated frequently. The Jacobian s role in sensitivity analysis is quite different. Here it is a coefficient in the definition of the sensitivity equations, as is 3f/9a matrix. Thus accurate computation of the sensitivity coefficients depends on accurate evaluation of these coefficient matrices. In general, for chemically reacting flow problems, it is usually difficult and often impractical to derive and program analytic expressions for the derivative matrices. However, advances in automatic-differentiation software are proving valuable for this task [36]. 

The effect is to keep the iteration matrix (Jacobian) banded, which considerably improves the efficiency of the Newton iteration that is used to solve the discrete problem. This procedure is equivalent to solving a simple first-order differential equation,... 

The Jacobian of the system is a square matrix, but importantly, because the residuals at any mesh point depend only on variables at the next-nearest-neighbor mesh point, the Jacobian is banded in a block-tridiagonal form. Figure 16.10 illustrates the structure of the Jacobian in the form used by the linear-equation solution at a step of the Newton iteration,... 

Because the T matrix is symmetric and has real elements, it turns out that the three roots A i, A.2, and A3 must be real numbers. For any arbitrary stress state the cubic equation is usually most conveniently solved by numerical iteration. For example, one could graph the determinant as a function of A, observing the approximate A values of the zeros. Then, taking the approximate values of A as starting iterates, a Newton iteration could be used to determine each of the exact roots. The three roots are the principal stresses. 

The bvp4c MATLAB code can deal with singular ODEs and we shall explain its use and the necessary preparations in Chapter 5. In fact there we show how to modify the inner workings of the built-in MATLAB BVP code bvp4c so that it does not stop when an intermediate Newton iteration encounters a singular or near singular Jacobian matrix, but rather continues with the least squares solution. The modifications to bvp4c will be explained when there is need in Chapter 5. 

If A is not a square matrix and we command A b in MATLAB, then the SVD is invoked and finds the least squares solution to the minimization problem min., Ax — b. A slight variant that uses only the QR factorization mentioned in subsection (F) for a singular but square system matrix A Rn,n is used inside our modified boundary value solver bvp4cf singhouseqr. m in Chapter 5 in order to deal successfully with singular Jacobian matrices inside its embedded Newton iteration. 

Iwork(18) indicates the number of damped Newton iterations allowed, with updated matrix G to) each time (see Eq. (B.1-2)), before changing the update interval to Iwork(18) iterations in the completion of the initial state. If the input value of Iwork(18) is zero, the code resets it to 5. 

The deferred updating of the Jacobian matrix PD evidently works very well 23 updates sufficed for the 231 integration steps. Here 563 evaluations of the vector f were needed besides those used in updating PD this corresponds to about two Newton iterations per integration step. 

The computational domain is the unit square in u and v, and this was divided into a 15 x 15 mesh i.e., 225 elements, and 16 x 16 = 256 nodes, so 256 basis functions and 256 residual equations. The Jacobian matrix was banded with a total bandwidth of 35. The first solution computed was the minimal surface, for which the initial estimate was an hyperbolic paraboloid. The nonlinear system of residual equations was solved by Newton iteration on a Cyber 124, each iteration using about 1 second cpu time. For nearly all the surfaces calculated, the mesh was an even mesh over the entire unit square. However, for the surfaces just near the close-packed spheres (CPS) limit, the nodes were evenly spaced in the u-direction but placed as follows in the i -direction i = 0,1/60,1/30,0.05,0.075,0.1,0.15,0.2,0.3,0.4,0.5,0.6,0.7,... 

This equation is known as the Gauss-Newton method [49]. For square matrices K, Eq. (62) can be reduced to Newton iterations using matrix identity (K C K) = K C(K ) This is why solving Eq. (2) with square K, as well as, Newton method also can be considered as a minimization of quadratic form [45]. 

The Newton iteration is well defined if the matrix of the system (KKT matrix) is nonsingular. It requires that the following two conditions are fulfilled ... 

For the solution of the entire system (2.1), (2.3), (2.4) it suggests itself to use the Newton iteration. But due to the structure of (2.4) the resulting Jacobian matrix would be non-sparse The oxygen concentration in the upper part of the reactor depends on the values of both the carbon concentration and the temperature in the layers lying underneath. Thus (2.1), (2.3) have been discretized by means of the finite element method and afterwards been solved each individually. The latter was realized through two fixed point iterations for (2.1), (2.4) with fixed temperature T and for (2.3) with fixed concentrations Cc, Coaj respectively. For this (2.1) (including boundary conditions) is written as... 

The better the Newton iteration matrix is approximated the smaller k. The compatibility condition can be transformed to... 

For performing the Gaufi-Newton iteration we have to integrate Eq. (7.3.2) and to compute the sensitivity matrix with respect to the initial values. In both subtasks we make use of the solution properties by applying coordinate projection. 

The lowest eigenvalue of the Hamiltonian matrix of Exercise 11.11.2 may also be determined by NewKai s method. Perform one Newton iteration with Cj = (1,0,0) as the starting vector. 

To start, we typically use a crude approximation ofthe Jacobian such as 5l l = /.In general, convergence is slower with Broyden s method than with Newton s method. Note, however, that the quadratic convergence of Newton s method is only obtained near the solution, and that this update formula does not require the N additional function evaluations of the finite difference approximation. This reduction in workload per Newton iteration makes the Broyden method a popular choice when we do not supply a routine to evaluate the Jacobian matrix analytically. For more information on qnasi-Newton methods and. Jacobian estimation, consult Nocedal Wright (1999). 

As X changes little from one iteration to the next and a i is fixed by Wh (/aed leading-coefficient BDF method), B varies slowly and we can save much CPU time through LU factorization. This algorithm marches forward in time similarly to an ODE system however, for the Newton iterations to be successfiil, the matrix must be nonsingular. This is unfortunately not always the case. We can identify the condition that must be met for 5 + to be invertible for the special case of (4.188),... 

However, if the direct sparse matrix solution method is used in the solution, then the calculation will take approximately twice the time without indicating that the equation is linear as the SPM method will require a complete matrix solution at the second Newton iteration. The approximate COE iterative solution method is used here to illustrate another possible solution approach with the pde2fe() function. 

The printed output, not shown here, also indicates that 3 to 4 Newton iterative cycles are required at each time point to solve the nonlinear set of PDE equations. The net result of all this is that the set of 2996 by 2996 matrix equations must be solved a total of over 500 times (or about 0.7 sec per solution). The PDE solutions of this and the previous chapter are some of flie longest running programs in this text. The reader is encouraged to re- execute the code in Listing 13.19 using the COE and SPM parameters and determine which approximate solution method is faster for this particular problem. 

---