gaussGrad.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-2019 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 "gaussGrad.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 template<class Type>
33 <
35  <
36  typename Foam::outerProduct<Foam::vector, Type>::type,
37  Foam::fvPatchField,
39  >
40 >
42 (
43  const GeometricField<Type, fvsPatchField, surfaceMesh>& ssf,
44  const word& name
45 )
46 {
47  typedef typename outerProduct<vector, Type>::type GradType;
48 
49  const fvMesh& mesh = ssf.mesh();
50 
51  tmp<GeometricField<GradType, fvPatchField, volMesh>> tgGrad
52  (
54  (
55  name,
56  mesh,
57  dimensioned<GradType>
58  (
59  "0",
60  ssf.dimensions()/dimLength,
61  Zero
62  ),
63  extrapolatedCalculatedFvPatchField<GradType>::typeName
64  )
65  );
66  GeometricField<GradType, fvPatchField, volMesh>& gGrad = tgGrad.ref();
67 
68  const labelUList& owner = mesh.owner();
69  const labelUList& neighbour = mesh.neighbour();
70  const vectorField& Sf = mesh.Sf();
71 
72  Field<GradType>& igGrad = gGrad;
73  const Field<Type>& issf = ssf;
74 
75  forAll(owner, facei)
76  {
77  GradType Sfssf = Sf[facei]*issf[facei];
78 
79  igGrad[owner[facei]] += Sfssf;
80  igGrad[neighbour[facei]] -= Sfssf;
81  }
82 
83  forAll(mesh.boundary(), patchi)
84  {
85  const labelUList& pFaceCells =
86  mesh.boundary()[patchi].faceCells();
87 
88  const vectorField& pSf = mesh.Sf().boundaryField()[patchi];
89 
90  const fvsPatchField<Type>& pssf = ssf.boundaryField()[patchi];
91 
92  forAll(mesh.boundary()[patchi], facei)
93  {
94  igGrad[pFaceCells[facei]] += pSf[facei]*pssf[facei];
95  }
96  }
97 
98  igGrad /= mesh.V();
99 
100  gGrad.correctBoundaryConditions();
101 
102  return tgGrad;
103 }
104 
105 
106 template<class Type>
107 Foam::tmp
108 <
110  <
111  typename Foam::outerProduct<Foam::vector, Type>::type,
112  Foam::fvPatchField,
114  >
115 >
117 (
118  const GeometricField<Type, fvPatchField, volMesh>& vsf,
119  const word& name
120 ) const
121 {
122  typedef typename outerProduct<vector, Type>::type GradType;
123 
124  tmp<GeometricField<GradType, fvPatchField, volMesh>> tgGrad
125  (
126  gradf(tinterpScheme_().interpolate(vsf), name)
127  );
128  GeometricField<GradType, fvPatchField, volMesh>& gGrad = tgGrad.ref();
129 
130  correctBoundaryConditions(vsf, gGrad);
131 
132  return tgGrad;
133 }
134 
135 
136 template<class Type>
138 (
141  <
143  >& gGrad
144 )
145 {
146  typename GeometricField
147  <
149  >::Boundary& gGradbf = gGrad.boundaryFieldRef();
150 
151  forAll(vsf.boundaryField(), patchi)
152  {
153  if (!vsf.boundaryField()[patchi].coupled())
154  {
155  const vectorField n
156  (
157  vsf.mesh().Sf().boundaryField()[patchi]
158  / vsf.mesh().magSf().boundaryField()[patchi]
159  );
160 
161  gGradbf[patchi] += n *
162  (
163  vsf.boundaryField()[patchi].snGrad()
164  - (n & gGradbf[patchi])
165  );
166  }
167  }
168 }
169 
170 
171 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
rDeltaT correctBoundaryConditions()
const Boundary & boundaryField() const
Return const-reference to the boundary field.
static tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > gradf(const GeometricField< Type, fvsPatchField, surfaceMesh > &, const word &name)
Return the gradient of the given field.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
volVectorField vectorField(fieldObject, mesh)
typeOfRank< typename pTraits< arg1 >::cmptType, direction(pTraits< arg1 >::rank)+direction(pTraits< arg2 >::rank) >::type type
Definition: products.H:90
Generic GeometricField class.
UList< label > labelUList
Definition: UList.H:65
virtual tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > calcGrad(const GeometricField< Type, fvPatchField, volMesh > &vsf, const word &name) const
Return the gradient of the given field to the gradScheme::grad.
Mesh data needed to do the Finite Volume discretisation.
Definition: volMesh.H:53
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
static const zero Zero
Definition: zero.H:97
static void correctBoundaryConditions(const GeometricField< Type, fvPatchField, volMesh > &, GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > &)
Correct the boundary values of the gradient using the patchField.
Definition: gaussGrad.C:138
const Mesh & mesh() const
Return mesh.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
label patchi
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
label n
A class for managing temporary objects.
Definition: PtrList.H:53