Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab (811443), страница 100
Текст из файла (страница 100)
For compressible flow, the conditions to be imposed are dictated by theMach number values. For an inviscid flow, the mathematical type of the equationschanges from elliptic to hyperbolic as the flow changes from subsonic to supersonic. Details regarding their implementation are given next.16.8Boundary Conditions66916.8.1 Inlet Boundary ConditionsAt the inlet to a domain the flow may be subsonic or supersonic necessitatingdifferent treatment since the flow equations may be either of the elliptic or thehyperbolic type.16.8.1.1 Subsonic Flow at InletAt subsonic speeds several conditions at inlet to a domain can be imposed. Thisinclude specified velocity, specified static pressure and velocity direction, orspecified stagnation pressure and velocity direction.
The last type should be usedwhen transition to supersonic speed occurs within the domain.Specified Velocityðpb ¼ ?; m_ b ¼ ?; vb specified ÞUnlike incompressible flow, since for compressible flow the density depends onpressure,the mass flux remains unknown even with a specified velocity at inleti:e, m_ 0b ¼ q0b vb Sb 6¼ 0 . At an inlet boundary the coefficient multiplying thepressure correction p0b is given by0apb ¼ Cq;bm_ bqbð16:35ÞFor implementation in the pressure correction equation, p0b is expressed in termsof internal nodes and the coefficients at these nodes modified accordingly. For theconstant profile case (i.e., pb ¼ pC ), the aC coefficient is obtained as!XXVCCm_ bCqq;f_þaCp ¼;0qDþCþmfq;bffDtqfqbf nbðC Þf nbðC Þ|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} boundary|fflfflffl{zfflfflffl}face contribution0interior faces contributionð16:36ÞThe value of the pressure pb is again obtained by extrapolation from the interioras explained in the incompressible flow chapter.Specified Static Pressure and Velocity Direction pb ¼ pspecified ; ev specified;m_ b ¼ ?; vb ¼ ?ÞIn the case of a specified static pressure at inlet, pb is known and thus p0b is set tozero and consequently q0b is also zero.
Therefore the implementation is similar to theincompressible case with the inlet treated as a Dirichlet boundary condition.Knowing the velocity direction, its magnitude is computed as in the incompressible67016Fluid Flow Computation: Compressible Flowscase using Eq. (16.33) leading to an equation similar to Eq. (15.137). The coefficient of the pressure-correction equation becomes0aCp!XXVC CqCq; f þ¼qf Df þqb DCm_ f ; 0 þDtq|fflffl{zfflffl}ff nbðCÞf nbðCÞ|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} boundary face contributioninterior faces contributionð16:37ÞSpecified Total Pressure and Velocity Direction ðpo; b ¼ po; specified ; ev specified;m_ b ¼ ?; vb ¼ ?ÞFor this case, the magnitude of the velocity and the pressure at the boundary areunknown while related through the total pressure equation given byc 1 2 c=ðc1ÞMbpo; b ¼ pb 1 þ2ð16:38Þwhere b refers to the boundary, po; b is the total pressure, pb the static pressure, c theratio of specific heats, and Mb the Mach number which is equivalent toMb ¼rffiffiffiffiffiffiffiffiffiffiffiffiffivb vbcRTbð16:39ÞEquation (16.37) can be rearranged to give the static pressure in terms of thetotal pressure aspb ¼ po; b!c=ðc1Þ 2m_ bc11þ 2 q 2 ðev nSb Þ2 cRTbð16:40Þbwhere ev is the unit vector in the direction of the velocity vector.
DifferentiatingEq. (16.40) with respect to m_ b givesÞ! ðð2c1 2c1Þm_ bdpbcm_ b po; bc1¼ 21þð16:41Þ 2 q 2 ðev nSb Þ2 cRTbd m_ bq ðev nSb Þ2 cRTbbbSubstituting Eq. (16.41) into Eq. (15.163) an equation for pressure correctionfunction of the mass flux correction is obtained as16.8Boundary Conditions671cm_ b po;bp_ 0b ¼ 2qb ðev nSb Þ2 cRTbÞ! ðð2c1 2c1Þm_ bc1m_ 0b1þ 2 q 2 ðev nSb Þ2 cRTbb¼ cb m_ 0bð16:42ÞReplacing p0b in Eq. (16.34) by its equivalent expression given by Eq. (16.42),the mass flux correction becomesm_ 0b ¼qb Db p0cm_ b1 þ qb Db cb q Cq;b cbð16:43ÞbThe modified boundary cell coefficient is obtained by substituting m_ 0b fromEq. (16.43) in the expanded continuity equation and is given by0aCp!X Cq; f XVC Cqqb DC þ¼qf Df þm_ f ; 0 þm_Dtqf1 þ qb DC cb qb Cq;b cbf nbðCÞf nbðC Þb|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}interior faces contributionboundary face contributionð16:44ÞIt should be mentioned that the boundary condition for the energy equation at asubsonic inlet is usually either a specified static temperature Tb or a specifiedstagnation temperature To;b .
If the static temperature is specified, then this is similarto a Dirichlet condition. If the stagnation temperature is specified then at eachiteration the value of the static temperature is extracted from the stagnation temperature equation usingTo;b ¼ Tb þvb vb2cpð16:45Þand the obtained value treated as known. Thus a Dirichlet-type boundary conditionis also applied.16.8.1.2 Supersonic Flow at InletSpecified static pressure, velocity, and temperaturepb ¼ pspecified ; vb ¼ vspecified ; T ¼ TspecifiedAt a supersonic inlet, values for all variables must be specified (pressure, velocity,and temperature).
This is equivalent to a Dirichlet-type condition implying that0m_ 0b ¼ p0b ¼ 0. Therefore the aCp coefficient of the boundary cell is found to be672160apCFluid Flow Computation: Compressible Flows!XXVC CqCq; f þ¼qf Dfm_ f ; 0 þDtqff nbðC Þf nbðC Þ|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}ð16:46Þinterior faces contribution16.8.2 Outlet Boundary Conditions16.8.2.1 Subsonic Flow at OutletSpecified Pressure pb ¼ pspecified ; m_ b ?; vb ¼ ?At a subsonic outlet, the pressure is usually prescribed.
Therefore the pressurecorrection p0b is set to zero while the mass flow rate correction m_ 0b is computed asm_ 0b¼qb DCp0bp0C m_ bþCq;b p0b ¼ qb DC p0Cqbð16:47ÞSince the velocity vb is not known, it is customary to assume its direction to bethat of the upwind velocity vC . The expression of the ac coefficient in thepressure-correction equation may be written as0aCp!XXVC CqCq; f þ¼qf Df þqb Dcm_ f ; 0 þ|ffl{zffl}Dtqff nbðCÞf nbðCÞ|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} boundary face contributioninterior faces contributionð16:48ÞFor the energy equation a zero flux Neumann-type boundary condition isapplied.Specified Mass Flow Rate m_ b ¼ m_ specified ; pb ? vb ?For a specified mass flow rate at outlet, m_ 0b is zero and is simply dropped from thepressure correction equation with no modifications required for the coefficients ofthe boundary elements. By setting m_ 0b to zero in Eq.
(16.34), an expression for thepressure correction at the boundary as a function of the pressure correction at theboundary cell centroid is obtained as16.8Boundary Conditions673p0b ¼qb DC p0Cm_qb DC qb Cq;bð16:49Þballowing the boundary pressure and density to be computed.
For the energyequation a zero flux Neumann-type boundary condition is applied.16.8.2.2 Supersonic Flow at OutletAt a supersonic outlet none of the variables should be specified and the values ofpressure, velocity, density, and temperature are extrapolated from the interior of thedomain. Thus both m_ b and pb are extrapolated from interior cells. This is equivalentto applying a Neumann boundary condition on pressure-correction, leading to thefollowing modified aC coefficient0apC! XXVC C qCq;f m_ b _ f ; 0 þþ¼q f Df þð16:50ÞCq;b mDtqqbff nbðC Þf nbðCÞ|fflfflfflfflfflffl{zfflfflfflfflfflffl}|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} boundary face contributioninterior faces contribution16.9 Computational Pointers16.9.1 uFVMFor compressible flows a major modification to the algorithm arises in the assemblyof the pressure equation through the inclusion of the convection like term.
This isadded to the cfdAssembleMdotTerm shown in Listing 16.1.% assemble terms X (for compressible flow)%(mdot_f/density_f)*( rho/ p) P'% where mdot_f is the newly computed mdot_f%local_mdot_f = local_FLUXCf_1*(pressure[iElement1]+ Pref) +local_FLUXCf_2*(pressure[iElement2]+ Pref) + local_FLUXVf;%local_FLUXCf1 = local_FLUXCf1 + max(local_mdot_f/density_f(iFace),0.0)*drhodp_f(iFace);local_FLUXCf2 = local_FLUXCf2 - max(-local_mdot_f/density_f(iFace),0.0)*drhodp_f(iFace);local_FLUXVf = local_FLUXVf -(max(local_mdot_f/density_f(iFace),0.0)*drhodp_f(iFace)*(pressure(iElement1)+ Pref)- max(-local_mdot_f/density_f(iFace),0.0)*drhodp_f(iFace)*(pressure(iElement2)+ Pref));local_FLUXVf+= -local_mdot_f;Listing 16.1 Script used to assemble the additional terms for the compressible pressure correctionequation67416Fluid Flow Computation: Compressible FlowsThe other important change is in the treatment of boundary conditions that nownecessitates accounting for a variable density field.
For example the supersonicinlet condition is implemented for the pressure correction equation as shown inListing 16.2.% assemble terms X (for compressible flow)%(mdot_f/density_f)*(?rho/?p) P'% where mdot_f is the newly computed mdot_f%local_mdot_f = local_FLUXCf_1*(pressure[iElement1]+ Pref) +local_FLUXCf_2*(pressure[iElement2]+ Pref) + local_FLUXVf;%local_FLUXCf1 = local_FLUXCf1 + max(local_mdot_f/density_f(iFace),0.0)*drhodp_f(iFace);local_FLUXCf2 = local_FLUXCf2 - max(-local_mdot_f/density_f(iFace),0.0)*drhodp_f(iFace);local_FLUXVf = local_FLUXVf -(max(local_mdot_f/density_f(iFace),0.0)*drhodp_f(iFace)*(pressure(iElement1)+ Pref)- max(-local_mdot_f/density_f(iFace),0.0)*drhodp_f(iFace)*(pressure(iElement2)+ Pref));local_FLUXVf+= -local_mdot_f;Listing 16.2 Implementation of a supersonic inlet condition16.9.2 OpenFOAM®In this section simpleFoam is extended to handle compressible fluid flow at allspeeds.
This entails the following modifications to simpleFoam: (i) the addition ofthe energy equation to be solved simultaneously with the continuity and momentumequations, (ii) the use of an equation of state relating density to temperature andpressure, (iii) and the introduction of the necessary modifications to the pressurecorrection equation and to a number of boundary conditions.
The resulting code isdenoted simpleFoamCompressible with many of the extensions added in the formof supplemental include files as shown in Listing 16.3.16.9Computational Pointers675The upwind.H and gaussConvectionScheme.H are used to force an upwinddiscretization on the convective term of the pressure correction equation.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
