54class fieldInterpolator
84 selectedFields_(selectedFields),
87 interpolator_(interpolator),
88 timeNames_(timeNames),
92 template<
class GeoFieldType>
97template<
class GeoFieldType>
98void fieldInterpolator::interpolate()
100 const word& clsName = GeoFieldType::typeName;
104 selectedFields_.empty()
105 ? objects_.sortedNames(clsName)
106 : objects_.sortedNames(clsName, selectedFields_)
109 if (fieldNames.
size())
111 Info<<
" " << clsName <<
's';
114 for (
const word& fieldName : fieldNames)
116 Info<<
' ' << fieldName <<
'(';
118 const scalar deltaT = (ti1_.value() - ti_.value())/(divisions_ + 1);
120 for (
int j=0; j<divisions_; j++)
128 if (j < divisions_-1)
136 interpolator_.valueWeights
154 ReadFields<GeoFieldType>
163 uniformInterpolate<GeoFieldType>
169 objects_[fieldName]->db(),
192int main(
int argc,
char *argv[])
196 "Interpolate fields between time-steps. Eg, for animation."
199 timeSelector::addOptions();
205 "The fields (or field) to be interpolated."
206 " Eg, '(U T p \"Y.*\")' or a single field 'U'"
212 "Specify number of temporal sub-divisions to create (default = 1)."
218 "The type of interpolation (linear or spline)"
221 argList::noFunctionObjects();
229 if (selectedFields.empty())
231 Info<<
"Interpolating all fields" <<
nl <<
endl;
241 Info<<
"Using " << divisions <<
" per time interval" <<
nl <<
endl;
244 const word interpolationType =
247 Info<<
"Using interpolation " << interpolationType <<
nl <<
endl;
256 timeVals[i] = timeDirs[i].value();
257 timeNames[i] = timeDirs[i].name();
261 interpolationWeights::New
271 Info<<
"Interpolating fields for times:" <<
endl;
273 for (label timei = 0; timei < timeDirs.
size() - 1; timei++)
275 runTime.setTime(timeDirs[timei], timei);
280 fieldInterpolator interpolator
Field reading functions for post-processing utilities.
List of IOobjects with searching and retrieving facilities.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const T & name() const noexcept
The name/key (const access)
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A List with indirect addressing. Like IndirectList but does not store addressing.
void size(const label n)
Older name for setAddressableSize.
List< T > getList(const label index) const
Get a List of values from the argument at index.
T getOrDefault(const word &optName, const T &deflt) const
Get a value from the named option if present, or return default.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Mesh data needed to do the Finite Volume discretisation.
An instant of time. Contains the time value and name. Uses Foam::Time when formatting the name.
Abstract base class for interpolating in 1D.
A wordRe is a Foam::word, but can contain a regular expression for matching words or strings.
A List of wordRe with additional matching capabilities.
A class for handling words, derived from Foam::string.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
constexpr char nl
The newline '\n' character (0x0a)
Foam::argList args(argc, argv)
#define forAll(list, i)
Loop across all elements in list.