ReadFields.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) 2011-2012 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 
26 #include "ReadFields.H"
27 #include "HashSet.H"
28 #include "Pstream.H"
29 #include "IOobjectList.H"
30 
31 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
32 
33 // Read all fields of type. Returns names of fields read. Guarantees all
34 // processors to read fields in same order.
35 template<class GeoField, class Mesh>
37 (
38  const Mesh& mesh,
39  const IOobjectList& objects,
40  PtrList<GeoField>& fields,
41  const bool syncPar
42 )
43 {
44  // Search list of objects for wanted type
45  IOobjectList fieldObjects(objects.lookupClass(GeoField::typeName));
46 
47  wordList masterNames(fieldObjects.names());
48 
49  if (syncPar && Pstream::parRun())
50  {
51  // Check that I have the same fields as the master
52  const wordList localNames(masterNames);
53  Pstream::scatter(masterNames);
54 
55  HashSet<word> localNamesSet(localNames);
56 
57  forAll(masterNames, i)
58  {
59  const word& masterFld = masterNames[i];
60 
61  HashSet<word>::iterator iter = localNamesSet.find(masterFld);
62 
63  if (iter == localNamesSet.end())
64  {
66  (
67  "ReadFields<class GeoField, class Mesh>"
68  "(const Mesh&, const IOobjectList&, PtrList<GeoField>&"
69  ", const bool)"
70  ) << "Fields not synchronised across processors." << endl
71  << "Master has fields " << masterNames
72  << " processor " << Pstream::myProcNo()
73  << " has fields " << localNames << exit(FatalError);
74  }
75  else
76  {
77  localNamesSet.erase(iter);
78  }
79  }
80 
81  forAllConstIter(HashSet<word>, localNamesSet, iter)
82  {
84  (
85  "ReadFields<class GeoField, class Mesh>"
86  "(const Mesh&, const IOobjectList&, PtrList<GeoField>&"
87  ", const bool)"
88  ) << "Fields not synchronised across processors." << endl
89  << "Master has fields " << masterNames
90  << " processor " << Pstream::myProcNo()
91  << " has fields " << localNames << exit(FatalError);
92  }
93  }
94 
95 
96  fields.setSize(masterNames.size());
97 
98  // Make sure to read in masterNames order.
99 
100  forAll(masterNames, i)
101  {
102  Info<< "Reading " << GeoField::typeName << ' ' << masterNames[i]
103  << endl;
104 
105  const IOobject& io = *fieldObjects[masterNames[i]];
106 
107  fields.set
108  (
109  i,
110  new GeoField
111  (
112  IOobject
113  (
114  io.name(),
115  io.instance(),
116  io.local(),
117  io.db(),
118  IOobject::MUST_READ,
119  IOobject::AUTO_WRITE,
120  io.registerObject()
121  ),
122  mesh
123  )
124  );
125  }
126  return masterNames;
127 }
128 
129 
130 template<class GeoField>
131 void Foam::ReadFields
132 (
133  const word& fieldName,
134  const typename GeoField::Mesh& mesh,
135  const wordList& timeNames,
136  objectRegistry& fieldsCache
137 )
138 {
139  // Collect all times that are no longer used
140  {
141  HashSet<word> usedTimes(timeNames);
142 
143  DynamicList<word> unusedTimes(fieldsCache.size());
144 
145  forAllIter(objectRegistry, fieldsCache, timeIter)
146  {
147  const word& tm = timeIter.key();
148  if (!usedTimes.found(tm))
149  {
150  unusedTimes.append(tm);
151  }
152  }
153 
154  //Info<< "Unloading times " << unusedTimes << endl;
155 
156  forAll(unusedTimes, i)
157  {
158  objectRegistry& timeCache = const_cast<objectRegistry&>
159  (
160  fieldsCache.lookupObject<objectRegistry>(unusedTimes[i])
161  );
162  fieldsCache.checkOut(timeCache);
163  }
164  }
165 
166 
167  // Load any new fields
168  forAll(timeNames, i)
169  {
170  const word& tm = timeNames[i];
171 
172  // Create if not found
173  if (!fieldsCache.found(tm))
174  {
175  //Info<< "Creating registry for time " << tm << endl;
176 
177  // Create objectRegistry if not found
178  objectRegistry* timeCachePtr = new objectRegistry
179  (
180  IOobject
181  (
182  tm,
183  tm,
184  fieldsCache,
185  IOobject::NO_READ,
186  IOobject::NO_WRITE
187  )
188  );
189  timeCachePtr->store();
190  }
191 
192  // Obtain cache for current time
193  const objectRegistry& timeCache =
194  fieldsCache.lookupObject<objectRegistry>
195  (
196  tm
197  );
198 
199  // Store field if not found
200  if (!timeCache.found(fieldName))
201  {
202  //Info<< "Loading field " << fieldName
203  // << " for time " << tm << endl;
204 
205  GeoField loadedFld
206  (
207  IOobject
208  (
209  fieldName,
210  tm,
211  mesh.thisDb(),
212  IOobject::MUST_READ,
213  IOobject::NO_WRITE,
214  false
215  ),
216  mesh
217  );
218 
219  // Transfer to timeCache (new objectRegistry and store flag)
220  GeoField* fldPtr = new GeoField
221  (
222  IOobject
223  (
224  fieldName,
225  tm,
226  timeCache,
227  IOobject::NO_READ,
228  IOobject::NO_WRITE
229  ),
230  loadedFld
231  );
232  fldPtr->store();
233  }
234  }
235 }
236 
237 
238 template<class GeoField>
239 void Foam::ReadFields
240 (
241  const word& fieldName,
242  const typename GeoField::Mesh& mesh,
243  const wordList& timeNames,
244  const word& registryName
245 )
246 {
247  ReadFields<GeoField>
248  (
249  fieldName,
250  mesh,
251  timeNames,
252  const_cast<objectRegistry&>
253  (
254  mesh.thisDb().subRegistry(registryName, true)
255  )
256  );
257 }
258 
259 
260 // ************************************************************************* //
#define forAllIter(Container, container, iter)
Definition: UList.H:440
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
messageStream Info
dynamicFvMesh & mesh
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
wordList ReadFields(const Mesh &mesh, const IOobjectList &objects, PtrList< GeoField > &fields, const bool syncPar=true)
Helper routine to read fields.
#define forAll(list, i)
Definition: UList.H:421
Helper routine to read fields.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
List< word > wordList
A List of words.
Definition: fileName.H:54
error FatalError