Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab (Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab.pdf), страница 102

The boundaryconditions are defined in the directory “derivedFvPatchFields” and are next presented. For a better understanding, it may be beneficial to read Chap. 18 prior togoing over the implementation process presented below.• totalPressureComp: This implements the total pressure condition at subsonicinlet. For that purpose the updateCoeffs() function is modified as shown inListing 16.16, which indicates that after gathering the necessary data from thesolver, the boundary static pressure is computed using Eq. (16.40) and theobtained values are stored in the defined newp variable.const fvPatchScalarField& TB =patch().lookupPatchField<volScalarField, scalar>("T");const fvPatchField<scalar>& RB =patch().lookupPatchField<volScalarField, scalar>("RGas");const fvPatchField<scalar>& gammaB =patch().lookupPatchField<volScalarField, scalar>("gammaGas");const fvsPatchField<scalar>& sphi =patch().lookupPatchField<surfaceScalarField, scalar>("mDot");const fvPatchField<scalar>& srho =patch().lookupPatchField<volScalarField, scalar>("rho");const fvPatchField<vector>& UF =patch().lookupPatchField<volVectorField, vector>("U");scalarField newp = max(min(p0_/pow((scalar(1.0) + (gammaB scalar(1.0) )*(sqr(sphi)/(sqr(srho) * sqr(patch().magSf())))/(scalar(2.0) * gammaB * RB * TB )), gammaB/(gammaB - scalar(1.0))), p0_),SMALL);operator==(newp);Listing 16.16 Modified updateCoeffs() function for implementing the total pressure boundarycondition at inlet• totalPressureCorrectorComp: This is also needed with a total inlet pressureboundary condition.

Its role is to deal with the mass flow rate correction at inletaffecting the diagonal coefficient of the boundary element, as expressed by68416Fluid Flow Computation: Compressible FlowsEqs. (16.43) and (16.44). In the implementation, the contributions to the coefficients are reset to zero, as can be inferred from Listing 16.17, except for thegradientInternalCoeffs() that is required to return the correct values. This ensuresthat for both the divergence and laplacian operator one boundary condition isapplied.tmp<Field<scalar> >totalPressureCorrectorCompFvPatchScalarField::valueInternalCoeffs(const tmp<scalarField>&) const{return tmp<Field<scalar> >(new Field<scalar>(this->size(), pTraits<scalar>::zero));}tmp<Field<scalar> >totalPressureCorrectorCompFvPatchScalarField::gradientBoundaryCoeffs()const{return tmp<Field<scalar> >(new Field<scalar>(this->size(), pTraits<scalar>::zero));}tmp<Field<scalar> >totalPressureCorrectorCompFvPatchScalarField::valueBoundaryCoeffs(const tmp<scalarField>&) const{return tmp<Field<scalar> >(new Field<scalar>(this->size(), pTraits<scalar>::zero));}tmp<Field<scalar> >totalPressureCorrectorCompFvPatchScalarField::gradientInternalCoeffs()const{const fvsPatchField<scalar>& srUA =patch().lookupPatchField<surfaceScalarField, scalar>("srUA");const fvPatchField<scalar>& srho =patch().lookupPatchField<volScalarField, scalar>("rho");return( deltaM()/(-srUA*patch().magSf()*srho) ); //to remove thelaplacian operator (see gaussLaplacianScheme)}Listing 16.17 Script used to reset the diagonal coefficients and to return only the diffusioncontribution or a total pressure boundary condition at inlet16.9Computational Pointers685The modified diagonal coefficient of the pressure correction equation for aboundary element is computed in Listing 16.18.

Here the function deltaM()implements Eq. (16.44) as is, while the diffusion term “(-srUA*patch().magSf()*srho)” is removed from the coefficient of the pressure correction equation, asdefined by the laplacian operator (see Chap. 8, computational pointers).Field<scalar> totalPressureCorrectorCompFvPatchScalarField::deltaM()const{const fvPatchField<scalar>& T =patch().lookupPatchField<volScalarField, scalar>("T");const fvPatchField<scalar>& srho =patch().lookupPatchField<volScalarField, scalar>("rho");const fvPatchField<scalar>& RB =patch().lookupPatchField<volScalarField, scalar>("RGas");const fvsPatchField<scalar>& srUA =patch().lookupPatchField<surfaceScalarField, scalar>("srUA");const fvsPatchField<scalar>& sphi =patch().lookupPatchField<surfaceScalarField, scalar>("mDot");scalarField Dp = patch().magSf()*patch().deltaCoeffs()*srUA;scalarField coeff = srho*Dp/(scalar(1.0) + srho*Dp*Cu() - (sphi/srho)*Cu()/(RB*T));return (coeff);}Listing 16.18 Script used to modify the diagonal coefficient of the boundary element for a totalpressure boundary condition at inlet• totalTemp: This function implements for the energy equation the total temperature boundary condition at a subsonic inlet.

The idea is to impose a statictemperature as a Dirichlet boundary condition using Eq. (16.45). For that purpose the updateCoeffs() function is modified as in Listing 16.19.68616Fluid Flow Computation: Compressible Flowsconst fvPatchField<vector>& Up =patch().lookupPatchField<volVectorField, vector>("U");const fvPatchField<scalar>& gammaB =patch().lookupPatchField<volScalarField, scalar>("gammaGas");scalarField gM1ByG = (gammaB - 1.0)/gammaB;const fvPatchScalarField& TB =patch().lookupPatchField<volScalarField, scalar>("T");const fvPatchField<scalar>& RB =patch().lookupPatchField<volScalarField, scalar>("RGas");scalarField psip = scalar(1.0)/(RB * TB);operator==(T0_/(1.0 + 0.5*psip*gM1ByG*magSqr(Up)));Listing 16.19 Modified updateCoeffs() function for implementing the total temperature boundarycondition at a subsonic inletAt each iteration the temperature value is updated based on the inlet totaltemperature and the boundary velocity.• totalVelocity: This function, described in Listing 16.20, implements the updatesto the velocity field required with total conditions applied at a subsonic inlet.

Forthat purpose a Dirichlet boundary condition is used, with velocity values iteratively computed based on the calculated fluxes “mDot” at the boundary itself.The algorithm is based on grabbing the flux field “mDot” (updated with newvalues after solving the pressure correction equation by invoking the “flux()”function) and dividing the flux by the face area and the corresponding density.const fvsPatchField<scalar>& sphi =patch().lookupPatchField<surfaceScalarField, scalar>("mDot");const fvPatchField<scalar>& rhop =patch().lookupPatchField<volScalarField, scalar>("rho");vectorField n = patch().nf();scalarField ndmagS = (n & inletDir())*patch().magSf();scalarField clip = neg(sphi);scalarField newvel = (sphi*clip)/(rhop*ndmagS);operator==(inletDir()*newvel););Listing 16.20 Updating the velocity at a subsonic inlet for a total pressure boundary condition16.9Computational Pointers687Additionally, a clipping variable is introduced in order to prevent any “backflow” at the inlet.

