All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
fourthGrad.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "fourthGrad.H"
27 #include "leastSquaresGrad.H"
28 #include "gaussGrad.H"
29 #include "fvMesh.H"
30 #include "volMesh.H"
31 #include "surfaceMesh.H"
32 #include "GeometricField.H"
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 template<class Type>
39 <
41 >
43 (
44  const VolField<Type>& vsf,
45  const word& name
46 ) const
47 {
48  // The fourth-order gradient is calculated in two passes. First,
49  // the standard least-square gradient is assembled. Then, the
50  // fourth-order correction is added to the second-order accurate
51  // gradient to complete the accuracy.
52 
53  typedef typename outerProduct<vector, Type>::type GradType;
54 
55  const fvMesh& mesh = vsf.mesh();
56 
57  // Assemble the second-order least-square gradient
58  // Calculate the second-order least-square gradient
59  tmp<VolField<GradType>> tsecondfGrad
60  = leastSquaresGrad<Type>(mesh).grad
61  (
62  vsf,
63  "leastSquaresGrad(" + vsf.name() + ")"
64  );
65  const VolField<GradType>& secondfGrad =
66  tsecondfGrad();
67 
68  tmp<VolField<GradType>> tfGrad
69  (
71  (
72  name,
73  secondfGrad
74  )
75  );
76  VolField<GradType>& fGrad = tfGrad.ref();
77 
78  const vectorField& C = mesh.C();
79 
80  const surfaceScalarField& lambda = mesh.weights();
81 
82  // Get reference to least square vectors
83  const leastSquaresVectors& lsv = leastSquaresVectors::New(mesh);
84  const surfaceVectorField& ownLs = lsv.pVectors();
85  const surfaceVectorField& neiLs = lsv.nVectors();
86 
87  // owner/neighbour addressing
88  const labelUList& own = mesh.owner();
89  const labelUList& nei = mesh.neighbour();
90 
91  // Assemble the fourth-order gradient
92 
93  // Internal faces
94  forAll(own, facei)
95  {
96  Type dDotGradDelta = 0.5*
97  (
98  (C[nei[facei]] - C[own[facei]])
99  & (secondfGrad[nei[facei]] - secondfGrad[own[facei]])
100  );
101 
102  fGrad[own[facei]] -= lambda[facei]*ownLs[facei]*dDotGradDelta;
103  fGrad[nei[facei]] -= (1.0 - lambda[facei])*neiLs[facei]*dDotGradDelta;
104  }
105 
106  // Boundary faces
107  forAll(vsf.boundaryField(), patchi)
108  {
109  if (secondfGrad.boundaryField()[patchi].coupled())
110  {
111  const fvsPatchVectorField& patchOwnLs =
112  ownLs.boundaryField()[patchi];
113 
114  const scalarField& lambdap = lambda.boundaryField()[patchi];
115 
116  const fvPatch& p = vsf.boundaryField()[patchi].patch();
117 
118  const labelUList& faceCells = p.faceCells();
119 
120  // Build the d-vectors
121  vectorField pd(p.delta());
122 
123  const Field<GradType> neighbourSecondfGrad
124  (
125  secondfGrad.boundaryField()[patchi].patchNeighbourField()
126  );
127 
128  forAll(faceCells, patchFacei)
129  {
130  fGrad[faceCells[patchFacei]] -=
131  0.5*lambdap[patchFacei]*patchOwnLs[patchFacei]
132  *(
133  pd[patchFacei]
134  & (
135  neighbourSecondfGrad[patchFacei]
136  - secondfGrad[faceCells[patchFacei]]
137  )
138  );
139  }
140  }
141  }
142 
143  fGrad.correctBoundaryConditions();
145 
146  return tfGrad;
147 }
148 
149 
150 // ************************************************************************* //
static const Foam::dimensionedScalar C("C", Foam::dimTemperature, 234.5)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Generic GeometricField class.
virtual tmp< VolField< typename outerProduct< vector, Type >::type > > calcGrad(const VolField< Type > &vsf, const word &name) const
Return the gradient of the given field to the gradScheme::grad.
A class for managing temporary objects.
Definition: tmp.H:55
volScalarField scalarField(fieldObject, mesh)
volVectorField vectorField(fieldObject, mesh)
label patchi
U correctBoundaryConditions()
dimensionedScalar lambda(viscosity->lookup("lambda"))
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
SurfaceField< scalar > surfaceScalarField
SurfaceField< vector > surfaceVectorField
UList< label > labelUList
Definition: UList.H:65
fvsPatchField< vector > fvsPatchVectorField
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
volScalarField & p