Go to the documentation of this file.
38 const dictionary&
dict,
44 return dict.get<Type>(
"average");
58 const word& fieldName,
59 const bool setAverage,
61 const word& interpolationScheme
65 patchField_(patchField),
66 fieldName_(fieldName),
67 setAverage_(setAverage),
69 interpolationScheme_(interpolationScheme)
82 patchField_(patchField),
85 dict.template getOrDefault<word>
88 patchField_.internalField().name()
91 setAverage_(
dict.getOrDefault(
"setAverage",
false)),
92 average_(getAverage(
dict, setAverage_)),
95 if (mapper_.mode() == mappedPatchBase::NEARESTCELL)
97 dict.readEntry(
"interpolationScheme", interpolationScheme_);
110 patchField_(patchField),
111 fieldName_(patchField_.internalField().name()),
142 patchField_(patchField),
158 if (mapper_.sameRegion())
160 if (fieldName_ == patchField_.internalField().name())
164 dynamic_cast<const fieldType&
>
166 patchField_.internalField()
172 return thisMesh.template lookupObject<fieldType>(fieldName_);
176 const fvMesh& nbrMesh = refCast<const fvMesh>(mapper_.sampleMesh());
177 return nbrMesh.template lookupObject<fieldType>(fieldName_);
189 int oldTag = UPstream::msgType();
190 UPstream::msgType() = oldTag + 1;
193 const fvMesh& nbrMesh = refCast<const fvMesh>(mapper_.sampleMesh());
197 auto& newValues = tnewValues.ref();
199 switch (mapper_.mode())
201 case mappedPatchBase::NEARESTCELL:
223 interpolationScheme_,
227 const auto& interp = *interpolator;
234 newValues[celli] = interp.interpolate
244 newValues = sampleField();
251 case mappedPatchBase::NEARESTPATCHFACE:
252 case mappedPatchBase::NEARESTPATCHFACEAMI:
254 const label nbrPatchID =
260 <<
"Unable to find sample patch " << mapper_.samplePatch()
261 <<
" in region " << mapper_.sampleRegion()
262 <<
" for patch " << patchField_.patch().name() <<
nl
266 const fieldType& nbrField = sampleField();
268 newValues = nbrField.boundaryField()[nbrPatchID];
269 mapper_.distribute(newValues);
273 case mappedPatchBase::NEARESTFACE:
277 const fieldType& nbrField = sampleField();
281 label faceStart = pf.patch().start();
285 allValues[faceStart++] = pf[facei];
289 mapper_.distribute(allValues);
290 newValues.transfer(allValues);
297 <<
"Unknown sampling mode: " << mapper_.mode() <<
nl
305 gSum(patchField_.patch().magSf()*newValues)
306 /
gSum(patchField_.patch().magSf());
308 if (
mag(averagePsi)/
mag(average_) > 0.5)
310 newValues *=
mag(average_)/
mag(averagePsi);
314 newValues += (average_ - averagePsi);
319 UPstream::msgType() = oldTag;
336 os.
writeEntry(
"interpolationScheme", interpolationScheme_);
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
A class for handling words, derived from Foam::string.
mappedPatchFieldBase(const mappedPatchBase &mapper, const fvPatchField< Type > &patchField, const word &fieldName, const bool setAverage, const Type average, const word &interpolationScheme)
Construct from components.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const fvPatchField< Type > & patchField_
Underlying patch field.
label nFaces() const
Number of mesh faces.
virtual tmp< Field< Type > > mappedField() const
Map sampleField onto *this patch.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Type gSum(const FieldField< Field, Type > &f)
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
#define forAll(list, i)
Loop across all elements in list.
Uses the cell value for any point in the cell.
word fieldName_
Name of field to sample.
label nCells() const
Number of mesh cells.
tmp< GeometricField< Type, faPatchField, areaMesh > > average(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
Area-weighted average a edgeField creating a areaField.
Class containing processor-to-processor mapping information.
scalarField samples(nIntervals, Zero)
const polyMesh & mesh() const
Return the mesh reference.
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
void distribute(List< T > &fld, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Distribute data using default commsType.
Functionality for sampling fields using mappedPatchBase. Every call to mappedField() returns a sample...
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Mesh data needed to do the Finite Volume discretisation.
errorManip< error > abort(error &err)
word interpolationScheme_
Interpolation scheme to use for nearestcell mode.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
virtual void write(Ostream &os) const
Write.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const mappedPatchBase & mapper_
Mapping engine.
const GeometricField< Type, fvPatchField, volMesh > & sampleField() const
Field to sample. Either on my or nbr mesh.
Traits class for primitives.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
void reverseDistribute(const label constructSize, List< T > &, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Reverse distribute data using default commsType.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
const bool setAverage_
If true adjust the mapped field to maintain average value average_.