FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "GeometricField.H"
27 #include <OpenFOAM/Time.H>
29 #include <OpenFOAM/dictionary.H>
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 // check mesh for two fields
34 
35 #define checkField(gf1, gf2, op) \
36 if ((gf1).mesh() != (gf2).mesh()) \
37 { \
38  FatalErrorIn("checkField(gf1, gf2, op)") \
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 <
52  GeometricBoundaryField
53 >
55 (
56  const dictionary& fieldDict
57 )
58 {
59  DimensionedField<Type, GeoMesh>::readField(fieldDict, "internalField");
60 
61  tmp<GeometricBoundaryField> tboundaryField
62  (
63  new GeometricBoundaryField
64  (
65  this->mesh().boundary(),
66  *this,
67  fieldDict.subDict("boundaryField")
68  )
69  );
70 
71  if (fieldDict.found("referenceLevel"))
72  {
73  Type fieldAverage(pTraits<Type>(fieldDict.lookup("referenceLevel")));
74 
75  Field<Type>::operator+=(fieldAverage);
76 
77  GeometricBoundaryField& boundaryField = tboundaryField();
78 
79  forAll(boundaryField, patchi)
80  {
81  boundaryField[patchi] == boundaryField[patchi] + fieldAverage;
82  }
83  }
84 
85  return tboundaryField;
86 }
87 
88 
89 template<class Type, template<class> class PatchField, class GeoMesh>
91 <
93  GeometricBoundaryField
94 >
96 {
97  if (is.version() < 2.0)
98  {
100  (
101  "GeometricField<Type, PatchField, GeoMesh>::readField(Istream&)",
102  is
103  ) << "IO versions < 2.0 are not supported."
104  << exit(FatalIOError);
105  }
106 
107  return readField(dictionary(is));
108 }
109 
110 
111 template<class Type, template<class> class PatchField, class GeoMesh>
113 {
114  if (this->readOpt() == IOobject::MUST_READ)
115  {
116  WarningIn
117  (
118  "GeometricField<Type, PatchField, GeoMesh>::readIfPresent()"
119  ) << "read option IOobject::MUST_READ "
120  << "suggests that a read constructor for field " << this->name()
121  << " would be more appropriate." << endl;
122  }
123  else if (this->readOpt() == IOobject::READ_IF_PRESENT && this->headerOk())
124  {
125  boundaryField_.transfer(readField(this->readStream(typeName))());
126  this->close();
127 
128  // Check compatibility between field and mesh
129  if (this->size() != GeoMesh::size(this->mesh()))
130  {
132  (
133  "GeometricField<Type, PatchField, GeoMesh>::"
134  "readIfPresent()",
135  this->readStream(typeName)
136  ) << " number of field elements = " << this->size()
137  << " number of mesh elements = "
138  << GeoMesh::size(this->mesh())
139  << exit(FatalIOError);
140  }
141 
142  readOldTimeIfPresent();
143 
144  return true;
145  }
146 
147  return false;
148 }
149 
150 
151 template<class Type, template<class> class PatchField, class GeoMesh>
153 {
154  // Read the old time field if present
155  IOobject field0
156  (
157  this->name() + "_0",
158  this->time().timeName(),
159  this->db(),
160  IOobject::READ_IF_PRESENT,
161  IOobject::AUTO_WRITE
162  );
163 
164  if (field0.headerOk())
165  {
166  if (debug)
167  {
168  Info<< "Reading old time level for field"
169  << endl << this->info() << endl;
170  }
171 
172  field0Ptr_ = new GeometricField<Type, PatchField, GeoMesh>
173  (
174  field0,
175  this->mesh()
176  );
177 
178  field0Ptr_->timeIndex_ = timeIndex_ - 1;
179 
180  if (!field0Ptr_->readOldTimeIfPresent())
181  {
182  field0Ptr_->oldTime();
183  }
184 
185  return true;
186  }
187 
188  return false;
189 }
190 
191 
192 // * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * //
193 
194 // Constructor given a GeometricField and dimensionSet
195 // This allocates storage for the field but not values.
196 // Note : This constructor should only be used to
197 // construct TEMPORARY variables
198 template<class Type, template<class> class PatchField, class GeoMesh>
200 (
201  const IOobject& io,
202  const Mesh& mesh,
203  const dimensionSet& ds,
204  const word& patchFieldType
205 )
206 :
208  timeIndex_(this->time().timeIndex()),
209  field0Ptr_(NULL),
210  fieldPrevIterPtr_(NULL),
211  boundaryField_(mesh.boundary(), *this, patchFieldType)
212 {
213  if (debug)
214  {
215  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
216  "creating temporary"
217  << endl << this->info() << endl;
218  }
219 
220  readIfPresent();
221 }
222 
223 
224 // Constructor given a GeometricField and dimensionSet
225 // This allocates storage for the field but not values.
226 // Note : This constructor should only be used to
227 // construct TEMPORARY variables
228 template<class Type, template<class> class PatchField, class GeoMesh>
230 (
231  const IOobject& io,
232  const Mesh& mesh,
233  const dimensionSet& ds,
234  const wordList& patchFieldTypes
235 )
236 :
238  timeIndex_(this->time().timeIndex()),
239  field0Ptr_(NULL),
240  fieldPrevIterPtr_(NULL),
241  boundaryField_(mesh.boundary(), *this, patchFieldTypes)
242 {
243  if (debug)
244  {
245  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
246  "creating temporary"
247  << endl << this->info() << endl;
248  }
249 
250  readIfPresent();
251 }
252 
253 
254 // Constructor given a GeometricField and dimensioned<Type>
255 template<class Type, template<class> class PatchField, class GeoMesh>
257 (
258  const IOobject& io,
259  const Mesh& mesh,
260  const dimensioned<Type>& dt,
261  const word& patchFieldType
262 )
263 :
265  timeIndex_(this->time().timeIndex()),
266  field0Ptr_(NULL),
267  fieldPrevIterPtr_(NULL),
268  boundaryField_(mesh.boundary(), *this, patchFieldType)
269 {
270  if (debug)
271  {
272  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
273  "creating temporary"
274  << endl << this->info() << endl;
275  }
276 
277  boundaryField_ == dt.value();
278 
279  readIfPresent();
280 }
281 
282 
283 // Constructor given a GeometricField and dimensioned<Type>
284 template<class Type, template<class> class PatchField, class GeoMesh>
286 (
287  const IOobject& io,
288  const Mesh& mesh,
289  const dimensioned<Type>& dt,
290  const wordList& patchFieldTypes
291 )
292 :
294  timeIndex_(this->time().timeIndex()),
295  field0Ptr_(NULL),
296  fieldPrevIterPtr_(NULL),
297  boundaryField_(mesh.boundary(), *this, patchFieldTypes)
298 {
299  if (debug)
300  {
301  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
302  "creating temporary"
303  << endl << this->info() << endl;
304  }
305 
306  boundaryField_ == dt.value();
307 
308  readIfPresent();
309 }
310 
311 
312 // construct from components
313 template<class Type, template<class> class PatchField, class GeoMesh>
315 (
316  const IOobject& io,
317  const Mesh& mesh,
318  const dimensionSet& ds,
319  const Field<Type>& iField,
320  const PtrList<PatchField<Type> >& ptfl
321 )
322 :
323  DimensionedField<Type, GeoMesh>(io, mesh, ds, iField),
324  timeIndex_(this->time().timeIndex()),
325  field0Ptr_(NULL),
326  fieldPrevIterPtr_(NULL),
327  boundaryField_(mesh.boundary(), *this, ptfl)
328 {
329  if (debug)
330  {
331  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
332  "constructing from components"
333  << endl << this->info() << endl;
334  }
335 
336  readIfPresent();
337 }
338 
339 
340 template<class Type, template<class> class PatchField, class GeoMesh>
342 (
343  const IOobject& io,
344  const Mesh& mesh
345 )
346 :
348  timeIndex_(this->time().timeIndex()),
349  field0Ptr_(NULL),
350  fieldPrevIterPtr_(NULL),
351  boundaryField_(*this, readField(this->readStream(typeName)))
352 {
353  this->close();
354 
355  // Check compatibility between field and mesh
356 
357  if (this->size() != GeoMesh::size(this->mesh()))
358  {
360  (
361  "GeometricField<Type, PatchField, GeoMesh>::GeometricField"
362  "(const IOobject&, const Mesh&)",
363  this->readStream(typeName)
364  ) << " number of field elements = " << this->size()
365  << " number of mesh elements = " << GeoMesh::size(this->mesh())
366  << exit(FatalIOError);
367  }
368 
369  readOldTimeIfPresent();
370 
371  if (debug)
372  {
373  Info<< "Finishing read-construct of "
374  "GeometricField<Type, PatchField, GeoMesh>"
375  << endl << this->info() << endl;
376  }
377 }
378 
379 
380 template<class Type, template<class> class PatchField, class GeoMesh>
382 (
383  const IOobject& io,
384  const Mesh& mesh,
385  Istream& is
386 )
387 :
389  timeIndex_(this->time().timeIndex()),
390  field0Ptr_(NULL),
391  fieldPrevIterPtr_(NULL),
392  boundaryField_(*this, readField(is))
393 {
394  // Check compatibility between field and mesh
395 
396  if (this->size() != GeoMesh::size(this->mesh()))
397  {
399  (
400  "GeometricField<Type, PatchField, GeoMesh>::GeometricField"
401  "(const IOobject&, const Mesh&, Istream&)",
402  is
403  ) << " number of field elements = " << this->size()
404  << " number of mesh elements = " << GeoMesh::size(this->mesh())
405  << exit(FatalIOError);
406  }
407 
408  readOldTimeIfPresent();
409 
410  if (debug)
411  {
412  Info<< "Finishing read-construct of "
413  "GeometricField<Type, PatchField, GeoMesh>"
414  << endl << this->info() << endl;
415  }
416 }
417 
418 
419 template<class Type, template<class> class PatchField, class GeoMesh>
421 (
422  const IOobject& io,
423  const Mesh& mesh,
424  const dictionary& dict
425 )
426 :
428  timeIndex_(this->time().timeIndex()),
429  field0Ptr_(NULL),
430  fieldPrevIterPtr_(NULL),
431  boundaryField_(*this, readField(dict))
432 {
433  // Check compatibility between field and mesh
434 
435  if (this->size() != GeoMesh::size(this->mesh()))
436  {
438  (
439  "GeometricField<Type, PatchField, GeoMesh>::GeometricField"
440  "(const IOobject&, const Mesh&, const dictionary&)"
441  ) << " number of field elements = " << this->size()
442  << " number of mesh elements = " << GeoMesh::size(this->mesh())
443  << exit(FatalIOError);
444  }
445 
446  readOldTimeIfPresent();
447 
448  if (debug)
449  {
450  Info<< "Finishing dictionary-construct of "
451  "GeometricField<Type, PatchField, GeoMesh>"
452  << endl << this->info() << endl;
453  }
454 }
455 
456 
457 // construct as copy
458 template<class Type, template<class> class PatchField, class GeoMesh>
460 (
462 )
463 :
465  timeIndex_(gf.timeIndex()),
466  field0Ptr_(NULL),
467  fieldPrevIterPtr_(NULL),
468  boundaryField_(*this, gf.boundaryField_)
469 {
470  if (debug)
471  {
472  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
473  "constructing as copy"
474  << endl << this->info() << endl;
475  }
476 
477  if (gf.field0Ptr_)
478  {
480  (
481  *gf.field0Ptr_
482  );
483  }
484 
485  this->writeOpt() = IOobject::NO_WRITE;
486 }
487 
488 // construct as copy of tmp<GeometricField> deleting argument
489 #ifdef ConstructFromTmp
490 template<class Type, template<class> class PatchField, class GeoMesh>
492 (
494 )
495 :
497  (
498  const_cast<GeometricField<Type, PatchField, GeoMesh>&>(tgf()),
499  tgf.isTmp()
500  ),
501  timeIndex_(tgf().timeIndex()),
502  field0Ptr_(NULL),
503  fieldPrevIterPtr_(NULL),
504  boundaryField_(*this, tgf().boundaryField_)
505 {
506  if (debug)
507  {
508  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
509  "constructing as copy"
510  << endl << this->info() << endl;
511  }
512 
513  this->writeOpt() = IOobject::NO_WRITE;
514 
515  tgf.clear();
516 }
517 #endif
518 
519 
520 // construct as copy resetting IO parameters
521 template<class Type, template<class> class PatchField, class GeoMesh>
522 Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
523 (
524  const IOobject& io,
526 )
527 :
529  timeIndex_(gf.timeIndex()),
530  field0Ptr_(NULL),
531  fieldPrevIterPtr_(NULL),
532  boundaryField_(*this, gf.boundaryField_)
533 {
534  if (debug)
535  {
536  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
537  "constructing as copy resetting IO params"
538  << endl << this->info() << endl;
539  }
540 
541  if (!readIfPresent() && gf.field0Ptr_)
542  {
544  (
545  io.name() + "_0",
546  *gf.field0Ptr_
547  );
548  }
549 }
550 
551 
552 // construct as copy resetting name
553 template<class Type, template<class> class PatchField, class GeoMesh>
554 Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
555 (
556  const word& newName,
558 )
559 :
560  DimensionedField<Type, GeoMesh>(newName, gf),
561  timeIndex_(gf.timeIndex()),
562  field0Ptr_(NULL),
563  fieldPrevIterPtr_(NULL),
564  boundaryField_(*this, gf.boundaryField_)
565 {
566  if (debug)
567  {
568  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
569  "constructing as copy resetting name"
570  << endl << this->info() << endl;
571  }
572 
573  if (!readIfPresent() && gf.field0Ptr_)
574  {
576  (
577  newName + "_0",
578  *gf.field0Ptr_
579  );
580  }
581 }
582 
583 
584 // construct as copy resetting name
585 #ifdef ConstructFromTmp
586 template<class Type, template<class> class PatchField, class GeoMesh>
587 Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
588 (
589  const word& newName,
591 )
592 :
594  (
595  newName,
596  const_cast<GeometricField<Type, PatchField, GeoMesh>&>(tgf()),
597  tgf.isTmp()
598  ),
599  timeIndex_(tgf().timeIndex()),
600  field0Ptr_(NULL),
601  fieldPrevIterPtr_(NULL),
602  boundaryField_(*this, tgf().boundaryField_)
603 {
604  if (debug)
605  {
606  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
607  "constructing from tmp resetting name"
608  << endl << this->info() << endl;
609  }
610 
611  tgf.clear();
612 }
613 #endif
614 
615 // construct as copy resetting IO parameters and patch type
616 template<class Type, template<class> class PatchField, class GeoMesh>
617 Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
618 (
619  const IOobject& io,
621  const word& patchFieldType
622 )
623 :
625  timeIndex_(gf.timeIndex()),
626  field0Ptr_(NULL),
627  fieldPrevIterPtr_(NULL),
628  boundaryField_(this->mesh().boundary(), *this, patchFieldType)
629 {
630  if (debug)
631  {
632  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
633  "constructing as copy resetting IO params"
634  << endl << this->info() << endl;
635  }
636 
637  boundaryField_ == gf.boundaryField_;
638 
639  if (!readIfPresent() && gf.field0Ptr_)
640  {
642  (
643  io.name() + "_0",
644  *gf.field0Ptr_
645  );
646  }
647 }
648 
649 
650 // construct as copy resetting IO parameters and boundary types
651 template<class Type, template<class> class PatchField, class GeoMesh>
652 Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
653 (
654  const IOobject& io,
656  const wordList& patchFieldTypes
657 )
658 :
660  timeIndex_(gf.timeIndex()),
661  field0Ptr_(NULL),
662  fieldPrevIterPtr_(NULL),
663  boundaryField_(this->mesh().boundary(), *this, patchFieldTypes)
664 {
665  if (debug)
666  {
667  Info<< "GeometricField<Type, PatchField, GeoMesh>::GeometricField : "
668  "constructing as copy resetting IO params"
669  << endl << this->info() << endl;
670  }
671 
672  boundaryField_ == gf.boundaryField_;
673 
674  if (!readIfPresent() && gf.field0Ptr_)
675  {
677  (
678  io.name() + "_0",
679  *gf.field0Ptr_
680  );
681  }
682 }
683 
684 
685 // * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * * //
686 
687 template<class Type, template<class> class PatchField, class GeoMesh>
689 {
690  deleteDemandDrivenData(field0Ptr_);
691  deleteDemandDrivenData(fieldPrevIterPtr_);
692 }
693 
694 
695 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
696 
697 template<class Type, template<class> class PatchField, class GeoMesh>
698 typename
701 {
702  this->setUpToDate();
703  storeOldTimes();
704  return *this;
705 }
706 
707 
708 template<class Type, template<class> class PatchField, class GeoMesh>
709 typename
712 {
713  this->setUpToDate();
714  storeOldTimes();
715  return *this;
716 }
717 
718 
719 // Return reference to GeometricBoundaryField
720 template<class Type, template<class> class PatchField, class GeoMesh>
721 typename
724 {
725  this->setUpToDate();
726  storeOldTimes();
727  return boundaryField_;
728 }
729 
730 
731 // Store old-time field
732 template<class Type, template<class> class PatchField, class GeoMesh>
734 {
735  if
736  (
737  field0Ptr_
738  && timeIndex_ != this->time().timeIndex()
739  && !(
740  this->name().size() > 2
741  && this->name()(this->name().size()-2, 2) == "_0"
742  )
743  )
744  {
745  storeOldTime();
746  }
747 
748  // Correct time index
749  timeIndex_ = this->time().timeIndex();
750 }
751 
752 // Store old-time field
753 template<class Type, template<class> class PatchField, class GeoMesh>
755 {
756  if (field0Ptr_)
757  {
758  field0Ptr_->storeOldTime();
759 
760  if (debug)
761  {
762  Info<< "Storing old time field for field" << endl
763  << this->info() << endl;
764  }
765 
766  *field0Ptr_ == *this;
767  field0Ptr_->timeIndex_ = timeIndex_;
768 
769  if (field0Ptr_->field0Ptr_)
770  {
771  field0Ptr_->writeOpt() = this->writeOpt();
772  }
773  }
774 }
775 
776 // Return the number of old time fields stored
777 template<class Type, template<class> class PatchField, class GeoMesh>
779 {
780  if (field0Ptr_)
781  {
782  return field0Ptr_->nOldTimes() + 1;
783  }
784  else
785  {
786  return 0;
787  }
788 }
789 
790 // Return old time internal field
791 template<class Type, template<class> class PatchField, class GeoMesh>
794 {
795  if (!field0Ptr_)
796  {
798  (
799  IOobject
800  (
801  this->name() + "_0",
802  this->time().timeName(),
803  this->db()
804  ),
805  *this
806  );
807  }
808  else
809  {
810  storeOldTimes();
811  }
812 
813  return *field0Ptr_;
814 }
815 
816 // Return old time internal field
817 template<class Type, template<class> class PatchField, class GeoMesh>
820 {
821  static_cast<const GeometricField<Type, PatchField, GeoMesh>&>(*this)
822  .oldTime();
823 
824  return *field0Ptr_;
825 }
826 
827 
828 // Store previous iteration field
829 template<class Type, template<class> class PatchField, class GeoMesh>
831 {
832  if (!fieldPrevIterPtr_)
833  {
834  if (debug)
835  {
836  Info<< "Allocating previous iteration field" << endl
837  << this->info() << endl;
838  }
839 
840  fieldPrevIterPtr_ = new GeometricField<Type, PatchField, GeoMesh>
841  (
842  this->name() + "PrevIter",
843  *this
844  );
845  }
846  else
847  {
848  *fieldPrevIterPtr_ == *this;
849  }
850 }
851 
852 
853 // Return previous iteration field
854 template<class Type, template<class> class PatchField, class GeoMesh>
857 {
858  if (!fieldPrevIterPtr_)
859  {
861  (
862  "GeometricField<Type, PatchField, GeoMesh>::prevIter() const"
863  ) << "previous iteration field" << endl << this->info() << endl
864  << " not stored."
865  << " Use field.storePrevIter() at start of iteration."
866  << abort(FatalError);
867  }
868 
869  return *fieldPrevIterPtr_;
870 }
871 
872 
873 // Correct the boundary conditions
874 template<class Type, template<class> class PatchField, class GeoMesh>
877 {
878  this->setUpToDate();
879  storeOldTimes();
880  boundaryField_.evaluate();
881 }
882 
883 
884 // Does the field need a reference level for solution
885 template<class Type, template<class> class PatchField, class GeoMesh>
887 {
888  // Search all boundary conditions, if any are
889  // fixed-value or mixed (Robin) do not set reference level for solution.
890 
891  bool needRef = true;
892 
893  forAll(boundaryField_, patchi)
894  {
895  if (boundaryField_[patchi].fixesValue())
896  {
897  needRef = false;
898  break;
899  }
900  }
901 
902  reduce(needRef, andOp<bool>());
903 
904  return needRef;
905 }
906 
907 
908 template<class Type, template<class> class PatchField, class GeoMesh>
910 {
911  operator==(prevIter() + alpha*(*this - prevIter()));
912 }
913 
914 
915 template<class Type, template<class> class PatchField, class GeoMesh>
917 {
918  scalar alpha = 0;
919 
920  if (this->mesh().relax(this->name()))
921  {
922  alpha = this->mesh().relaxationFactor(this->name());
923  }
924 
925  if (alpha > 0)
926  {
927  relax(alpha);
928  }
929 }
930 
931 
932 template<class Type, template<class> class PatchField, class GeoMesh>
934 (
935  bool final
936 ) const
937 {
938  if (final)
939  {
940  return this->name() + "Final";
941  }
942  else
943  {
944  return this->name();
945  }
946 }
947 
948 
949 // writeData member function required by regIOobject
950 template<class Type, template<class> class PatchField, class GeoMesh>
952 writeData(Ostream& os) const
953 {
954  os << *this;
955  return os.good();
956 }
957 
958 
959 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
960 
961 template<class Type, template<class> class PatchField, class GeoMesh>
964 {
966  (
968  (
969  IOobject
970  (
971  this->name() + ".T()",
972  this->instance(),
973  this->db()
974  ),
975  this->mesh(),
976  this->dimensions()
977  )
978  );
979 
980  Foam::T(result().internalField(), internalField());
981  Foam::T(result().boundaryField(), boundaryField());
982 
983  return result;
984 }
985 
986 
987 template<class Type, template<class> class PatchField, class GeoMesh>
988 Foam::tmp
989 <
991  <
993  PatchField,
994  GeoMesh
995  >
996 >
998 (
999  const direction d
1000 ) const
1001 {
1003  (
1005  (
1006  IOobject
1007  (
1008  this->name() + ".component(" + Foam::name(d) + ')',
1009  this->instance(),
1010  this->db()
1011  ),
1012  this->mesh(),
1013  this->dimensions()
1014  )
1015  );
1016 
1017  Foam::component(Component().internalField(), internalField(), d);
1018  Foam::component(Component().boundaryField(), boundaryField(), d);
1019 
1020  return Component;
1021 }
1022 
1023 
1024 template<class Type, template<class> class PatchField, class GeoMesh>
1026 (
1027  const direction d,
1028  const GeometricField
1029  <
1030  typename GeometricField<Type, PatchField, GeoMesh>::cmptType,
1031  PatchField,
1032  GeoMesh
1033  >& gcf
1034 )
1035 {
1036  internalField().replace(d, gcf.internalField());
1037  boundaryField().replace(d, gcf.boundaryField());
1038 }
1039 
1040 
1041 template<class Type, template<class> class PatchField, class GeoMesh>
1044  const direction d,
1045  const dimensioned<cmptType>& ds
1046 )
1047 {
1048  internalField().replace(d, ds.value());
1049  boundaryField().replace(d, ds.value());
1050 }
1051 
1052 
1053 template<class Type, template<class> class PatchField, class GeoMesh>
1056  const dimensioned<Type>& dt
1057 )
1058 {
1059  Foam::max(internalField(), internalField(), dt.value());
1061 }
1062 
1063 
1064 template<class Type, template<class> class PatchField, class GeoMesh>
1067  const dimensioned<Type>& dt
1068 )
1069 {
1070  Foam::min(internalField(), internalField(), dt.value());
1072 }
1073 
1074 
1075 template<class Type, template<class> class PatchField, class GeoMesh>
1077 {
1078  internalField().negate();
1079  boundaryField().negate();
1080 }
1081 
1082 
1083 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
1084 
1085 template<class Type, template<class> class PatchField, class GeoMesh>
1089 )
1090 {
1091  if (this == &gf)
1092  {
1093  FatalErrorIn
1094  (
1095  "GeometricField<Type, PatchField, GeoMesh>::operator="
1096  "(const GeometricField<Type, PatchField, GeoMesh>&)"
1097  ) << "attempted assignment to self"
1098  << abort(FatalError);
1099  }
1100 
1101  checkField(*this, gf, "=");
1102 
1103  // only equate field contents not ID
1104 
1105  dimensionedInternalField() = gf.dimensionedInternalField();
1106  boundaryField() = gf.boundaryField();
1107 }
1108 
1109 
1110 template<class Type, template<class> class PatchField, class GeoMesh>
1111 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator=
1114 )
1115 {
1116  if (this == &(tgf()))
1117  {
1118  FatalErrorIn
1119  (
1120  "GeometricField<Type, PatchField, GeoMesh>::operator="
1121  "(const tmp<GeometricField<Type, PatchField, GeoMesh> >&)"
1122  ) << "attempted assignment to self"
1123  << abort(FatalError);
1124  }
1125 
1127 
1128  checkField(*this, gf, "=");
1129 
1130  // only equate field contents not ID
1131 
1132  this->dimensions() = gf.dimensions();
1133 
1134  // This is dodgy stuff, don't try it at home.
1135  internalField().transfer
1136  (
1137  const_cast<Field<Type>&>(gf.internalField())
1138  );
1139 
1140  boundaryField() = gf.boundaryField();
1141 
1142  tgf.clear();
1143 }
1144 
1145 
1146 template<class Type, template<class> class PatchField, class GeoMesh>
1147 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator=
1149  const dimensioned<Type>& dt
1150 )
1151 {
1152  dimensionedInternalField() = dt;
1153  boundaryField() = dt.value();
1154 }
1155 
1156 
1157 template<class Type, template<class> class PatchField, class GeoMesh>
1158 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator==
1161 )
1162 {
1164 
1165  checkField(*this, gf, "==");
1166 
1167  // only equate field contents not ID
1168 
1169  dimensionedInternalField() = gf.dimensionedInternalField();
1170  boundaryField() == gf.boundaryField();
1171 
1172  tgf.clear();
1173 }
1174 
1175 
1176 template<class Type, template<class> class PatchField, class GeoMesh>
1177 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator==
1179  const dimensioned<Type>& dt
1180 )
1181 {
1182  dimensionedInternalField() = dt;
1183  boundaryField() == dt.value();
1184 }
1185 
1186 
1187 #define COMPUTED_ASSIGNMENT(TYPE, op) \
1188  \
1189 template<class Type, template<class> class PatchField, class GeoMesh> \
1190 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \
1191 ( \
1192  const GeometricField<TYPE, PatchField, GeoMesh>& gf \
1193 ) \
1194 { \
1195  checkField(*this, gf, #op); \
1196  \
1197  dimensionedInternalField() op gf.dimensionedInternalField(); \
1198  boundaryField() op gf.boundaryField(); \
1199 } \
1200  \
1201 template<class Type, template<class> class PatchField, class GeoMesh> \
1202 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \
1203 ( \
1204  const tmp<GeometricField<TYPE, PatchField, GeoMesh> >& tgf \
1205 ) \
1206 { \
1207  operator op(tgf()); \
1208  tgf.clear(); \
1209 } \
1210  \
1211 template<class Type, template<class> class PatchField, class GeoMesh> \
1212 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \
1213 ( \
1214  const dimensioned<TYPE>& dt \
1215 ) \
1216 { \
1217  dimensionedInternalField() op dt; \
1218  boundaryField() op dt.value(); \
1219 }
1220 
1225 
1226 #undef COMPUTED_ASSIGNMENT
1227 
1228 
1229 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
1230 
1231 template<class Type, template<class> class PatchField, class GeoMesh>
1232 Foam::Ostream& Foam::operator<<
1234  Ostream& os,
1236 )
1237 {
1238  gf.dimensionedInternalField().writeData(os, "internalField");
1239  os << nl;
1240  gf.boundaryField().writeEntry("boundaryField", os);
1241 
1242  // Check state of IOstream
1243  os.check
1244  (
1245  "Ostream& operator<<(Ostream&, "
1246  "const GeometricField<Type, PatchField, GeoMesh>&)"
1247  );
1248 
1249  return (os);
1250 }
1251 
1252 
1253 template<class Type, template<class> class PatchField, class GeoMesh>
1254 Foam::Ostream& Foam::operator<<
1256  Ostream& os,
1258 )
1259 {
1260  os << tgf();
1261  tgf.clear();
1262  return os;
1263 }
1264 
1265 
1266 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1267 
1268 #undef checkField
1269 
1270 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1271 
1272 #include "GeometricBoundaryField.C"
1273 #include "GeometricFieldFunctions.C"
1274 
1275 // ************************ vim: set sw=4 sts=4 et: ************************ //