cubic.H
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 Class
25  Foam::cubic
26 
27 Description
28  Cubic interpolation scheme class derived from linear and returns
29  linear weighting factors but also applies an explicit correction.
30 
31 SourceFiles
32  cubic.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef cubic_H
37 #define cubic_H
38 
39 #include "linear.H"
40 #include "gaussGrad.H"
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 namespace Foam
45 {
46 
47 /*---------------------------------------------------------------------------*\
48  Class cubic Declaration
49 \*---------------------------------------------------------------------------*/
50 
51 template<class Type>
52 class cubic
53 :
54  public linear<Type>
55 {
56 
57 public:
58 
59  //- Runtime type information
60  TypeName("cubic");
61 
62 
63  // Constructors
64 
65  //- Construct from mesh
66  cubic(const fvMesh& mesh)
67  :
68  linear<Type>(mesh)
69  {}
70 
71  //- Construct from mesh and Istream
72  cubic
73  (
74  const fvMesh& mesh,
75  Istream&
76  )
77  :
78  linear<Type>(mesh)
79  {}
80 
81  //- Construct from mesh, faceFlux and Istream
82  cubic
83  (
84  const fvMesh& mesh,
85  const surfaceScalarField&,
86  Istream&
87  )
88  :
89  linear<Type>(mesh)
90  {}
91 
92  //- Disallow default bitwise copy construction
93  cubic(const cubic&) = delete;
94 
95 
96  // Member Functions
97 
98  //- Return true if this scheme uses an explicit correction
99  virtual bool corrected() const
100  {
101  return true;
102  }
103 
104  //- Return the explicit correction to the face-interpolate
106  correction
107  (
108  const VolField<Type>& vf
109  ) const
110  {
111  const fvMesh& mesh = this->mesh();
112 
113  // calculate the appropriate interpolation factors
115 
116  const surfaceScalarField kSc
117  (
118  lambda*(scalar(1) - lambda*(scalar(3) - scalar(2)*lambda))
119  );
120 
121  const surfaceScalarField kVecP(sqr(scalar(1) - lambda)*lambda);
122  const surfaceScalarField kVecN(sqr(lambda)*(lambda - scalar(1)));
123 
124  tmp<SurfaceField<Type>> tsfCorr
125  (
127  (
128  "cubic::correction(" + vf.name() +')',
130  )
131  );
132  SurfaceField<Type>& sfCorr =
133  tsfCorr.ref();
134 
135  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
136  {
137  sfCorr.replace
138  (
139  cmpt,
140  sfCorr.component(cmpt)
141  + (
143  <
144  typename outerProduct
145  <
146  vector,
147  typename pTraits<Type>::cmptType
148  >::type
149  >::interpolate
150  (
152  <typename pTraits<Type>::cmptType>(mesh)
153  .grad(vf.component(cmpt)),
154  kVecP,
155  kVecN
156  ) & mesh.Sf()
157  )/mesh.magSf()/mesh.surfaceInterpolation::deltaCoeffs()
158  );
159  }
160 
161  typename SurfaceField<Type>::
162  Boundary& sfCorrbf = sfCorr.boundaryFieldRef();
163 
164  forAll(sfCorrbf, pi)
165  {
166  if (!sfCorrbf[pi].coupled())
167  {
168  sfCorrbf[pi] = Zero;
169  }
170  }
171 
172  return tsfCorr;
173  }
174 
175 
176  // Member Operators
177 
178  //- Disallow default bitwise assignment
179  void operator=(const cubic&) = delete;
180 };
181 
182 
183 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
184 
185 } // End namespace Foam
186 
187 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
188 
189 #endif
190 
191 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Generic GeometricField class.
void replace(const direction, const GeometricField< cmptType, PatchField, GeoMesh > &)
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const word & name() const
Return name.
Definition: IOobject.H:310
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
Cubic interpolation scheme class derived from linear and returns linear weighting factors but also ap...
Definition: cubic.H:54
virtual bool corrected() const
Return true if this scheme uses an explicit correction.
Definition: cubic.H:98
TypeName("cubic")
Runtime type information.
void operator=(const cubic &)=delete
Disallow default bitwise assignment.
cubic(const fvMesh &mesh)
Construct from mesh.
Definition: cubic.H:65
virtual tmp< SurfaceField< Type > > correction(const VolField< Type > &vf) const
Return the explicit correction to the face-interpolate.
Definition: cubic.H:106
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const surfaceVectorField & Sf() const
Return cell face area vectors.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
Basic second-order gradient scheme using face-interpolation and Gauss' theorem.
Definition: gaussGrad.H:60
Centred interpolation interpolation scheme class.
Definition: linear.H:53
Traits class for primitives.
Definition: pTraits.H:53
Abstract base class for surface interpolation schemes.
const fvMesh & mesh() const
Return mesh reference.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &, const tmp< surfaceScalarField > &, const tmp< surfaceScalarField > &)
Return the face-interpolate of the given cell field.
const surfaceScalarField & weights() const
Return reference to linear difference weighting factors.
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
dimensionedScalar lambda(viscosity->lookup("lambda"))
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
dimensionedSymmTensor sqr(const dimensionedVector &dv)
uint8_t direction
Definition: direction.H:45
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488