Loading...
Searching...
No Matches
uniformFixedValuePointPatchField.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2019-2024 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
30#include "SubField.H"
31#include "polyPatch.H"
32
33// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
34
35// Alternative
36// {
37// if (const auto* fpp = isA<facePointPatch>(p).patch())
38// {
39// return fpp->patch();
40// }
41// else
42// {
43// return nullptr;
44// }
45// }
46
47template<class Type>
48const Foam::polyPatch*
49Foam::uniformFixedValuePointPatchField<Type>::getPolyPatch(const pointPatch& p)
50{
51 const polyMesh& mesh = p.boundaryMesh().mesh()();
52 return mesh.boundaryMesh().cfindPatch(p.name());
53}
54
55
56// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * //
57
58template<class Type>
61(
62 const pointPatch& p,
65:
67{}
68
69
70template<class Type>
73(
74 const pointPatch& p,
76 const dictionary& dict
77)
78:
80{
81 if (const polyPatch* pp = this->getPolyPatch(this->patch()))
82 {
83 refValueFunc_ = PatchFunction1<Type>::New
84 (
85 *pp,
86 "uniformValue",
87 dict,
88 false // point values (faceValues = false)
89 );
90 }
91 // Fallback
92 refPointValueFunc_ = Function1<Type>::New
93 (
94 "uniformValue",
95 dict,
96 &this->internalField().db()
97 );
98
99 if (!this->readValueEntry(dict))
100 {
101 // Ensure field has reasonable initial values
102 this->extrapolateInternal();
103
104 // Evaluate to assign a value
105 this->evaluate();
106 }
107}
108
109
110template<class Type>
113(
114 const uniformFixedValuePointPatchField<Type>& ptf,
115 const pointPatch& p,
116 const DimensionedField<Type, pointMesh>& iF,
117 const pointPatchFieldMapper& mapper
118)
119:
120 fixedValuePointPatchField<Type>(ptf, p, iF, mapper)
121{
122 if (const polyPatch* pp = this->getPolyPatch(this->patch()))
123 {
124 refValueFunc_ = ptf.refValueFunc_.clone(*pp);
125 }
126 // Fallback
127 refPointValueFunc_ = ptf.refPointValueFunc_.clone();
128
129 if (mapper.direct() && !mapper.hasUnmapped())
130 {
131 // Use mapping instead of re-evaluation
132 this->map(ptf, mapper);
133 }
134 else
135 {
136 // Evaluate since value not mapped
137 this->evaluate();
138 }
139}
140
141
142template<class Type>
145(
148)
149:
150 fixedValuePointPatchField<Type>(pfld, iF)
151{
152 if (const polyPatch* pp = this->getPolyPatch(this->patch()))
153 {
154 refValueFunc_ = pfld.refValueFunc_.clone(*pp);
155 }
156 // Fallback
157 refPointValueFunc_ = pfld.refPointValueFunc_.clone();
158}
159
160
161// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
162
163template<class Type>
165(
166 const pointPatchFieldMapper& mapper
167)
168{
170
171 bool canEvaluate(false);
172
173 if (refValueFunc_)
174 {
175 refValueFunc_().autoMap(mapper);
176
177 // If mapper is not dependent on time we're ok to evaluate
178 if (refValueFunc_->constant())
179 {
180 canEvaluate = true;
181 }
182 }
183 if (refPointValueFunc_)
184 {
185 // If mapper is not dependent on time we're ok to evaluate
186 if (refPointValueFunc_->constant())
187 {
188 canEvaluate = true;
189 }
190 }
191
192 if (canEvaluate)
194 this->evaluate();
195 }
196}
197
198
199template<class Type>
201(
202 const pointPatchField<Type>& ptf,
203 const labelList& addr
204)
205{
207
209
210 if (refValueFunc_ && tiptf.refValueFunc_)
212 refValueFunc_().rmap(tiptf.refValueFunc_(), addr);
213 }
214}
215
216
217template<class Type>
219{
220 if (this->updated())
221 {
222 return;
223 }
224 const scalar t = this->db().time().timeOutputValue();
225
226 if (refValueFunc_)
227 {
228 valuePointPatchField<Type>::operator=(refValueFunc_->value(t));
229 }
230 else
231 {
232 valuePointPatchField<Type>::operator=(refPointValueFunc_->value(t));
233 }
235}
236
237
238template<class Type>
240write(Ostream& os) const
241{
242 // Note: write value
244 if (refValueFunc_)
245 {
246 refValueFunc_->writeData(os);
247 }
248 else if (refPointValueFunc_)
249 {
250 refPointValueFunc_->writeData(os);
251 }
252}
253
254
255// ************************************************************************* //
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual bool direct() const =0
Is it a direct (non-interpolating) mapper?
virtual bool hasUnmapped() const =0
Any unmapped values?
void autoMap(const FieldMapper &map, const bool applyFlip=true)
Map from self.
Definition Field.C:465
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition Field.C:528
A simple container of IOobject preferences. Can also be used for general handling of read/no-read/rea...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition Ostream.H:59
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition dictionary.H:133
A FixedValue boundary condition for pointField.
fixedValuePointPatchField(const pointPatch &, const DimensionedField< Type, pointMesh > &)
Construct from patch and internal field.
const pointPatch & patch() const noexcept
Return the patch.
const objectRegistry & db() const
The associated objectRegistry.
bool updated() const noexcept
True if the boundary condition has already been updated.
Foam::pointPatchFieldMapper.
Abstract base class for point-mesh patch fields.
const DimensionedField< Type, pointMesh > & internalField() const noexcept
Return const-reference to the dimensioned internal field.
Basic pointPatch represents a set of points from the mesh.
Definition pointPatch.H:67
A patch is a list of labels that address the faces in the global face list.
Definition polyPatch.H:73
Enables the specification of a uniform fixed value condition.
uniformFixedValuePointPatchField(const pointPatch &, const DimensionedField< Type, pointMesh > &)
Construct from patch and internal field.
virtual void autoMap(const pointPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void rmap(const pointPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void write(Ostream &) const
Write.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void operator=(const valuePointPatchField< Type > &)
bool readValueEntry(const dictionary &dict, IOobjectOption::readOption readOpt=IOobjectOption::LAZY_READ)
Read the "value" entry into *this.
void extrapolateInternal()
Assign the patch field from the internal field.
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::commsTypes::buffered)
Evaluate the patch field.
volScalarField & p
dynamicFvMesh & mesh
OBJstream os(runTime.globalPath()/outputName)
Type & refCast(U &obj)
A dynamic_cast (for references) to Type reference.
Definition typeInfo.H:172
List< label > labelList
A List of labels.
Definition List.H:62
dictionary dict