GeometricField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 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 "GeometricField.H"
27 #include "Time.H"
28 #include "demandDrivenData.H"
29 #include "dictionary.H"
30 #include "localIOdictionary.H"
31 #include "data.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 #define checkField(gf1, gf2, op) \
36 if ((gf1).mesh() != (gf2).mesh()) \
37 { \
38  FatalErrorInFunction \
39  << "different mesh for fields " \
40  << (gf1).name() << " and " << (gf2).name() \
41  << " during operatrion " << op \
42  << abort(FatalError); \
43 }
44 
45 
46 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
47 
48 template<class Type, template<class> class PatchField, class GeoMesh>
50 (
51  const dictionary& dict
52 )
53 {
54  Internal::readField(dict, "internalField");
55 
56  boundaryField_.readField(*this, dict.subDict("boundaryField"));
57 
58  if (dict.found("referenceLevel"))
59  {
60  Type fieldAverage(pTraits<Type>(dict.lookup("referenceLevel")));
61 
62  Field<Type>::operator+=(fieldAverage);
63 
64  forAll(boundaryField_, patchi)
65  {
66  boundaryField_[patchi] == boundaryField_[patchi] + fieldAverage;
67  }
68  }
69 }
70 
71 
72 template<class Type, template<class> class PatchField, class GeoMesh>
74 {
75  const localIOdictionary dict
76  (
77  IOobject
78  (
79  this->name(),
80  this->instance(),
81  this->local(),
82  this->db(),
83  IOobject::MUST_READ,
84  IOobject::NO_WRITE,
85  false
86  ),
87  typeName
88  );
89 
90  this->close();
91 
92  readFields(dict);
93 }
94 
95 
96 template<class Type, template<class> class PatchField, class GeoMesh>
98 {
99  if
100  (
101  this->readOpt() == IOobject::MUST_READ
102  || this->readOpt() == IOobject::MUST_READ_IF_MODIFIED
103  )
104  {
106  << "read option IOobject::MUST_READ or MUST_READ_IF_MODIFIED"
107  << " suggests that a read constructor for field " << this->name()
108  << " would be more appropriate." << endl;
109  }
110  else if
111  (
112  this->readOpt() == IOobject::READ_IF_PRESENT
113  && this->template typeHeaderOk<GeometricField<Type, PatchField, GeoMesh>>
114  (
115  true
116  )
117  )
118  {
119  readFields();
120 
121  // Check compatibility between field and mesh
122  if (this->size() != GeoMesh::size(this->mesh()))
123  {
124  FatalIOErrorInFunction(this->readStream(typeName))
125  << " number of field elements = " << this->size()
126  << " number of mesh elements = "
127  << GeoMesh::size(this->mesh())
128  << exit(FatalIOError);
129  }
130 
131  readOldTimeIfPresent();
132 
133  return true;
134  }
135 
136  return false;
137 }
138 
139 
140 template<class Type, template<class> class PatchField, class GeoMesh>
142 {
143  // Read the old time field if present
144  IOobject field0
145  (
146  this->name() + "_0",
147  this->time().timeName(),
148  this->db(),
149  IOobject::READ_IF_PRESENT,
150  IOobject::AUTO_WRITE,
151  this->registerObject()
152  );
153 
154  if
155  (
156  field0.template typeHeaderOk<GeometricField<Type, PatchField, GeoMesh>>
157  (
158  true
159  )
160  )
161  {
162  if (debug)
163  {
164  InfoInFunction << "Reading old time level for field"
165  << endl << this->info() << endl;
166  }
167 
168  field0Ptr_ = new GeometricField<Type, PatchField, GeoMesh>
169  (
170  field0,
171  this->mesh()
172  );
173 
174  field0Ptr_->timeIndex_ = timeIndex_ - 1;
175 
176  if (!field0Ptr_->readOldTimeIfPresent())
177  {
178  field0Ptr_->oldTime();
179  }
180 
181  return true;
182  }
183 
184  return false;
185 }
186 
187 
188 // * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * //
189 
190 template<class Type, template<class> class PatchField, class GeoMesh>
192 (
193  const IOobject& io,
194  const Mesh& mesh,
195  const dimensionSet& ds,
196  const word& patchFieldType
197 )
198 :
199  Internal(io, mesh, ds, false),
200  timeIndex_(this->time().timeIndex()),
201  field0Ptr_(nullptr),
202  fieldPrevIterPtr_(nullptr),
203  boundaryField_(mesh.boundary(), *this, patchFieldType)
204 {
205  if (debug)
206  {
207  InfoInFunction << "Creating temporary" << endl << this->info() << endl;
208  }
209 
210  readIfPresent();
211 }
212 
213 
214 template<class Type, template<class> class PatchField, class GeoMesh>
216 (
217  const IOobject& io,
218  const Mesh& mesh,
219  const dimensionSet& ds,
220  const wordList& patchFieldTypes,
221  const wordList& actualPatchTypes
222 )
223 :
224  Internal(io, mesh, ds, false),
225  timeIndex_(this->time().timeIndex()),
226  field0Ptr_(nullptr),
227  fieldPrevIterPtr_(nullptr),
228  boundaryField_(mesh.boundary(), *this, patchFieldTypes, actualPatchTypes)
229 {
230  if (debug)
231  {
232  InfoInFunction << "Creating temporary" << endl << this->info() << endl;
233  }
234 
235  readIfPresent();
236 }
237 
238 
239 template<class Type, template<class> class PatchField, class GeoMesh>
241 (
242  const IOobject& io,
243  const Mesh& mesh,
244  const dimensioned<Type>& dt,
245  const word& patchFieldType
246 )
247 :
248  Internal(io, mesh, dt, false),
249  timeIndex_(this->time().timeIndex()),
250  field0Ptr_(nullptr),
251  fieldPrevIterPtr_(nullptr),
252  boundaryField_(mesh.boundary(), *this, patchFieldType)
253 {
254  if (debug)
255  {
256  InfoInFunction << "Creating temporary" << endl << this->info() << endl;
257  }
258 
259  boundaryField_ == dt.value();
260 
261  readIfPresent();
262 }
263 
264 
265 template<class Type, template<class> class PatchField, class GeoMesh>
267 (
268  const IOobject& io,
269  const Mesh& mesh,
270  const dimensioned<Type>& dt,
271  const wordList& patchFieldTypes,
272  const wordList& actualPatchTypes
273 )
274 :
275  Internal(io, mesh, dt, false),
276  timeIndex_(this->time().timeIndex()),
277  field0Ptr_(nullptr),
278  fieldPrevIterPtr_(nullptr),
279  boundaryField_(mesh.boundary(), *this, patchFieldTypes, actualPatchTypes)
280 {
281  if (debug)
282  {
283  InfoInFunction << "Creating temporary" << endl << this->info() << endl;
284  }
285 
286  boundaryField_ == dt.value();
287 
288  readIfPresent();
289 }
290 
291 
292 template<class Type, template<class> class PatchField, class GeoMesh>
294 (
295  const IOobject& io,
296  const Internal& diField,
297  const PtrList<PatchField<Type>>& ptfl
298 )
299 :
300  Internal(io, diField),
301  timeIndex_(this->time().timeIndex()),
302  field0Ptr_(nullptr),
303  fieldPrevIterPtr_(nullptr),
304  boundaryField_(this->mesh().boundary(), *this, ptfl)
305 {
306  if (debug)
307  {
309  << "Constructing from components" << endl << this->info() << endl;
310  }
311 
312  readIfPresent();
313 }
314 
315 
316 template<class Type, template<class> class PatchField, class GeoMesh>
318 (
319  const IOobject& io,
320  const Mesh& mesh,
321  const dimensionSet& ds,
322  const Field<Type>& iField,
323  const PtrList<PatchField<Type>>& ptfl
324 )
325 :
326  Internal(io, mesh, ds, iField),
327  timeIndex_(this->time().timeIndex()),
328  field0Ptr_(nullptr),
329  fieldPrevIterPtr_(nullptr),
330  boundaryField_(mesh.boundary(), *this, ptfl)
331 {
332  if (debug)
333  {
335  << "Constructing from components" << endl << this->info() << endl;
336  }
337 
338  readIfPresent();
339 }
340 
341 
342 template<class Type, template<class> class PatchField, class GeoMesh>
344 (
345  const IOobject& io,
346  const Mesh& mesh
347 )
348 :
349  Internal(io, mesh, dimless, false),
350  timeIndex_(this->time().timeIndex()),
351  field0Ptr_(nullptr),
352  fieldPrevIterPtr_(nullptr),
353  boundaryField_(mesh.boundary())
354 {
355  readFields();
356 
357  // Check compatibility between field and mesh
358 
359  if (this->size() != GeoMesh::size(this->mesh()))
360  {
361  FatalIOErrorInFunction(this->readStream(typeName))
362  << " number of field elements = " << this->size()
363  << " number of mesh elements = " << GeoMesh::size(this->mesh())
364  << exit(FatalIOError);
365  }
366 
367  readOldTimeIfPresent();
368 
369  if (debug)
370  {
372  << "Finishing read-construction of" << endl << this->info() << endl;
373  }
374 }
375 
376 
377 template<class Type, template<class> class PatchField, class GeoMesh>
379 (
380  const IOobject& io,
381  const Mesh& mesh,
382  const dictionary& dict
383 )
384 :
385  Internal(io, mesh, dimless, false),
386  timeIndex_(this->time().timeIndex()),
387  field0Ptr_(nullptr),
388  fieldPrevIterPtr_(nullptr),
389  boundaryField_(mesh.boundary())
390 {
391  readFields(dict);
392 
393  // Check compatibility between field and mesh
394 
395  if (this->size() != GeoMesh::size(this->mesh()))
396  {
398  << " number of field elements = " << this->size()
399  << " number of mesh elements = " << GeoMesh::size(this->mesh())
400  << exit(FatalIOError);
401  }
402 
403  readOldTimeIfPresent();
404 
405  if (debug)
406  {
408  << "Finishing dictionary-construct of "
409  << endl << this->info() << endl;
410  }
411 }
412 
413 
414 template<class Type, template<class> class PatchField, class GeoMesh>
416 (
418 )
419 :
420  Internal(gf),
421  timeIndex_(gf.timeIndex()),
422  field0Ptr_(nullptr),
423  fieldPrevIterPtr_(nullptr),
424  boundaryField_(*this, gf.boundaryField_)
425 {
426  if (debug)
427  {
429  << "Constructing as copy" << endl << this->info() << endl;
430  }
431 
432  if (gf.field0Ptr_)
433  {
435  (
436  *gf.field0Ptr_
437  );
438  }
439 
440  this->writeOpt() = IOobject::NO_WRITE;
441 }
442 
443 
444 #ifndef NoConstructFromTmp
445 template<class Type, template<class> class PatchField, class GeoMesh>
447 (
449 )
450 :
451  Internal
452  (
453  const_cast<GeometricField<Type, PatchField, GeoMesh>&>(tgf()),
454  tgf.isTmp()
455  ),
456  timeIndex_(tgf().timeIndex()),
457  field0Ptr_(nullptr),
458  fieldPrevIterPtr_(nullptr),
459  boundaryField_(*this, tgf().boundaryField_)
460 {
461  if (debug)
462  {
464  << "Constructing from tmp" << endl << this->info() << endl;
465  }
466 
467  this->writeOpt() = IOobject::NO_WRITE;
468 
469  tgf.clear();
470 }
471 #endif
472 
473 
474 template<class Type, template<class> class PatchField, class GeoMesh>
476 (
477  const IOobject& io,
479 )
480 :
481  Internal(io, gf),
482  timeIndex_(gf.timeIndex()),
483  field0Ptr_(nullptr),
484  fieldPrevIterPtr_(nullptr),
485  boundaryField_(*this, gf.boundaryField_)
486 {
487  if (debug)
488  {
490  << "Constructing as copy resetting IO params"
491  << endl << this->info() << endl;
492  }
493 
494  if (!readIfPresent() && gf.field0Ptr_)
495  {
497  (
498  io.name() + "_0",
499  *gf.field0Ptr_
500  );
501  }
502 }
503 
504 
505 #ifndef NoConstructFromTmp
506 template<class Type, template<class> class PatchField, class GeoMesh>
508 (
509  const IOobject& io,
511 )
512 :
513  Internal
514  (
515  io,
516  const_cast<GeometricField<Type, PatchField, GeoMesh>&>(tgf()),
517  tgf.isTmp()
518  ),
519  timeIndex_(tgf().timeIndex()),
520  field0Ptr_(nullptr),
521  fieldPrevIterPtr_(nullptr),
522  boundaryField_(*this, tgf().boundaryField_)
523 {
524  if (debug)
525  {
527  << "Constructing from tmp resetting IO params"
528  << endl << this->info() << endl;
529  }
530 
531  tgf.clear();
532 
533  readIfPresent();
534 }
535 #endif
536 
537 
538 template<class Type, template<class> class PatchField, class GeoMesh>
540 (
541  const word& newName,
543 )
544 :
545  Internal(newName, gf),
546  timeIndex_(gf.timeIndex()),
547  field0Ptr_(nullptr),
548  fieldPrevIterPtr_(nullptr),
549  boundaryField_(*this, gf.boundaryField_)
550 {
551  if (debug)
552  {
554  << "Constructing as copy resetting name"
555  << endl << this->info() << endl;
556  }
557 
558  if (!readIfPresent() && gf.field0Ptr_)
559  {
561  (
562  newName + "_0",
563  *gf.field0Ptr_
564  );
565  }
566 }
567 
568 
569 #ifndef NoConstructFromTmp
570 template<class Type, template<class> class PatchField, class GeoMesh>
572 (
573  const word& newName,
575 )
576 :
577  Internal
578  (
579  newName,
580  const_cast<GeometricField<Type, PatchField, GeoMesh>&>(tgf()),
581  tgf.isTmp()
582  ),
583  timeIndex_(tgf().timeIndex()),
584  field0Ptr_(nullptr),
585  fieldPrevIterPtr_(nullptr),
586  boundaryField_(*this, tgf().boundaryField_)
587 {
588  if (debug)
589  {
591  << "Constructing from tmp resetting name"
592  << endl << this->info() << endl;
593  }
594 
595  tgf.clear();
596 }
597 #endif
598 
599 
600 template<class Type, template<class> class PatchField, class GeoMesh>
602 (
603  const IOobject& io,
605  const word& patchFieldType
606 )
607 :
608  Internal(io, gf),
609  timeIndex_(gf.timeIndex()),
610  field0Ptr_(nullptr),
611  fieldPrevIterPtr_(nullptr),
612  boundaryField_(this->mesh().boundary(), *this, patchFieldType)
613 {
614  if (debug)
615  {
617  << "Constructing as copy resetting IO params"
618  << endl << this->info() << endl;
619  }
620 
621  boundaryField_ == gf.boundaryField_;
622 
623  if (!readIfPresent() && gf.field0Ptr_)
624  {
626  (
627  io.name() + "_0",
628  *gf.field0Ptr_
629  );
630  }
631 }
632 
633 
634 template<class Type, template<class> class PatchField, class GeoMesh>
636 (
637  const IOobject& io,
639  const wordList& patchFieldTypes,
640  const wordList& actualPatchTypes
641 
642 )
643 :
644  Internal(io, gf),
645  timeIndex_(gf.timeIndex()),
646  field0Ptr_(nullptr),
647  fieldPrevIterPtr_(nullptr),
648  boundaryField_
649  (
650  this->mesh().boundary(),
651  *this,
652  patchFieldTypes,
653  actualPatchTypes
654  )
655 {
656  if (debug)
657  {
659  << "Constructing as copy resetting IO params and patch types"
660  << endl << this->info() << endl;
661  }
662 
663  boundaryField_ == gf.boundaryField_;
664 
665  if (!readIfPresent() && gf.field0Ptr_)
666  {
668  (
669  io.name() + "_0",
670  *gf.field0Ptr_
671  );
672  }
673 }
674 
675 
676 #ifndef NoConstructFromTmp
677 template<class Type, template<class> class PatchField, class GeoMesh>
679 (
680  const IOobject& io,
682  const wordList& patchFieldTypes,
683  const wordList& actualPatchTypes
684 )
685 :
686  Internal
687  (
688  io,
689  const_cast<GeometricField<Type, PatchField, GeoMesh>&>(tgf()),
690  tgf.isTmp()
691  ),
692  timeIndex_(tgf().timeIndex()),
693  field0Ptr_(nullptr),
694  fieldPrevIterPtr_(nullptr),
695  boundaryField_
696  (
697  this->mesh().boundary(),
698  *this,
699  patchFieldTypes,
700  actualPatchTypes
701  )
702 {
703  if (debug)
704  {
706  << "Constructing from tmp resetting IO params and patch types"
707  << endl << this->info() << endl;
708  }
709 
710  boundaryField_ == tgf().boundaryField_;
711 
712  tgf.clear();
713 }
714 #endif
715 
716 
717 template<class Type, template<class> class PatchField, class GeoMesh>
720 {
722  (
724  );
725 }
726 
727 
728 // * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * * //
729 
730 template<class Type, template<class> class PatchField, class GeoMesh>
732 {
733  deleteDemandDrivenData(field0Ptr_);
734  deleteDemandDrivenData(fieldPrevIterPtr_);
735 }
736 
737 
738 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
739 
740 template<class Type, template<class> class PatchField, class GeoMesh>
741 typename
744 {
745  this->setUpToDate();
746  storeOldTimes();
747  return *this;
748 }
749 
750 
751 template<class Type, template<class> class PatchField, class GeoMesh>
752 typename
755 {
756  this->setUpToDate();
757  storeOldTimes();
758  return *this;
759 }
760 
761 
762 template<class Type, template<class> class PatchField, class GeoMesh>
763 typename
766 {
767  this->setUpToDate();
768  storeOldTimes();
769  return boundaryField_;
770 }
771 
772 
773 template<class Type, template<class> class PatchField, class GeoMesh>
775 {
776  if
777  (
778  field0Ptr_
779  && timeIndex_ != this->time().timeIndex()
780  && !(
781  this->name().size() > 2
782  && this->name()(this->name().size()-2, 2) == "_0"
783  )
784  )
785  {
786  storeOldTime();
787  }
788 
789  // Correct time index
790  timeIndex_ = this->time().timeIndex();
791 }
792 
793 
794 template<class Type, template<class> class PatchField, class GeoMesh>
796 {
797  if (field0Ptr_)
798  {
799  field0Ptr_->storeOldTime();
800 
801  if (debug)
802  {
804  << "Storing old time field for field" << endl
805  << this->info() << endl;
806  }
807 
808  *field0Ptr_ == *this;
809  field0Ptr_->timeIndex_ = timeIndex_;
810 
811  if (field0Ptr_->field0Ptr_)
812  {
813  field0Ptr_->writeOpt() = this->writeOpt();
814  }
815  }
816 }
817 
818 
819 template<class Type, template<class> class PatchField, class GeoMesh>
821 {
822  if (field0Ptr_)
823  {
824  return field0Ptr_->nOldTimes() + 1;
825  }
826  else
827  {
828  return 0;
829  }
830 }
831 
832 
833 template<class Type, template<class> class PatchField, class GeoMesh>
836 {
837  if (!field0Ptr_)
838  {
840  (
841  IOobject
842  (
843  this->name() + "_0",
844  this->time().timeName(),
845  this->db(),
846  IOobject::NO_READ,
847  IOobject::NO_WRITE,
848  this->registerObject()
849  ),
850  *this
851  );
852  }
853  else
854  {
855  storeOldTimes();
856  }
857 
858  return *field0Ptr_;
859 }
860 
861 
862 template<class Type, template<class> class PatchField, class GeoMesh>
865 {
866  static_cast<const GeometricField<Type, PatchField, GeoMesh>&>(*this)
867  .oldTime();
868 
869  return *field0Ptr_;
870 }
871 
872 
873 template<class Type, template<class> class PatchField, class GeoMesh>
875 {
876  if (!fieldPrevIterPtr_)
877  {
878  if (debug)
879  {
881  << "Allocating previous iteration field" << endl
882  << this->info() << endl;
883  }
884 
885  fieldPrevIterPtr_ = new GeometricField<Type, PatchField, GeoMesh>
886  (
887  this->name() + "PrevIter",
888  *this
889  );
890  }
891  else
892  {
893  *fieldPrevIterPtr_ == *this;
894  }
895 }
896 
897 
898 template<class Type, template<class> class PatchField, class GeoMesh>
901 {
902  if (!fieldPrevIterPtr_)
903  {
905  << "previous iteration field" << endl << this->info() << endl
906  << " not stored."
907  << " Use field.storePrevIter() at start of iteration."
908  << abort(FatalError);
909  }
910 
911  return *fieldPrevIterPtr_;
912 }
913 
914 
915 template<class Type, template<class> class PatchField, class GeoMesh>
918 {
919  this->setUpToDate();
920  storeOldTimes();
921  boundaryField_.evaluate();
922 }
923 
924 
925 template<class Type, template<class> class PatchField, class GeoMesh>
927 {
928  // Search all boundary conditions, if any are
929  // fixed-value or mixed (Robin) do not set reference level for solution.
930 
931  bool needRef = true;
932 
933  forAll(boundaryField_, patchi)
934  {
935  if (boundaryField_[patchi].fixesValue())
936  {
937  needRef = false;
938  break;
939  }
940  }
941 
942  reduce(needRef, andOp<bool>());
943 
944  return needRef;
945 }
946 
947 
948 template<class Type, template<class> class PatchField, class GeoMesh>
950 {
951  if (debug)
952  {
954  << "Relaxing" << endl << this->info() << " by " << alpha << endl;
955  }
956 
957  operator==(prevIter() + alpha*(*this - prevIter()));
958 }
959 
960 
961 template<class Type, template<class> class PatchField, class GeoMesh>
963 {
964  word name = this->name();
965 
966  if
967  (
968  this->mesh().data::template lookupOrDefault<bool>
969  (
970  "finalIteration",
971  false
972  )
973  )
974  {
975  name += "Final";
976  }
977 
978  if (this->mesh().relaxField(name))
979  {
980  relax(this->mesh().fieldRelaxationFactor(name));
981  }
982 }
983 
984 
985 template<class Type, template<class> class PatchField, class GeoMesh>
987 (
988  bool final
989 ) const
990 {
991  if (final)
992  {
993  return this->name() + "Final";
994  }
995  else
996  {
997  return this->name();
998  }
999 }
1000 
1001 
1002 template<class Type, template<class> class PatchField, class GeoMesh>
1005  Ostream& os
1006 ) const
1007 {
1008  os << "min/max(" << this->name() << ") = "
1009  << Foam::min(*this).value() << ", "
1010  << Foam::max(*this).value()
1011  << endl;
1012 }
1013 
1014 
1015 template<class Type, template<class> class PatchField, class GeoMesh>
1017 writeData(Ostream& os) const
1018 {
1019  os << *this;
1020  return os.good();
1021 }
1022 
1023 
1024 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
1025 
1026 template<class Type, template<class> class PatchField, class GeoMesh>
1029 {
1031  (
1033  (
1034  IOobject
1035  (
1036  this->name() + ".T()",
1037  this->instance(),
1038  this->db()
1039  ),
1040  this->mesh(),
1041  this->dimensions()
1042  )
1043  );
1044 
1045  Foam::T(result.ref().primitiveFieldRef(), primitiveField());
1046  Foam::T(result.ref().boundaryFieldRef(), boundaryField());
1047 
1048  return result;
1049 }
1050 
1051 
1052 template<class Type, template<class> class PatchField, class GeoMesh>
1053 Foam::tmp
1054 <
1056  <
1057  typename Foam::GeometricField<Type, PatchField, GeoMesh>::cmptType,
1058  PatchField,
1059  GeoMesh
1060  >
1061 >
1063 (
1064  const direction d
1065 ) const
1066 {
1068  (
1070  (
1071  IOobject
1072  (
1073  this->name() + ".component(" + Foam::name(d) + ')',
1074  this->instance(),
1075  this->db()
1076  ),
1077  this->mesh(),
1078  this->dimensions()
1079  )
1080  );
1081 
1082  Foam::component(Component.ref().primitiveFieldRef(), primitiveField(), d);
1083  Foam::component(Component.ref().boundaryFieldRef(), boundaryField(), d);
1084 
1085  return Component;
1086 }
1087 
1088 
1089 template<class Type, template<class> class PatchField, class GeoMesh>
1091 (
1092  const direction d,
1093  const GeometricField
1094  <
1096  PatchField,
1097  GeoMesh
1098  >& gcf
1099 )
1100 {
1101  primitiveFieldRef().replace(d, gcf.primitiveField());
1102  boundaryFieldRef().replace(d, gcf.boundaryField());
1103 }
1104 
1105 
1106 template<class Type, template<class> class PatchField, class GeoMesh>
1109  const direction d,
1110  const dimensioned<cmptType>& ds
1111 )
1112 {
1113  primitiveFieldRef().replace(d, ds.value());
1114  boundaryFieldRef().replace(d, ds.value());
1115 }
1116 
1117 
1118 template<class Type, template<class> class PatchField, class GeoMesh>
1121  const dimensioned<Type>& dt
1122 )
1123 {
1124  Foam::max(primitiveFieldRef(), primitiveField(), dt.value());
1125  Foam::max(boundaryFieldRef(), boundaryField(), dt.value());
1126 }
1127 
1128 
1129 template<class Type, template<class> class PatchField, class GeoMesh>
1132  const dimensioned<Type>& dt
1133 )
1134 {
1135  Foam::min(primitiveFieldRef(), primitiveField(), dt.value());
1136  Foam::min(boundaryFieldRef(), boundaryField(), dt.value());
1137 }
1138 
1139 
1140 template<class Type, template<class> class PatchField, class GeoMesh>
1143  const dimensioned<Type>& minDt,
1144  const dimensioned<Type>& maxDt
1145 )
1146 {
1147  Foam::max(primitiveFieldRef(), primitiveField(), minDt.value());
1148  Foam::max(boundaryFieldRef(), boundaryField(), minDt.value());
1149  Foam::min(primitiveFieldRef(), primitiveField(), maxDt.value());
1150  Foam::min(boundaryFieldRef(), boundaryField(), maxDt.value());
1151 }
1152 
1153 
1154 template<class Type, template<class> class PatchField, class GeoMesh>
1156 {
1157  primitiveFieldRef().negate();
1158  boundaryFieldRef().negate();
1159 }
1160 
1161 
1162 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
1163 
1164 template<class Type, template<class> class PatchField, class GeoMesh>
1165 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator=
1168 )
1169 {
1170  if (this == &gf)
1171  {
1173  << "attempted assignment to self"
1174  << abort(FatalError);
1175  }
1176 
1177  checkField(*this, gf, "=");
1178 
1179  // Only assign field contents not ID
1180 
1181  ref() = gf();
1182  boundaryFieldRef() = gf.boundaryField();
1183 }
1184 
1185 
1186 template<class Type, template<class> class PatchField, class GeoMesh>
1187 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator=
1190 )
1191 {
1192  if (this == &(tgf()))
1193  {
1195  << "attempted assignment to self"
1196  << abort(FatalError);
1197  }
1198 
1200 
1201  checkField(*this, gf, "=");
1202 
1203  // Only assign field contents not ID
1204 
1205  this->dimensions() = gf.dimensions();
1206 
1207  if (tgf.isTmp())
1208  {
1209  // Transfer the storage from the tmp
1210  primitiveFieldRef().transfer
1211  (
1212  const_cast<Field<Type>&>(gf.primitiveField())
1213  );
1214  }
1215  else
1216  {
1218  }
1219 
1220  boundaryFieldRef() = gf.boundaryField();
1221 
1222  tgf.clear();
1223 }
1224 
1225 
1226 template<class Type, template<class> class PatchField, class GeoMesh>
1227 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator=
1229  const dimensioned<Type>& dt
1230 )
1231 {
1232  ref() = dt;
1233  boundaryFieldRef() = dt.value();
1234 }
1235 
1236 
1237 template<class Type, template<class> class PatchField, class GeoMesh>
1238 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator==
1241 )
1242 {
1244 
1245  checkField(*this, gf, "==");
1246 
1247  // Only assign field contents not ID
1248 
1249  ref() = gf();
1250  boundaryFieldRef() == gf.boundaryField();
1251 
1252  tgf.clear();
1253 }
1254 
1255 
1256 template<class Type, template<class> class PatchField, class GeoMesh>
1257 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator==
1259  const dimensioned<Type>& dt
1260 )
1261 {
1262  ref() = dt;
1263  boundaryFieldRef() == dt.value();
1264 }
1265 
1266 
1267 #define COMPUTED_ASSIGNMENT(TYPE, op) \
1268  \
1269 template<class Type, template<class> class PatchField, class GeoMesh> \
1270 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \
1271 ( \
1272  const GeometricField<TYPE, PatchField, GeoMesh>& gf \
1273 ) \
1274 { \
1275  checkField(*this, gf, #op); \
1276  \
1277  ref() op gf(); \
1278  boundaryFieldRef() op gf.boundaryField(); \
1279 } \
1280  \
1281 template<class Type, template<class> class PatchField, class GeoMesh> \
1282 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \
1283 ( \
1284  const tmp<GeometricField<TYPE, PatchField, GeoMesh>>& tgf \
1285 ) \
1286 { \
1287  operator op(tgf()); \
1288  tgf.clear(); \
1289 } \
1290  \
1291 template<class Type, template<class> class PatchField, class GeoMesh> \
1292 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \
1293 ( \
1294  const dimensioned<TYPE>& dt \
1295 ) \
1296 { \
1297  ref() op dt; \
1298  boundaryFieldRef() op dt.value(); \
1299 }
1300 
1305 
1306 #undef COMPUTED_ASSIGNMENT
1307 
1308 
1309 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
1310 
1311 template<class Type, template<class> class PatchField, class GeoMesh>
1312 Foam::Ostream& Foam::operator<<
1313 (
1314  Ostream& os,
1316 )
1317 {
1318  gf().writeData(os, "internalField");
1319  os << nl;
1320  gf.boundaryField().writeEntry("boundaryField", os);
1321 
1322  // Check state of IOstream
1323  os.check
1324  (
1325  "Ostream& operator<<(Ostream&, "
1326  "const GeometricField<Type, PatchField, GeoMesh>&)"
1327  );
1328 
1329  return (os);
1330 }
1331 
1332 
1333 template<class Type, template<class> class PatchField, class GeoMesh>
1334 Foam::Ostream& Foam::operator<<
1335 (
1336  Ostream& os,
1338 )
1339 {
1340  os << tgf();
1341  tgf.clear();
1342  return os;
1343 }
1344 
1345 
1346 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1347 
1348 #undef checkField
1349 
1350 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1351 
1352 #include "GeometricBoundaryField.C"
1353 #include "GeometricFieldFunctions.C"
1354 
1355 // ************************************************************************* //
void replace(const direction, const GeometricField< cmptType, PatchField, GeoMesh > &)
#define COMPUTED_ASSIGNMENT(TYPE, op)
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
#define checkField(gf1, gf2, op)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const word & name() const
Return name.
Definition: IOobject.H:297
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
UEqn relax()
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject *> &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const Boundary & boundaryField() const
Return const-reference to the boundary field.
bool writeData(Ostream &) const
WriteData member function required by regIOobject.
Field< Type >::cmptType cmptType
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
uint8_t direction
Definition: direction.H:45
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
Generic dimensioned Type class.
label nOldTimes() const
Return the number of old time fields stored.
void maxMin(const dimensioned< Type > &minDt, const dimensioned< Type > &maxDt)
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:333
word select(bool final) const
Select the final iteration parameters if `final&#39; is true.
conserve primitiveFieldRef()+
bool needReference() const
Does the field need a reference level for solution.
const dimensionSet & dimensions() const
Return dimensions.
Dimension set for the base types.
Definition: dimensionSet.H:120
dynamicFvMesh & mesh
void storeOldTimes() const
Store the old-time fields.
Pre-declare SubField and related Field type.
Definition: Field.H:57
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
void storeOldTime() const
Store the old-time field.
const Type & value() const
Return const reference to value.
void writeMinMax(Ostream &os) const
Helper function to write the min and max to an Ostream.
word timeName
Definition: getTimeIndex.H:3
void min(const dimensioned< Type > &)
faceListList boundary(nPatches)
GeoMesh::Mesh Mesh
Type of mesh on which this DimensionedField is instantiated.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
GeometricField(const IOobject &, const Mesh &, const dimensionSet &, const word &patchFieldType=PatchField< Type >::calculatedType())
Constructor given IOobject, mesh, dimensions and patch type.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
static const char nl
Definition: Ostream.H:265
const GeometricField< Type, PatchField, GeoMesh > & prevIter() const
Return previous iteration field.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
virtual ~GeometricField()
Destructor.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label timeIndex() const
Return the time index of the field.
Internal & ref()
Return a reference to the dimensioned internal field.
Template functions to aid in the implementation of demand driven data.
tmp< GeometricField< Type, PatchField, GeoMesh > > clone() const
Clone.
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
#define WarningInFunction
Report a warning using Foam::Warning.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
void max(const dimensioned< Type > &)
label timeIndex
Definition: getTimeIndex.H:4
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void correctBoundaryConditions()
Correct boundary field.
void storePrevIter() const
Store the field as the previous iteration value.
rDeltaT ref()
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
void deleteDemandDrivenData(DataPtr &dataPtr)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
void relax()
Relax field (for steady-state solution).
IOerror FatalIOError
#define InfoInFunction
Report an information message using Foam::Info.