Here the “clip” variable may assume either a value of zero or one.In fact the “neg” function returns 1 and 0 for negative and positive values,respectively, preventing outward velocities to be accepted. On the other hand, the“inletDir” variable gives the velocity direction at inlet, as defined by the user.16.10 ClosureIn this chapter the incompressible segregated pressure based approach developed inthe previous chapter was extended to handle compressible fluid flow at all speeds.This involved modifying the pressure correction equation to include a convectionlike term that changes its type from elliptic to hyperbolic.

It also required alterationsto the momentum equations, the solution of the energy equation, as well as theaddition of an equation of state. Just as critical, are the special boundary conditionsneeded in the simulation of compressible flows. A number of boundary conditionswere presented as well as some implementation details.The needed modifications to the base incompressible code represent a relativelysmall change to the bulk of any code and yet allow a drastic extension of itscapabilities. The next chapter will present the additional techniques needed fordealing with the time averaged Navier-Stokes equations required for solving turbulent flow problems.16.11 ExercisesExercise 1 A portion of a gas-supply system is shown in Fig.

16.4. The mass flowrate m_ in a pipe section is given bym_ ¼ qCDpwhere Dp is the pressure drop over the length of the pipe section, q is the gasdensity, and C is the gas conductance. The following data is known:p1 ¼ 400; p2 ¼ 350m_ F ¼ 25CA ¼ CC ¼ 1:2; CB ¼ 1:4; CD ¼ 1:6; CE ¼ 1:8with the density related to the pressure via p ¼ qR0 with R0 ¼ 2000.If the direction of the flow is as shown in the figure, find p3 ; p4 ; p5 ; m_ A ; m_ B ;m_ C ; m_ D and m_ E using the following procedure:688••••16Fluid Flow Computation: Compressible FlowsStart with a guess for p3 ; p4 ; and p5 .Compute m_ values based on the guessed pressures and densities.Construct the pressure-correction equations and solve for p03 ; p04 and p05 .Update the pressures and the m_ valuesDo you need to iterate? Why?Fig.

16.4 A portion of a gassupply systemp2 = 350mBp3mEmF = 25mCp5mDp1 = 400p4mAExercise 2 Consider the flow of an ideal gas in a converging<b>Текст обрезан, так как является слишком большим</b>.