Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab (811443), страница 89
Текст из файла (страница 89)
(15.107) with Eq. (15.94), Df is found to beDf ¼Duf Sxf2 2 2þ Dvf Syf þ Dwf SzfyzxdCFDuf Sxf þ dCFDvf Syf þ dCFDwf Szfð15:108ÞAny of the above approaches can be adopted to calculate the value of Df :15.5.4 The Collocated SIMPLE AlgorithmThe sequence of events in the collocated SIMPLE algorithm is displayed inFig. 15.16 and can be summarized as follows:1. To compute the solution at iteration n + 1, start with the solution at time t forpressure, velocity, and mass flow rate fields, i.e., p(n), u(n), and m_ ðnÞ , respectively,as the initial guess.2. Solve the momentum equation given by Eq. (15.70) to obtain a new velocityfield v*.3. Update the mass flow rate at the element faces using the Rhie-Chow interpolation (Eq.
15.100) to compute a momentum satisfying mass flow rate field m_ .4. Using the new mass flow rates assemble the pressure correction equation(Eq. 15.98) and solve it to obtain a pressure correction field p0 .5. With the pressure correction field update the pressure and velocity fields at theelement centroids and the mass flow rate at the element faces to obtaincontinuity-satisfying fields using Eq. (15.101).
These resulting fields aredenoted by u**, m_ , and p*, respectively.6. Treat the obtained solution as a new initial guess, return to step 2 and repeatuntil convergence.7. Set the solution at time n + 1 (i.e., t + Δt) to be equal to the converged solutionand set the current time n + 1 (i.e., t + Δt) to be n (i.e., t).8. Advance to the next time step.9. Return to step 1 and repeat until the last time step is reached.59815 Fluid Flow Computation: Incompressible FlowsFig.
15.16 A flow chart ofthe SIMPLE algorithmset initial guess m( ) ,v ( ) ,n( n)n,and p(n)at timet + t toconverged values at timetassemble and solve momentumequation for v *compute m*f using theRhie Chow interpolationassemble and solve pressurecorrectionequation for prepeatcorrect m*f ,v * ,repeat( n)obtain m**f ,v ** ,set m(f ) = m**f ,v ( ) = v ** ,nnno*,and p( ) ton,and p*( n)=*,and p( ) = p*nconverged ?yesset solution at time t + t tobeequal tothe converged solutionadvanceintimeset t = t + tnotimeexceeded ?yesstopExample 3In the two dimensional problem shown below, the following quantities aregiven pW = 100, pN = 20 and pE = 50 and an inlet condition at face s withvs = 20 and zero pressure gradient.The flow is steady state and the density is uniform of value 1.
Themomentum equations for ue and vn are15.5General Derivation599uC ¼ dx ðpe pw ÞandvC ¼ dy ðpn ps Þwhere the constants dx = 1 and dy = 0.25. The element shown hasΔx = Δy = 1. Use the collocated SIMPLE algorithm to derive the pressurecorrection equation and solve it to find the pressure for element C. TakeðnÞpC ¼ 70 as an initial guess for pressure at C (Fig.
15.17).pN = 20pW = 100 NWpE = 50CEvs = 20Fig. 15.17 The two dimensional domain used in Example 3SolutionAt the inlet the zero gradient condition can be used to computeðnÞps ¼ pC ¼ pC ¼ 70ðnÞðnÞnow compute uC ; vC usinguC ¼ dx peðnÞ pwðnÞ¼ 1ð60 85Þ¼ 25andvC ¼ dy pnðnÞ psðnÞ¼ 0:25ð45 70Þ¼ 6:25ðnÞsince peðnÞ ¼ 0:5 pEðnÞ þ pCðnÞ ; pnðnÞ ¼ 0:5 pEðnÞ þ pNðnÞ and pwðnÞ ¼ 0:5 pEðnÞ þ pWNow the pressure correction equation is constructed by substituting themass flux and its correction into the mass conservation equation. Using theRhie-Chow interpolation the face fluxes are computed as60015 Fluid Flow Computation: Incompressible Flows32766 776 ðnÞðnÞðnÞðnÞðnÞm_ e ¼ ue Dy ¼ ue dx 6 pE pC 0:5 pðnÞppþp7eweee76 |fflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflffl} |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}54pressure differenceacross faceaverage pressure differenceacross face32766h i 776 ðnÞðnÞðnÞðnÞðnÞðnÞðnÞðnÞðnÞ¼ 0:5 dx pðnÞpppdpp0:5pppþdþp76xxECeweewweweee76|fflfflfflfflfflfflfflffl|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}ffl{zfflfflfflfflfflfflfflfflffl}|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}54ueðnÞðnÞ¼ dx pE pCpressure differenceacross faceaverage pressure differenceacross face¼ 1ð50 70Þ¼ 20similarly for the n and w face32m_ n¼vn Dx¼vn766 776 ðnÞðnÞðnÞðnÞðnÞðnÞ dy 6 pN pC 0:5 pn ps þ pnn pn76 |fflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflffl} |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}754pressure differenceacross faceðnÞðnÞ¼ dy pN pCaverage pressure differenceacross face¼ 0:25ð20 70Þ¼ 12:532m_ w¼uw Dy766 776ðnÞðnÞðnÞðnÞðnÞðnÞ¼ uw þ dx 6 pC pW 0:5 pe pw þ pw pww 76 |fflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflffl} |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}754ðnÞðnÞ¼ dx pC pWpressure differenceacross faceaverage pressure differenceacross face¼ 1ð70 100Þ¼ 30with m_ s ¼ 20Dx ¼ 20 the pressure correction equation is constructedfrom m_ e þ m_ 0e þ m_ n þ m_ 0n þ m_ w þ m_ 0w þ m_ s ¼ 0and15.5General Derivation601m_ 0e þ m_ 0n þ m_ 0w ¼ m_ e þ m_ n þ m_ w þ m_ s¼ ð20 þ 12:5 30 20Þ¼ 17:5nowm_ 0e ¼ dx p0E p0Cm_ 0n ¼ dy p0N p0Cm_ 0w ¼ dx p0C p0Wthus we have p0E p0C 0:25 p0N p0C þ p0C p0W ¼ 17:5or2:25p0C p0E 0:25p0N p0W ¼ 17:5if the E, N and W values were exact then we would have2:25p0C ¼ 17:5orp0C ¼ 7:78then we would proceed with correcting the pressure and velocity componentsat C to yieldðnÞpC ¼ pC þ p0C ¼ 77:78u0C ¼ 0 ) u¼ 25 Cv0C ¼ dy p0n p0s¼ 0:25 0:5p0C p0C¼ 0:25ð0:5 7:78 7:78Þ¼ 0:9725sovC ¼ 6:25 þ 0:9725¼ 7:222560215 Fluid Flow Computation: Incompressible Flows15.6 Boundary ConditionsA boundary element has at least one face located at a boundary patch, which isdenoted as a boundary face (Fig.
15.18). The treatment of boundary conditions at aboundary face is critical to the accuracy and robustness of any CFD code. Thus forany pressure-based code to succeed, it is imperative for the boundary conditions ofboth momentum and pressure-correction equations to be properly implemented.This section describes in detail the implementation of a variety of boundaryconditions.A first note of interest is the expression of the Rhie-Chow interpolation at aboundary face, which has to be modified since the averaging cannot be performedin a fashion similar to an interior face.
The average at a boundary face is written interms of the element value ashb ¼ hCð15:109Þwhere b refers to the boundary face centroid and C to the element centroid.Adopting this practice, the averages in the Rhie-Chow interpolation, the velocity,and the mass flow rate at a boundary face becomevb ¼ vCðnÞðnÞð15:110Þrpb ¼ rpCDvb ¼ DvCFig. 15.18 An example of aboundary elementboundaryelementebbSbnboundaryfaceC15.6Boundary Conditions603ðnÞðnÞ¼ vb DvC rpb rpb|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}boundary facestandard RhieChowðnÞðnÞ¼ vC DvC rpb rpC|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}vb|{z}ð15:111Þboundary RhieChowm_ b ¼ qb vb SbhiðnÞðnÞ¼ qb vC DvC rpb rpC SbðnÞðnÞðnÞðnÞ¼ qb vC Sb qb DC pb pC qb DvC rpb Tb DvC rpC Sbð15:112ÞIn what follows the implementation of boundary conditions is first presented forthe momentum equation, followed by the implementation of the boundary conditions for the pressure (correction) equation.
For the cases when the boundaryconditions for the momentum and pressure correction equations are co-dependent,their full treatment is detailed in the pressure correction equation section.15.6.1 Boundary Conditions for the Momentum EquationThe semi-discretized form of the momentum equation can be expressed asX X X ðqvÞC ðqvÞCVC þm_ f vf ¼ pf Sf þsf Sf þ |{z}BDt|fflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflffl} f nbðCÞf nbðCÞf nbðCÞelement|fflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflffl}|fflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflffl} |fflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflffl} discretizationelement discretizationface discretizationface discretizationface discretizationð15:113Þwhere the discretization type of the various terms is explicitly stated. Terms that areevaluated at the element faces should be modified along a boundary face inaccordance with the type of boundary condition used. Therefore, for boundaryelements, the terms evaluated at the element faces are written asX Xð15:114Þm_ f vf ¼m_ f vf þ m_ b vb|ffl{zffl}f nbðCÞf interior nbðCÞboundary face|fflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflffl}face discretization60415 Fluid Flow Computation: Incompressible FlowsX sf Sf ¼f nbðCÞXf interior nbðCÞ|fflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflffl}face discretization¼Xsf Sf þ sb Sb|fflffl{zfflffl}boundary facesf Sf þf interior nbðCÞX pf Sf ¼f nbðCÞ|fflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflffl}Xð15:115ÞFb|{z}boundaryfacepf S f þf interior nbðCÞpb Sb|ffl{zffl}ð15:116Þboundary faceface discretizationwhere subscript b refers to a value at the boundary.
As previously stated, the pressureterm may be discretized following either an element or a face discretization.PIn eithercase the expanded form is the same since VC ðrpÞC is calculated aspf S ff nbðCÞimplying that the value at the boundary is always required.
To show the wayboundary pressure affects solution, the expanded form of the pressure gradient (i.e.,face discretization) is adopted in the implementation of boundary conditions.15.6.1.1 Wall Boundary ConditionsGenerally a no-slip or a slip boundary condition may be applied to a moving or astationary wall. The implementation involves computing and linearizing the shearstress at the wall. This is different from the Dirichlet condition, though for acartesian grid the two conditions may be shown to be identical.No-Slip Wall Boundary ðpb ¼ ?; m_ b ¼ 0; vb ¼ vwall ÞA no slip boundary condition implies that the velocity of the fluid at the wall vb isequal to the velocity of the wall vwall.
For a stationary wall, the boundary velocityvb is zero. The known velocity at the wall could be wrongly viewed as a Dirichletboundary condition. However this is not really the case, since what is needed is tohave a zero normal boundary flux while also accounting for the shear stress. Thiscannot be satisfied by simply setting vb = vwall. Figure 15.19 shows that this can beguaranteed by ensuring a shear stress that is tangential to the wall in addition to aboundary velocity equation that is equal to the wall velocity. The force Fb exertedby the wall on the fluid can be written asFb ¼ F? þ Fkð15:117Þ15.6Boundary Conditions605vCebvCvebbwallSbvbwallnSbnFig.
15.19 A schematic of a no-slip wall boundary conditionwhere Fk is in the tangential direction to the wall and F? is in the normal direction,which is supposed to be zero. ThusFb ¼ Fk ¼ swall Sbð15:118Þwhere swall is the shear stress exerted by the wall on the fluid given byswall ¼ l@vk:@d?ð15:119ÞIn Eq. (15.119) vk is the velocity vector in the direction parallel to the wall andd? is the normal distance from the centroid of the boundary element to the wallcomputed asSb¼ nx i þ ny j þ nz kSbdCb Sbd? ¼ dCb n ¼Sbn¼ð15:120Þand n the wall normal unit vector. The velocity vector vk can be written as thedifference between v and its normal component as8 9< u unx þ vny þ wnz nx =vk ¼ v ðv nÞn ¼ v unx þ vny þ wnz ny:;w unx þ vny þ wnz nzð15:121Þ60615 Fluid Flow Computation: Incompressible FlowsFrom Eq.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
