Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab (811443), страница 62
Текст из файла (страница 62)
(11.75) resulting in the following RHSconvection term:1111RHSconvection ¼ /C /W þ /E km_ e ; 0k þ /E /EE þ /C km_ e ; 0k44441111 /C /E þ /W km_ w ; 0k þ /W /WW þ /C km_ w ; 0k4444ð11:116Þ39811 Discretization of the Convection TermThe rate of change of this term with respect to /C is found as@ ðRHSconvection Þ331¼ km_ e ; 0k km_ w ; 0k ðm_ e þ m_ w Þ@/C444ð11:117Þwhich, for a constant velocity field, is always negative indicating a stable scheme,but not for the general case when the velocity is varying.11.5.15Comparison of the Various SchemesBy comparing the stability of the various schemes presented so far, it is clear thatthe most negative (with a coefficient of −3/2) is the one associated with the SOUscheme.
This is followed by the upwind (with a coefficient of −1), then FROMM(with a coefficient of −3/4), after that QUICK(with a coefficient of −3/8), andfinally the central difference scheme (with a coefficient of zero, i.e., a neutralsensitivity to changes in /C ). The self corrective action discussed above is the sumof both convection and diffusion contributions. The false diffusion produced by theupwind scheme adds to its stability and even though its coefficient is −1, is the moststable scheme. This is demonstrated by the profiles presented in Fig.
11.14, whichrepresent the solutions using the various numerical schemes of theconvection-diffusion equation over a domain of length L = 1 and subject to theDirichlet conditions of / ¼ 1 at x = 0 and / ¼ 0 at x = 1. While all solutions arestable at low Péclet number (i.e., Pe = 1, Fig. 11.14a), the CD, FROMM, andQUICK scheme solutions shown in Fig. 11.14b are seen to be wiggly at Pe = 10.At low Pe, the accuracy of the CD and QUICK schemes is comparable and theirsolutions are very close to the exact solution. The least accurate is the solutionproduced by the first order upwind scheme.
The solution of the FROMM scheme is(a)(b)Fig. 11.14 Stability of solutions generated using several convective schemes at two values of cellPéclet number of a 1 and b 1011.5Higher Order Upwind Schemes399more accurate than the SOU scheme, which, in turn, is more accurate than theupwind scheme solution, while both the FROMM and SOU scheme solutions areless accurate than the QUICK scheme solution.At high Pe, the behavior of the schemes changes. The only stable solutions arethe ones obtained by the upwind and SOU schemes with their profiles being almostidentical indicating that the SOU scheme is still highly diffusive.
The CD scheme iswiggly over most of the domain. The FROMM and QUICK schemes, on the otherhand, show wiggles but of smaller amplitudes. The reason for these wiggles is theimposed boundary condition at exit from the domain. As the phenomenon isconvection dominated, it is affected by upstream values. At exit from the domain,the solution faces an imposed value that it has to satisfy, resulting in a largeunexpected gradient causing the over and under shoots to occur. The diffusiveupwind and SOU schemes being based on upstream values only, are smooth and donot show any sign of under/overshoots for all Pe values as the solution is independent of the imposed value at exit.
However the SOU scheme is expected to giverise to oscillations (over/under shoots) in the presence of a high gradient in thedomain like in the presence of a shock wave.11.5.16Functional Relationships for Uniformand Non-uniform GridsThe various interpolation profiles presented so far have been derived on a onedimensional Cartesian grid. Along a curvilinear coordinate axis, the same functionalrelationships can be used by replacing the Cartesian x-axis by a curvilinear f-axis,as shown in Fig.
11.15. For uniform grid, the functional relationships remainexactly the same, independent of whether a Cartesian or a curvilinear grid system isused. For non-uniform grid, the independent variable x appearing in the functionalrelationships should be replaced by the more general independent variable f, whichrepresents distance along the coordinate axis.
If O is the origin (Fig. 11.15), then fUfor example is calculated asOUUUvfCfFig. 11.15 Notation on a curvilinear coordinate axis fDDD40011 Discretization of the Convection TermTable 11.1 Functional relationships for several convection schemes on uniform and non-uniformgridsSchemeUniform gridNon-uniform gridUpwind/f ¼ /C/f ¼ /CDownwind/f ¼ /D/f ¼ /DCD/f ¼ 0:5ð/C þ /D ÞSOU31/f ¼ /C /U22313/f ¼ /C /U þ /D488QUICKFROMM/D /C ff fCfD fC/C /U ff fC/f ¼ /C þfC fUff fU ff fC/f ¼ /U þð/ /U ÞðfD fU ÞðfD fC Þ Dff fU ff fDþð/ /U ÞðfC fU ÞðfC fD Þ Cff fC/f ¼ /C þð/ /U ÞfD fU D/f ¼ /C þ1/f ¼ /C þ ð/D /U Þ4fU ¼ fUU þ ðfU fUU ÞqffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffifUU ¼ ðxUU xO Þ2 þðyUU yO Þ2 þðzUU zO Þ2qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffifU fUU ¼ ðxU xUU Þ2 þðyU yUU Þ2 þðzU zUU Þ2ð11:118ÞTherefore distances are calculated by subdividing the curvilinear line into a numberof straight lines and in general the following applies:qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffif1 f2 ¼ ðx1 x2 Þ2 þðy1 y2 Þ2 þðz1 z2 Þ2ð11:119ÞAdopting a general coordinate system, the functional relationships for the variousschemes presented so far on uniform and non-uniform grids can be derived to be asshown in Table 11.1.11.6Steady Two Dimensional AdvectionThe steady two dimensional advection equation is given byr ðqv/Þ ¼ 0:ð11:120ÞIntegrating over the two-dimensional element of volume VC shown in Fig.
11.16,using the divergence theorem, and replacing the surface integral by a summationover the element faces, Eq. (11.120) becomes11.6Steady Two Dimensional Advection401( x )w( x )eNENNWSn( y )nnSeSwwWeC( y )s( y )CEsSsySWSSE( x )CxFig. 11.16 Notation for a two dimensional Cartesian grid system01X BZC@ J/;C dSA ¼ 0f nbðC Þð11:121ÞfUsing a single Gaussian point for the face integral, the left hand side of Eq. (11.121)is transformed to10X BZX /;CXCJf Sf ¼ðqv/Þf Sf@ J/;C dSA ¼f nbðC Þf nbðC Þfð11:122Þf nbðC ÞSubstitution of Eq. (11.122) in Eq.
(11.121) yieldsXðqv/Þf Sf ¼ 0ð11:123Þf nbðC ÞThe full discretized form of Eq. (11.120) over a Cartesian grid is obtained as01 0 m_ w 1 0 m_ n 1 0 m_ s 1m_ ezfflfflffl}|fflfflffl{zfflfflffl}|fflfflffl{zfflfflffl}|fflfflffl{zfflfflffl}|fflfflffl{@quDy/ A @quDy/ A þ@qvDx/ A @qvDx/A ¼ 0ewnsð11:124Þ40211 Discretization of the Convection TermAdopting an upwind scheme in each coordinate direction by treating the flow aslocally one dimensional, Eq. (11.124) becomesaC /C þ aE /E þ aW /W þ aN /N þ aS /S ¼ 0ð11:125Þwith the coefficients given byaE ¼ FLuxFe ¼ km_ e ; 0kaW ¼ FLuxFw ¼ km_ w ; 0kaN ¼ FLuxFn ¼ km_ n ; 0kaS ¼ FLuxFs ¼ km_ s ; 0kXaC ¼FluxCf ¼ km_ e ; 0k þ km_ w ; 0k þ km_ n ; 0k þ km_ n ; 0kf nbðC Þð11:126Þ¼ ðaE þ aW þ aN þ aS Þ þm_ e þ m_ w þ m_ n þ m_ s|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflPffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}aFF NBðC ÞIf the QUICK scheme is used, the discretized equation is changed intoaC /C þ aE /E þ aW /W þ aEE /EE þ aWW /WWþ aN /N þ aS /S þ aNN /NN þ aSS /SS ¼ 0ð11:127Þwith its coefficients computed as331aE ¼ FluxFe ¼ km_ e ; 0k þ km_ e ; 0k km_ w ; 0k488331aW ¼ FluxFw ¼ km_ w ; 0k þ km_ w ; 0k km_ e ; 0k48811aEE ¼ FluxFee ¼ km_ e ; 0k aWW ¼ FluxFww ¼ km_ w ; 0k88331aN ¼ FluxFn ¼ km_ n ; 0k þ km_ n ; 0k km_ s ; 0k488331aS ¼ FluxFs ¼ km_ s ; 0k þ km_ s ; 0k km_ n ; 0k48811aNN ¼ FluxFnn ¼ km_ n ; 0k aSS ¼ FluxFss ¼ km_ s ; 0k88XXaC ¼FluxCf ¼ aF þ ðm_ e þ m_ w þ m_ n þ m_ s Þf nbðC Þð11:128ÞF NBðC ÞBoth discretized equations are used to solve the pure advection of a step profile inan oblique velocity field problem.
The physical domain is schematically depicted inFig. 11.17a. It represents a square domain with a property / being convected in afield with a velocity v(1,1). The value of / is 1 on the left side of the domain and 011.6Steady Two Dimensional Advection(a)403(b)L1=1vExact solutionUpwind solution0.8QUICK solution=0=1vv0.60.40.2=0v000.20.4y0.60.81Fig. 11.17 a Physical domain and b / profiles along the vertical centerline of the domain for pureadvection of a step profile in an oblique velocity fieldon the bottom. In the absence of any diffusion, the exact solution is / ¼ 1 above thediagonal shown in Fig. 11.17a and / ¼ 0 below it. Figure 11.17b compares theexact profile at x = 0.5 with similar ones obtained numerically using the upwind andQUICK schemes.In comparison with the exact solution, the profile generated by the upwindscheme is smeared and highly inaccurate but very smooth.
This inaccuracy is due toa new type of error known as cross-stream diffusion, which is caused by theone-dimensional interpolation profile used, i.e., it is due to treating the flow aslocally one dimensional. The origin of cross-flow diffusion was identified byPatankar [4] and Stubley [13] as being a multi-dimensional phenomenon. It occursonly when the velocity field is not aligned with the grid.
An approximate expressionfor the cross flow diffusion has been given by de Vahl Davis and Mallinson [14] fortwo dimensional flows asqjvjDxDy sinð2hÞC/false ¼ 4 Dy sin3 ðhÞ þ Dx cos3 ðhÞð11:129Þwhere jvj is the velocity magnitude and h the angle made by the velocity vectorwith the x coordinate axis. This error can be reduced by using higher order interpolation schemes as demonstrated by the profile generated with the QUICKscheme.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
