regionModelTemplates.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-2016 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 template<class Type>
29 (
30  const regionModel& nbrRegion,
31  const label regionPatchi,
32  const label nbrPatchi,
33  const Field<Type>& nbrField,
34  const bool flip
35 ) const
36 {
37  int oldTag = UPstream::msgType();
38  UPstream::msgType() = oldTag + 1;
39 
40  const AMIPatchToPatchInterpolation& ami =
41  interRegionAMI(nbrRegion, regionPatchi, nbrPatchi, flip);
42 
43  tmp<Field<Type>> tresult(ami.interpolateToSource(nbrField));
44 
45  UPstream::msgType() = oldTag;
46 
47  return tresult;
48 }
49 
50 
51 template<class Type>
54 (
55  const regionModel& nbrRegion,
56  const word& fieldName,
57  const label regionPatchi,
58  const bool flip
59 ) const
60 {
62 
63  const fvMesh& nbrRegionMesh = nbrRegion.regionMesh();
64 
65  if (nbrRegionMesh.foundObject<fieldType>(fieldName))
66  {
67  const label nbrPatchi = nbrCoupledPatchID(nbrRegion, regionPatchi);
68 
69  int oldTag = UPstream::msgType();
70  UPstream::msgType() = oldTag + 1;
71 
72  const AMIPatchToPatchInterpolation& ami =
73  interRegionAMI(nbrRegion, regionPatchi, nbrPatchi, flip);
74 
75  const fieldType& nbrField =
76  nbrRegionMesh.lookupObject<fieldType>(fieldName);
77 
78  const Field<Type>& nbrFieldp = nbrField.boundaryField()[nbrPatchi];
79 
80  tmp<Field<Type>> tresult(ami.interpolateToSource(nbrFieldp));
81 
82  UPstream::msgType() = oldTag;
83 
84  return tresult;
85  }
86  else
87  {
88  const polyPatch& p = regionMesh().boundaryMesh()[regionPatchi];
89 
90  return
92  (
93  new Field<Type>
94  (
95  p.size(),
96  Zero
97  )
98  );
99  }
100 }
101 
102 
103 template<class Type>
106 (
107  const regionModel& nbrRegion,
108  const word& fieldName,
109  const label regionPatchi,
110  const bool flip
111 ) const
112 {
114 
115  const fvMesh& nbrRegionMesh = nbrRegion.regionMesh();
116 
117  if (nbrRegionMesh.foundObject<fieldType>(fieldName))
118  {
119  const label nbrPatchi = nbrCoupledPatchID(nbrRegion, regionPatchi);
120 
121  int oldTag = UPstream::msgType();
122  UPstream::msgType() = oldTag + 1;
123 
124  const AMIPatchToPatchInterpolation& ami =
125  interRegionAMI(nbrRegion, regionPatchi, nbrPatchi, flip);
126 
127  const fieldType& nbrField =
128  nbrRegionMesh.lookupObject<fieldType>(fieldName);
129 
130  const fvPatchField<Type>& nbrFieldp =
131  nbrField.boundaryField()[nbrPatchi];
132 
133  tmp<Field<Type>> tresult
134  (
135  ami.interpolateToSource(nbrFieldp.patchInternalField())
136  );
137 
138  UPstream::msgType() = oldTag;
139 
140  return tresult;
141  }
142  else
143  {
144  const polyPatch& p = regionMesh().boundaryMesh()[regionPatchi];
145 
146  return
148  (
149  new Field<Type>
150  (
151  p.size(),
152  Zero
153  )
154  );
155  }
156 }
157 
158 
159 template<class Type>
161 (
162  const label regionPatchi,
163  List<Type>& regionField
164 ) const
165 {
166  forAll(intCoupledPatchIDs_, i)
167  {
168  if (intCoupledPatchIDs_[i] == regionPatchi)
169  {
170  const mappedPatchBase& mpb =
171  refCast<const mappedPatchBase>
172  (
173  regionMesh().boundaryMesh()[regionPatchi]
174  );
175  mpb.reverseDistribute(regionField);
176  return;
177  }
178  }
179 
181  << "Region patch ID " << regionPatchi << " not found in region mesh"
182  << abort(FatalError);
183 }
184 
185 
186 template<class Type>
188 (
189  const label regionPatchi,
190  List<Type>& primaryField
191 ) const
192 {
193  forAll(intCoupledPatchIDs_, i)
194  {
195  if (intCoupledPatchIDs_[i] == regionPatchi)
196  {
197  const mappedPatchBase& mpb =
198  refCast<const mappedPatchBase>
199  (
200  regionMesh().boundaryMesh()[regionPatchi]
201  );
202  mpb.distribute(primaryField);
203  return;
204  }
205  }
206 
208  << "Region patch ID " << regionPatchi << " not found in region mesh"
209  << abort(FatalError);
210 }
211 
212 
213 template<class Type, class CombineOp>
215 (
216  const label regionPatchi,
217  List<Type>& regionField,
218  const CombineOp& cop
219 ) const
220 {
221  forAll(intCoupledPatchIDs_, i)
222  {
223  if (intCoupledPatchIDs_[i] == regionPatchi)
224  {
225  const mappedPatchBase& mpb =
226  refCast<const mappedPatchBase>
227  (
228  regionMesh().boundaryMesh()[regionPatchi]
229  );
230  mpb.reverseDistribute(regionField, cop);
231  return;
232  }
233  }
234 
236  << "Region patch ID " << regionPatchi << " not found in region mesh"
237  << abort(FatalError);
238 }
239 
240 
241 template<class Type, class CombineOp>
243 (
244  const label regionPatchi,
245  List<Type>& primaryField,
246  const CombineOp& cop
247 ) const
248 {
249  forAll(intCoupledPatchIDs_, i)
250  {
251  if (intCoupledPatchIDs_[i] == regionPatchi)
252  {
253  const mappedPatchBase& mpb =
254  refCast<const mappedPatchBase>
255  (
256  regionMesh().boundaryMesh()[regionPatchi]
257  );
258  mpb.distribute(primaryField, cop);
259  return;
260  }
261  }
262 
264  << "Region patch ID " << regionPatchi << " not found in region mesh"
265  << abort(FatalError);
266 }
267 
268 
269 // ************************************************************************* //
#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
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
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.
Pre-declare SubField and related Field type.
Definition: Field.H:57
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:61
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:91
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.
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
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.
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