mappedFixedInternalValueFvPatchField.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-2018 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 "UIndirectList.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class Type>
34 (
35  const fvPatch& p,
37 )
38 :
40 {}
41 
42 
43 template<class Type>
46 (
48  const fvPatch& p,
50  const fvPatchFieldMapper& mapper
51 )
52 :
53  mappedFixedValueFvPatchField<Type>(ptf, p, iF, mapper)
54 {}
55 
56 
57 template<class Type>
60 (
61  const fvPatch& p,
63  const dictionary& dict
64 )
65 :
67 {}
68 
69 
70 template<class Type>
73 (
75 )
76 :
78 {}
79 
80 
81 template<class Type>
84 (
87 )
88 :
90 {}
91 
92 
93 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
94 
95 template<class Type>
97 {
99 
100  if (this->updated())
101  {
102  return;
103  }
104 
105  // Since we're inside initEvaluate/evaluate there might be processor
106  // comms underway. Change the tag we use.
107  int oldTag = UPstream::msgType();
108  UPstream::msgType() = oldTag + 1;
109 
110  // Retrieve the neighbour values and assign to this patch boundary field
112 
113  // Get the coupling information from the mappedPatchBase
114  const mappedPatchBase& mpp =
115  refCast<const mappedPatchBase>(this->patch().patch());
116  const fvMesh& nbrMesh = refCast<const fvMesh>(mpp.sampleMesh());
117 
118  Field<Type> nbrIntFld;
119 
120  switch (mpp.mode())
121  {
122  case mappedPatchBase::NEARESTCELL:
123  {
125  << "Cannot apply "
126  << mappedPatchBase::sampleModeNames_
127  [
128  mappedPatchBase::NEARESTCELL
129  ]
130  << " mapping mode for patch " << this->patch().name()
131  << exit(FatalError);
132 
133  break;
134  }
135  case mappedPatchBase::NEARESTPATCHFACE:
136  case mappedPatchBase::NEARESTPATCHFACEAMI:
137  {
138  const label samplePatchi = mpp.samplePolyPatch().index();
139  const fvPatchField<Type>& nbrPatchField =
140  this->sampleField().boundaryField()[samplePatchi];
141  nbrIntFld = nbrPatchField.patchInternalField();
142  mpp.distribute(nbrIntFld);
143 
144  break;
145  }
146  case mappedPatchBase::NEARESTFACE:
147  {
148  Field<Type> allValues(nbrMesh.nFaces(), Zero);
149 
150  const FieldType& nbrField = this->sampleField();
151 
152  forAll(nbrField.boundaryField(), patchi)
153  {
154  const fvPatchField<Type>& pf = nbrField.boundaryField()[patchi];
155  const Field<Type> pif(pf.patchInternalField());
156 
157  label faceStart = pf.patch().start();
158 
159  forAll(pf, facei)
160  {
161  allValues[faceStart++] = pif[facei];
162  }
163  }
164 
165  mpp.distribute(allValues);
166  nbrIntFld.transfer(allValues);
167 
168  break;
169  }
170  default:
171  {
173  << "Unknown sampling mode: " << mpp.mode()
174  << abort(FatalError);
175  }
176  }
177 
178  // Restore tag
179  UPstream::msgType() = oldTag;
180 
181  // Assign to (this) patch internal field its neighbour values
182  Field<Type>& intFld = const_cast<Field<Type>&>(this->primitiveField());
183  UIndirectList<Type>(intFld, this->patch().faceCells()) = nbrIntFld;
184 }
185 
186 
187 template<class Type>
189 (
190  Ostream& os
191 ) const
192 {
194 }
195 
196 
197 // ************************************************************************* //
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const polyMesh & sampleMesh() const
Get the region mesh.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
This boundary condition maps the value at a set of cells or patch faces back to *this.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
Generic GeometricField class.
const sampleMode & mode() const
What to sample.
void write(Ostream &, const label, const dictionary &)
Write with dictionary lookup.
Pre-declare SubField and related Field type.
Definition: Field.H:57
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:340
mappedFixedInternalValueFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:224
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
label patchi
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A List with indirect addressing.
Definition: fvMatrix.H:106
label index() const
Return the index of this patch in the boundaryMesh.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
volScalarField & p
const polyPatch & samplePolyPatch() const
Get the patch on the region.
This boundary condition maps the boundary and internal values of a neighbour patch field to the bound...