uniformInterpolatedDisplacementPointPatchVectorField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2012-2013 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 \*---------------------------------------------------------------------------*/
25 
27 #include "pointFields.H"
29 #include "Time.H"
30 #include "polyMesh.H"
31 #include "interpolationWeights.H"
32 #include "uniformInterpolate.H"
33 #include "ReadFields.H"
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
44 (
45  const pointPatch& p,
47 )
48 :
50 {}
51 
52 
55 (
56  const pointPatch& p,
58  const dictionary& dict
59 )
60 :
62  fieldName_(dict.lookup("fieldName")),
63  interpolationScheme_(dict.lookup("interpolationScheme"))
64 {
65  const pointMesh& pMesh = this->dimensionedInternalField().mesh();
66 
67  // Read time values
68  instantList allTimes = Time::findTimes(pMesh().time().path());
69 
70  // Only keep those that contain the field
71  DynamicList<word> names(allTimes.size());
72  DynamicList<scalar> values(allTimes.size());
73 
74  forAll(allTimes, i)
75  {
76  IOobject io
77  (
78  fieldName_,
79  allTimes[i].name(),
80  pMesh(),
83  false
84  );
85  if (io.headerOk())
86  {
87  names.append(allTimes[i].name());
88  values.append(allTimes[i].value());
89  }
90  }
91  timeNames_.transfer(names);
92  timeVals_.transfer(values);
93 
94  Info<< type() << " : found " << fieldName_ << " for times "
95  << timeNames_ << endl;
96 
97  if (timeNames_.size() < 1)
98  {
100  (
101  "uniformInterpolatedDisplacementPointPatchVectorField::\n"
102  "uniformInterpolatedDisplacementPointPatchVectorField\n"
103  "(\n"
104  " const pointPatch&,\n"
105  " const DimensionedField<vector, pointMesh>&,\n"
106  " const dictionary&\n"
107  ")\n"
108  ) << "Did not find any times with " << fieldName_
109  << exit(FatalError);
110  }
111 
112 
113  if (!dict.found("value"))
114  {
115  updateCoeffs();
116  }
117 }
118 
119 
122 (
124  const pointPatch& p,
126  const pointPatchFieldMapper& mapper
127 )
128 :
129  fixedValuePointPatchField<vector>(ptf, p, iF, mapper),
130  fieldName_(ptf.fieldName_),
131  interpolationScheme_(ptf.interpolationScheme_),
132  timeNames_(ptf.timeNames_),
133  timeVals_(ptf.timeVals_),
134  interpolatorPtr_(ptf.interpolatorPtr_)
135 {}
136 
137 
140 (
143 )
144 :
146  fieldName_(ptf.fieldName_),
147  interpolationScheme_(ptf.interpolationScheme_),
148  timeNames_(ptf.timeNames_),
149  timeVals_(ptf.timeVals_),
150  interpolatorPtr_(ptf.interpolatorPtr_)
151 {}
152 
153 
154 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
155 
157 {
158  if (this->updated())
159  {
160  return;
161  }
162 
163  if (!interpolatorPtr_.valid())
164  {
165  interpolatorPtr_ = interpolationWeights::New
166  (
167  interpolationScheme_,
168  timeVals_
169  );
170  }
171 
172  const pointMesh& pMesh = this->dimensionedInternalField().mesh();
173  const Time& t = pMesh().time();
174 
175  // Update indices of times and weights
176  bool timesChanged = interpolatorPtr_->valueWeights
177  (
178  t.timeOutputValue(),
179  currentIndices_,
180  currentWeights_
181  );
182 
183  const wordList currentTimeNames
184  (
185  UIndirectList<word>(timeNames_, currentIndices_)
186  );
187 
188 
189  // Load if necessary fields for this interpolation
190  if (timesChanged)
191  {
192  objectRegistry& fieldsCache = const_cast<objectRegistry&>
193  (
194  pMesh.thisDb().subRegistry("fieldsCache", true)
195  );
196  // Save old times so we now which ones have been loaded and need
197  // 'correctBoundaryConditions'. Bit messy.
198  HashSet<word> oldTimes(fieldsCache.toc());
199 
200  ReadFields<pointVectorField>
201  (
202  fieldName_,
203  pMesh,
204  currentTimeNames
205  );
206 
207  forAllConstIter(objectRegistry, fieldsCache, fieldsCacheIter)
208  {
209  if (!oldTimes.found(fieldsCacheIter.key()))
210  {
211  // Newly loaded fields. Make sure the internal
212  // values are consistent with the boundary conditions.
213  // This is quite often not the case since these
214  // fields typically are constructed 'by hand'
215 
216  const objectRegistry& timeCache = dynamic_cast
217  <
218  const objectRegistry&
219  >(*fieldsCacheIter());
220 
221 
222  pointVectorField& d = const_cast<pointVectorField&>
223  (
224  timeCache.lookupObject<pointVectorField>
225  (
226  fieldName_
227  )
228  );
230  }
231  }
232  }
233 
234 
235  // Interpolate the whole field
236  pointVectorField result
237  (
238  uniformInterpolate<pointVectorField>
239  (
240  IOobject
241  (
242  word("uniformInterpolate(")
243  + this->dimensionedInternalField().name()
244  + ')',
245  pMesh.time().timeName(),
246  pMesh.thisDb(),
249  ),
250  fieldName_,
251  currentTimeNames,
252  currentWeights_
253  )
254  );
255 
256 
257  // Extract back from the internal field
258  this->operator==
259  (
261  );
262 
264 }
265 
266 
268 const
269 {
271  os.writeKeyword("fieldName")
272  << fieldName_ << token::END_STATEMENT << nl;
273  os.writeKeyword("interpolationScheme")
274  << interpolationScheme_ << token::END_STATEMENT << nl;
275  writeEntry("value", os);
276 }
277 
278 
279 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
280 
282 (
285 );
286 
287 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
288 
289 } // End namespace Foam
290 
291 // ************************************************************************* //
const objectRegistry & thisDb() const
Return database. For now is its polyMesh.
Definition: pointMesh.H:118
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:56
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const objectRegistry & subRegistry(const word &name, const bool forceCreate=false) const
Lookup and return a const sub-objectRegistry. Optionally create.
A class for handling words, derived from string.
Definition: word.H:59
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
DimensionedInternalField & dimensionedInternalField()
Return dimensioned internal field.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
messageStream Info
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
scalar timeOutputValue() const
Return current time value.
Definition: TimeState.C:67
Namespace for OpenFOAM.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
makePointPatchTypeField(pointPatchVectorField, solidBodyMotionDisplacementPointPatchVectorField)
dictionary dict
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
A List with indirect addressing.
Definition: fvMatrix.H:106
volScalarField & p
Definition: createFields.H:51
const DimensionedField< vector, pointMesh > & dimensionedInternalField() const
Return dimensioned internal field reference.
#define forAll(list, i)
Definition: UList.H:421
void writeEntry(const word &keyword, Ostream &os) const
Write the field as a dictionary entry.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Macros for easy insertion into run-time selection tables.
List< Key > toc() const
Return the table of contents.
Definition: HashTable.C:201
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
Foam::pointPatchFieldMapper.
Helper routine to read fields.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
rDeltaT dimensionedInternalField()
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
tmp< Field< vector > > patchInternalField() const
Return field created from appropriate internal field values.
error FatalError
static autoPtr< interpolationWeights > New(const word &type, const scalarField &samples)
Return a reference to the selected interpolationWeights.
void correctBoundaryConditions()
Correct boundary field.
Registry of regIOobjects.
uniformInterpolatedDisplacementPointPatchVectorField(const pointPatch &, const DimensionedField< vector, pointMesh > &)
Construct from patch and internal field.
A HashTable with keys but without contents.
Definition: HashSet.H:59
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
const Time & time() const
Return time.
Definition: IOobject.C:251
bool headerOk()
Read and check header info.
Definition: IOobject.C:424
bool updated() const
Return true if the boundary condition has already been updated.
static instantList findTimes(const fileName &, const word &constantName="constant")
Search a given directory for valid time directories.
Definition: findTimes.C:38
virtual void write(Ostream &) const
Write.