temporalInterpolate.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-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 Description
25  Interpolate fields between time-steps e.g. for animation.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "argList.H"
30 #include "timeSelector.H"
31 
32 #include "fvMesh.H"
33 #include "Time.H"
34 #include "volMesh.H"
35 #include "surfaceMesh.H"
36 #include "volFields.H"
37 #include "surfaceFields.H"
38 #include "pointFields.H"
39 #include "ReadFields.H"
40 #include "interpolationWeights.H"
41 #include "uniformInterpolate.H"
42 
43 using namespace Foam;
44 
45 class fieldInterpolator
46 {
47  Time& runTime_;
48  const fvMesh& mesh_;
49  const IOobjectList& objects_;
50  const HashSet<word>& selectedFields_;
51  instant ti_;
52  instant ti1_;
53  const interpolationWeights& interpolator_;
54  const wordList& timeNames_;
55  int divisions_;
56 
57 public:
58 
59  fieldInterpolator
60  (
61  Time& runTime,
62  const fvMesh& mesh,
63  const IOobjectList& objects,
64  const HashSet<word>& selectedFields,
65  const instant& ti,
66  const instant& ti1,
67  const interpolationWeights& interpolator,
68  const wordList& timeNames,
69  int divisions
70  )
71  :
72  runTime_(runTime),
73  mesh_(mesh),
74  objects_(objects),
75  selectedFields_(selectedFields),
76  ti_(ti),
77  ti1_(ti1),
78  interpolator_(interpolator),
79  timeNames_(timeNames),
80  divisions_(divisions)
81  {}
82 
83  template<class GeoFieldType>
84  void interpolate();
85 };
86 
87 
88 template<class GeoFieldType>
90 {
91  const word& fieldClassName = GeoFieldType::typeName;
92 
93  IOobjectList fields = objects_.lookupClass(fieldClassName);
94 
95  if (fields.size())
96  {
97  Info<< " " << fieldClassName << "s:";
98 
99  forAllConstIter(IOobjectList, fields, fieldIter)
100  {
101  if
102  (
103  selectedFields_.empty()
104  || selectedFields_.found(fieldIter()->name())
105  )
106  {
107  Info<< " " << fieldIter()->name() << '(';
108 
109  scalar deltaT = (ti1_.value() - ti_.value())/(divisions_ + 1);
110 
111  for (int j=0; j<divisions_; j++)
112  {
113  instant timej = instant(ti_.value() + (j + 1)*deltaT);
114 
115  runTime_.setTime(instant(timej.name()), 0);
116 
117  Info<< timej.name();
118 
119  if (j < divisions_-1)
120  {
121  Info<< " ";
122  }
123 
124  // Calculate times to read and weights
125  labelList indices;
126  scalarField weights;
127  interpolator_.valueWeights
128  (
129  runTime_.value(),
130  indices,
131  weights
132  );
133 
134  const wordList selectedTimeNames
135  (
136  UIndirectList<word>(timeNames_, indices)()
137  );
138 
139  // Info<< "For time " << runTime_.value()
140  // << " need times " << selectedTimeNames
141  // << " need weights " << weights << endl;
142 
143 
144  // Read on the objectRegistry all the required fields
145  ReadFields<GeoFieldType>
146  (
147  fieldIter()->name(),
148  mesh_,
149  selectedTimeNames
150  );
151 
152  GeoFieldType fieldj
153  (
154  uniformInterpolate<GeoFieldType>
155  (
156  IOobject
157  (
158  fieldIter()->name(),
159  runTime_.timeName(),
160  fieldIter()->db(),
163  false
164  ),
165  fieldIter()->name(),
166  selectedTimeNames,
167  weights
168  )
169  );
170 
171  fieldj.write();
172  }
173 
174  Info<< ')';
175  }
176  }
177 
178  Info<< endl;
179  }
180 }
181 
182 
183 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
184 
185 int main(int argc, char *argv[])
186 {
188  #include "addRegionOption.H"
190  (
191  "fields",
192  "list",
193  "specify a list of fields to be interpolated. Eg, '(U T p)' - "
194  "regular expressions not currently supported"
195  );
197  (
198  "divisions",
199  "integer",
200  "specify number of temporal sub-divisions to create (default = 1)."
201  );
203  (
204  "interpolationType",
205  "word",
206  "specify type of interpolation (linear or spline)"
207  );
208 
209  #include "setRootCase.H"
210  #include "createTime.H"
212 
213  HashSet<word> selectedFields;
214  if (args.optionFound("fields"))
215  {
216  args.optionLookup("fields")() >> selectedFields;
217  }
218  if (selectedFields.size())
219  {
220  Info<< "Interpolating fields " << selectedFields << nl << endl;
221  }
222  else
223  {
224  Info<< "Interpolating all fields" << nl << endl;
225  }
226 
227 
228  int divisions = 1;
229  if (args.optionFound("divisions"))
230  {
231  args.optionLookup("divisions")() >> divisions;
232  }
233  Info<< "Using " << divisions << " per time interval" << nl << endl;
234 
235 
236  const word interpolationType = args.optionLookupOrDefault<word>
237  (
238  "interpolationType",
239  "linear"
240  );
241  Info<< "Using interpolation " << interpolationType << nl << endl;
242 
243 
245 
246  scalarField timeVals(timeDirs.size());
247  wordList timeNames(timeDirs.size());
248  forAll(timeDirs, i)
249  {
250  timeVals[i] = timeDirs[i].value();
251  timeNames[i] = timeDirs[i].name();
252  }
253  autoPtr<interpolationWeights> interpolatorPtr
254  (
256  (
257  interpolationType,
258  timeVals
259  )
260  );
261 
262 
263  #include "createNamedMesh.H"
264 
265  Info<< "Interpolating fields for times:" << endl;
266 
267  for (label timei = 0; timei < timeDirs.size() - 1; timei++)
268  {
269  runTime.setTime(timeDirs[timei], timei);
270 
271  // Read objects in time directory
272  IOobjectList objects(mesh, runTime.timeName());
273 
274  fieldInterpolator interpolator
275  (
276  runTime,
277  mesh,
278  objects,
279  selectedFields,
280  timeDirs[timei],
281  timeDirs[timei+1],
282  interpolatorPtr(),
283  timeNames,
284  divisions
285  );
286 
287  // Interpolate vol fields
288  interpolator.interpolate<volScalarField>();
289  interpolator.interpolate<volVectorField>();
290  interpolator.interpolate<volSphericalTensorField>();
291  interpolator.interpolate<volSymmTensorField>();
292  interpolator.interpolate<volTensorField>();
293 
294  // Interpolate surface fields
295  interpolator.interpolate<surfaceScalarField>();
296  interpolator.interpolate<surfaceVectorField>();
297  interpolator.interpolate<surfaceSphericalTensorField>();
298  interpolator.interpolate<surfaceSymmTensorField>();
299  interpolator.interpolate<surfaceTensorField>();
300  }
301 
302  Info<< "End\n" << endl;
303 
304  return 0;
305 }
306 
307 
308 // ************************************************************************* //
Foam::surfaceFields.
static autoPtr< interpolationWeights > New(const word &type, const scalarField &samples)
Return a reference to the selected interpolationWeights.
A HashTable with keys but without contents.
Definition: HashSet.H:59
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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 word & name() const
Return name.
Definition: IOobject.H:297
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
void off()
Switch the function objects off.
const word & name() const
Name (const access)
Definition: instant.H:126
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
Field reading functions for post-processing utilities.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
Abstract base class for interpolating in 1D.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
T optionLookupOrDefault(const word &opt, const T &deflt) const
Read a value from the named option if present.
Definition: argListI.H:237
dynamicFvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
static void addOption(const word &opt, const string &param="", const string &usage="")
Add to an option to validOptions with usage information.
Definition: argList.C:120
virtual void setTime(const Time &)
Reset the time and time-index to those of the given time.
Definition: Time.C:873
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
static const char nl
Definition: Ostream.H:265
static instantList select0(Time &runTime, const argList &args)
Return the set of times selected based on the argList options.
Definition: timeSelector.C:257
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
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.
An instant of time. Contains the time value and name.
Definition: instant.H:66
IOobjectList lookupClass(const word &className) const
Return the list for all IOobjects of a given class.
Definition: IOobjectList.C:192
const functionObjectList & functionObjects() const
Return the list of function objects.
Definition: Time.H:410
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A List with indirect addressing.
Definition: fvMatrix.H:106
messageStream Info
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Foam::argList args(argc, argv)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
static void addOptions(const bool constant=true, const bool withZero=false)
Add the options handled by timeSelector to argList::validOptions.
Definition: timeSelector.C:114
Namespace for OpenFOAM.
IStringStream optionLookup(const word &opt) const
Return an IStringStream from the named option.
Definition: argListI.H:114