Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab (811443), страница 19
Текст из файла (страница 19)
(4.11) is sought. The process starts by integrating Eq. (4.1) over element C thatenables recovering its integral balance form, which was described in Chap. 3, asZZVC34ZZ_qdVr ðkrT ÞdV ¼ð4:3ÞVCF25f2814915Fig. 4.11 Discretization stencil1016f3F3Cf1f4F4F14.1 The Discretization Process95Then, using the divergence theorem, the volume integral is transformed into asurface integral yieldingZð4:4Þ ðkrT Þ dS ¼ q_ C VCSCThis equation is actually a heat balance over element C. It is basically the integralform of the original partial differential equation and involves no approximation.Replacing the surface integral by a summation over the control volume faces,Eq. (4.4) becomesXðkrT Þf Sf ¼ q_ C VCð4:5Þf nbðC Þwhere f represents the integration point at the centroid of the bounding face.
Thistransformation is the first approximation introduced. Therefore the integral inEq. (4.4) is numerically approximated by the fluxes at the centroids of the faces.This is a second order approximation as will be demonstrated in a later chapter.Expanding the summation, Eq. (4.5) can be written asðkrT Þf1 Sf1 ðkrT Þf2 Sf2 ðkrT Þf3 Sf3 ðkrT Þf4 Sf4 ¼ q_ C VCð4:6ÞConsidering face f1 shown in Fig.
(4.12), the surface vector and temperaturegradient in Eq. (4.6) are given bySf1 ¼ Dyf1 idxf1 ¼ xF1 xC @T@TrTf1 ¼iþj@x f1@y f1wherexCxF1Dyf1Sf1rTf1isisisisisthethethethetheð4:7Þx-coordinate of the centroid of element C.x-coordinate of the centroid of element F1.area of face f1.surface vector of face f1 directed out of element C.gradient of T at the centroid of face f1.and by substitution, the first term in Eq. (4.6) is converted to@T@Tiþj@x@y @TDyf¼@x f1 1rTf1 Sf1 ¼Dyf1 if1ð4:8Þ964 The Discretization ProcessxF2S f2f2yF3f3S f1CS f3F1f1yS f 4 f4F4xFig.
4.12 Finite volume notationTo proceed further, a profile approximating the variation of T between C and F1 isneeded. Assuming linear variation of T, the x-component of the gradient at the facef1 can be written as@T@x¼f1TF1 TCdxf1ð4:9Þthus Eq. (4.8) can be approximated asrTf1 Sf1 ¼TF1 TCDyf1dxf1ð4:10Þor more generally asðkrT Þf1 Sf1 ¼ aF1 ðTF1 TC Þð4:11ÞwhereaF1 ¼ kDyf1dxf1ð4:12ÞRepeating for each of the remaining faces, the following coefficients are obtained:4.1 The Discretization Process97Dxf2dyf2Dyf¼ k 3dxf3Dxf¼ k 4dyf4aF2 ¼ kaF3aF4ð4:13Þwhich when substituted into Eq. (4.6) yieldsXXðkrT Þf Sf ¼f nbðC ÞaF ðTF TC ÞF NBðC Þ¼ ðaF1 þ aF2 þ aF3 þ aF4 ÞTC þ aF1 TF1 þ aF2 TF2 þ aF3 TF3 þ aF4 TF4¼ q_ C VCð4:14Þor more compactlyaC TC þXaF TF ¼ bCð4:15ÞF NBðC ÞwhereaC ¼ XaF ¼ ðaF1 þ aF2 þ aF3 þ aF4 ÞF NBðCÞð4:16ÞbC ¼ q_ C VCEquations similar to Eq.
(4.15) may be derived for all cells in the domain, yielding aset of algebraic equations, which can be solved using a variety of direct or iterativemethods. Focusing on element C in Fig. (4.13), Eq. (4.15) implies a relationF24S f2f2F3f389Sf3f1Sf 415F4Fig. 4.13 Element notationF1S f1Cf410984 The Discretization Process4891015.........9=........Fig. 4.14 System of equationsbetween TC and the temperatures at its four neighbors namely TF1 , TF2 , TF3 , and TF4 ,which in global assembly would be T9 , and T10 , T4 , T8 and T15 .Similar equations are also derived for boundary elements and their collectionyields the set of equations illustrated in Fig.
(4.14), which can be represented inmatrix form as given in Eq. (4.2), where A is the matrix of coefficients, [T] thesolution vector, and b is a vector composed of terms that cannot be included in A.The methodology to solve Eq. (4.2) is presented in the next section.Finally it should be stated that the properties of the finite volume method asrelated to accuracy, robustness, and other characteristics will be reviewed in laterchapters. This includes examining in more details the finite volume discretization ofthe diffusion term, which was presented above for a rectangular Cartesian grid.4.1.5Step IV: Solution of the Discretized EquationsThe discretization of the differential equation results in a set of discrete algebraicequations, which must be solved to obtain the discrete values of T.
The coefficientsof these equations may be independent of T (i.e., linear) or dependent on T (i.e.non-linear). The techniques to solve this algebraic system of equations are independent of the discretization method, and represent the various trajectories that canbe followed to obtain a solution. For the linear algebraic sets encountered in thisbook, the uniqueness of the solution is guaranteed. Therefore if the adopted solutionmethod gives a solution, it will be the desired solution.
All solution methods (i.e.,all paths to solution) which arrive at a solution will give the same solution for thesame set of discrete equations.4.1 The Discretization Process99The solution methods for solving systems of algebraic equations may be broadlyclassified as direct or iterative and are briefly reviewed below.4.1.5.1Direct MethodsIn a direct method the solution to the system of equations [e.g., Eq. (4.2)] isobtained by applying a relatively complex algorithm, in comparison with an iterative method, only once to obtain the solution for a given set of coefficients. Anexample of a direct method is matrix inversion whereby the solution is obtained as½T ¼ A1 bð4:17ÞTherefore a solution for [T] is guaranteed if A−1 can be found. However, theoperation count for the inversion of an N × N matrix is O(N3), which is computationally expensive.
Consequently, inversion is almost never employed in practicalproblems. More efficient methods for linear systems are available. For the discretization methods of interest here, A is sparse, and for structured meshes it isbanded. For certain types of equations (e.g., pure diffusion), the matrix is symmetric. Matrix manipulation can take into account the special structure of A indevising efficient solution techniques. Such methods will be reviewed in Chap. 10.In general, direct methods are rarely used in computational fluid dynamicsbecause of their large computational and storage requirements. Most industrial CFDproblems today involve hundreds of thousands of cells, with 5–10 unknowns percell even for simple problems. Thus the matrix A is usually very large, and mostdirect methods become impractical for these large problems.
Furthermore, thematrix A is usually non-linear, so that the direct method must be embedded withinan iterative loop to update nonlinearities in A. Thus, the direct method is appliedover and over again, making it all the more time-consuming.4.1.5.2Iterative MethodsIterative methods follow a guess-and-correct procedure to gradually refine theestimated solution by repeatedly solving the discrete system of equations. Let usconsider an extremely simple Gauss-Seidel iterative method.
The overall solutionloop for this method may be written as follows:(a) Guess the discrete values of T at all grid elements in the domain.(b) Visit each grid element in turn. Update T usingTC ¼PF NBðC ÞaF T F þ bCaCð4:18Þ1004 The Discretization ProcessThe neighboring values are required for the update of TC. These are assumedknown at prevailing values. Thus, grid elements which have already been visitedwill have updated values of T and those that have not will have old values.(c) Sweep the domain until all grid elements are covered.
This completes oneiteration.(d) Check if an appropriate convergence criterion is met. The requirement, forexample, could be that the maximum change in the grid-point values of T beless than 1 %. If the criterion is met, stop. Else, go back to step b and repeat.The iteration procedure described here is not guaranteed to converge to asolution for arbitrary combinations of aC and aNB. Convergence of the process isguaranteed for linear problems if the Scarborough criterion is satisfied. TheScarborough criterion requires that aC and aNB should satisfyPaF F NBðC Þ 1 for all grid pointsð4:19Þ\1 for at least one pointaCMatrices which satisfy the Scarborough criterion have diagonal dominance.The Gauss-Seidel scheme can be implemented with very little storage. All that isrequired is storage for the discrete values of T at the grid elements.
The coefficientscan be computed on the fly if desired, since the entire coefficient matrix for thedomain is not required when updating the value of T at any grid point. Also, theiterative nature of the scheme makes it particularly suitable for non-linear problems.If the coefficients depend on T, they may be updated using prevailing values of T asiterations proceed. Nevertheless, the Gauss-Seidel scheme is rarely used in practicefor solving the systems encountered in CFD. The rate of convergence of the schemedecreases to unacceptably low levels if the system of equations is large.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
