gaussLaplacianScheme.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 "gaussLaplacianScheme.H"
27 #include "surfaceInterpolate.H"
28 #include "fvcDiv.H"
29 #include "fvcGrad.H"
30 #include "fvMatrices.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 namespace fv
40 {
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 template<class Type, class GType>
45 tmp<fvMatrix<Type>>
47 (
48  const surfaceScalarField& gammaMagSf,
49  const surfaceScalarField& deltaCoeffs,
50  const VolField<Type>& vf
51 )
52 {
53  tmp<fvMatrix<Type>> tfvm
54  (
55  new fvMatrix<Type>
56  (
57  vf,
58  deltaCoeffs.dimensions()*gammaMagSf.dimensions()*vf.dimensions()
59  )
60  );
61  fvMatrix<Type>& fvm = tfvm.ref();
62 
63  fvm.upper() = deltaCoeffs.primitiveField()*gammaMagSf.primitiveField();
64  fvm.negSumDiag();
65 
67  {
68  const fvPatchField<Type>& pvf = vf.boundaryField()[patchi];
69  const fvsPatchScalarField& pGamma = gammaMagSf.boundaryField()[patchi];
70  const fvsPatchScalarField& pDeltaCoeffs =
71  deltaCoeffs.boundaryField()[patchi];
72 
73  if (pvf.coupled())
74  {
75  fvm.internalCoeffs()[patchi] =
76  pGamma*pvf.gradientInternalCoeffs(pDeltaCoeffs);
77  fvm.boundaryCoeffs()[patchi] =
78  -pGamma*pvf.gradientBoundaryCoeffs(pDeltaCoeffs);
79  }
80  else
81  {
82  fvm.internalCoeffs()[patchi] = pGamma*pvf.gradientInternalCoeffs();
83  fvm.boundaryCoeffs()[patchi] = -pGamma*pvf.gradientBoundaryCoeffs();
84  }
85  }
86 
87  return tfvm;
88 }
89 
90 
91 template<class Type, class GType>
94 (
95  const surfaceVectorField& SfGammaCorr,
96  const VolField<Type>& vf
97 )
98 {
99  const fvMesh& mesh = this->mesh();
100 
101  tmp<SurfaceField<Type>> tgammaSnGradCorr
102  (
104  (
105  "gammaSnGradCorr("+vf.name()+')',
106  mesh,
107  SfGammaCorr.dimensions()
108  *vf.dimensions()*mesh.deltaCoeffs().dimensions()
109  )
110  );
111 
112  for (direction cmpt = 0; cmpt < pTraits<Type>::nComponents; cmpt++)
113  {
114  tgammaSnGradCorr.ref().replace
115  (
116  cmpt,
117  fvc::dotInterpolate(SfGammaCorr, fvc::grad(vf.component(cmpt)))
118  );
119  }
120 
121  return tgammaSnGradCorr;
122 }
123 
124 
125 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
126 
127 template<class Type, class GType>
128 tmp<VolField<Type>>
130 (
131  const VolField<Type>& vf
132 )
133 {
134  const fvMesh& mesh = this->mesh();
135 
136  tmp<VolField<Type>> tLaplacian
137  (
138  fvc::div(this->tsnGradScheme_().snGrad(vf)*mesh.magSf())
139  );
140 
141  tLaplacian.ref().rename("laplacian(" + vf.name() + ')');
142 
143  return tLaplacian;
144 }
145 
146 
147 template<class Type, class GType>
150 (
151  const SurfaceField<GType>& gamma,
152  const VolField<Type>& vf
153 )
154 {
155  const fvMesh& mesh = this->mesh();
156 
157  const surfaceVectorField Sn(mesh.Sf()/mesh.magSf());
158 
159  const surfaceVectorField SfGamma(mesh.Sf() & gamma);
160  const SurfaceField<scalar> SfGammaSn
161  (
162  SfGamma & Sn
163  );
164  const surfaceVectorField SfGammaCorr(SfGamma - SfGammaSn*Sn);
165 
166  tmp<fvMatrix<Type>> tfvm = fvmLaplacianUncorrected
167  (
168  SfGammaSn,
169  this->tsnGradScheme_().deltaCoeffs(vf),
170  vf
171  );
172  fvMatrix<Type>& fvm = tfvm.ref();
173 
174  tmp<SurfaceField<Type>> tfaceFluxCorrection
175  = gammaSnGradCorr(SfGammaCorr, vf);
176 
177  if (this->tsnGradScheme_().corrected())
178  {
179  tfaceFluxCorrection.ref() +=
180  SfGammaSn*this->tsnGradScheme_().correction(vf);
181  }
182 
183  fvm.source() -= mesh.V()*fvc::div(tfaceFluxCorrection())().primitiveField();
184 
185  if (mesh.schemes().fluxRequired(vf.name()))
186  {
187  fvm.faceFluxCorrectionPtr() = tfaceFluxCorrection.ptr();
188  }
189 
190  return tfvm;
191 }
192 
193 
194 template<class Type, class GType>
197 (
198  const SurfaceField<GType>& gamma,
199  const VolField<Type>& vf
200 )
201 {
202  const fvMesh& mesh = this->mesh();
203 
204  const surfaceVectorField Sn(mesh.Sf()/mesh.magSf());
205  const surfaceVectorField SfGamma(mesh.Sf() & gamma);
206  const SurfaceField<scalar> SfGammaSn
207  (
208  SfGamma & Sn
209  );
210  const surfaceVectorField SfGammaCorr(SfGamma - SfGammaSn*Sn);
211 
212  tmp<VolField<Type>> tLaplacian
213  (
214  fvc::div
215  (
216  SfGammaSn*this->tsnGradScheme_().snGrad(vf)
217  + gammaSnGradCorr(SfGammaCorr, vf)
218  )
219  );
220 
221  tLaplacian.ref().rename
222  (
223  "laplacian(" + gamma.name() + ',' + vf.name() + ')'
224  );
225 
226  return tLaplacian;
227 }
228 
229 
230 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
231 
232 } // End namespace fv
233 
234 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
235 
236 } // End namespace Foam
237 
238 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const dimensionSet & dimensions() const
Return dimensions.
Generic GeometricField class.
const Internal::FieldType & primitiveField() const
Return a const-reference to the primitive field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
const word & name() const
Return name.
Definition: IOobject.H:310
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
SurfaceField< Type > *& faceFluxCorrectionPtr()
Return pointer to face-flux non-orthogonal correction field.
Definition: fvMatrix.H:332
Field< Type > & source()
Definition: fvMatrix.H:307
FieldField< Field, Type > & internalCoeffs()
fvBoundary scalar field containing pseudo-matrix coeffs
Definition: fvMatrix.H:319
FieldField< Field, Type > & boundaryCoeffs()
fvBoundary scalar field containing pseudo-matrix coeffs
Definition: fvMatrix.H:326
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
const fvSchemes & schemes() const
Return the fvSchemes.
Definition: fvMesh.C:1751
const surfaceVectorField & Sf() const
Return cell face area vectors.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:88
virtual tmp< Field< Type > > gradientBoundaryCoeffs() const
Return the matrix source coefficients corresponding to the.
Definition: fvPatchField.H:502
virtual bool coupled() const
Return true if this patch field is coupled.
Definition: fvPatchField.H:328
virtual tmp< Field< Type > > gradientInternalCoeffs() const
Return the matrix diagonal coefficients corresponding to the.
Definition: fvPatchField.H:482
bool fluxRequired(const word &name) const
Definition: fvSchemes.C:514
Basic second-order laplacian using face-gradients and Gauss' theorem.
tmp< fvMatrix< Type > > fvmLaplacian(const SurfaceField< GType > &, const VolField< Type > &)
static tmp< fvMatrix< Type > > fvmLaplacianUncorrected(const surfaceScalarField &gammaMagSf, const surfaceScalarField &deltaCoeffs, const VolField< Type > &)
tmp< VolField< Type > > fvcLaplacian(const VolField< Type > &)
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:83
scalarField & upper()
Definition: lduMatrix.C:197
const surfaceScalarField & deltaCoeffs() const
Return reference to cell-centre difference coefficients.
A class for managing temporary objects.
Definition: tmp.H:55
T * ptr() const
Return tmp pointer for reuse.
Definition: tmpI.H:205
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A special matrix type and solver, designed for finite volume solutions of scalar equations.
Calculate the divergence of the given field.
Calculate the gradient of the given field.
label patchi
static tmp< SurfaceField< typename innerProduct< vector, Type >::type > > dotInterpolate(const surfaceVectorField &Sf, const VolField< Type > &tvf)
Interpolate field onto faces.
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
Namespace for OpenFOAM.
static scalar Sn(const scalar a, const scalar x)
Definition: invIncGamma.C:83
uint8_t direction
Definition: direction.H:45
labelList fv(nPoints)