mappedVelocityFluxFixedValueFvPatchField.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 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 
27 #include "fvPatchFieldMapper.H"
28 #include "mappedPatchBase.H"
29 #include "volFields.H"
30 #include "surfaceFields.H"
32 
33 
34 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
35 
38 (
39  const fvPatch& p,
41 )
42 :
43  fixedValueFvPatchVectorField(p, iF),
44  phiName_("phi")
45 {}
46 
47 
50 (
52  const fvPatch& p,
54  const fvPatchFieldMapper& mapper
55 )
56 :
57  fixedValueFvPatchVectorField(ptf, p, iF, mapper),
58  phiName_(ptf.phiName_)
59 {
60  if (!isA<mappedPatchBase>(this->patch().patch()))
61  {
63  (
64  "mappedVelocityFluxFixedValueFvPatchField::"
65  "mappedVelocityFluxFixedValueFvPatchField"
66  "("
67  "const mappedVelocityFluxFixedValueFvPatchField&, "
68  "const fvPatch&, "
69  "const DimensionedField<vector, volMesh>&, "
70  "const fvPatchFieldMapper&"
71  ")"
72  ) << "Patch type '" << p.type()
73  << "' not type '" << mappedPatchBase::typeName << "'"
74  << " for patch " << p.name()
75  << " of field " << dimensionedInternalField().name()
76  << " in file " << dimensionedInternalField().objectPath()
77  << exit(FatalError);
78  }
79 }
80 
81 
84 (
85  const fvPatch& p,
87  const dictionary& dict
88 )
89 :
90  fixedValueFvPatchVectorField(p, iF, dict),
91  phiName_(dict.lookupOrDefault<word>("phi", "phi"))
92 {
93  if (!isA<mappedPatchBase>(this->patch().patch()))
94  {
96  (
97  "mappedVelocityFluxFixedValueFvPatchField::"
98  "mappedVelocityFluxFixedValueFvPatchField"
99  "("
100  "const fvPatch&, "
101  "const DimensionedField<vector, volMesh>&, "
102  "const dictionary&"
103  ")"
104  ) << "Patch type '" << p.type()
105  << "' not type '" << mappedPatchBase::typeName << "'"
106  << " for patch " << p.name()
107  << " of field " << dimensionedInternalField().name()
108  << " in file " << dimensionedInternalField().objectPath()
109  << exit(FatalError);
110  }
111 
112  const mappedPatchBase& mpp = refCast<const mappedPatchBase>
113  (
114  this->patch().patch()
115  );
116  if (mpp.mode() == mappedPolyPatch::NEARESTCELL)
117  {
119  (
120  "mappedVelocityFluxFixedValueFvPatchField::"
121  "mappedVelocityFluxFixedValueFvPatchField"
122  "("
123  "const fvPatch&, "
124  "const DimensionedField<vector, volMesh>&, "
125  "const dictionary&"
126  ")"
127  ) << "Patch " << p.name()
128  << " of type '" << p.type()
129  << "' can not be used in 'nearestCell' mode"
130  << " of field " << dimensionedInternalField().name()
131  << " in file " << dimensionedInternalField().objectPath()
132  << exit(FatalError);
133  }
134 }
135 
136 
139 (
141 )
142 :
143  fixedValueFvPatchVectorField(ptf),
144  phiName_(ptf.phiName_)
145 {}
146 
147 
150 (
153 )
154 :
155  fixedValueFvPatchVectorField(ptf, iF),
156  phiName_(ptf.phiName_)
157 {}
158 
159 
160 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
161 
163 {
164  if (updated())
165  {
166  return;
167  }
168 
169  // Since we're inside initEvaluate/evaluate there might be processor
170  // comms underway. Change the tag we use.
171  int oldTag = UPstream::msgType();
172  UPstream::msgType() = oldTag+1;
173 
174  // Get the mappedPatchBase
175  const mappedPatchBase& mpp = refCast<const mappedPatchBase>
176  (
177  mappedVelocityFluxFixedValueFvPatchField::patch().patch()
178  );
179  const fvMesh& nbrMesh = refCast<const fvMesh>(mpp.sampleMesh());
180  const word& fieldName = dimensionedInternalField().name();
181  const volVectorField& UField =
182  nbrMesh.lookupObject<volVectorField>(fieldName);
183 
184  surfaceScalarField& phiField = const_cast<surfaceScalarField&>
185  (
186  nbrMesh.lookupObject<surfaceScalarField>(phiName_)
187  );
188 
189  vectorField newUValues;
190  scalarField newPhiValues;
191 
192  switch (mpp.mode())
193  {
195  {
196  vectorField allUValues(nbrMesh.nFaces(), vector::zero);
197  scalarField allPhiValues(nbrMesh.nFaces(), 0.0);
198 
199  forAll(UField.boundaryField(), patchI)
200  {
201  const fvPatchVectorField& Upf = UField.boundaryField()[patchI];
202  const scalarField& phipf = phiField.boundaryField()[patchI];
203 
204  label faceStart = Upf.patch().start();
205 
206  forAll(Upf, faceI)
207  {
208  allUValues[faceStart + faceI] = Upf[faceI];
209  allPhiValues[faceStart + faceI] = phipf[faceI];
210  }
211  }
212 
213  mpp.distribute(allUValues);
214  newUValues.transfer(allUValues);
215 
216  mpp.distribute(allPhiValues);
217  newPhiValues.transfer(allPhiValues);
218 
219  break;
220  }
223  {
224  const label nbrPatchID =
225  nbrMesh.boundaryMesh().findPatchID(mpp.samplePatch());
226 
227  newUValues = UField.boundaryField()[nbrPatchID];
228  mpp.distribute(newUValues);
229 
230  newPhiValues = phiField.boundaryField()[nbrPatchID];
231  mpp.distribute(newPhiValues);
232 
233  break;
234  }
235  default:
236  {
238  (
239  "mappedVelocityFluxFixedValueFvPatchField::"
240  "updateCoeffs()"
241  ) << "patch can only be used in NEARESTPATCHFACE, "
242  << "NEARESTPATCHFACEAMI or NEARESTFACE mode" << nl
243  << abort(FatalError);
244  }
245  }
246 
247  operator==(newUValues);
248  phiField.boundaryField()[patch().index()] == newPhiValues;
249 
250  // Restore tag
251  UPstream::msgType() = oldTag;
252 
253  fixedValueFvPatchVectorField::updateCoeffs();
254 }
255 
256 
258 (
259  Ostream& os
260 ) const
261 {
263  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
264  this->writeEntry("value", os);
265 }
266 
267 
268 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
269 
270 namespace Foam
271 {
273  (
276  );
277 }
278 
279 
280 // ************************************************************************* //
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Foam::surfaceFields.
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...
const polyMesh & sampleMesh() const
Get the region mesh.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:390
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
const word & samplePatch() const
Patch (only if NEARESTPATCHFACE)
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
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
Foam::fvPatchFieldMapper.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:354
Namespace for OpenFOAM.
const word & name() const
Return name.
Definition: fvPatch.H:149
static const char nl
Definition: Ostream.H:260
#define forAll(list, i)
Definition: UList.H:421
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:451
Macros for easy insertion into run-time selection tables.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
rDeltaT dimensionedInternalField()
This boundary condition maps the velocity and flux from a neighbour patch to this patch...
const sampleMode & mode() const
What to sample.
error FatalError
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
static const Vector zero
Definition: Vector.H:80
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
mappedVelocityFluxFixedValueFvPatchField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:300
label start() const
Return start label of this patch in the polyMesh face list.
Definition: fvPatch.H:155
virtual void updateCoeffs()
Update the coefficients associated with the patch field.