GeometricField.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-2015 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 "data.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 #define checkField(gf1, gf2, op) \
35 if ((gf1).mesh() != (gf2).mesh()) \
36 { \
37  FatalErrorIn("checkField(gf1, gf2, op)") \
38  << "different mesh for fields " \
39  << (gf1).name() << " and " << (gf2).name() \
40  << " during operatrion " << op \
41  << abort(FatalError); \
42 }
43 
44 
45 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
46 
47 template<class Type, template<class> class PatchField, class GeoMesh>
49 (
50  const dictionary& dict
51 )
52 {
53  DimensionedField<Type, GeoMesh>::readField(dict, "internalField");
54 
55  boundaryField_.readField(*this, dict.subDict("boundaryField"));
56 
57  if (dict.found("referenceLevel"))
58  {
59  Type fieldAverage(pTraits<Type>(dict.lookup("referenceLevel")));
60 
61  Field<Type>::operator+=(fieldAverage);
62 
63  forAll(boundaryField_, patchi)
64  {
65  boundaryField_[patchi] == boundaryField_[patchi] + fieldAverage;
66  }
67  }
68 }
69 
70 
71 template<class Type, template<class> class PatchField, class GeoMesh>
73 {
74  const IOdictionary dict
75  (
76  IOobject
77  (
78  this->name(),
79  this->time().timeName(),
80  this->db(),
81  IOobject::NO_READ,
82  IOobject::NO_WRITE,
83  false
84  ),
85  this->readStream(typeName)
86  );
87 
88  this->close();
89 
90  readFields(dict);
91 }
92 
93 
94 template<class Type, template<class> class PatchField, class GeoMesh>
96 {
97  if
98  (
99  this->readOpt() == IOobject::MUST_READ
100  || this->readOpt() == IOobject::MUST_READ_IF_MODIFIED
101  )
102  {
103  WarningIn
104  (
105  "GeometricField<Type, PatchField, GeoMesh>::readIfPresent()"
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 (this->readOpt() == IOobject::READ_IF_PRESENT && this->headerOk())
111  {
112  readFields();
113 
114  // Check compatibility between field and mesh
115  if (this->size() != GeoMesh::size(this->mesh()))
116  {
118  (
119  "GeometricField<Type, PatchField, GeoMesh>::"
120  "readIfPresent()",
121  this->readStream(typeName)
122  ) << " number of field elements = " << this->size()
123  << " number of mesh elements = "
124  << GeoMesh::size(this->mesh())
125  << exit(FatalIOError);
126  }
127 
128  readOldTimeIfPresent();
129 
130  return true;
131  }
132 
133  return false;
134 }
135 
136 
137 template<class Type, template<class> class PatchField, class GeoMesh>
139 {
140  // Read the old time field if present
141  IOobject field0
142  (
143  this->name() + "_0",
144  this->time().timeName(),
145  this->db(),
146  IOobject::READ_IF_PRESENT,
147  IOobject::AUTO_WRITE,
148  this->registerObject()
149  );
150 
151  if (field0.headerOk())
152  {
153  if (debug)
154  {
155  Info<< "Reading old time level for field"
156  << endl << this->info() << endl;
157  }
158 
159  field0Ptr_ = new GeometricField<Type, PatchField, GeoMesh>
160  (
161  field0,
162  this->mesh()
163  );
164 
165  field0Ptr_->timeIndex_ = timeIndex_ - 1;
166 
167  if (!field0Ptr_->readOldTimeIfPresent())
168  {
169  field0Ptr_->oldTime();
170  }
171 
172  return true;
173  }
174 
175  return false;
176 }
177 
178 
179 // * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * //
180 
181 template<class Type, template<class> class PatchField, class GeoMesh>
183 (
184  const IOobject& io,
185  const Mesh& mesh,
186  const dimensionSet& ds,
187  const word& patchFieldType
188 )
189 :
190  DimensionedField<Type, GeoMesh>(io, mesh, ds, false),
191  timeIndex_(this->time().timeIndex()),
192  field0Ptr_(NULL),
193  fieldPrevIterPtr_(NULL),
194  boundaryField_(mesh.boundary(), *this, patchFieldType)
195 {
196  if (debug)
197  {
198  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
199  "creating temporary"
200  << endl << this->info() << endl;
201  }
202 
203  readIfPresent();
204 }
205 
206 
207 template<class Type, template<class> class PatchField, class GeoMesh>
209 (
210  const IOobject& io,
211  const Mesh& mesh,
212  const dimensionSet& ds,
213  const wordList& patchFieldTypes,
214  const wordList& actualPatchTypes
215 )
216 :
217  DimensionedField<Type, GeoMesh>(io, mesh, ds, false),
218  timeIndex_(this->time().timeIndex()),
219  field0Ptr_(NULL),
220  fieldPrevIterPtr_(NULL),
221  boundaryField_(mesh.boundary(), *this, patchFieldTypes, actualPatchTypes)
222 {
223  if (debug)
224  {
225  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
226  "creating temporary"
227  << endl << this->info() << endl;
228  }
229 
230  readIfPresent();
231 }
232 
233 
234 template<class Type, template<class> class PatchField, class GeoMesh>
236 (
237  const IOobject& io,
238  const Mesh& mesh,
239  const dimensioned<Type>& dt,
240  const word& patchFieldType
241 )
242 :
243  DimensionedField<Type, GeoMesh>(io, mesh, dt, false),
244  timeIndex_(this->time().timeIndex()),
245  field0Ptr_(NULL),
246  fieldPrevIterPtr_(NULL),
247  boundaryField_(mesh.boundary(), *this, patchFieldType)
248 {
249  if (debug)
250  {
251  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
252  "creating temporary"
253  << endl << this->info() << endl;
254  }
255 
256  boundaryField_ == dt.value();
257 
258  readIfPresent();
259 }
260 
261 
262 template<class Type, template<class> class PatchField, class GeoMesh>
264 (
265  const IOobject& io,
266  const Mesh& mesh,
267  const dimensioned<Type>& dt,
268  const wordList& patchFieldTypes,
269  const wordList& actualPatchTypes
270 )
271 :
272  DimensionedField<Type, GeoMesh>(io, mesh, dt, false),
273  timeIndex_(this->time().timeIndex()),
274  field0Ptr_(NULL),
275  fieldPrevIterPtr_(NULL),
276  boundaryField_(mesh.boundary(), *this, patchFieldTypes, actualPatchTypes)
277 {
278  if (debug)
279  {
280  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
281  "creating temporary"
282  << endl << this->info() << endl;
283  }
284 
285  boundaryField_ == dt.value();
286 
287  readIfPresent();
288 }
289 
290 
291 template<class Type, template<class> class PatchField, class GeoMesh>
293 (
294  const IOobject& io,
295  const Mesh& mesh,
296  const dimensionSet& ds,
297  const Field<Type>& iField,
298  const PtrList<PatchField<Type> >& ptfl
299 )
300 :
301  DimensionedField<Type, GeoMesh>(io, mesh, ds, iField),
302  timeIndex_(this->time().timeIndex()),
303  field0Ptr_(NULL),
304  fieldPrevIterPtr_(NULL),
305  boundaryField_(mesh.boundary(), *this, ptfl)
306 {
307  if (debug)
308  {
309  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
310  "constructing from components"
311  << endl << this->info() << endl;
312  }
313 
314  readIfPresent();
315 }
316 
317 
318 template<class Type, template<class> class PatchField, class GeoMesh>
320 (
321  const IOobject& io,
322  const Mesh& mesh
323 )
324 :
326  timeIndex_(this->time().timeIndex()),
327  field0Ptr_(NULL),
328  fieldPrevIterPtr_(NULL),
329  boundaryField_(mesh.boundary())
330 {
331  readFields();
332 
333  // Check compatibility between field and mesh
334 
335  if (this->size() != GeoMesh::size(this->mesh()))
336  {
338  (
339  "GeometricField<Type, PatchField, GeoMesh>::GeometricField"
340  "(const IOobject&, const Mesh&)",
341  this->readStream(typeName)
342  ) << " number of field elements = " << this->size()
343  << " number of mesh elements = " << GeoMesh::size(this->mesh())
344  << exit(FatalIOError);
345  }
346 
347  readOldTimeIfPresent();
348 
349  if (debug)
350  {
351  Info<< "Finishing read-construct of "
352  "GeometricField<Type, PatchField, GeoMesh>"
353  << endl << this->info() << endl;
354  }
355 }
356 
357 
358 template<class Type, template<class> class PatchField, class GeoMesh>
360 (
361  const IOobject& io,
362  const Mesh& mesh,
363  const dictionary& dict
364 )
365 :
367  timeIndex_(this->time().timeIndex()),
368  field0Ptr_(NULL),
369  fieldPrevIterPtr_(NULL),
370  boundaryField_(mesh.boundary())
371 {
372  readFields(dict);
373 
374  // Check compatibility between field and mesh
375 
376  if (this->size() != GeoMesh::size(this->mesh()))
377  {
379  (
380  "GeometricField<Type, PatchField, GeoMesh>::GeometricField"
381  "(const IOobject&, const Mesh&, const dictionary&)"
382  ) << " number of field elements = " << this->size()
383  << " number of mesh elements = " << GeoMesh::size(this->mesh())
384  << exit(FatalIOError);
385  }
386 
387  readOldTimeIfPresent();
388 
389  if (debug)
390  {
391  Info<< "Finishing dictionary-construct of "
392  "GeometricField<Type, PatchField, GeoMesh>"
393  << endl << this->info() << endl;
394  }
395 }
396 
397 
398 template<class Type, template<class> class PatchField, class GeoMesh>
400 (
402 )
403 :
405  timeIndex_(gf.timeIndex()),
406  field0Ptr_(NULL),
407  fieldPrevIterPtr_(NULL),
408  boundaryField_(*this, gf.boundaryField_)
409 {
410  if (debug)
411  {
412  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
413  "constructing as copy"
414  << endl << this->info() << endl;
415  }
416 
417  if (gf.field0Ptr_)
418  {
420  (
421  *gf.field0Ptr_
422  );
423  }
424 
425  this->writeOpt() = IOobject::NO_WRITE;
426 }
427 
428 
429 #ifndef NoConstructFromTmp
430 template<class Type, template<class> class PatchField, class GeoMesh>
432 (
434 )
435 :
437  (
438  const_cast<GeometricField<Type, PatchField, GeoMesh>&>(tgf()),
439  tgf.isTmp()
440  ),
441  timeIndex_(tgf().timeIndex()),
442  field0Ptr_(NULL),
443  fieldPrevIterPtr_(NULL),
444  boundaryField_(*this, tgf().boundaryField_)
445 {
446  if (debug)
447  {
448  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
449  "constructing as copy"
450  << endl << this->info() << endl;
451  }
452 
453  this->writeOpt() = IOobject::NO_WRITE;
454 
455  tgf.clear();
456 }
457 #endif
458 
459 
460 template<class Type, template<class> class PatchField, class GeoMesh>
462 (
463  const IOobject& io,
465 )
466 :
468  timeIndex_(gf.timeIndex()),
469  field0Ptr_(NULL),
470  fieldPrevIterPtr_(NULL),
471  boundaryField_(*this, gf.boundaryField_)
472 {
473  if (debug)
474  {
475  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
476  "constructing as copy resetting IO params"
477  << endl << this->info() << endl;
478  }
479 
480  if (!readIfPresent() && gf.field0Ptr_)
481  {
483  (
484  io.name() + "_0",
485  *gf.field0Ptr_
486  );
487  }
488 }
489 
490 
491 #ifndef NoConstructFromTmp
492 template<class Type, template<class> class PatchField, class GeoMesh>
494 (
495  const IOobject& io,
497 )
498 :
500  (
501  io,
502  const_cast<GeometricField<Type, PatchField, GeoMesh>&>(tgf()),
503  tgf.isTmp()
504  ),
505  timeIndex_(tgf().timeIndex()),
506  field0Ptr_(NULL),
507  fieldPrevIterPtr_(NULL),
508  boundaryField_(*this, tgf().boundaryField_)
509 {
510  if (debug)
511  {
512  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
513  "constructing from tmp resetting IO params"
514  << endl << this->info() << endl;
515  }
516 
517  tgf.clear();
518 
519  readIfPresent();
520 }
521 #endif
522 
523 
524 template<class Type, template<class> class PatchField, class GeoMesh>
526 (
527  const word& newName,
529 )
530 :
531  DimensionedField<Type, GeoMesh>(newName, gf),
532  timeIndex_(gf.timeIndex()),
533  field0Ptr_(NULL),
534  fieldPrevIterPtr_(NULL),
535  boundaryField_(*this, gf.boundaryField_)
536 {
537  if (debug)
538  {
539  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
540  "constructing as copy resetting name"
541  << endl << this->info() << endl;
542  }
543 
544  if (!readIfPresent() && gf.field0Ptr_)
545  {
547  (
548  newName + "_0",
549  *gf.field0Ptr_
550  );
551  }
552 }
553 
554 
555 #ifndef NoConstructFromTmp
556 template<class Type, template<class> class PatchField, class GeoMesh>
558 (
559  const word& newName,
561 )
562 :
564  (
565  newName,
566  const_cast<GeometricField<Type, PatchField, GeoMesh>&>(tgf()),
567  tgf.isTmp()
568  ),
569  timeIndex_(tgf().timeIndex()),
570  field0Ptr_(NULL),
571  fieldPrevIterPtr_(NULL),
572  boundaryField_(*this, tgf().boundaryField_)
573 {
574  if (debug)
575  {
576  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
577  "constructing from tmp resetting name"
578  << endl << this->info() << endl;
579  }
580 
581  tgf.clear();
582 }
583 #endif
584 
585 
586 template<class Type, template<class> class PatchField, class GeoMesh>
588 (
589  const IOobject& io,
591  const word& patchFieldType
592 )
593 :
595  timeIndex_(gf.timeIndex()),
596  field0Ptr_(NULL),
597  fieldPrevIterPtr_(NULL),
598  boundaryField_(this->mesh().boundary(), *this, patchFieldType)
599 {
600  if (debug)
601  {
602  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
603  "constructing as copy resetting IO params"
604  << endl << this->info() << endl;
605  }
606 
607  boundaryField_ == gf.boundaryField_;
608 
609  if (!readIfPresent() && gf.field0Ptr_)
610  {
612  (
613  io.name() + "_0",
614  *gf.field0Ptr_
615  );
616  }
617 }
618 
619 
620 template<class Type, template<class> class PatchField, class GeoMesh>
622 (
623  const IOobject& io,
625  const wordList& patchFieldTypes,
626  const wordList& actualPatchTypes
627 
628 )
629 :
631  timeIndex_(gf.timeIndex()),
632  field0Ptr_(NULL),
633  fieldPrevIterPtr_(NULL),
634  boundaryField_
635  (
636  this->mesh().boundary(),
637  *this,
638  patchFieldTypes,
639  actualPatchTypes
640  )
641 {
642  if (debug)
643  {
644  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
645  "constructing as copy resetting IO params and patch types"
646  << endl << this->info() << endl;
647  }
648 
649  boundaryField_ == gf.boundaryField_;
650 
651  if (!readIfPresent() && gf.field0Ptr_)
652  {
654  (
655  io.name() + "_0",
656  *gf.field0Ptr_
657  );
658  }
659 }
660 
661 
662 // * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * * //
663 
664 template<class Type, template<class> class PatchField, class GeoMesh>
666 {
667  deleteDemandDrivenData(field0Ptr_);
668  deleteDemandDrivenData(fieldPrevIterPtr_);
669 }
670 
671 
672 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
673 
674 template<class Type, template<class> class PatchField, class GeoMesh>
675 typename
678 {
679  this->setUpToDate();
680  storeOldTimes();
681  return *this;
682 }
683 
684 
685 template<class Type, template<class> class PatchField, class GeoMesh>
686 typename
689 {
690  this->setUpToDate();
691  storeOldTimes();
692  return *this;
693 }
694 
695 
696 template<class Type, template<class> class PatchField, class GeoMesh>
697 typename
700 {
701  this->setUpToDate();
702  storeOldTimes();
703  return boundaryField_;
704 }
705 
706 
707 template<class Type, template<class> class PatchField, class GeoMesh>
709 {
710  if
711  (
712  field0Ptr_
713  && timeIndex_ != this->time().timeIndex()
714  && !(
715  this->name().size() > 2
716  && this->name()(this->name().size()-2, 2) == "_0"
717  )
718  )
719  {
720  storeOldTime();
721  }
722 
723  // Correct time index
724  timeIndex_ = this->time().timeIndex();
725 }
726 
727 
728 template<class Type, template<class> class PatchField, class GeoMesh>
730 {
731  if (field0Ptr_)
732  {
733  field0Ptr_->storeOldTime();
734 
735  if (debug)
736  {
737  Info<< "Storing old time field for field" << endl
738  << this->info() << endl;
739  }
740 
741  *field0Ptr_ == *this;
742  field0Ptr_->timeIndex_ = timeIndex_;
743 
744  if (field0Ptr_->field0Ptr_)
745  {
746  field0Ptr_->writeOpt() = this->writeOpt();
747  }
748  }
749 }
750 
751 
752 template<class Type, template<class> class PatchField, class GeoMesh>
754 {
755  if (field0Ptr_)
756  {
757  return field0Ptr_->nOldTimes() + 1;
758  }
759  else
760  {
761  return 0;
762  }
763 }
764 
765 
766 template<class Type, template<class> class PatchField, class GeoMesh>
769 {
770  if (!field0Ptr_)
771  {
773  (
774  IOobject
775  (
776  this->name() + "_0",
777  this->time().timeName(),
778  this->db(),
779  IOobject::NO_READ,
780  IOobject::NO_WRITE,
781  this->registerObject()
782  ),
783  *this
784  );
785  }
786  else
787  {
788  storeOldTimes();
789  }
790 
791  return *field0Ptr_;
792 }
793 
794 
795 template<class Type, template<class> class PatchField, class GeoMesh>
798 {
799  static_cast<const GeometricField<Type, PatchField, GeoMesh>&>(*this)
800  .oldTime();
801 
802  return *field0Ptr_;
803 }
804 
805 
806 template<class Type, template<class> class PatchField, class GeoMesh>
808 {
809  if (!fieldPrevIterPtr_)
810  {
811  if (debug)
812  {
813  Info<< "Allocating previous iteration field" << endl
814  << this->info() << endl;
815  }
816 
817  fieldPrevIterPtr_ = new GeometricField<Type, PatchField, GeoMesh>
818  (
819  this->name() + "PrevIter",
820  *this
821  );
822  }
823  else
824  {
825  *fieldPrevIterPtr_ == *this;
826  }
827 }
828 
829 
830 template<class Type, template<class> class PatchField, class GeoMesh>
833 {
834  if (!fieldPrevIterPtr_)
835  {
837  (
838  "GeometricField<Type, PatchField, GeoMesh>::prevIter() const"
839  ) << "previous iteration field" << endl << this->info() << endl
840  << " not stored."
841  << " Use field.storePrevIter() at start of iteration."
842  << abort(FatalError);
843  }
844 
845  return *fieldPrevIterPtr_;
846 }
847 
848 
849 template<class Type, template<class> class PatchField, class GeoMesh>
852 {
853  this->setUpToDate();
854  storeOldTimes();
855  boundaryField_.evaluate();
856 }
857 
858 
859 template<class Type, template<class> class PatchField, class GeoMesh>
861 {
862  // Search all boundary conditions, if any are
863  // fixed-value or mixed (Robin) do not set reference level for solution.
864 
865  bool needRef = true;
866 
867  forAll(boundaryField_, patchi)
868  {
869  if (boundaryField_[patchi].fixesValue())
870  {
871  needRef = false;
872  break;
873  }
874  }
875 
876  reduce(needRef, andOp<bool>());
877 
878  return needRef;
879 }
880 
881 
882 template<class Type, template<class> class PatchField, class GeoMesh>
884 {
885  if (debug)
886  {
887  InfoIn
888  (
889  "GeometricField<Type, PatchField, GeoMesh>::relax"
890  "(const scalar alpha)"
891  ) << "Relaxing" << endl << this->info() << " by " << alpha << endl;
892  }
893 
894  operator==(prevIter() + alpha*(*this - prevIter()));
895 }
896 
897 
898 template<class Type, template<class> class PatchField, class GeoMesh>
900 {
901  word name = this->name();
902 
903  if
904  (
905  this->mesh().data::template lookupOrDefault<bool>
906  (
907  "finalIteration",
908  false
909  )
910  )
911  {
912  name += "Final";
913  }
914 
915  if (this->mesh().relaxField(name))
916  {
917  relax(this->mesh().fieldRelaxationFactor(name));
918  }
919 }
920 
921 
922 template<class Type, template<class> class PatchField, class GeoMesh>
924 (
925  bool final
926 ) const
927 {
928  if (final)
929  {
930  return this->name() + "Final";
931  }
932  else
933  {
934  return this->name();
935  }
936 }
937 
938 
939 template<class Type, template<class> class PatchField, class GeoMesh>
941 (
942  Ostream& os
943 ) const
944 {
945  os << "min/max(" << this->name() << ") = "
946  << Foam::min(*this).value() << ", "
947  << Foam::max(*this).value()
948  << endl;
949 }
950 
951 
952 template<class Type, template<class> class PatchField, class GeoMesh>
954 writeData(Ostream& os) const
955 {
956  os << *this;
957  return os.good();
958 }
959 
960 
961 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
962 
963 template<class Type, template<class> class PatchField, class GeoMesh>
966 {
968  (
970  (
971  IOobject
972  (
973  this->name() + ".T()",
974  this->instance(),
975  this->db()
976  ),
977  this->mesh(),
978  this->dimensions()
979  )
980  );
981 
982  Foam::T(result().internalField(), internalField());
983  Foam::T(result().boundaryField(), boundaryField());
984 
985  return result;
986 }
987 
988 
989 template<class Type, template<class> class PatchField, class GeoMesh>
990 Foam::tmp
991 <
993  <
994  typename Foam::GeometricField<Type, PatchField, GeoMesh>::cmptType,
995  PatchField,
996  GeoMesh
997  >
998 >
1000 (
1001  const direction d
1002 ) const
1003 {
1005  (
1007  (
1008  IOobject
1009  (
1010  this->name() + ".component(" + Foam::name(d) + ')',
1011  this->instance(),
1012  this->db()
1013  ),
1014  this->mesh(),
1015  this->dimensions()
1016  )
1017  );
1018 
1019  Foam::component(Component().internalField(), internalField(), d);
1020  Foam::component(Component().boundaryField(), boundaryField(), d);
1021 
1022  return Component;
1023 }
1024 
1025 
1026 template<class Type, template<class> class PatchField, class GeoMesh>
1028 (
1029  const direction d,
1030  const GeometricField
1031  <
1033  PatchField,
1034  GeoMesh
1035  >& gcf
1036 )
1037 {
1038  internalField().replace(d, gcf.internalField());
1039  boundaryField().replace(d, gcf.boundaryField());
1040 }
1041 
1042 
1043 template<class Type, template<class> class PatchField, class GeoMesh>
1046  const direction d,
1047  const dimensioned<cmptType>& ds
1048 )
1049 {
1050  internalField().replace(d, ds.value());
1051  boundaryField().replace(d, ds.value());
1052 }
1053 
1054 
1055 template<class Type, template<class> class PatchField, class GeoMesh>
1058  const dimensioned<Type>& dt
1059 )
1060 {
1062  Foam::max(boundaryField(), boundaryField(), dt.value());
1063 }
1064 
1065 
1066 template<class Type, template<class> class PatchField, class GeoMesh>
1069  const dimensioned<Type>& dt
1070 )
1071 {
1073  Foam::min(boundaryField(), boundaryField(), dt.value());
1074 }
1075 
1076 
1077 template<class Type, template<class> class PatchField, class GeoMesh>
1079 {
1080  internalField().negate();
1081  boundaryField().negate();
1082 }
1083 
1084 
1085 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
1086 
1087 template<class Type, template<class> class PatchField, class GeoMesh>
1088 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator=
1091 )
1092 {
1093  if (this == &gf)
1094  {
1095  FatalErrorIn
1096  (
1097  "GeometricField<Type, PatchField, GeoMesh>::operator="
1098  "(const GeometricField<Type, PatchField, GeoMesh>&)"
1099  ) << "attempted assignment to self"
1100  << abort(FatalError);
1101  }
1102 
1103  checkField(*this, gf, "=");
1104 
1105  // only equate field contents not ID
1106 
1108  boundaryField() = gf.boundaryField();
1109 }
1110 
1111 
1112 template<class Type, template<class> class PatchField, class GeoMesh>
1113 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator=
1116 )
1117 {
1118  if (this == &(tgf()))
1119  {
1120  FatalErrorIn
1121  (
1122  "GeometricField<Type, PatchField, GeoMesh>::operator="
1123  "(const tmp<GeometricField<Type, PatchField, GeoMesh> >&)"
1124  ) << "attempted assignment to self"
1125  << abort(FatalError);
1126  }
1127 
1129 
1130  checkField(*this, gf, "=");
1131 
1132  // only equate field contents not ID
1133 
1134  this->dimensions() = gf.dimensions();
1135 
1136  // This is dodgy stuff, don't try it at home.
1137  internalField().transfer
1138  (
1139  const_cast<Field<Type>&>(gf.internalField())
1140  );
1141 
1142  boundaryField() = gf.boundaryField();
1143 
1144  tgf.clear();
1145 }
1146 
1147 
1148 template<class Type, template<class> class PatchField, class GeoMesh>
1149 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator=
1151  const dimensioned<Type>& dt
1152 )
1153 {
1154  dimensionedInternalField() = dt;
1155  boundaryField() = dt.value();
1156 }
1157 
1158 
1159 template<class Type, template<class> class PatchField, class GeoMesh>
1160 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator==
1163 )
1164 {
1166 
1167  checkField(*this, gf, "==");
1168 
1169  // only equate field contents not ID
1170 
1172  boundaryField() == gf.boundaryField();
1173 
1174  tgf.clear();
1175 }
1176 
1177 
1178 template<class Type, template<class> class PatchField, class GeoMesh>
1179 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator==
1181  const dimensioned<Type>& dt
1182 )
1183 {
1184  dimensionedInternalField() = dt;
1185  boundaryField() == dt.value();
1186 }
1187 
1188 
1189 #define COMPUTED_ASSIGNMENT(TYPE, op) \
1190  \
1191 template<class Type, template<class> class PatchField, class GeoMesh> \
1192 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \
1193 ( \
1194  const GeometricField<TYPE, PatchField, GeoMesh>& gf \
1195 ) \
1196 { \
1197  checkField(*this, gf, #op); \
1198  \
1199  dimensionedInternalField() op gf.dimensionedInternalField(); \
1200  boundaryField() op gf.boundaryField(); \
1201 } \
1202  \
1203 template<class Type, template<class> class PatchField, class GeoMesh> \
1204 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \
1205 ( \
1206  const tmp<GeometricField<TYPE, PatchField, GeoMesh> >& tgf \
1207 ) \
1208 { \
1209  operator op(tgf()); \
1210  tgf.clear(); \
1211 } \
1212  \
1213 template<class Type, template<class> class PatchField, class GeoMesh> \
1214 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \
1215 ( \
1216  const dimensioned<TYPE>& dt \
1217 ) \
1218 { \
1219  dimensionedInternalField() op dt; \
1220  boundaryField() op dt.value(); \
1221 }
1222 
1227 
1228 #undef COMPUTED_ASSIGNMENT
1229 
1230 
1231 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
1232 
1233 template<class Type, template<class> class PatchField, class GeoMesh>
1234 Foam::Ostream& Foam::operator<<
1235 (
1236  Ostream& os,
1238 )
1239 {
1240  gf.dimensionedInternalField().writeData(os, "internalField");
1241  os << nl;
1242  gf.boundaryField().writeEntry("boundaryField", os);
1243 
1244  // Check state of IOstream
1245  os.check
1246  (
1247  "Ostream& operator<<(Ostream&, "
1248  "const GeometricField<Type, PatchField, GeoMesh>&)"
1249  );
1250 
1251  return (os);
1252 }
1253 
1254 
1255 template<class Type, template<class> class PatchField, class GeoMesh>
1256 Foam::Ostream& Foam::operator<<
1257 (
1258  Ostream& os,
1260 )
1261 {
1262  os << tgf();
1263  tgf.clear();
1264  return os;
1265 }
1266 
1267 
1268 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1269 
1270 #undef checkField
1271 
1272 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1273 
1274 #include "GeometricBoundaryField.C"
1275 #include "GeometricFieldFunctions.C"
1276 
1277 // ************************************************************************* //
#define COMPUTED_ASSIGNMENT(TYPE, op)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
UEqn relax()
unsigned char direction
Definition: direction.H:43
faceListList boundary(nPatches)
Field< Type >::cmptType cmptType
virtual ~GeometricField()
Destructor.
bool needReference() const
Does the field need a reference level for solution.
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
void deleteDemandDrivenData(DataPtr &dataPtr)
label timeIndex
Definition: getTimeIndex.H:4
A class for handling words, derived from string.
Definition: word.H:59
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
InternalField & internalField()
Return internal field.
DimensionedInternalField & dimensionedInternalField()
Return dimensioned internal field.
messageStream Info
dynamicFvMesh & mesh
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
void writeMinMax(Ostream &os) const
Helper function to write the min and max to an Ostream.
#define checkField(gf1, gf2, op)
dictionary dict
Generic dimensioned Type class.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOerror FatalIOError
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define WarningIn(functionName)
Report a warning using Foam::Warning.
void readFields(const Mesh &mesh, const IOobjectList &objects, PtrList< GeoField > &fields)
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
GeometricField(const IOobject &, const Mesh &, const dimensionSet &, const word &patchFieldType=PatchField< Type >::calculatedType())
Constructor given IOobject, mesh, dimensions and patch type.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Dimension set for the base types.
Definition: dimensionSet.H:116
const GeometricField< Type, PatchField, GeoMesh > & prevIter() const
Return previous iteration field.
#define forAll(list, i)
Definition: UList.H:421
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:333
label patchi
label nOldTimes() const
Return the number of old time fields stored.
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
Template functions to aid in the implementation of demand driven data.
const dimensionSet & dimensions() const
Return dimensions.
#define InfoIn(functionName)
Report a information message using Foam::Info.
Pre-declare SubField and related Field type.
Definition: Field.H:57
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const word & name() const
Return name.
Definition: IOobject.H:260
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
rDeltaT dimensionedInternalField()
Generic GeometricField class.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
Definition: GeoMesh.H:46
error FatalError
void correctBoundaryConditions()
Correct boundary field.
void replace(const direction, const GeometricField< cmptType, PatchField, GeoMesh > &)
void max(const dimensioned< Type > &)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
label timeIndex() const
Return the time index of the field.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
void storeOldTimes() const
Store the old-time fields.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
bool writeData(Ostream &) const
WriteData member function required by regIOobject.
void storePrevIter() const
Store the field as the previous iteration value.
void relax()
Relax field (for steady-state solution).
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
void min(const dimensioned< Type > &)
A class for managing temporary objects.
Definition: PtrList.H:118
void storeOldTime() const
Store the old-time field.
conserve internalField()+
const Type & value() const
Return const reference to value.
word select(bool final) const
Select the final iteration parameters if `final&#39; is true.
word timeName
Definition: getTimeIndex.H:3