Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab (811443), страница 75
Текст из файла (страница 75)
They showed that in order for the solution of a differenceequation to converge to the solution of the partial differential equation the numericalscheme must use all the information contained in the initial data that influence thesolution. This requirement has become later known as the CFL condition.13.2The Finite Difference Approach495In reality the CFL condition can be interpreted simply as one of the basic rules thatshould be satisfied by the coefficients, namely the opposite signs rule extended toinclude the transient coefficients.
Thus just as /F is considered a ‘spatial’ neighbor of/C , /C is a ‘temporal’ neighbor of /C , and the opposite signs rule should equallyapply to both. Noting that the diagonal coefficient is now aC and the coefficient of its‘temporal’ neighbor is aC þ aC , the opposite signs requirement becomesaC þ aC 0:13.2.2.1ð13:13ÞStability of a Transient-Advection CaseFor the one dimensional pure advection problem with a flow moving right wiseshown in Fig.
13.5, the aC and aC coefficients in the discretized equation of elementC, using the upwind scheme for the interpolation of all variables at an element face,are given byaC ¼ m_ e ¼ qC uC DyCaC ¼ qC VCq DxC DyC¼ CDtDtð13:14ÞTherefore, the CFL condition requiresaC þ aC 0qC uC DyC )qC DxC DyC0Dtð13:15ÞorDt DxC:uCð13:16ÞFor convection dominated flows, defining a CFL number asCFLconv v Dt¼ CDxCð13:17ÞxCuwwWWuwueeueWECxwEExeFig.
13.5 A portion of the discretized domain for a one dimensional convection problemyC49613 Temporal Discretization: The Transient Termimplies that for numerical stability the CFL number should satisfyCFLconv 1:13.2.2.2ð13:18ÞStability of a Transient-Diffusion CaseFor pure diffusion problems, the expression for the CFL number is different. Forthat purpose, the one dimensional pure diffusion problem schematically depicted inFig. (13.6) is considered.The aC and aC coefficients in the discretized equation of element C using linearinterpolation profiles are given byaC ¼C/e DyC C/w DyCþdxedxwaC ¼ qC VCq DxC DyC¼ CDtDtð13:19ÞTherefore, the CFL condition requiresaC þ aC 0)C/e DyC C/w DyC qC DxC DyC0þdxedxwDtð13:20ÞorDt qC DxC:C/eC/wþdxe dxwð13:21ÞFor the case when the grid is uniform and the diffusion coefficient is constant,Eq.
(13.21) becomesDt qC ðDxC Þ22C/C:ð13:22ÞxCWWWECxwEExeFig. 13.6 A portion of the discretized domain for a one dimensional diffusion problemyC13.2The Finite Difference Approach497For diffusion dominated problems, a CFL number is defined asCFLdiff ¼C/C Dtð13:23ÞqC ðDxC Þ2implying that for stability the following condition should be satisfied:CFLdiff 13.2.2.31:2ð13:24ÞStability of a Transient-Convection-Diffusion CaseFor the case of a multi dimensional unsteady convection and diffusion problem(Fig. 13.7) and based on the derivations presented in Chap.
12, the coefficients inEq. (13.14) are given byqC VCDt P / EfCfþ m_ f ; 0aC ¼dCFf nbðC ÞaC ¼ ð13:25ÞSubstituting the expressions for the coefficients from Eq. (13.25) in Eq. (13.13), theCFL condition becomes q VX / EfC0ð13:26ÞCfþ m_ f ; 0 CdDtCFf nbðCÞleading to the following constraint on the time step:Dt qC VC :PEfC/fþ m_ f ; 0dCFf nbðC Þð13:27ÞFig. 13.7 A portion of thediscretized domain for a multidimensional convectionproblemSffCdCFEfF49813 Temporal Discretization: The Transient TermEquation (13.27) is the general requirement for stability of explicit transient schemes.In fact the conditions obtained earlier for pure convection and pure diffusion in onedimensional domains can be derived as special cases of Eq.
(13.27). For the case of aone dimensional diffusion problem with a uniform grid of cell size Dx, constantdensity q, and a uniform diffusion coefficient C/ , Eq. (13.27) reduces toDxC DyCDt qC0z}|{VC1 ) Dt ¼DyqC ðDxC Þ2Cz}|{P B / Ef CB C_þm;0Cff@ |{z}AdCFf nbðC Þ|{z}|fflfflfflffl{zfflfflfflffl}¼C/e þC/w2C/C:ð13:28Þ¼0DxCWhile for the case of a one dimensional advection problem discretized using theupwind scheme and with the flow moving from left to right, Eq.
(13.27) reduces to¼DxC DyCDt 0qCz}|{VC1 ) Dt P B / Ef CBC_þm;0 Cf@ fAf nbðC Þ |{z} dCF|fflfflfflffl{zfflfflfflffl}¼0DxC:uCð13:29Þm_ e ¼qC uC DyCThis stability constraint is stringent and very restrictive as it forces the use ofextremely small time steps when solving transient problems. That is, whereas thecomputational cost at each time step is small in comparison to what would berequired to solve a system of equations at that level, the imposed limitation by theCFL condition necessitates a larger number of steps to move the solution in time.Therefore the benefit of reducing the calculations per time step is lost by the muchlarger number of time steps required.
Moreover, Eq. (13.27) indicates also thatimproving the spatial accuracy by decreasing the grid size, decreases further themaximum time step size that can be used without causing instabilities.As shown next, such a constraint does not apply to implicit schemes for whichthe transient term has always the proper sign.13.2.3Backward Euler SchemeTo derive the backward Euler scheme, the value of the function T at time t Dt isexpressed using a Taylor series using the values of T and its derivates at time t as13.2The Finite Difference Approach499T ðt DtÞ ¼ T ðtÞ @T ðtÞ@ 2 T ðtÞ Dt2þ Dt þ@t@t2 2!ð13:30ÞManipulating Eq. (13.30), an equation for the first derivative is obtained as@T ðtÞ T ðtÞ T ðt DtÞ @ 2 T ðtÞ Dt¼þþ @tDt@t2 2!ð13:31ÞReplacing T by ðq/Þ in Eq.
(13.31) and substituting the resulting expression for thederivative in Eq. (13.3), the discretized equation becomes ðqC /C Þt ðqC /C ÞtDtVC þ L /tC ¼ 0:Dtð13:32ÞThen invoking the algebraic relation of the spatial operator and the suggestednotation, the complete algebraic form of the transient scalar equation is obtained asX aF /F ¼ bC þaC /Cð13:33ÞaC þ aC /C þF NBðC Þwith the coefficients given byqC VCDtq VCaC ¼ CDtaC ¼ð13:34ÞThe stencil for Eq. (13.33) is shown in Fig. 13.8.
It is clear that with the spatial operatorevaluated at the same time level as the new temporal coefficient, resolving the ϕ field atFig. 13.8 Stencil for thebackward Euler stencilCC50013 Temporal Discretization: The Transient Terma new time level requires solving a system of equations. This type of schemesrequiring the solution of a system of equations is denoted by implicit schemes [5–12].As can be inferred from Eq. (13.34) aC and aC are of opposite signs guaranteeingthat /C is bounded by the values of its spatial neighbors at the current time stept and by the value of its temporal neighbor at the previous time step t Dt. Thisimplies that the scheme is always stable independent of the time step used, allowingfor the solution to proceed rapidly by using large time steps. Nonetheless this is notthe ideal scheme as it is of low order and solutions obtained with this scheme are oflow accuracy unless small time steps are used, which puts its use in a quandary.Adopting large time steps for computational efficiency results in a solution of lowaccuracy and using small time steps for higher accuracy is associated with lowcomputational efficiency.13.2.4Crank-Nicolson SchemeIn the Crank-Nicolson scheme [2, 14] a more accurate representation of the transient term is derived by expressing the values of the function T at times t Dt andt þ Dt in terms of the values of T and its derivates at time t asT ðt þ DtÞ ¼ T ðtÞ þ@T ðtÞ@ 2 T ðtÞ Dt2 @ 3 T ðtÞ Dt3þþ Dt þ@t@t2 2!@t3 3!@T ðtÞ@ 2 T ðtÞ Dt2 @ 3 T ðtÞ Dt3þ T ðt DtÞ ¼ T ðtÞ Dt þ@t@t2 2!@t3 3!ð13:35ÞThen, subtracting T ðt þ DtÞ from T ðt DtÞ given in Eq.
(13.35), an equation forthe first derivative is obtained as @T ðtÞ T ðt þ DtÞ T ðt DtÞ¼þ O Dt2@t2Dtð13:36ÞNote that the order of accuracy of the derivative is now OðDt2 Þ since the secondorder derivative is completely eliminated.Substituting the time derivative given by Eq. (13.36) into Eq. (13.3) yields ðqC /C ÞtþDt ðqC /C ÞtDtVC þ L /tC ¼ 02Dtð13:37ÞThen invoking the algebraic relation of the spatial operator, and using the suggestednotation, the complete algebraic form of the transient scalar equation is obtained as0aC /C ¼ bC @aC /C þXF NBðC Þ1aF /F A aC /Cð13:38Þ13.2The Finite Difference Approach501with the coefficients given byqC V2Dt q V¼ C2DtaC ¼aCð13:39ÞThe stencil for Eq.
(13.38) is shown in Fig. 13.9. It is clear that the scheme is anexplicit type scheme, since the evaluation of ðq/ÞtþDt can be performed using onlyold values. However two old levels are now needed, with the spatial operator beingevaluated at one of these levels.An analysis of the stability of the CN scheme can be performed after slightlymodifying the original equation.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
