Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab (811443), страница 65
Текст из файла (страница 65)
11.21 UML Graph for the surfaceInterpolationScheme class, where a box with a black borderdenotes a documented struct or class for which not all inheritance/containment relations are shownwhere the upper and lower vectors are now defined using the interpolation weightsdenoted tweights, as described above. Further, Eq.
(11.163) indicates that theimplementation of the implicit operator in OpenFOAM® is based on the so called“downwind” discretization (as explained later) in which all schemes are discretizedin a fully implicit manner independently of the order of accuracy and without usinga deferred correction approach. This technique is also similar to the downwindweighing factor method that will be described in the next chapter.As described above for both explicit and implicit discretization methods, thedivergence operator for the calculation of the convection term is based on the faceinterpolation and the calculation of the weights. OpenFOAM®performs these tasksin a base class denoted by surfaceInterpolationScheme.
From this base class alarge number of interpolation schemes are derived and implemented, as shown inthe UML graph depicted in Fig. 11.21. For better understanding, additional detailsabout this class are given next.41811 Discretization of the Convection TermAs mentioned above, the fvc::div and fvm::div operators use thesurfaceInterpolationScheme class to perform the needed tasks for each discretization method. In this class, the functions used to calculate the values at the faces(vf) and the weights are displayed in Listings 11.9 and 11.10, respectively, andwhere tinterpScheme_ (Listing 11.10) is a surfaceInterpolationScheme object.template<class Type>tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >gaussConvectionScheme<Type>::interpolate(const surfaceScalarField&,const GeometricField<Type, fvPatchField, volMesh>& vf) const{return tinterpScheme_().interpolate(vf);}Listing 11.9 The interpolation function where the value at the face (vf) is computedandtemplate<class Type>tmp<fvMatrix<Type> >gaussConvectionScheme<Type>::fvmDiv(const surfaceScalarField& faceFlux,const GeometricField<Type, fvPatchField, volMesh>& vf) const{tmp<surfaceScalarField> tweights = tinterpScheme_().weights(vf);const surfaceScalarField& weights = tweights();Listing 11.10 The function used to calculate the weightsThe definition of the class can be found in FOAM_SRC/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/Surface InterpolationScheme.H where the two main functions (Listings 11.11 and 11.12),necessary for the calculation of the divergence operator are defined.//- Return the face-interpolate of the given cell field// with explicit correctionvirtual tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >interpolate(const GeometricField<Type, fvPatchField, volMesh>&) const;Listing 11.11 Definition of the main function for face interpolation11.9Computational Pointers419//- Return the interpolation weighting factors for the given fieldvirtual tmp<surfaceScalarField> weights(const GeometricField<Type, fvPatchField, volMesh>&) const = 0;Listing 11.12 Definition of the main function for weight factors calculationThe interpolate function in Listing 11.9 (here OpenFOAM® adopts the samefunction name but for a different implementation) is used with the explicit operatorfvc::div and is defined as a normal virtual function meaning that a derived classmay or may not adopt it.
For most of the derived classes this function is not selectedand the definition in the surfaceInterpolationScheme class is based on a modifiedform of Eq. (11.163) written as/f ¼ |{z}- /O þ ð1 -Þ /N ¼ /N þ -ð/O /N Þ|fflfflfflffl{zfflfflfflffl}ownercoefficientð11:164ÞneighbourcoefficientTo implement Eq. (11.164), OpenFOAM® uses auxiliary functions with namessimilar to the ones used earlier in the surfaceInterpolationScheme class. First aninterpolate class with a single argument is instantiated from the divergenceoperator. Then inside the interpolate function a new interpolate function with twoarguments is introduced along with a new weights function, using the script shownin Listing 11.13.template<class Type>tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >surfaceInterpolationScheme<Type>::interpolate(const GeometricField<Type, fvPatchField, volMesh>& vf) const{tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsf= interpolate(vf, weights(vf));Listing 11.13 The new interpolate function with two argumentsThe new interpolate function (the one with the two arguments) compute theface value according to Eq.
(11.164). The script used to perform this task is givenby (Listing 11.14)42011 Discretization of the Convection Term//- Return the face-interpolate of the given cell field// with the given weighting factorstemplate<class Type>tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >surfaceInterpolationScheme<Type>::interpolate(const GeometricField<Type, fvPatchField, volMesh>& vf,const tmp<surfaceScalarField>& tlambdas){const surfaceScalarField& lambdas = tlambdas();const Field<Type>& vfi = vf.internalField();const scalarField& lambda = lambdas.internalField();const fvMesh& mesh = vf.mesh();const labelUList& P = mesh.owner();const labelUList& N = mesh.neighbour();GeometricField<Type, fvsPatchField, surfaceMesh>& sf = tsf();Field<Type>& sfi = sf.internalField();for (register label fi=0; fi<P.size(); fi++){sfi[fi] = lambda[fi]*(vfi[P[fi]] - vfi[N[fi]]) + vfi[N[fi]];}Listing 11.14 Script used to compute face values according to Eq.
(11.164)where lambda represents the weight - in Eq. (11.164).As briefly introduced above, the weights function is now a pure virtual functionand it defines the weights of the interpolation. Contrary to the interpolate function,the weights function being pure virtual has to be redefined for every derived class(C++ requirement). Looking again at Fig. 11.20, all derived classes have to definethe weights function. The redefinition must be performed by constructing theproper weights corresponding to the scheme adopted.For the case of the central difference scheme the weights function is defined in aclass named linear<Type> with its script given by (Listing 11.15),11.9Computational Pointers421//Member Functions//- Return the interpolation weighting factorsvirtual tmp<surfaceScalarField> weights(const GeometricField<Type, fvPatchField, volMesh>&) const{return this->mesh().surfaceInterpolation::weights();}Listing 11.15 Script for calculating the weights of the central difference schemewhere surfaceInterpolation::weights() returns the distance weights based on themesh, thus defining the central difference discretization.Another example is the downwind scheme defined by Eq.
11.44, and for whichthe weights function is given by (Listing 11.16)//- Return the interpolation weighting factorsvirtual tmp<surfaceScalarField> weights(const GeometricField<Type, fvPatchField, volMesh>&) const{return neg(faceFlux_);}Listing 11.16 Script for calculating the weights of the downwind schemewhere the neg function returns 0 for positive fluxes and 1 for negative ones, whichis the opposite of the pos function used with the upwind scheme.11.10ClosureThe chapter dealt with the discretization of the convection term.
The difficultiesassociated with the use of a linear symmetrical profile were discussed and remediessuggested. The upwind scheme, although stable and generates physically plausibleresults, is highly diffusive smearing sharp gradients and producing results that arefirst order accurate. Upwind-biased HO convection schemes were shown toimprove accuracy but to produce a new type of error known as the dispersion error,which manifest itself in the solution through wiggles, oscillations, and over/undershoots across sharp gradients. No attempt was made to address this error in thischapter.
Among other issues, this will form the subject of the next chapter thatfurther expands on the discretization of the convection term.42211 Discretization of the Convection Term11.11ExercisesExercise 1(a) Find third order accurate expressions, i.e., OðDx3 Þ, of /U , /C , and /D byusing a Taylor series expansion about the face f for the one dimensionaluniform grid shown in Fig.