** EXACT-FACTORIZATION, COMPACT METHODS FOR
MULTIDIMENSIONAL HEAT TRANSFER: APPLICATION TO COMBUSTION**

J. I. Ramos

Departamento de Lenguajes y Ciencias de la Computación

E. T. S. Ingenieros Industriales, Universidad de Málaga

Plaza El Ejido, s/n, 29013-Málaga, Spain

It is well known that the use of e-methods for the solution of multidimensional advection-diffusion- reaction equations such as those occurring in heat and mass transfer, combustion, fluid dynamics, etc., provides a large system of non-linearly coupled algebraic equations whose solution requires iterative techniques. Iterations may, however, be eliminated by employing linearized e-methods which are based on the discretization of the time derivative and the time-linearization of the nonlinear terms; however, a large system of linearly coupled algebraic equations still results upon the discretization of the spatial derivatives. In order to eliminate the coupling between different spatial directions in linearized e-methods, approximate factorization techniques which reduce a multi-dimensional problem to the solution of a sequence of one-dimensional ones have been proposed. The linear one- dimensional operators may be easily and inexpensively solved by means of a block tridiagonal solver. Unfortunately, these approximate factorization techniques introduce a factorization error due to the factorization of the three-dimensional operator into one-dimensional ones, and the magnitude of this error may be very large unless small time and space steps are employed in the calculations.

In this paper, an iterative, exact-factorization method based on linearized e-techniques is proposed. Since it is an exact factorization method, it does not have approximate factorization errors. Furthermore, since a linearized e-technique is employed, the factorized, one-dimensional operators in each spatial direction are linear and may be easily and inexpensively solved. However, the exact factorization introduces higher-order derivatives and a coupling between the different spatial directions. This coupling has been accounted for by means of a predictor-corrector strategy as follows. In the predictor step, an approximate factorization method which introduces factorization errors is employed, whereas the factorization errors are accounted for in an iterative manner in the corrector step. Moreover, the mixed fourth-order spatial derivatives have been eliminated by using the one-dimensional operators so that only second-order derivatives appear in the one-dimensional operators which are solved in an iterative manner. Furthermore, in order to accurately resolve the steep fronts or boundary layers that may occur, the spatial derivatives have been discretized by means of three-point, compact difference expressions which are fourth-order accurate in space, so that the resulting method is second-order accurate in time and fourth-order accurate in space, and the difference equations are formulated in delta form so that, for problems which have steady state solutions, the asymptotic numerical solution is independent of the time step employed in the calculations if approximate factorization techniques are employed to obtain the numerical solution. Since the approximate factorization errors depend on the curvature of the solution, they are expected to be largest where the curvature of the dependent variables is largest. In addition, since these errors depend on the time step, the steady state solution, if any, depends on the second power of the time step.

The use of three-point, compact difference expressions for the discretization of the spatial derivatives requires that, in addition to the unknown dependent variables, their first- and second-order spatial derivatives be considered as unknowns. This implies that, if U is the vector of dependent variables and n denotes the number of spatial directions, the blocks of the one-dimensional finite difference operators which result upon the factorization of compact, linearized e-methods are 3Nx3N for advection-diffusion-reaction equations compared with NxN for the case that the spatial derivatives are discretized with second-order accurate formula, where N denotes the number of dependent variables. The large increase in the dimensionality of the matrices is compensated by both the ease and inexpensiveness with which one-dimensional operators can be solved and the values of the first- and second-order spatial derivatives that they provide. These derivatives may be used to concentrate the grid points were steep gradients occur. Moreover, the use of the one-dimensional operators at three successive grid points together with the fourth-order accurate finite difference expressions which relate the first- and second-order derivatives to the discrete dependent variables provides seven equations for nine unknowns, i.e., the dependent variables and its first- and second-order derivatives; therefore, the spatial derivatives may be eliminated to obtain a tridiagonal matrix for the dependent variables.

The exact-factorization, compact method proposed here has been used to study the combustion of a mixture which is governed by a one-step irreversible chemical reaction in two-dimensions, i.e., n = 2, by solving the energy and mass conservation equations, i.e., N = 2, subject to homogeneous Dirichlet boundary conditions at the four boundaries as a function of time starting with a uniform mixture composition and a temperature spike at the center of the computational domain which serves as an ignition source. Calculations have been performed with second- and fourth-order accurate discretizations of the spatial derivatives, and with exact and approximate factorization methods in order to determine the influence of the spatial discretization and the factorization errors on the numerical solution. Some sample results are shown in Figures 1 and 2 which show the temperature (u) and fuel mass fraction (v) distributions at different times (t) and the errors that result from the approximate factorization of the two-dimensional operator into a sequence of two one-dimensional ones. Although the approximate factorization errors are second-order in time, their largest magnitude occurs at the steep fronts of both u and v.

Figure 1. u (top left), v (top right), and errors in u (bottom left) and v (bottom right) at t =10. (The
top figures were obtained with the implicit, linearized compact method presented here Ot=0.04,
Ox=Dy=0.5 and e= 0.5. The bottom figures correspond to the difference between the results of a
second-order accurate approximate factorization method in both space and time and the compact,
linearized technique presented here.)

Numerical experiments indicate that the exact factorization technique proposed here is able to accurately resolve the steep gradients which occur in both the temperature and fuel mass fraction as the combustion front propagates through the combustible mixture and approaches the boundaries of the domain. The differences between the results obtained with second- and fourth-order accurate discretizations of the spatial derivatives are largest at the steep fronts; in fact, fourth-order accurate discretizations provide steeper fronts than second-order ones. Moreover, the iterative, exact- factorization method presented here converges in only two or three iterations and is very efficient because the one-dimensional finite difference operators provide block tridiagonal matrices and their solution was obtained by means of LU decomposition.

Numerical experiments have also been performed in order to determine the effects of 8, i.e., the implicitness parameter, and the allocation of the reaction terms to the one-dimensional operators which result from the factorization. These experiments indicate that the factorization errors of the fully-implicit method, i.e., 8 =1, are second-order in time which is the same order as that of the Crank-Nicolson method, i.e., e = 1/2; however, the discretization errors of the former, i.e., first- order accuracy in time, are larger than those of the latter, i.e., second-order, and the fully-implicit method predicts a slightly faster flame speed than the Crank-Nicoloson technique when the combustion front is sufficiently far from the boundaries.

The allocation of the reaction terms to the one-dimensional operators results in errors which are also second-order in time and the differences between different allocations and the one which splits the reaction terms in equal amounts amongst the different spatial operators are largest at the steep combustion front. However, equal allocation of the reaction terms to the one-dimensional operators is computationally more demanding that those which assign the source terms to only one-dimensional operator because the Jacobian matrix of the reaction terms has to be evaluated for the two one- dimensional operators. Furthermore, for symmetric problems such as the one considered here, the unequal allocation of the reaction terms to the one-dimensional operators may result in solutions which are not symmetric because the non-symmetric distribution of the reaction terms increases the splitting errors associated with the factorization of the two-dimensional reaction-diffusion operator into one-dimensional ones, although, if the source terms are allocated to only one-dimensional operator, the other one is of the diffusion type.

Figure 2. u (top left), v (top right), and errors in u (bottom left) and v (bottom right) at t = 20. (The
top figures were obtained with the implicit, linearized compact method presented here Ot=0.04,
Ox=Dy=0.5 and e= 0.5. The bottom figures correspond to the difference between the results of a
second-order accurate approximate factorization method in both space and time and the compact,
linearized technique presented here.)

**ACKNOWLEDGEMENTS**

This work has been supported by Project No. PB94-1494 from the DGICYT of Spain.