DAUtility.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2 
3  DAFoam : Discrete Adjoint with OpenFOAM
4  Version : v4
5 
6 \*---------------------------------------------------------------------------*/
7 
8 #include "DAUtility.H"
9 
10 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
11 
12 namespace Foam
13 {
14 
15 // Constructors
17 {
18 }
19 
21 {
22 }
23 
25  PyObject* pyDict,
26  dictionary& ofDict)
27 {
28  /*
29  Description:
30  Parse a Python dictionary to an OpenFOAM dictionary
31  Only support a certain number of data types
32 
33  Input:
34  pyDict: Pytion dictionary
35 
36  Output:
37  ofDict: OpenFoam dictionary
38 
39  Example:
40  We support two types of pyDict, one is to set the key value
41  as a list and have the value type as the first element.
42  This is needed in pyDAFoam.py to reminder users the date type
43  of a key value
44 
45  pyDict = {
46  "solverName": [str, "DASimpleFoam"],
47  "patches": [list, [1.0, 2.0]]
48  }
49 
50  The second type of pyDict is just like a usual dict:
51 
52  pyDict = {
53  "solverName": "DASimpleFoam",
54  "patches": [1.0, 2.0]
55  }
56 
57  Note: one can have multiple levels of sub dictionaries
58  */
59 
60  //PyObject_Print(pyDict,stdout,0);Info<<endl;
61 
62  // size of the pyOpions keys
63  Py_ssize_t dictSize = PyDict_Size(pyDict);
64  // all the keys
65  PyObject* keys = PyDict_Keys(pyDict);
66  // loop over all the keys in pyDict and assign their values
67  // to ofDict
68  for (label i = 0; i < dictSize; i++)
69  {
70  // the ith key
71  PyObject* keyI = PyList_GetItem(keys, i);
72  // convert it to UTF8 such that we can use it in C++
73  const char* keyUTF8 = PyUnicode_AsUTF8(keyI);
74 
75  //std::cout << "Key is "<<keyUTF8<<std::endl;
76  // the actual value of this key, NOTE: it is a list
77  // the 1st one is its type and the 2nd one is its value
78  // this is because we need to make sure users prescribe
79  // the right data type in pyDAFoam.py, so we set the
80  // value like this in pyDAFoam.py:
81  // defOptions = {
82  // "solverName": [str, "DASimpleFOam"] }
83  // however, the above has one exception when using subDict
84  // in this case, we only have one element
85  PyObject* value = PyDict_GetItem(pyDict, keyI);
86  const char* valueTypeTmp = Py_TYPE(value)->tp_name;
87  PyObject* value1;
88  if (word(valueTypeTmp) == "list")
89  {
90  // there are two possibilities for this case
91  // for nonSubDict case, we need to have pyDict format
92  // like this: { "solverName": [str, "DASimpleDAFoam"] }
93  // however, for subDicts, we dont specify the type property
94  PyObject* value0 = PyList_GetItem(value, 0);
95  const char* valueTypeTmp0 = Py_TYPE(value0)->tp_name;
96  if (word(valueTypeTmp0) == "type")
97  {
98  // this is nonSubdict case
99  // the second element of value is the actual value
100  value1 = PyList_GetItem(value, 1);
101  }
102  else
103  {
104  // it is for subdicts, so value1=value
105  value1 = value;
106  }
107  }
108  else
109  {
110  // if we only have one element, set value1 = value
111  value1 = value;
112  }
113 
114  //PyObject_Print(value1,stdout,0);Info<<endl;
115  // get the type of value1
116  const char* valueType = Py_TYPE(value1)->tp_name;
117  bool overwrite = true;
118  if (word(valueType) == "str")
119  {
120  const char* valSet = PyUnicode_AsUTF8(value1);
121  ofDict.add(keyUTF8, word(valSet), overwrite);
122  }
123  else if (word(valueType) == "int")
124  {
125  long valSet = PyLong_AsLong(value1);
126  ofDict.add(keyUTF8, label(valSet), overwrite);
127  }
128  else if (word(valueType) == "bool")
129  {
130  label valSet = PyObject_IsTrue(value1);
131  ofDict.add(keyUTF8, valSet, overwrite);
132  }
133  else if (word(valueType) == "float")
134  {
135  scalar valSet = PyFloat_AS_DOUBLE(value1);
136  ofDict.add(keyUTF8, valSet, overwrite);
137  }
138  else if (word(valueType) == "list")
139  {
140  // size of the list
141  Py_ssize_t listSize = PyList_Size(value1);
142 
143  //Info<<listSize<<endl;
144  // create OpenFOAM lists to hold this list
145  // we create all the possible types
146  scalarList valSetScalar;
147  valSetScalar.setSize(label(listSize));
148  labelList valSetLabel;
149  valSetLabel.setSize(label(listSize));
150  List<word> valSetWord;
151  valSetWord.setSize(label(listSize));
152  // we also support scalarListList or labelListList
153  List<List<scalar>> valSetScalarList;
154  valSetScalarList.setSize(label(listSize));
155  List<List<label>> valSetLabelList;
156  valSetLabelList.setSize(label(listSize));
157 
158  label isListListScalar = 0;
159  label isListListLabel = 0;
160 
161  // need to check what type of list this is
162  // by checking its first element
163  PyObject* tmp = PyList_GetItem(value1, 0);
164  const char* tmpType = Py_TYPE(tmp)->tp_name;
165  word tmpTypeWord = word(tmpType);
166  // assign value to the OpenFOAM list
167  for (label j = 0; j < listSize; j++)
168  {
169  PyObject* valueListJ = PyList_GetItem(value1, j);
170  if (tmpTypeWord == "str")
171  {
172  const char* valSet = PyUnicode_AsUTF8(valueListJ);
173  valSetWord[j] = word(valSet);
174  }
175  else if (tmpTypeWord == "int")
176  {
177  long valSet = PyLong_AsLong(valueListJ);
178  valSetLabel[j] = label(valSet);
179  }
180  else if (tmpTypeWord == "float")
181  {
182  scalar valSet = PyFloat_AS_DOUBLE(valueListJ);
183  valSetScalar[j] = valSet;
184  }
185  else if (tmpTypeWord == "bool")
186  {
187  label valSet = PyObject_IsTrue(valueListJ);
188  valSetLabel[j] = valSet;
189  }
190  else if (tmpTypeWord == "list")
191  {
192  // it is a ListList, we need to add values to
193  // valSetScalarList or valSetLabelList
194  // NOTE: we support only labelListList or scalarListList
195  // size of the list
196  Py_ssize_t listSize1 = PyList_Size(valueListJ);
197  PyObject* tmp1 = PyList_GetItem(valueListJ, 0);
198  const char* tmpType1 = Py_TYPE(tmp1)->tp_name;
199  word tmpTypeWord1 = word(tmpType1);
200 
201  valSetLabelList[j].setSize(listSize1);
202  valSetScalarList[j].setSize(listSize1);
203 
204  // assign value to the OpenFOAM listList
205  for (label k = 0; k < listSize1; k++)
206  {
207  PyObject* valueListK = PyList_GetItem(valueListJ, k);
208  if (tmpTypeWord1 == "int")
209  {
210  long valSet = PyLong_AsLong(valueListK);
211  valSetLabelList[j][k] = label(valSet);
212  isListListLabel = 1;
213  }
214  else if (tmpTypeWord1 == "float")
215  {
216  scalar valSet = PyFloat_AS_DOUBLE(valueListK);
217  valSetScalarList[j][k] = valSet;
218  isListListScalar = 1;
219  }
220  else
221  {
222  FatalErrorIn("pyDict2OFDict") << "listList only supports double and int"
223  << abort(FatalError);
224  }
225  }
226  }
227  else
228  {
229  FatalErrorIn("pyDict2OFDict") << "Type: <" << tmpTypeWord << "> for " << keyUTF8
230  << " list is not supported! Options are:"
231  << " str, int, bool, float, and list!"
232  << abort(FatalError);
233  }
234  }
235 
236  // add the list to the ofDict dict
237  if (tmpTypeWord == "str")
238  {
239  ofDict.add(keyUTF8, valSetWord, overwrite);
240  }
241  else if (tmpTypeWord == "int")
242  {
243  ofDict.add(keyUTF8, valSetLabel, overwrite);
244  }
245  else if (tmpTypeWord == "float")
246  {
247  ofDict.add(keyUTF8, valSetScalar, overwrite);
248  }
249  else if (tmpTypeWord == "bool")
250  {
251  ofDict.add(keyUTF8, valSetLabel, overwrite);
252  }
253 
254  // if there are listLists, add them separately
255  if (isListListLabel)
256  {
257  ofDict.add(keyUTF8, valSetLabelList, overwrite);
258  }
259 
260  if (isListListScalar)
261  {
262  ofDict.add(keyUTF8, valSetScalarList, overwrite);
263  }
264  }
265  else if (word(valueType) == "dict")
266  {
267  // if its a subdict, recursely call this function
268  dictionary subDict;
269  DAUtility::pyDict2OFDict(value1, subDict);
270  ofDict.add(keyUTF8, subDict, overwrite);
271  }
272  else
273  {
274  FatalErrorIn("pyDict2OFDict") << "Type: " << valueType << " for " << keyUTF8
275  << " is not supported! Options are: str, int, float, bool, list, and dict!"
276  << abort(FatalError);
277  }
278  //std::cout << "My type is " << valueType << std::endl;
279  }
280 }
281 
283  Vec vecIn,
284  const word prefix)
285 {
286  /*
287  Description:
288  Read a vector in binary form
289 
290  Input:
291  vecIn: a Petsc vector to read values into (also output)
292  prefix: Name of the Petsc vector from disk
293 
294  Example:
295  If the vector storing in the disk reads: dFdWVector.bin
296  Then read the vector using:
297  Vec dFdW;
298  codes to initialize vector dFdW....
299  readVectorBinary(dFdW,"dFdwVector");
300  NOTE: the prefix does not include ".bin"
301  */
302 
303  std::ostringstream fileNameStream("");
304  fileNameStream << prefix << ".bin";
305  word fileName = fileNameStream.str();
306 
307  PetscViewer viewer;
308  PetscViewerBinaryOpen(PETSC_COMM_WORLD, fileName.c_str(), FILE_MODE_READ, &viewer);
309  VecLoad(vecIn, viewer);
310  PetscViewerDestroy(&viewer);
311 
312  return;
313 }
314 
316  const Vec vecIn,
317  const word prefix)
318 {
319  /*
320  Description:
321  Write a vector in binary form
322 
323  Input:
324  vecIn: a Petsc vector to write to disk
325  prefix: Name of the Petsc vector to write to disk
326 
327  Example:
328  Vec dFdW;
329  codes to initialize vector dFdW....
330  writeVectorBinary(dFdW,"dFdwVector");
331  This will write the dFdW vector to disk with name "dFdWVector.bin"
332  NOTE: the prefix does not include ".bin"
333  */
334 
335  std::ostringstream fileNameStream("");
336  fileNameStream << prefix << ".bin";
337  word fileName = fileNameStream.str();
338 
339  PetscViewer viewer;
340  PetscViewerBinaryOpen(PETSC_COMM_WORLD, fileName.c_str(), FILE_MODE_WRITE, &viewer);
341  VecView(vecIn, viewer);
342  PetscViewerDestroy(&viewer);
343 
344  return;
345 }
346 
348  const Vec vecIn,
349  const word prefix)
350 {
351  /*
352  Description:
353  Write a vector in ASCII form
354  Input:
355  vecIn: a Petsc vector to write to disk
356  prefix: Name of the Petsc vector to write to disk
357  Example:
358  Vec dFdW;
359  codes to initialize vector dFdW....
360  writeVectorASCII(dFdW,"dFdwVector");
361  This will write the dFdW vector to disk with name "dFdWVector.dat"
362  NOTE: the prefix does not include ".dat"
363  */
364 
365  std::ostringstream fileNameStream("");
366  fileNameStream << prefix << ".dat";
367  word fileName = fileNameStream.str();
368 
369  PetscViewer viewer;
370  PetscViewerASCIIOpen(PETSC_COMM_WORLD, fileName.c_str(), &viewer);
371  PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_MATLAB); // write all the digits
372  VecView(vecIn, viewer);
373  PetscViewerDestroy(&viewer);
374 
375  return;
376 }
377 
379  Mat matIn,
380  const word prefix)
381 {
382  /*
383  Description:
384  Read a matrix in binary form
385 
386  Input:
387  matIn: a Petsc matrix to read values into (also output)
388  prefix: Name of the Petsc matrix from disk
389 
390  Example:
391  If the matrix storing in the disk reads: dRdWMat.bin
392  Then read the matrix using:
393  Mat dRdW;
394  codes to initialize matrix dRdW....
395  readMatrixBinary(dRdW,"dRdwVector");
396  NOTE: the prefix does not include ".bin"
397  */
398 
399  std::ostringstream fileNameStream("");
400  fileNameStream << prefix << ".bin";
401  word fileName = fileNameStream.str();
402 
403  PetscViewer viewer;
404  PetscViewerBinaryOpen(PETSC_COMM_WORLD, fileName.c_str(), FILE_MODE_READ, &viewer);
405  MatLoad(matIn, viewer);
406  PetscViewerDestroy(&viewer);
407 
408  return;
409 }
410 
412  const Mat matIn,
413  const word prefix)
414 {
415  /*
416  Description:
417  Write a matrix in binary form
418 
419  Input:
420  matIn: a Petsc matrix to write to disk
421  prefix: Name of the Petsc matrix to write to disk
422 
423  Example:
424  Mat dRdW;
425  codes to initialize matrix dRdW....
426  writeMatrixBinary(dRdW,"dRdWMat");
427  This will write the dRdW matrix to disk with name "dRdWMat.bin"
428  NOTE: the prefix does not include ".bin"
429  */
430 
431  std::ostringstream fileNameStream("");
432  fileNameStream << prefix << ".bin";
433  word fileName = fileNameStream.str();
434 
435  PetscViewer viewer;
436  PetscViewerBinaryOpen(PETSC_COMM_WORLD, fileName.c_str(), FILE_MODE_WRITE, &viewer);
437  MatView(matIn, viewer);
438  PetscViewerDestroy(&viewer);
439 
440  return;
441 }
442 
444  const Mat matIn,
445  const word prefix)
446 {
447  /*
448  Description:
449  Write a matrix in ASCII form
450 
451  Input:
452  matIn: a Petsc matrix to write to disk
453  prefix: Name of the Petsc matrix to write to disk
454 
455  Example:
456  Mat dRdW;
457  codes to initialize matrix dRdW....
458  writeMatrixASCII(dRdW,"dRdWMat");
459  This will write the dRdW matrix to disk with name "dRdWMat.dat"
460  NOTE: the prefix does not include ".dat"
461  */
462 
463  std::ostringstream fileNameStream("");
464  fileNameStream << prefix << ".dat";
465  word fileName = fileNameStream.str();
466 
467  PetscViewer viewer;
468  PetscViewerASCIIOpen(PETSC_COMM_WORLD, fileName.c_str(), &viewer);
469  MatView(matIn, viewer);
470  PetscViewerDestroy(&viewer);
471 
472  return;
473 }
474 
476  const dictionary& allOptions,
477  volScalarField& var,
478  const label printToScreen)
479 {
480  /*
481  Description:
482  Bound a field variable according to the bounds defined
483  in the allOptions dict. This is a overload function for volScalarField
484 
485  Input:
486  allOptions: a dictionary that has upper and lower bound values
487  We need to give specific name for the bounds, i.e.,
488 
489  variable name + Max -> setting upper bound
490 
491  variable name + Min -> setting lower bound
492 
493  Input & Output:
494  var: an OpenFOAM field variable to be bounded
495 
496  Example:
497  dictionary allOptions;
498  allOptions.set("pMax", 120000.0);
499  allOptions.set("pMin", 80000.0);
500  DAUtility daUtil;
501  volScalarField p = .... // initialize p
502  daUtil.boundVar(allOptions, p);
503  */
504 
505  const scalar vGreat = 1e200;
506  label useUpperBound = 0, useLowerBound = 0;
507 
508  dictionary varBoundsDict = allOptions.subDict("primalVarBounds");
509 
510  word lowerBoundName = var.name() + "Min";
511  word upperBoundName = var.name() + "Max";
512 
513  scalar varMin = varBoundsDict.lookupOrDefault<scalar>(lowerBoundName, -vGreat);
514  scalar varMax = varBoundsDict.lookupOrDefault<scalar>(upperBoundName, vGreat);
515 
516  forAll(var, cellI)
517  {
518  if (var[cellI] <= varMin)
519  {
520  var[cellI] = varMin;
521  useLowerBound = 1;
522  }
523  if (var[cellI] >= varMax)
524  {
525  var[cellI] = varMax;
526  useUpperBound = 1;
527  }
528  }
529 
530  forAll(var.boundaryField(), patchI)
531  {
532  forAll(var.boundaryField()[patchI], faceI)
533  {
534  if (var.boundaryFieldRef()[patchI][faceI] <= varMin)
535  {
536  var.boundaryFieldRef()[patchI][faceI] = varMin;
537  useLowerBound = 1;
538  }
539  if (var.boundaryFieldRef()[patchI][faceI] >= varMax)
540  {
541  var.boundaryFieldRef()[patchI][faceI] = varMax;
542  useUpperBound = 1;
543  }
544  }
545  }
546 
547  if (printToScreen)
548  {
549  if (useUpperBound)
550  {
551  Info << "Bounding " << var.name() << "<" << varMax << endl;
552  }
553  if (useLowerBound)
554  {
555  Info << "Bounding " << var.name() << ">" << varMin << endl;
556  }
557  }
558 
559  return;
560 }
561 
563  const dictionary& allOptions,
564  volVectorField& var,
565  const label printToScreen)
566 {
567  /*
568  Description:
569  Bound a field variable according to the bounds defined
570  in the allOptions dict. This is a overload function for volVectorField
571 
572  Input:
573  allOptions: a dictionary that has upper and lower bound values
574  We need to give specific name for the bounds, i.e.,
575 
576  variable name + Max -> setting upper bound
577 
578  variable name + Min -> setting lower bound
579 
580  Input & Output:
581  var: an OpenFOAM field variable to be bounded
582 
583  Example:
584  dictionary allOptions;
585  allOptions.set("pMax", 120000.0);
586  allOptions.set("pMin", 80000.0);
587  DAUtility daUtil;
588  volScalarField p = .... // initialize p
589  daUtil.boundVar(allOptions, p);
590  */
591 
592  const scalar vGreat = 1e200;
593  label useUpperBound = 0, useLowerBound = 0;
594 
595  dictionary varBoundsDict = allOptions.subDict("primalVarBounds");
596 
597  word lowerBoundName = var.name() + "Min";
598  word upperBoundName = var.name() + "Max";
599 
600  scalar varMin = varBoundsDict.lookupOrDefault<scalar>(lowerBoundName, -vGreat);
601  scalar varMax = varBoundsDict.lookupOrDefault<scalar>(upperBoundName, vGreat);
602 
603  forAll(var, cellI)
604  {
605  for (label i = 0; i < 3; i++)
606  {
607  if (var[cellI][i] <= varMin)
608  {
609  var[cellI][i] = varMin;
610  useLowerBound = 1;
611  }
612  if (var[cellI][i] >= varMax)
613  {
614  var[cellI][i] = varMax;
615  useUpperBound = 1;
616  }
617  }
618  }
619 
620  forAll(var.boundaryField(), patchI)
621  {
622  forAll(var.boundaryField()[patchI], faceI)
623  {
624  for (label i = 0; i < 3; i++)
625  {
626  if (var.boundaryFieldRef()[patchI][faceI][i] <= varMin)
627  {
628  var.boundaryFieldRef()[patchI][faceI][i] = varMin;
629  useLowerBound = 1;
630  }
631  if (var.boundaryFieldRef()[patchI][faceI][i] >= varMax)
632  {
633  var.boundaryFieldRef()[patchI][faceI][i] = varMax;
634  useUpperBound = 1;
635  }
636  }
637  }
638  }
639 
640  if (printToScreen)
641  {
642  if (useUpperBound)
643  {
644  Info << "Bounding " << var.name() << "<" << varMax << endl;
645  }
646  if (useLowerBound)
647  {
648  Info << "Bounding " << var.name() << ">" << varMin << endl;
649  }
650  }
651 
652  return;
653 }
654 
655 globalIndex DAUtility::genGlobalIndex(const label localIndexSize)
656 {
657  /*
658  Generate a glocal index system based on the local index size
659  such that we can use it to map a local index to a global one
660 
661  Input:
662  -----
663  localIndexSize: the SIZE of local index
664 
665  Output:
666  ------
667  globalIndex object: the global index object to map a local index
668  to a global index
669 
670  Example:
671  --------
672  If the local index reads:
673 
674  On processor 0:
675  labelList sampleList = {0, 1, 2};
676  globalIndex glbSample = genGlobalIndex(sampleList.size());
677 
678  On processor 1:
679  labelList sampleList = {0, 1};
680  globalIndex glbSample = genGlobalIndex(sampleList.size());
681 
682  After each processor calls genGlobalIndex and get the glbSample
683  object, we can use it to map a local index to a global one,
684  e.g., on processor 0, if we call:
685 
686  label glxIdx = glbSample.toGlobal(1);
687 
688  it will return glbIdx = 1;
689  However, on processor 1, if we call
690 
691  label glxIdx = glbSample.toGlobal(1);
692 
693  it will return glbIdx = 4;
694 
695  The date storage structure is illustrated as follows
696 
697  global index -> 0 1 2 3 4
698  local index -> 0 1 2 0 1
699  ----- ===
700  proc0 proc1
701 
702  */
703  globalIndex result(localIndexSize);
704  return result;
705 }
706 
708  const scalar val,
709  const scalar refVal,
710  const scalar tol)
711 {
712  /*
713  Description:
714  check whether a value is close to a reference value by a tolerance
715 
716  Input:
717  val: value to check
718  refVal: reference value
719  tol: tolerance
720 
721  Output:
722  return: 1 means fabs(val-refVal) < tol
723  */
724  if (fabs(val - refVal) < tol)
725  {
726  return 1;
727  }
728  else
729  {
730  return 0;
731  }
732 }
733 
735  const SolverPerformance<scalar>& solverP,
736  const label printToScreen,
737  const word varName,
738  scalar& primalMaxRes)
739 {
740  /*
741  Description:
742  Setup maximal residual control and print the residual as needed
743  */
744 
745  // calculate the initial residual mag and set it to primalResidualNorms_
746 
747  scalar initRes = solverP.initialResidual();
748 
749  if (initRes > primalMaxRes)
750  {
751  primalMaxRes = initRes;
752  }
753 
754  if (printToScreen)
755  {
756  Info << varName << " initRes: " << solverP.initialResidual()
757  << " finalRes: " << solverP.finalResidual()
758  << " nIters: " << solverP.nIterations() << endl;
759  }
760 }
761 
763  const SolverPerformance<vector>& solverP,
764  const label printToScreen,
765  const word varName,
766  scalar& primalMaxRes)
767 {
768  /*
769  Description:
770  Setup maximal residual control and print the residual as needed
771  */
772 
773  // calculate the initial residual mag and set it to primalResidualNorms_
774 
775  // for vectors, we need to use the median value for the residual
776  // this is because we often need to run 2D simulations with symmetry
777  // BC, so one component of the residual vector, which is related to the symmetry BC,
778  // may be high while the other two components' residuals are low.
779  // In this case, we can use the median value for the residual vector, which better
780  // represents the convergence.
781 
782  vector initRes = solverP.initialResidual();
783  scalarList initResList = {initRes[0], initRes[1], initRes[2]};
784  sort(initResList);
785 
786  if (initResList[1] > primalMaxRes)
787  {
788  primalMaxRes = initResList[1];
789  }
790 
791  if (printToScreen)
792  {
793  for (label i = 0; i < 3; i++)
794  {
795  Info << varName << i
796  << " initRes: " << solverP.initialResidual()[i]
797  << " finalRes: " << solverP.finalResidual()[i]
798  << " nIters: " << solverP.nIterations()[i] << endl;
799  }
800  }
801 }
802 
804  const primitiveMesh& mesh,
805  const point& point)
806 {
807  /*
808  A self-defined findCell function. We will need to cast fvMesh to primitiveMesh
809  and then call primitiveMesh's findCell. For some reasons, the fvMesh's findCell
810  did not work correctly in ADR mode...
811  */
812 
813  label cellI = mesh.findCell(point);
814  return cellI;
815 }
816 
818  const fvMesh& mesh,
819  const word fieldName,
820  const word fieldType)
821 {
822  /*
823  Whether a field is readable from the disk. We will just check if the
824  field exists on the disk (typically in the 0 folder)
825  */
826 
827  label readable = 0;
828 
829  IOobject headerFile(
830  fieldName,
831  mesh.time().timeName(),
832  mesh,
833  IOobject::NO_READ);
834 
835  if (fieldType == "volScalarField")
836  {
837  if (headerFile.typeHeaderOk<volScalarField>(true))
838  {
839  readable = 1;
840  }
841  }
842  else if (fieldType == "volVectorField")
843  {
844  if (headerFile.typeHeaderOk<volVectorField>(true))
845  {
846  readable = 1;
847  }
848  }
849  else
850  {
851  FatalErrorIn("DAUtility::isFieldReadable")
852  << "fieldType not supported! Options are volScalarField or volVectorField"
853  << abort(FatalError);
854  }
855 
856  return readable;
857 }
858 
859 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
860 
861 } // End namespace Foam
862 
863 // ************************************************************************* //
Foam::DAUtility::readMatrixBinary
static void readMatrixBinary(Mat matIn, const word prefix)
read petsc matrix in binary format
Definition: DAUtility.C:378
Foam::DAUtility::readVectorBinary
static void readVectorBinary(Vec vecIn, const word prefix)
read petsc vector in binary format
Definition: DAUtility.C:282
Foam::DAUtility::DAUtility
DAUtility()
Constructors.
Definition: DAUtility.C:16
Foam::DAUtility::~DAUtility
virtual ~DAUtility()
Destructor.
Definition: DAUtility.C:20
DAUtility.H
Foam::DAUtility::isValueCloseToRef
static label isValueCloseToRef(const scalar val, const scalar refVal, const scalar tol=1.0e-6)
check whether a value is close to a reference value by a tolerance
Definition: DAUtility.C:707
mesh
fvMesh & mesh
Definition: createRefsHeatTransfer.H:4
Foam::DAUtility::pyDict2OFDict
static void pyDict2OFDict(PyObject *pyDict, dictionary &ofDict)
convert a python dictionary object to OpenFoam dictionary
Definition: DAUtility.C:24
Foam::DAUtility::boundVar
static void boundVar(const dictionary &allOptions, volScalarField &var, const label printToScreen)
bound a volScalar variable based on parametes defined in DAOption::allOptions_
Definition: DAUtility.C:475
Foam::DAUtility::writeVectorBinary
static void writeVectorBinary(const Vec vecIn, const word prefix)
write petsc vector in binary format
Definition: DAUtility.C:315
Foam::DAUtility::primalResidualControl
static void primalResidualControl(const SolverPerformance< scalar > &solverP, const label printToScreen, const word varName, scalar &primalMaxRes)
control when to print the residual and also compute the maxInitRes
Definition: DAUtility.C:734
k
dimensionedScalar & k
Definition: createRefsHeatTransfer.H:6
Foam
Definition: checkGeometry.C:32
allOptions
const dictionary & allOptions
Definition: createRefsPimpleDyM.H:14
forAll
forAll(nuTilda1, cellI)
Definition: nuTilda1EqnIrkPimple.H:19
Foam::DAUtility::writeMatrixBinary
static void writeMatrixBinary(const Mat matIn, const word prefix)
write petsc matrix in binary format
Definition: DAUtility.C:411
Foam::DAUtility::myFindCell
static label myFindCell(const primitiveMesh &mesh, const point &point)
Definition: DAUtility.C:803
Foam::DAUtility::writeVectorASCII
static void writeVectorASCII(const Vec vecIn, const word prefix)
write petsc vector in ascii format
Definition: DAUtility.C:347
Foam::DAUtility::writeMatrixASCII
static void writeMatrixASCII(const Mat matIn, const word prefix)
write petsc matrix in ascii format
Definition: DAUtility.C:443
Foam::DAUtility::isFieldReadable
static label isFieldReadable(const fvMesh &mesh, const word fieldName, const word fieldType)
Definition: DAUtility.C:817
Foam::DAUtility::genGlobalIndex
static globalIndex genGlobalIndex(const label localIndexSize)
generate global index numbering for local-global index transferring
Definition: DAUtility.C:655