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-2023 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 "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 protected:
58 
59  // Protected data
60 
61  //- The vol field to interpolate
62  const VolField<Type>& psi_;
63 
64  //- Reference to the mesh
65  const polyMesh& mesh_;
66 
67 
68 public:
69 
70  //- Runtime type information
71  TypeName("interpolation");
72 
73 
74  // Declare run-time constructor selection table
75 
77  (
78  autoPtr,
80  dictionary,
81  (
82  const VolField<Type>& psi
83  ),
84  (psi)
85  );
86 
87 
88  // Selectors
89 
90  //- Return a reference to the specified interpolation scheme
92  (
93  const word& interpolationType,
94  const VolField<Type>& psi
95  );
96 
97  //- Return a reference to the selected interpolation scheme
99  (
100  const dictionary& interpolationSchemes,
101  const VolField<Type>& psi
102  );
103 
104 
105  // Constructors
106 
107  //- Construct from components
109  (
110  const VolField<Type>& psi
111  );
112 
113  //- Copy constructor
115 
116  //- Clone
117  virtual autoPtr<interpolation<Type>> clone() const = 0;
118 
119 
120 
121  //- Destructor
122  virtual ~interpolation()
123  {}
124 
125 
126  // Member Functions
127 
128  //- Return the field to be interpolated
129  const VolField<Type>& psi() const
130  {
131  return psi_;
132  }
133 
134  //- Interpolate field to the given point in the given cell
135  virtual Type interpolate
136  (
137  const vector& position,
138  const label celli,
139  const label facei = -1
140  ) const = 0;
141 
142  //- As above, but for a field
144  (
145  const vectorField& position,
146  const labelList& celli,
147  const labelList& facei = NullObjectRef<labelList>()
148  ) const = 0;
149 
150  //- Interpolate field to the given coordinates in the tetrahedron
151  // defined by the given indices. Calls interpolate function
152  // above here except where overridden by derived
153  // interpolation types.
154  virtual Type interpolate
155  (
156  const barycentric& coordinates,
157  const tetIndices& tetIs,
158  const label facei = -1
159  ) const;
160 
161  //- As above, but for a field
163  (
164  const Field<barycentric>& coordinates,
165  const labelList& celli,
166  const labelList& tetFacei,
167  const labelList& tetPti,
168  const labelList& facei = NullObjectRef<labelList>()
169  ) const = 0;
170 };
171 
172 
173 /*---------------------------------------------------------------------------*\
174  Class fieldInterpolation Declaration
175 \*---------------------------------------------------------------------------*/
176 
177 template<class Type, class InterpolationType>
178 class fieldInterpolation
179 :
180  public interpolation<Type>
181 {
182 public:
183 
184  // Constructors
185 
186  //- Inherit base class constructors
188 
189 
190  // Member Functions
191 
192  //- Interpolate field to the given points in the given cells
194  (
195  const vectorField& position,
196  const labelList& celli,
197  const labelList& facei = NullObjectRef<labelList>()
198  ) const;
199 
200  //- Interpolate field to the given coordinates in the given tetrahedra
202  (
203  const Field<barycentric>& coordinates,
204  const labelList& celli,
205  const labelList& tetFacei,
206  const labelList& tetPti,
207  const labelList& facei = NullObjectRef<labelList>()
208  ) const;
209 };
210 
211 
212 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
213 
214 } // End namespace Foam
215 
216 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
217 
218 #define defineInterpolation(Type, InterpolationType) \
219  defineNamedTemplateTypeNameAndDebug(InterpolationType<Type>, 0);
220 
221 #define makeInterpolation(Type, InterpolationType) \
222  addTemplatedToRunTimeSelectionTable \
223  ( \
224  interpolation, \
225  InterpolationType, \
226  Type, \
227  dictionary \
228  );
229 
230 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
231 
232 #ifdef NoRepository
233  #include "interpolation.C"
234 #endif
235 
236 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
237 
238 #endif
239 
240 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
virtual tmp< Field< Type > > interpolate(const vectorField &position, const labelList &celli, const labelList &facei=NullObjectRef< labelList >()) const
Interpolate field to the given points in the given cells.
Abstract base class for interpolation.
Definition: interpolation.H:55
static autoPtr< interpolation< Type > > New(const word &interpolationType, const VolField< Type > &psi)
Return a reference to the specified interpolation scheme.
Definition: interpolation.C:56
interpolation(const VolField< Type > &psi)
Construct from components.
Definition: interpolation.C:35
virtual autoPtr< interpolation< Type > > clone() const =0
Clone.
const VolField< Type > & psi() const
Return the field to be interpolated.
const VolField< Type > & psi_
The vol field to interpolate.
Definition: interpolation.H:61
declareRunTimeSelectionTable(autoPtr, interpolation, dictionary,(const VolField< Type > &psi),(psi))
const polyMesh & mesh_
Reference to the mesh.
Definition: interpolation.H:64
TypeName("interpolation")
Runtime type information.
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.
virtual ~interpolation()
Destructor.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
Storage and named access for the indices of a tet which is part of the decomposition of a cell.
Definition: tetIndices.H:82
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
Namespace for OpenFOAM.
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
Basic run-time type information using word as the type's name. Used to enhance the standard RTTI to c...