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-2021 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 (
47  const fvPatch& p,
49  const dictionary& dict
50 )
51 :
53 {}
54 
55 
56 template<class Type>
59 (
61  const fvPatch& p,
63  const fvPatchFieldMapper& mapper
64 )
65 :
66  mappedFixedValueFvPatchField<Type>(ptf, p, iF, mapper)
67 {}
68 
69 
70 template<class Type>
73 (
76 )
77 :
79 {}
80 
81 
82 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
83 
84 template<class Type>
86 {
88 
89  if (this->updated())
90  {
91  return;
92  }
93 
94  // Since we're inside initEvaluate/evaluate there might be processor
95  // comms underway. Change the tag we use.
96  int oldTag = UPstream::msgType();
97  UPstream::msgType() = oldTag + 1;
98 
99  // Retrieve the neighbour values and assign to this patch boundary field
101 
102  // Get the coupling information from the mappedPatchBase
103  const mappedPatchBase& mpp =
104  refCast<const mappedPatchBase>(this->patch().patch());
105  const fvMesh& nbrMesh = refCast<const fvMesh>(mpp.sampleMesh());
106 
107  Field<Type> nbrIntFld;
108 
109  switch (mpp.mode())
110  {
111  case mappedPatchBase::NEARESTCELL:
112  {
114  << "Cannot apply "
115  << mappedPatchBase::sampleModeNames_
116  [
117  mappedPatchBase::NEARESTCELL
118  ]
119  << " mapping mode for patch " << this->patch().name()
120  << exit(FatalError);
121 
122  break;
123  }
124  case mappedPatchBase::NEARESTPATCHFACE:
125  case mappedPatchBase::NEARESTPATCHFACEAMI:
126  {
127  const label samplePatchi = mpp.samplePolyPatch().index();
128  const fvPatchField<Type>& nbrPatchField =
129  this->sampleField().boundaryField()[samplePatchi];
130  nbrIntFld = nbrPatchField.patchInternalField();
131  mpp.distribute(nbrIntFld);
132 
133  break;
134  }
135  case mappedPatchBase::NEARESTFACE:
136  {
137  Field<Type> allValues(nbrMesh.nFaces(), Zero);
138 
139  const FieldType& nbrField = this->sampleField();
140 
141  forAll(nbrField.boundaryField(), patchi)
142  {
143  const fvPatchField<Type>& pf = nbrField.boundaryField()[patchi];
144  const Field<Type> pif(pf.patchInternalField());
145 
146  label faceStart = pf.patch().start();
147 
148  forAll(pf, facei)
149  {
150  allValues[faceStart++] = pif[facei];
151  }
152  }
153 
154  mpp.distribute(allValues);
155  nbrIntFld.transfer(allValues);
156 
157  break;
158  }
159  default:
160  {
162  << "Unknown sampling mode: " << mpp.mode()
163  << abort(FatalError);
164  }
165  }
166 
167  // Restore tag
168  UPstream::msgType() = oldTag;
169 
170  // Assign to (this) patch internal field its neighbour values
171  Field<Type>& intFld = const_cast<Field<Type>&>(this->primitiveField());
172  UIndirectList<Type>(intFld, this->patch().faceCells()) = nbrIntFld;
173 }
174 
175 
176 template<class Type>
178 (
179  Ostream& os
180 ) const
181 {
183 }
184 
185 
186 // ************************************************************************* //
dictionary dict
#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:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
const polyMesh & sampleMesh() const
Get the region mesh.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
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.
Pre-declare SubField and related Field type.
Definition: Field.H:56
Foam::fvPatchFieldMapper.
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
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:337
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:54
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:174
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...