regionModelTemplates.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 
26 #include "mappedPatchFieldBase.H"
27 
28 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
29 
30 template<class Type>
33 (
34  const regionModel& nbrRegion,
35  const label regionPatchi,
36  const label nbrPatchi,
37  const Field<Type>& nbrField,
38  const bool flip
39 ) const
40 {
41  int oldTag = UPstream::msgType();
42  UPstream::msgType() = oldTag + 1;
43 
44  const AMIInterpolation& ami =
45  interRegionAMI(nbrRegion, regionPatchi, nbrPatchi, flip);
46 
47  tmp<Field<Type>> tresult(ami.interpolateToSource(nbrField));
48 
49  UPstream::msgType() = oldTag;
50 
51  return tresult;
52 }
53 
54 
55 template<class Type>
58 (
59  const regionModel& nbrRegion,
60  const word& fieldName,
61  const label regionPatchi,
62  const bool flip
63 ) const
64 {
66 
67  const fvMesh& nbrRegionMesh = nbrRegion.regionMesh();
68 
69  if (nbrRegionMesh.foundObject<fieldType>(fieldName))
70  {
71  const label nbrPatchi = nbrCoupledPatchID(nbrRegion, regionPatchi);
72 
73  int oldTag = UPstream::msgType();
74  UPstream::msgType() = oldTag + 1;
75 
76  const AMIInterpolation& ami =
77  interRegionAMI(nbrRegion, regionPatchi, nbrPatchi, flip);
78 
79  const fieldType& nbrField =
80  nbrRegionMesh.lookupObject<fieldType>(fieldName);
81 
82  const Field<Type>& nbrFieldp = nbrField.boundaryField()[nbrPatchi];
83 
84  tmp<Field<Type>> tresult(ami.interpolateToSource(nbrFieldp));
85 
86  UPstream::msgType() = oldTag;
87 
88  return tresult;
89  }
90  else
91  {
92  const polyPatch& p = regionMesh().boundaryMesh()[regionPatchi];
93 
94  return
96  (
97  new Field<Type>
98  (
99  p.size(),
100  Zero
101  )
102  );
103  }
104 }
105 
106 
107 template<class Type>
110 (
111  const regionModel& nbrRegion,
112  const word& fieldName,
113  const label regionPatchi,
114  const bool flip
115 ) const
116 {
118 
119  const fvMesh& nbrRegionMesh = nbrRegion.regionMesh();
120 
121  if (nbrRegionMesh.foundObject<fieldType>(fieldName))
122  {
123  const label nbrPatchi = nbrCoupledPatchID(nbrRegion, regionPatchi);
124 
125  int oldTag = UPstream::msgType();
126  UPstream::msgType() = oldTag + 1;
127 
128  const AMIInterpolation& ami =
129  interRegionAMI(nbrRegion, regionPatchi, nbrPatchi, flip);
130 
131  const fieldType& nbrField =
132  nbrRegionMesh.lookupObject<fieldType>(fieldName);
133 
134  const fvPatchField<Type>& nbrFieldp =
135  nbrField.boundaryField()[nbrPatchi];
136 
137  tmp<Field<Type>> tresult
138  (
139  ami.interpolateToSource(nbrFieldp.patchInternalField())
140  );
141 
142  UPstream::msgType() = oldTag;
143 
144  return tresult;
145  }
146  else
147  {
148  const polyPatch& p = regionMesh().boundaryMesh()[regionPatchi];
149 
150  return
152  (
153  new Field<Type>
154  (
155  p.size(),
156  Zero
157  )
158  );
159  }
160 }
161 
162 
163 template<class Type>
165 (
166  const label regionPatchi,
167  List<Type>& regionField
168 ) const
169 {
170  forAll(intCoupledPatchIDs_, i)
171  {
172  if (intCoupledPatchIDs_[i] == regionPatchi)
173  {
174  const mappedPatchBase& mpb =
175  refCast<const mappedPatchBase>
176  (
177  regionMesh().boundaryMesh()[regionPatchi]
178  );
179  mpb.reverseDistribute(regionField);
180  return;
181  }
182  }
183 
185  << "Region patch ID " << regionPatchi << " not found in region mesh"
186  << abort(FatalError);
187 }
188 
189 
190 template<class Type, class CombineOp>
192 (
193  const label regionPatchi,
194  List<Type>& regionField,
195  const CombineOp& cop
196 ) const
197 {
198  forAll(intCoupledPatchIDs_, i)
199  {
200  if (intCoupledPatchIDs_[i] == regionPatchi)
201  {
202  const mappedPatchBase& mpb =
203  refCast<const mappedPatchBase>
204  (
205  regionMesh().boundaryMesh()[regionPatchi]
206  );
207  mpb.reverseDistribute(regionField, cop);
208  return;
209  }
210  }
211 
213  << "Region patch ID " << regionPatchi << " not found in region mesh"
214  << abort(FatalError);
215 }
216 
217 
218 template<class Type>
220 (
221  const label regionPatchi,
222  List<Type>& primaryField
223 ) const
224 {
225  forAll(intCoupledPatchIDs_, i)
226  {
227  if (intCoupledPatchIDs_[i] == regionPatchi)
228  {
229  const mappedPatchBase& mpb =
230  refCast<const mappedPatchBase>
231  (
232  regionMesh().boundaryMesh()[regionPatchi]
233  );
234  mpb.distribute(primaryField);
235  return;
236  }
237  }
238 
240  << "Region patch ID " << regionPatchi << " not found in region mesh"
241  << abort(FatalError);
242 }
243 
244 
245 template<class Type>
247 (
248  Field<Type>& regionField,
249  const label regionPatchi,
250  const fvPatchField<Type>& primaryPatchField
251 ) const
252 {
253  const polyPatch& regionPatch = regionMesh().boundaryMesh()[regionPatchi];
254  const mappedPatchBase& mpb = refCast<const mappedPatchBase>(regionPatch);
255 
256  mappedPatchFieldBase<Type> mpf(mpb, primaryPatchField);
257 
258  UIndirectList<Type>(regionField, regionPatch.faceCells())
259  = mpf.mappedField();
260 }
261 
262 
263 template<class Type>
265 (
266  Field<Type>& rf,
268 ) const
269 {
270  forAll(intCoupledPatchIDs_, i)
271  {
272  const label regionPatchi = intCoupledPatchIDs_[i];
273  const label primaryPatchi = primaryPatchIDs_[i];
274 
275  toRegion(rf, regionPatchi, pBf[primaryPatchi]);
276  }
277 }
278 
279 
280 // ************************************************************************* //
#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
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:291
bool foundObject(const word &name) const
Is the named Type found?
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 Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
void toPrimary(const label regionPatchi, List< Type > &regionField) const
Convert a local region field to the primary region.
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:334
Pre-declare SubField and related Field type.
Definition: Field.H:56
A class for handling words, derived from string.
Definition: word.H:59
tmp< Field< Type > > mapRegionPatchInternalField(const regionModel &nbrRegion, const word &fieldName, const label regionPatchi, const bool flip=false) const
Map patch internal field from another region model to local.
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:55
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
void interpolateToSource(const UList< Type > &fld, const CombineOp &cop, List< Type > &result, const UList< Type > &defaultValues=UList< Type >::null()) const
Interpolate from target to source with supplied op.
Functionality for sampling fields using mappedPatchBase. Every call to mappedField() returns a sample...
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A List with indirect addressing.
Definition: fvMatrix.H:106
tmp< Foam::Field< Type > > mapRegionPatchField(const regionModel &nbrRegion, const label regionPatchi, const label nbrPatchi, const Field< Type > &nbrField, const bool flip=false) const
Map patch field from another region model to local patch.
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
void toRegion(const label regionPatchi, List< Type > &primaryFieldField) const
Convert a primary region field to the local region.
Base class for region models.
Definition: regionModel.H:55
void reverseDistribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66