40 #define checkField(fld1, fld2, op) \ 41 if (&(fld1).mesh() != &(fld2).mesh()) \ 43 FatalErrorInFunction \ 44 << "Different mesh for fields " \ 45 << (fld1).name() << " and " << (fld2).name() \ 46 << " during operation " << op \ 47 << abort(FatalError); \ 53 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
56 const dictionary&
dict 59 Internal::readField(
dict,
"internalField");
61 boundaryField_.readField(*
this,
dict.subDict(
"boundaryField"));
65 if (
dict.readIfPresent(
"referenceLevel", refLevel))
67 Field<Type>::operator+=(refLevel);
69 forAll(boundaryField_, patchi)
71 boundaryField_[patchi] == boundaryField_[patchi] + refLevel;
77 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
80 const localIOdictionary
dict 88 IOobjectOption::MUST_READ,
89 IOobjectOption::NO_WRITE,
90 IOobjectOption::NO_REGISTER
101 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
104 if (this->isReadRequired())
107 <<
"The readOption MUST_READ or READ_MODIFIED" 108 <<
" suggests that a read constructor for field " << this->
name()
109 <<
" would be more appropriate." <<
endl;
113 this->isReadOptional()
114 && this->
template typeHeaderOk<GeometricField<Type, PatchField, GeoMesh>>
121 readOldTimeIfPresent();
130 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
139 IOobjectOption::LAZY_READ,
140 IOobjectOption::AUTO_WRITE,
141 this->registerObject()
146 field0.template typeHeaderOk<GeometricField<Type, PatchField, GeoMesh>>
153 <<
"Reading old time level for field" <<
nl << this->info() <<
endl;
155 field0Ptr_ =
new GeometricField<Type, PatchField, GeoMesh>
164 field0Ptr_->oriented() = this->oriented();
166 field0Ptr_->timeIndex_ = timeIndex_ - 1;
168 if (!field0Ptr_->readOldTimeIfPresent())
170 field0Ptr_->oldTime();
182 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
188 const word& patchFieldType
191 Internal(
io,
mesh, dims, false),
194 fieldPrevIterPtr_(nullptr),
198 <<
"Creating" <<
nl << this->info() <<
endl;
204 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
214 Internal(
io,
mesh, dims, false),
217 fieldPrevIterPtr_(nullptr),
218 boundaryField_(
mesh.
boundary(), *this, patchFieldTypes, actualPatchTypes)
221 <<
"Creating" <<
nl << this->info() <<
endl;
227 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
234 const word& patchFieldType
237 Internal(
io,
mesh, value, dims, false),
240 fieldPrevIterPtr_(nullptr),
244 <<
"Creating" <<
nl << this->info() <<
endl;
246 boundaryField_ == value;
252 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
263 Internal(
io,
mesh, value, dims, false),
266 fieldPrevIterPtr_(nullptr),
267 boundaryField_(
mesh.
boundary(), *this, patchFieldTypes, actualPatchTypes)
270 <<
"Creating" <<
nl << this->info() <<
endl;
272 boundaryField_ == value;
278 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
284 const word& patchFieldType
298 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
320 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
324 const Internal& diField,
325 const PtrList<PatchField<Type>>& ptfl
328 Internal(
io, diField),
331 fieldPrevIterPtr_(nullptr),
335 <<
"Copy construct from components" <<
nl << this->info() <<
endl;
341 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
346 const PtrList<PatchField<Type>>& ptfl
349 Internal(
io,
std::move(diField)),
352 fieldPrevIterPtr_(nullptr),
356 <<
"Move construct from components" <<
nl << this->info() <<
endl;
362 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
367 const PtrList<PatchField<Type>>& ptfl
370 Internal(
io, tfield),
373 fieldPrevIterPtr_(nullptr),
377 <<
"Construct from tmp internalField" <<
nl << this->info() <<
endl;
383 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
386 const Internal& diField,
387 const PtrList<PatchField<Type>>& ptfl
393 fieldPrevIterPtr_(nullptr),
397 <<
"Copy construct from components" <<
nl << this->info() <<
endl;
403 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
407 const PtrList<PatchField<Type>>& ptfl
410 Internal(
std::move(diField)),
413 fieldPrevIterPtr_(nullptr),
417 <<
"Move construct from components" <<
nl << this->info() <<
endl;
423 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
430 const word& patchFieldType
433 Internal(
io,
mesh, dims, iField),
436 fieldPrevIterPtr_(nullptr),
440 <<
"Copy construct from internal field" <<
nl << this->info() <<
endl;
446 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
453 const word& patchFieldType
456 Internal(
io,
mesh, dims,
std::move(iField)),
459 fieldPrevIterPtr_(nullptr),
463 <<
"Move construct from internal field" <<
nl << this->info() <<
endl;
469 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
476 const PtrList<PatchField<Type>>& ptfl
479 Internal(
io,
mesh, dims, iField),
482 fieldPrevIterPtr_(nullptr),
486 <<
"Copy construct from components" <<
nl << this->info() <<
endl;
492 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
499 const PtrList<PatchField<Type>>& ptfl
502 Internal(
io,
mesh, dims,
std::move(iField)),
505 fieldPrevIterPtr_(nullptr),
509 <<
"Move construct from components" <<
nl << this->info() <<
endl;
515 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
522 const PtrList<PatchField<Type>>& ptfl
525 Internal(
io,
mesh, dims, tfield),
528 fieldPrevIterPtr_(nullptr),
532 <<
"Construct from tmp internalField" <<
nl << this->info() <<
endl;
538 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
543 const bool readOldTime
549 fieldPrevIterPtr_(nullptr),
553 <<
"Read construct" <<
nl << this->info() <<
endl;
555 if (!this->isAnyRead())
560 <<
"Had readOption NO_READ for field " 561 << this->
name() <<
", but constructor always reads field!" 569 readOldTimeIfPresent();
573 <<
"Finishing read-construction" <<
nl << this->info() <<
endl;
577 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
588 fieldPrevIterPtr_(nullptr),
594 <<
"Finishing dictionary-construct" <<
nl << this->info() <<
endl;
598 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
607 fieldPrevIterPtr_(nullptr),
608 boundaryField_(*this, gf.boundaryField_)
611 <<
"Copy construct" <<
nl << this->info() <<
endl;
621 this->writeOpt(IOobjectOption::NO_WRITE);
625 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
631 Internal(tgf.constCast(), tgf.movable()),
634 fieldPrevIterPtr_(nullptr),
635 boundaryField_(*this, tgf().boundaryField_)
638 <<
"Constructing from tmp" <<
nl << this->info() <<
endl;
640 this->writeOpt(IOobjectOption::NO_WRITE);
646 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
656 fieldPrevIterPtr_(nullptr),
657 boundaryField_(*this, gf.boundaryField_)
660 <<
"Copy construct, resetting IO params" <<
nl 661 << this->info() <<
endl;
674 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
678 const tmp<GeometricField<Type, PatchField, GeoMesh>>& tgf
681 Internal(
io, tgf.constCast(), tgf.movable()),
684 fieldPrevIterPtr_(nullptr),
685 boundaryField_(*this, tgf().boundaryField_)
688 <<
"Constructing from tmp resetting IO params" <<
nl 689 << this->info() <<
endl;
697 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
704 Internal(newName, gf),
707 fieldPrevIterPtr_(nullptr),
708 boundaryField_(*this, gf.boundaryField_)
711 <<
"Copy construct, resetting name" <<
nl 712 << this->info() <<
endl;
725 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
729 const tmp<GeometricField<Type, PatchField, GeoMesh>>& tgf
732 Internal(newName, tgf.constCast(), tgf.movable()),
735 fieldPrevIterPtr_(nullptr),
736 boundaryField_(*this, tgf().boundaryField_)
739 <<
"Constructing from tmp resetting name" <<
nl 740 << this->info() <<
endl;
746 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
751 const word& patchFieldType
757 fieldPrevIterPtr_(nullptr),
758 boundaryField_(this->
mesh().
boundary(), *this, patchFieldType)
761 <<
"Copy construct, resetting IO params" <<
nl 762 << this->info() <<
endl;
764 boundaryField_ == gf.boundaryField_;
777 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
781 const GeometricField<Type, PatchField, GeoMesh>& gf,
789 fieldPrevIterPtr_(nullptr),
799 <<
"Copy construct, resetting IO params and patch types" <<
nl 800 << this->info() <<
endl;
802 boundaryField_ == gf.boundaryField_;
806 field0Ptr_ =
new GeometricField<Type, PatchField, GeoMesh>
815 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
819 const GeometricField<Type, PatchField, GeoMesh>& gf,
821 const word& patchFieldType
827 fieldPrevIterPtr_(nullptr),
828 boundaryField_(*this, gf.boundaryField_,
patchIDs, patchFieldType)
831 <<
"Copy construct, resetting IO params and setting patchFieldType " 832 <<
"for patchIDs" <<
nl 833 << this->info() <<
endl;
837 field0Ptr_ =
new GeometricField<Type, PatchField, GeoMesh>
846 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
850 const tmp<GeometricField<Type, PatchField, GeoMesh>>& tgf,
855 Internal(
io, tgf.constCast(), tgf.movable()),
858 fieldPrevIterPtr_(nullptr),
868 <<
"Constructing from tmp resetting IO params and patch types" <<
nl 869 << this->info() <<
endl;
871 boundaryField_ == tgf().boundaryField_;
877 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
887 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
912 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
917 const bool updateAccessTime
920 if (updateAccessTime)
929 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
934 const bool updateAccessTime
937 if (updateAccessTime)
946 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
951 const bool updateAccessTime
954 if (updateAccessTime)
959 return boundaryField_;
963 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
967 return (field0Ptr_ ? (field0Ptr_->nOldTimes() + 1) : 0);
971 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
977 && timeIndex_ != this->time().
timeIndex()
978 && !this->
name().ends_with(
"_0")
982 timeIndex_ = this->time().timeIndex();
990 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
995 field0Ptr_->storeOldTime();
998 <<
"Storing old time field for field" <<
nl << this->info() <<
endl;
1000 *field0Ptr_ == *
this;
1001 field0Ptr_->timeIndex_ = timeIndex_;
1003 if (field0Ptr_->field0Ptr_)
1005 field0Ptr_->writeOpt(this->writeOpt());
1011 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1017 field0Ptr_ =
new GeometricField<Type, PatchField, GeoMesh>
1021 this->
name() +
"_0",
1024 IOobjectOption::NO_READ,
1025 IOobjectOption::NO_WRITE,
1026 this->registerObject()
1034 <<
"created old time field " << field0Ptr_->info() <<
endl;
1038 error::printStack(
Info);
1051 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1062 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1065 if (!fieldPrevIterPtr_)
1068 <<
"Allocating previous iteration field" <<
nl 1069 << this->info() <<
endl;
1073 this->
name() +
"PrevIter",
1079 *fieldPrevIterPtr_ == *
this;
1084 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1088 if (!fieldPrevIterPtr_)
1091 <<
"previous iteration field" <<
endl << this->info() <<
endl 1093 <<
" Use field.storePrevIter() at start of iteration." 1097 return *fieldPrevIterPtr_;
1101 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1109 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1113 this->setUpToDate();
1115 boundaryField_.evaluate();
1119 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1123 this->setUpToDate();
1125 boundaryField_.evaluateLocal();
1129 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1135 bool needRef =
true;
1137 for (
const auto& pf : boundaryField_)
1139 if (pf.fixesValue())
1150 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1160 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1165 if (this->
mesh().data().isFinalIteration())
1170 scalar relaxCoeff = 1;
1172 if (this->
mesh().relaxField(
name, relaxCoeff))
1179 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1187 return this->
name() +
"Final";
1190 return this->
name();
1194 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1207 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1218 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1224 this->
name() +
".T()",
1229 Foam::T(tresult.ref().primitiveFieldRef(), primitiveField());
1230 Foam::T(tresult.ref().boundaryFieldRef(), boundaryField());
1236 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1258 Foam::component(tresult.ref().primitiveFieldRef(), primitiveField(), d);
1259 Foam::component(tresult.ref().boundaryFieldRef(), boundaryField(), d);
1265 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1269 const GeometricField
1271 typename GeometricField<Type, PatchField, GeoMesh>::cmptType,
1282 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1294 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1305 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1316 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1326 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1332 this->clamp_max(
upper.value());
1336 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1348 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1359 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1366 this->clamp_range(
lower.value(),
upper.value());
1370 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1376 this->clamp_range(
range.value());
1380 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1388 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1398 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1418 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1424 const auto& gf = tgf();
1454 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1460 internalFieldRef() = dt;
1465 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1471 const auto& gf = tgf();
1484 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1490 internalFieldRef() = dt;
1495 #define COMPUTED_ASSIGNMENT(TYPE, op) \ 1497 template<class Type, template<class> class PatchField, class GeoMesh> \ 1498 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \ 1500 const GeometricField<TYPE, PatchField, GeoMesh>& gf \ 1503 checkField(*this, gf, #op); \ 1505 internalFieldRef() op gf.internalField(); \ 1506 boundaryFieldRef() op gf.boundaryField(); \ 1509 template<class Type, template<class> class PatchField, class GeoMesh> \ 1510 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \ 1512 const tmp<GeometricField<TYPE, PatchField, GeoMesh>>& tgf \ 1515 operator op(tgf()); \ 1519 template<class Type, template<class> class PatchField, class GeoMesh> \ 1520 void Foam::GeometricField<Type, PatchField, GeoMesh>::operator op \ 1522 const dimensioned<TYPE>& dt \ 1525 internalFieldRef() op dt; \ 1526 boundaryFieldRef() op dt.value(); \ 1534 #undef COMPUTED_ASSIGNMENT 1539 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
1555 template<
class Type,
template<
class>
class PatchField,
class GeoMesh>
void clamp_min(const Type &lower)
Impose lower (floor) clamp on the field values (in-place)
const labelList patchIDs(pbm.indices(polyPatchNames, true))
void clamp_range(const dimensioned< MinMax< Type >> &range)
Clamp field values (in-place) to the specified range.
#define COMPUTED_ASSIGNMENT(TYPE, op)
const Type & value() const noexcept
Return const reference to value.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
const Internal::FieldType & primitiveField() const noexcept
Return a const-reference to the internal field values.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Y [inertIndex] clamp_min(0)
void clearOldTimes()
Remove old-time and prev-iter fields.
Field< Type >::cmptType cmptType
Component type of the field elements.
const Internal & internalField() const noexcept
Return a const-reference to the dimensioned internal field.
constexpr char nl
The newline '\n' character (0x0a)
A min/max value pair with additional methods. In addition to conveniently storing values...
GeometricField(const IOobject &io, const Mesh &mesh, const dimensionSet &dims, const word &patchFieldType=PatchField< Type >::calculatedType())
Construct given IOobject, mesh, dimensions and patch type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
label nOldTimes() const noexcept
The number of old time fields stored.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
orientedType oriented() const noexcept
Return oriented type.
string upper(const std::string &s)
Return string copy transformed with std::toupper on each character.
Generic GeometricField class.
Generic dimensioned Type class.
const dimensionSet dimless
Dimensionless.
bool writeData(Ostream &os) const
The writeData function (required by regIOobject), calls operator<<.
word select(bool final) const
Select the final iteration parameters if final is true by returning the field name + "Final" otherwis...
conserve primitiveFieldRef()+
MinMax< label > minMax(const labelHashSet &set)
Find the min/max values of labelHashSet.
#define forAll(list, i)
Loop across all elements in list.
bool needReference() const
Does the field need a reference level for solution.
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
const T & min() const noexcept
The min value (first)
bool returnReduceAnd(const bool value, const label comm=UPstream::worldComm)
Perform logical (and) MPI Allreduce on a copy. Uses UPstream::reduceAnd.
void storeOldTimes() const
Store the old-time fields.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
Generic templated field type.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const NameMatchPredicate &selectedFields, DynamicList< regIOobject *> &storedObjects)
Read the selected GeometricFields of the templated type and store on the objectRegistry.
A class for handling words, derived from Foam::string.
#define DebugInFunction
Report an information message using Foam::Info.
void storeOldTime() const
Store the old-time field.
void writeMinMax(Ostream &os) const
Helper function to write the min and max to an Ostream.
errorManip< error > abort(error &err)
Us boundaryFieldRef().evaluateCoupled< coupledFaPatch >()
void negate()
Negate the field inplace. See notes in Field.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Internal & internalFieldRef(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
int debug
Static debugging option.
void normalise()
Normalise the field inplace. See notes in Field.
const GeometricField< Type, PatchField, GeoMesh > & prevIter() const
Return previous iteration field.
OBJstream os(runTime.globalPath()/outputName)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual ~GeometricField()
Destructor.
Generic GeometricBoundaryField class.
void clamp_max(const Type &upper)
Impose upper (ceiling) clamp on the field values (in-place)
List< word > wordList
List of word.
Template functions to aid in the implementation of demand driven data.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field values.
tmp< GeometricField< Type, PatchField, GeoMesh > > clone() const
Clone.
string lower(const std::string &s)
Return string copy transformed with std::tolower on each character.
#define WarningInFunction
Report a warning using Foam::Warning.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
meshDefDict readIfPresent("polyMeshPatches", polyPatchNames)
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
void correctBoundaryConditions()
Correct boundary field.
messageStream Info
Information stream (stdout output on master, null elsewhere)
void storePrevIter() const
Store the field as the previous iteration value.
Generic mesh wrapper used by volMesh, surfaceMesh, pointMesh etc.
List< label > labelList
A List of labels.
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, IOobject::NO_REGISTER)
A class for managing temporary objects.
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
void deleteDemandDrivenData(DataPtr &dataPtr)
Defines the attributes of an object for which implicit objectRegistry management is supported...
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
const Boundary & boundaryField() const noexcept
Return const-reference to the boundary field.
void relax()
Relax field (for steady-state solution).
void correctLocalBoundaryConditions()
Correct boundary conditions after a purely local operation. Is.
const dimensionSet & dimensions() const noexcept
Return dimensions.
#define checkField(fld1, fld2, op)
#define InfoInFunction
Report an information message using Foam::Info.
void replace(const direction d, const GeometricField< cmptType, PatchField, GeoMesh > &gcf)
Replace specified field component with content from another field.