FreeFOAM The Cross-Platform CFD Toolkit
mapFields.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) 1991-2010 OpenCFD Ltd.
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 Application
25  mapFields
26 
27 Description
28  Maps volume fields from one mesh to another.
29 
30  Maps volume fields from one mesh to another, reading and
31  interpolating all fields present in the time directory of both cases.
32  Parallel and non-parallel cases are handled without the need to reconstruct
33  them first.
34 
35 Usage
36 
37  - mapFields [OPTIONS] <source dir>
38 
39  @param <source dir> \n
40  @todo Detailed description of argument.
41 
42  @param -consistent \n
43  Meshes are consistent.
44 
45  @param -parallelSource \n
46  The source case is decomposed.
47 
48  @param -sourceTime <scalar>\n
49  Time of the source.
50 
51  @param -parallelTarget \n
52  The target case is decomposed.
53 
54  @param -case <dir>\n
55  Case directory.
56 
57  @param -help \n
58  Display help message.
59 
60  @param -doc \n
61  Display Doxygen API documentation page for this application.
62 
63  @param -srcDoc \n
64  Display Doxygen source documentation page for this application.
65 
66 \*---------------------------------------------------------------------------*/
67 
68 #include <finiteVolume/fvCFD.H>
69 #include <sampling/meshToMesh.H>
70 #include "MapVolFields.H"
71 #include "MapConsistentVolFields.H"
72 #include "UnMapped.H"
74 #include "mapLagrangian.H"
75 
76 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
77 
78 void mapConsistentMesh
79 (
80  const fvMesh& meshSource,
81  const fvMesh& meshTarget
82 )
83 {
84  // Create the interpolation scheme
85  meshToMesh meshToMeshInterp(meshSource, meshTarget);
86 
87  Info<< nl
88  << "Consistently creating and mapping fields for time "
89  << meshSource.time().timeName() << nl << endl;
90 
91  {
92  // Search for list of objects for this time
93  IOobjectList objects(meshSource, meshSource.time().timeName());
94 
95  // Map volFields
96  // ~~~~~~~~~~~~~
97  MapConsistentVolFields<scalar>(objects, meshToMeshInterp);
98  MapConsistentVolFields<vector>(objects, meshToMeshInterp);
99  MapConsistentVolFields<sphericalTensor>(objects, meshToMeshInterp);
100  MapConsistentVolFields<symmTensor>(objects, meshToMeshInterp);
101  MapConsistentVolFields<tensor>(objects, meshToMeshInterp);
102  }
103 
104  {
105  // Search for list of target objects for this time
106  IOobjectList objects(meshTarget, meshTarget.time().timeName());
107 
108  // Mark surfaceFields as unmapped
109  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
110  UnMapped<surfaceScalarField>(objects);
111  UnMapped<surfaceVectorField>(objects);
112  UnMapped<surfaceSphericalTensorField>(objects);
113  UnMapped<surfaceSymmTensorField>(objects);
114  UnMapped<surfaceTensorField>(objects);
115 
116  // Mark pointFields as unmapped
117  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
118  UnMapped<pointScalarField>(objects);
119  UnMapped<pointVectorField>(objects);
120  UnMapped<pointSphericalTensorField>(objects);
121  UnMapped<pointSymmTensorField>(objects);
122  UnMapped<pointTensorField>(objects);
123  }
124 
125  mapLagrangian(meshToMeshInterp);
126 }
127 
128 
129 void mapSubMesh
130 (
131  const fvMesh& meshSource,
132  const fvMesh& meshTarget,
133  const HashTable<word>& patchMap,
134  const wordList& cuttingPatches
135 )
136 {
137  // Create the interpolation scheme
138  meshToMesh meshToMeshInterp
139  (
140  meshSource,
141  meshTarget,
142  patchMap,
143  cuttingPatches
144  );
145 
146  Info<< nl
147  << "Mapping fields for time " << meshSource.time().timeName()
148  << nl << endl;
149 
150  {
151  // Search for list of source objects for this time
152  IOobjectList objects(meshSource, meshSource.time().timeName());
153 
154  // Map volFields
155  // ~~~~~~~~~~~~~
156  MapVolFields<scalar>(objects, meshToMeshInterp);
157  MapVolFields<vector>(objects, meshToMeshInterp);
158  MapVolFields<sphericalTensor>(objects, meshToMeshInterp);
159  MapVolFields<symmTensor>(objects, meshToMeshInterp);
160  MapVolFields<tensor>(objects, meshToMeshInterp);
161  }
162 
163  {
164  // Search for list of target objects for this time
165  IOobjectList objects(meshTarget, meshTarget.time().timeName());
166 
167  // Mark surfaceFields as unmapped
168  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
169  UnMapped<surfaceScalarField>(objects);
170  UnMapped<surfaceVectorField>(objects);
171  UnMapped<surfaceSphericalTensorField>(objects);
172  UnMapped<surfaceSymmTensorField>(objects);
173  UnMapped<surfaceTensorField>(objects);
174 
175  // Mark pointFields as unmapped
176  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
177  UnMapped<pointScalarField>(objects);
178  UnMapped<pointVectorField>(objects);
179  UnMapped<pointSphericalTensorField>(objects);
180  UnMapped<pointSymmTensorField>(objects);
181  UnMapped<pointTensorField>(objects);
182  }
183 
184  mapLagrangian(meshToMeshInterp);
185 }
186 
187 
188 void mapConsistentSubMesh
189 (
190  const fvMesh& meshSource,
191  const fvMesh& meshTarget
192 )
193 {
194  HashTable<word> patchMap;
195  HashTable<label> cuttingPatchTable;
196 
197  forAll(meshTarget.boundary(), patchi)
198  {
199  if (!isA<processorFvPatch>(meshTarget.boundary()[patchi]))
200  {
201  patchMap.insert
202  (
203  meshTarget.boundary()[patchi].name(),
204  meshTarget.boundary()[patchi].name()
205  );
206  }
207  else
208  {
209  cuttingPatchTable.insert
210  (
211  meshTarget.boundaryMesh()[patchi].name(),
212  -1
213  );
214  }
215  }
216 
217  mapSubMesh(meshSource, meshTarget, patchMap, cuttingPatchTable.toc());
218 }
219 
220 
221 wordList addProcessorPatches
222 (
223  const fvMesh& meshTarget,
224  const wordList& cuttingPatches
225 )
226 {
227  // Add the processor patches to the cutting list
228  HashTable<label> cuttingPatchTable;
229  forAll (cuttingPatches, i)
230  {
231  cuttingPatchTable.insert(cuttingPatches[i], i);
232  }
233 
234  forAll (meshTarget.boundary(), patchi)
235  {
236  if (isA<processorFvPatch>(meshTarget.boundary()[patchi]))
237  {
238  if
239  (
240  !cuttingPatchTable.found
241  (
242  meshTarget.boundaryMesh()[patchi].name()
243  )
244  )
245  {
246  cuttingPatchTable.insert
247  (
248  meshTarget.boundaryMesh()[patchi].name(),
249  -1
250  );
251  }
252  }
253  }
254 
255  return cuttingPatchTable.toc();
256 }
257 
258 
259 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
260 
261 int main(int argc, char *argv[])
262 {
263  #include "setRoots.H"
264  #include "createTimes.H"
265 
266  HashTable<word> patchMap;
267  wordList cuttingPatches;
268 
269  if (!consistent)
270  {
271  IOdictionary mapFieldsDict
272  (
273  IOobject
274  (
275  "mapFieldsDict",
276  runTimeTarget.system(),
277  runTimeTarget,
278  IOobject::MUST_READ,
279  IOobject::NO_WRITE,
280  false
281  )
282  );
283 
284  mapFieldsDict.lookup("patchMap") >> patchMap;
285 
286  mapFieldsDict.lookup("cuttingPatches") >> cuttingPatches;
287  }
288 
289  if (parallelSource && !parallelTarget)
290  {
291  IOdictionary decompositionDict
292  (
293  IOobject
294  (
295  "decomposeParDict",
296  runTimeSource.system(),
297  runTimeSource,
298  IOobject::MUST_READ,
299  IOobject::NO_WRITE
300  )
301  );
302 
303  int nProcs(readInt(decompositionDict.lookup("numberOfSubdomains")));
304 
305  Info<< "Create target mesh\n" << endl;
306 
307  fvMesh meshTarget
308  (
309  IOobject
310  (
311  fvMesh::defaultRegion,
312  runTimeTarget.timeName(),
313  runTimeTarget
314  )
315  );
316 
317  Info<< "Target mesh size: " << meshTarget.nCells() << endl;
318 
319  for (int procI=0; procI<nProcs; procI++)
320  {
321  Info<< nl << "Source processor " << procI << endl;
322 
323  Time runTimeSource
324  (
325  Time::controlDictName,
327  caseDirSource/fileName(word("processor") + name(procI))
328  );
329 
330  #include "setTimeIndex.H"
331 
332  fvMesh meshSource
333  (
334  IOobject
335  (
336  fvMesh::defaultRegion,
337  runTimeSource.timeName(),
338  runTimeSource
339  )
340  );
341 
342  Info<< "mesh size: " << meshSource.nCells() << endl;
343 
344  if (consistent)
345  {
346  mapConsistentSubMesh(meshSource, meshTarget);
347  }
348  else
349  {
350  mapSubMesh(meshSource, meshTarget, patchMap, cuttingPatches);
351  }
352  }
353  }
354  else if (!parallelSource && parallelTarget)
355  {
356  IOdictionary decompositionDict
357  (
358  IOobject
359  (
360  "decomposeParDict",
361  runTimeTarget.system(),
362  runTimeTarget,
363  IOobject::MUST_READ,
364  IOobject::NO_WRITE
365  )
366  );
367 
368  int nProcs(readInt(decompositionDict.lookup("numberOfSubdomains")));
369 
370  Info<< "Create source mesh\n" << endl;
371 
372  #include "setTimeIndex.H"
373 
374  fvMesh meshSource
375  (
376  IOobject
377  (
378  fvMesh::defaultRegion,
379  runTimeSource.timeName(),
380  runTimeSource
381  )
382  );
383 
384  Info<< "Source mesh size: " << meshSource.nCells() << endl;
385 
386  for (int procI=0; procI<nProcs; procI++)
387  {
388  Info<< nl << "Target processor " << procI << endl;
389 
390  Time runTimeTarget
391  (
392  Time::controlDictName,
394  caseDirTarget/fileName(word("processor") + name(procI))
395  );
396 
397  fvMesh meshTarget
398  (
399  IOobject
400  (
401  fvMesh::defaultRegion,
402  runTimeTarget.timeName(),
403  runTimeTarget
404  )
405  );
406 
407  Info<< "mesh size: " << meshTarget.nCells() << endl;
408 
409  if (consistent)
410  {
411  mapConsistentSubMesh(meshSource, meshTarget);
412  }
413  else
414  {
415  mapSubMesh
416  (
417  meshSource,
418  meshTarget,
419  patchMap,
420  addProcessorPatches(meshTarget, cuttingPatches)
421  );
422  }
423  }
424  }
425  else if (parallelSource && parallelTarget)
426  {
427  IOdictionary decompositionDictSource
428  (
429  IOobject
430  (
431  "decomposeParDict",
432  runTimeSource.system(),
433  runTimeSource,
434  IOobject::MUST_READ,
435  IOobject::NO_WRITE
436  )
437  );
438 
439  int nProcsSource
440  (
441  readInt(decompositionDictSource.lookup("numberOfSubdomains"))
442  );
443 
444 
445  IOdictionary decompositionDictTarget
446  (
447  IOobject
448  (
449  "decomposeParDict",
450  runTimeTarget.system(),
451  runTimeTarget,
452  IOobject::MUST_READ,
453  IOobject::NO_WRITE
454  )
455  );
456 
457  int nProcsTarget
458  (
459  readInt(decompositionDictTarget.lookup("numberOfSubdomains"))
460  );
461 
462  List<boundBox> bbsTarget(nProcsTarget);
463  List<bool> bbsTargetSet(nProcsTarget, false);
464 
465  for (int procISource=0; procISource<nProcsSource; procISource++)
466  {
467  Info<< nl << "Source processor " << procISource << endl;
468 
469  Time runTimeSource
470  (
471  Time::controlDictName,
473  caseDirSource/fileName(word("processor") + name(procISource))
474  );
475 
476  #include "setTimeIndex.H"
477 
478  fvMesh meshSource
479  (
480  IOobject
481  (
482  fvMesh::defaultRegion,
483  runTimeSource.timeName(),
484  runTimeSource
485  )
486  );
487 
488  Info<< "mesh size: " << meshSource.nCells() << endl;
489 
490  boundBox bbSource(meshSource.bounds());
491 
492  for (int procITarget=0; procITarget<nProcsTarget; procITarget++)
493  {
494  if
495  (
496  !bbsTargetSet[procITarget]
497  || (
498  bbsTargetSet[procITarget]
499  && bbsTarget[procITarget].overlaps(bbSource)
500  )
501  )
502  {
503  Info<< nl << "Target processor " << procITarget << endl;
504 
505  Time runTimeTarget
506  (
507  Time::controlDictName,
509  caseDirTarget/fileName(word("processor")
510  + name(procITarget))
511  );
512 
513  fvMesh meshTarget
514  (
515  IOobject
516  (
517  fvMesh::defaultRegion,
518  runTimeTarget.timeName(),
519  runTimeTarget
520  )
521  );
522 
523  Info<< "mesh size: " << meshTarget.nCells() << endl;
524 
525  bbsTarget[procITarget] = meshTarget.bounds();
526  bbsTargetSet[procITarget] = true;
527 
528  if (bbsTarget[procITarget].overlaps(bbSource))
529  {
530  if (consistent)
531  {
532  mapConsistentSubMesh(meshSource, meshTarget);
533  }
534  else
535  {
536  mapSubMesh
537  (
538  meshSource,
539  meshTarget,
540  patchMap,
541  addProcessorPatches(meshTarget, cuttingPatches)
542  );
543  }
544  }
545  }
546  }
547  }
548  }
549  else
550  {
551  #include "setTimeIndex.H"
552 
553  Info<< "Create meshes\n" << endl;
554 
555  fvMesh meshSource
556  (
557  IOobject
558  (
559  fvMesh::defaultRegion,
560  runTimeSource.timeName(),
561  runTimeSource
562  )
563  );
564 
565  fvMesh meshTarget
566  (
567  IOobject
568  (
569  fvMesh::defaultRegion,
570  runTimeTarget.timeName(),
571  runTimeTarget
572  )
573  );
574 
575  Info<< "Source mesh size: " << meshSource.nCells() << tab
576  << "Target mesh size: " << meshTarget.nCells() << nl << endl;
577 
578  if (consistent)
579  {
580  mapConsistentMesh(meshSource, meshTarget);
581  }
582  else
583  {
584  mapSubMesh(meshSource, meshTarget, patchMap, cuttingPatches);
585  }
586  }
587 
588  Info<< "\nEnd\n" << endl;
589 
590  return 0;
591 }
592 
593 
594 // ************************ vim: set sw=4 sts=4 et: ************************ //