mappedVelocityFluxFixedValueFvPatchField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 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  << "Patch type '" << p.type()
64  << "' not type '" << mappedPatchBase::typeName << "'"
65  << " for patch " << p.name()
66  << " of field " << internalField().name()
67  << " in file " << internalField().objectPath()
68  << exit(FatalError);
69  }
70 }
71 
72 
75 (
76  const fvPatch& p,
78  const dictionary& dict
79 )
80 :
81  fixedValueFvPatchVectorField(p, iF, dict),
82  phiName_(dict.lookupOrDefault<word>("phi", "phi"))
83 {
84  if (!isA<mappedPatchBase>(this->patch().patch()))
85  {
87  << "Patch type '" << p.type()
88  << "' not type '" << mappedPatchBase::typeName << "'"
89  << " for patch " << p.name()
90  << " of field " << internalField().name()
91  << " in file " << internalField().objectPath()
92  << exit(FatalError);
93  }
94 
95  const mappedPatchBase& mpp = refCast<const mappedPatchBase>
96  (
97  this->patch().patch()
98  );
99  if (mpp.mode() == mappedPolyPatch::NEARESTCELL)
100  {
102  << "Patch " << p.name()
103  << " of type '" << p.type()
104  << "' can not be used in 'nearestCell' mode"
105  << " of field " << internalField().name()
106  << " in file " << internalField().objectPath()
107  << exit(FatalError);
108  }
109 }
110 
111 
114 (
116 )
117 :
118  fixedValueFvPatchVectorField(ptf),
119  phiName_(ptf.phiName_)
120 {}
121 
122 
125 (
128 )
129 :
130  fixedValueFvPatchVectorField(ptf, iF),
131  phiName_(ptf.phiName_)
132 {}
133 
134 
135 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
136 
138 {
139  if (updated())
140  {
141  return;
142  }
143 
144  // Since we're inside initEvaluate/evaluate there might be processor
145  // comms underway. Change the tag we use.
146  int oldTag = UPstream::msgType();
147  UPstream::msgType() = oldTag+1;
148 
149  // Get the mappedPatchBase
150  const mappedPatchBase& mpp = refCast<const mappedPatchBase>
151  (
152  mappedVelocityFluxFixedValueFvPatchField::patch().patch()
153  );
154  const fvMesh& nbrMesh = refCast<const fvMesh>(mpp.sampleMesh());
155  const word& fieldName = internalField().name();
156  const volVectorField& UField =
157  nbrMesh.lookupObject<volVectorField>(fieldName);
158 
159  surfaceScalarField& phiField =
160  nbrMesh.lookupObjectRef<surfaceScalarField>(phiName_);
161 
162  vectorField newUValues;
163  scalarField newPhiValues;
164 
165  switch (mpp.mode())
166  {
168  {
169  vectorField allUValues(nbrMesh.nFaces(), Zero);
170  scalarField allPhiValues(nbrMesh.nFaces(), 0.0);
171 
172  forAll(UField.boundaryField(), patchi)
173  {
174  const fvPatchVectorField& Upf = UField.boundaryField()[patchi];
175  const scalarField& phipf = phiField.boundaryField()[patchi];
176 
177  label faceStart = Upf.patch().start();
178 
179  forAll(Upf, facei)
180  {
181  allUValues[faceStart + facei] = Upf[facei];
182  allPhiValues[faceStart + facei] = phipf[facei];
183  }
184  }
185 
186  mpp.distribute(allUValues);
187  newUValues.transfer(allUValues);
188 
189  mpp.distribute(allPhiValues);
190  newPhiValues.transfer(allPhiValues);
191 
192  break;
193  }
196  {
197  const label nbrPatchID =
198  nbrMesh.boundaryMesh().findPatchID(mpp.samplePatch());
199 
200  newUValues = UField.boundaryField()[nbrPatchID];
201  mpp.distribute(newUValues);
202 
203  newPhiValues = phiField.boundaryField()[nbrPatchID];
204  mpp.distribute(newPhiValues);
205 
206  break;
207  }
208  default:
209  {
211  << "patch can only be used in NEARESTPATCHFACE, "
212  << "NEARESTPATCHFACEAMI or NEARESTFACE mode" << nl
213  << abort(FatalError);
214  }
215  }
216 
217  operator==(newUValues);
218  phiField.boundaryFieldRef()[patch().index()] == newPhiValues;
219 
220  // Restore tag
221  UPstream::msgType() = oldTag;
222 
223  fixedValueFvPatchVectorField::updateCoeffs();
224 }
225 
226 
228 (
229  Ostream& os
230 ) const
231 {
233  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
234  writeEntry(os, "value", *this);
235 }
236 
237 
238 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
239 
240 namespace Foam
241 {
243  (
246  );
247 }
248 
249 
250 // ************************************************************************* //
Foam::surfaceFields.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual label start() const
Return start label of this patch in the polyMesh face list.
Definition: fvPatch.H:149
const polyMesh & sampleMesh() const
Get the region mesh.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
mappedVelocityFluxFixedValueFvPatchField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
This boundary condition maps the velocity and flux from a neighbour patch to this patch...
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:476
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:280
Macros for easy insertion into run-time selection tables.
const sampleMode & mode() const
What to sample.
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
static const zero Zero
Definition: zero.H:97
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:325
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:265
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
const word & samplePatch() const
Patch (only if NEARESTPATCHFACE)
virtual const word & name() const
Return name.
Definition: fvPatch.H:143
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:342
Namespace for OpenFOAM.