Moukalled F., Mangani L., Darwish M. The finite volume method in computational fluid dynamics. An advanced introduction with OpenFOAM and Matlab (811443), страница 29
Текст из файла (страница 29)
The process starts by computing the location ofthe geometric centre of the polyhedron element and decomposing it into a numberof polygonal pyramids. As shown in Fig. 6.21, each polygonal pyramid is formedof the geometric centre as the apex and a polygonal face of the element as the base,with its side faces being triangles.For a polygonal pyramid, the volume and centroid are readily computed. Thevolume is calculated as ð1=3Þ Base Height. The base is basically one of thesurfaces of the element, while the sub-element pyramid centroid, measured fromthe centroid of the base, is situated at 1=4 of the line joining the centroid of the baseto the apex of the pyramid.
The volume of the polyhedron element is the sum of thevolumes of the polygonal pyramids. As for the centroid it is computed as theFig. 6.21 A sub-elementpyramidGGd Gfd CffSfSf6.5 Geometric Quantities159volume-weighted average of the centroids of the pyramids. Mathematically this iscomputed from the following relations:xG ¼k1Xxik i¼1ðxCE Þpyramid ¼ 0:75ðxCE Þf þ 0:25ðxG ÞpyramidVpyramid ¼VC ¼dGf Sf3ð6:26ÞXVpyramid SubpyramidsðC ÞPðxCE ÞC ¼ xC ¼ SubpyramidsðC ÞðxCE Þpyramid VpyramidVCSince the centroid and volume of a polygonal pyramid are easily computed, thisprocedure allows accurate computations of both the volume and centroid of an element.6.5.2.3Face Weighting FactorConsider the one dimensional finite volume mesh system shown in Fig. 6.22. Thevalues of ϕ are known at the control volume centroids C and F, and are to be used tocompute the value of ϕ at the interface f.A simple linear interpolation will result in the following formula:/f ¼ gf /F þ 1 gf /Cð6:27Þwheregf ¼dCfdCf þ dfFð6:28ÞThe simplicity of this formula does not extend into multi-dimensional situations as intwo or three dimensions the circumstances become a bit more complicated.
In thiscase there is not a unique option for the definition of the geometric weighting factors.FfCfCFig. 6.22 One dimensional mesh systemF1606 The Finite Volume MeshOne choice would be to base the weighting factor on the respective volumes, suchthatgf ¼VCVC þ VFð6:29ÞThis however yields wrong results in certain cases, such as in the configurationshown in Fig. 6.23.Another difficulty arises when the points C, f, and F are not collinear as depictedin Fig.
6.24a.A better alternative for such cases, as displayed in Fig. 6.24b, is to base theinterpolation on the normal distances to the face, i.e., Cf 0 and Ff 00 . Thus theinterpolation factor is computed asgf ¼dCf efdCf ef þ dfF efð6:30Þwhere ef is the surface unit vector given bySfef ¼ Sf ð6:31ÞFfCFig. 6.23 Axisymmetric grid system(a)F(b)ffCfCFd fFdCf ffefSfFig. 6.24 Two dimensional control volume with the points C, f, and F being non collinear6.5 Geometric Quantities161Example 5Compute the weighing factor using Eqs. (6.26) and (6.27) for the two triangular elements shown in Fig. 6.25 (Table 6.5).3SffdCfd fFFfC214Fig.
6.25 Two neighboring polygonal elementsTable 6.5 Coordinates of the triangular elements for Example 4Node1234xy001.20.41120.1SolutionTo calculate the interpolation factor using Eq. (6.26) the volumes of the twoelements are needed and are computed as1V C ¼ ½ x1 ð y2 y3 Þ þ x2 ð y3 y1 Þ þ x3 ð y1 y2 Þ 21¼ ½0 þ 1:2ð1 0Þ þ 1ð0 0:4Þ ¼ 0:421V F ¼ ½ x2 ð y4 y3 Þ þ x4 ð y3 y2 Þ þ x3 ð y2 y4 Þ 21¼ ½1:2ð0:1 1Þ þ 2ð1 0:4Þ þ 1ð0:4 0:1Þ ¼ 0:422VC0:4¼ 0:4878gf ¼¼VC þ VF 0:4 þ 0:42The calculation of gf based on Eq. (6.27) is more involved.
The centroids ofthe two control volumes are required and are calculated as1xC ¼ ð0 þ 1 þ 1:2Þ ¼ 0:733331xF ¼ ð1:2 þ 1 þ 2Þ ¼ 1:431yC ¼ ð0 þ 1 þ 0:4Þ ¼ 0:466631yF ¼ ð0:4 þ 1 þ 0:1Þ ¼ 0:531626 The Finite Volume MeshThe face centroid is also required and is found to be1xf ¼ ð1 þ 1:2Þ ¼ 1:121yf ¼ ð1 þ 0:4Þ ¼ 0:72The surface vector is calculated asSf ¼ ðy3 y2 Þi ðx3 x2 Þj ¼ 0:6i þ 0:2jThe unit vector normal to the face becomesef ¼Sf0:6i þ 0:2j¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ 0:949i þ 0:316jSf0:62 þ 0:22The distance from the centroids of the control volumes to the face centroid aregiven by dCf ¼ xf xC i þ yf yC j ¼ 0:3667i þ 0:2334jdfF ¼ xF xf i þ yF yf j ¼ 0:3i 0:2jThe interpolation factor using Eq.
(6.27) is found asgf ¼dCf efð0:3667i þ 0:2334jÞ ð0:949i þ 0:316jÞ¼ 0:6556¼dCf ef þ dfF ef ð0:6667i þ 0:0334jÞ ð0:949i þ 0:316jÞthe difference in values is obvious.6.66.6.1Computational PointersuFVMIn ufvm, the processing of all geometric and topological data is performed in oneroutine denoted by cfdProcessOpenFoamMesh, which is executed right afterreading the OpenFOAM® mesh with cfdReadOpenFoamMesh. Reading theOpenFOAM® mesh in its native form, requires reading the various files defining anOpenFOAM® mesh, namely and in that order, points, faces, owners, neighbours,and boundaries files.In cfdProcessOpenFoamMesh the face geometry is initially processed (centroid,area and surface vector), as shown in Listing 6.1.6.6 Computational Pointers163%% Process basic Face Geometry%numberOfFaces = theMesh.numberOfFaces;for iFace=1:numberOfFacesiNodes = theMesh.faces(iFace).iNodes;numberOfiNodes = length(iNodes);%% Compute a rough centre of the face%centre = [0 0 0]';for iNode=iNodescentre = centre + theMesh.nodes(iNode).centroid;endcentre = centre/numberOfiNodes;centroid = [0 0 0]';Sf = [0 0 0]';area = 0;%% using the center compute the area and centroid% of virtual triangles based on the centre and the% face nodes%for iTriangle=1:numberOfiNodespoint1 = centre;point2 = theMesh.nodes(iNodes(iTriangle)).centroid;if(iTriangle<numberOfiNodes)point3 = theMesh.nodes(iNodes(iTriangle+1)).centroid;elsepoint3 = theMesh.nodes(iNodes(1)).centroid;endlocal_centroid = (point1+point2+point3)/3;local_Sf = 0.5*cross(point2-point1,point3-point1);local_area = cfdMagnitude(local_Sf);centroid = centroid + local_area*local_centroid;Sf = Sf + local_Sf;area = area + local_area;endcentroid = centroid/area;%theMesh.faces(iFace).centroid = centroid;theMesh.faces(iFace).Sf = Sf;theMesh.faces(iFace).area = area;endListing 6.1 Processing of basic face geometryThis is followed by processing the basic element geometry (Listing 6.2).
Theelement volume is computed by subdividing it into non-overlapping pyramids,whose geometric characteristics can be easily computed.1646 The Finite Volume Mesh%% compute volume and centroid of each element%numberOfElements = theMesh.numberOfElements;for iElement=1:numberOfElementsiFaces = theMesh.elements(iElement).iFaces;%% Compute a rough centre of the element%centre = [0 0 0]';for iFace=1:length(iFaces)centre = centre + theMesh.faces(iFace).centroid;endcentroid = [0 0 0]';Sf = [0 0 0]';centre = centre/length(iFaces);% using the centre, compute the area and centroid% of virtual triangles based on the centre and the% face nodes%localVolumeCentroidSum = [0 0 0]';localVolumeSum = 0;for iFace=1:length(iFaces)localFace = theMesh.faces(iFaces(iFace));localFaceSign = theMesh.elements(iElement).faceSign(iFace);Sf = localFace.Sf*localFaceSign;Cf = localFace.centroid - centre;% calculate face-pyramid volumelocalVolume = Sf'*Cf/3;% Calculate face-pyramid centrelocalCentroid0.25*centre;=0.75*localFace.centroid+%Accumulate volume-weighted face-pyramid centrelocalVolumeCentroidSum = localVolumeCentroidSum +localCentroid*localVolume;% Accumulate face-pyramid volumelocalVolumeSum = localVolumeSum + localVolume;endcentroid = localVolumeCentroidSum/localVolumeSum;volume = localVolumeSum;%theMesh.elements(iElement).volume = volume;theMesh.elements(iElement).centroid = centroid;endListing 6.2 Processing of basic element geometry6.6.2OpenFOAM®As OpenFOAM® [10] uses an unstructured grid platform, all its geometric entities,such as elements, volumes, areas, and centroids, as well as face weighting factors,have to be evaluated and stored.
This section provides an overview of the geometricrelations needed in the evaluation of the geometric quantities used in OpenFOAM®.6.6 Computational Pointers6.6.2.1165Area and Centroid of FacesIn order to evaluate areas and centers for a generic polygon face, the face isdecomposed into a series of triangular faces. The overall polygon face metrics arecalculated by summing up the properties of each triangular portion.
ThusEqs. (6.21)–(6.23) have to be applied for each computational cell.OpenFOAM® constructs centers of faces and calculates their areas in the file“$FOAM_SRC/OpenFOAM/meshes/primitiveMesh/primitiveMeshFaceCentresAndAreas.C” in which the following function is defined (Listing 6.3):void Foam::primitiveMesh::makeFaceCentresAndAreas(const pointField& p,vectorField& fCtrs,vectorField& fAreas) constListing 6.3 Function used for defining centers and areas of facesThe function has three arguments with the first, defined as const, representingdata read from files, while the second and third arguments designate returned lists ofobjects of dimensions equal to the number of faces in the domain, containing thecenters (fCtrs) and areas (fAreas) of the polygonal faces, respectively.The pointField data is a list of all mesh vertices with each vertex defined usingthree spatial coordinates. The second data needed is the face definition.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
