interpolation.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-2018 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::interpolation
26 
27 Description
28  Abstract base class for interpolation
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #ifndef interpolation_H
33 #define interpolation_H
34 
35 #include "faceList.H"
36 #include "volFieldsFwd.H"
37 #include "pointFields.H"
38 #include "typeInfo.H"
39 #include "autoPtr.H"
40 #include "runTimeSelectionTables.H"
41 #include "tetIndices.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 
48 class polyMesh;
49 
50 /*---------------------------------------------------------------------------*\
51  Class interpolation Declaration
52 \*---------------------------------------------------------------------------*/
53 
54 template<class Type>
55 class interpolation
56 {
57 
58 protected:
59 
60  // Protected data
61 
63 
64  const polyMesh& pMesh_;
69 
70 
71 public:
72 
73  //- Runtime type information
74  virtual const word& type() const = 0;
75 
76 
77  // Declare run-time constructor selection table
78 
80  (
81  autoPtr,
83  dictionary,
84  (
86  ),
87  (psi)
88  );
89 
90 
91  // Selectors
92 
93  //- Return a reference to the specified interpolation scheme
95  (
96  const word& interpolationType,
98  );
99 
100  //- Return a reference to the selected interpolation scheme
102  (
103  const dictionary& interpolationSchemes,
105  );
106 
107 
108  // Constructors
109 
110  //- Construct from components
112  (
114  );
115 
116 
117  //- Destructor
118  virtual ~interpolation()
119  {}
120 
121 
122  // Member Functions
123 
124  //- Return the field to be interpolated
126  {
127  return psi_;
128  }
129 
130  //- Interpolate field to the given point in the given cell
131  virtual Type interpolate
132  (
133  const vector& position,
134  const label celli,
135  const label facei = -1
136  ) const = 0;
137 
138  //- Interpolate field to the given coordinates in the tetrahedron
139  // defined by the given indices. Calls interpolate function
140  // above here execpt where overridden by derived
141  // interpolation types.
142  virtual Type interpolate
143  (
144  const barycentric& coordinates,
145  const tetIndices& tetIs,
146  const label facei = -1
147  ) const
148  {
149  return
151  (
152  tetIs.tet(pMesh_).barycentricToPoint(coordinates),
153  tetIs.cell(),
154  facei
155  );
156  }
157 };
158 
159 
160 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
161 
162 } // End namespace Foam
163 
164 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
166 #define makeInterpolationType(SS, Type) \
167  \
168 defineNamedTemplateTypeNameAndDebug(SS<Type>, 0); \
169  \
170 interpolation<Type>::adddictionaryConstructorToTable<SS<Type>> \
171  add##SS##Type##ConstructorToTable_;
172 
174 #define makeInterpolation(SS) \
175  \
176 makeInterpolationType(SS, scalar) \
177 makeInterpolationType(SS, vector) \
178 makeInterpolationType(SS, sphericalTensor) \
179 makeInterpolationType(SS, symmTensor) \
180 makeInterpolationType(SS, tensor)
181 
182 
183 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
184 
185 #ifdef NoRepository
186  #include "interpolation.C"
187 #endif
188 
189 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
190 
191 #endif
192 
193 // ************************************************************************* //
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const GeometricField< Type, fvPatchField, volMesh > & psi_
Definition: interpolation.H:61
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const polyMesh & pMesh_
Definition: interpolation.H:63
static autoPtr< interpolation< Type > > New(const word &interpolationType, const GeometricField< Type, fvPatchField, volMesh > &psi)
Return a reference to the specified interpolation scheme.
virtual ~interpolation()
Destructor.
virtual const word & type() const =0
Runtime type information.
Generic GeometricField class.
tetPointRef tet(const polyMesh &mesh) const
Return the geometry corresponding to this tet.
Definition: tetIndicesI.H:100
Point barycentricToPoint(const barycentric &bary) const
Calculate the point from the given barycentric coordinates.
Definition: tetrahedronI.H:254
virtual Type interpolate(const vector &position, const label celli, const label facei=-1) const =0
Interpolate field to the given point in the given cell.
const vectorField & pMeshFaceCentres_
Definition: interpolation.H:66
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Return the field to be interpolated.
const vectorField & pMeshPoints_
Definition: interpolation.H:64
A class for handling words, derived from string.
Definition: word.H:59
label cell() const
Return the cell.
Definition: tetIndicesI.H:28
const faceList & pMeshFaces_
Definition: interpolation.H:65
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
Definition: tetIndices.H:81
declareRunTimeSelectionTable(autoPtr, interpolation, dictionary,(const GeometricField< Type, fvPatchField, volMesh > &psi),(psi))
PtrList< coordinateSystem > coordinates(solidRegions.size())
Abstract base class for interpolation.
interpolation(const GeometricField< Type, fvPatchField, volMesh > &psi)
Construct from components.
Definition: interpolation.C:35
const vectorField & pMeshFaceAreas_
Definition: interpolation.H:67
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
Macros to ease declaration of run-time selection tables.
Namespace for OpenFOAM.