mapNearestAMI.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) 2013-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 #include "mapNearestAMI.H"
27 
28 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
29 
30 template<class SourcePatch, class TargetPatch>
32 (
33  const SourcePatch& srcPatch,
34  const TargetPatch& tgtPatch,
35  const label& srcFacei,
36  label& tgtFacei
37 ) const
38 {
39  const vectorField& srcCf = srcPatch.faceCentres();
40  const vectorField& tgtCf = tgtPatch.faceCentres();
41 
42  const vector srcP = srcCf[srcFacei];
43 
44  DynamicList<label> tgtFaces(10);
45  tgtFaces.append(tgtFacei);
46 
47  DynamicList<label> visitedFaces(10);
48 
49  scalar d = GREAT;
50 
51  do
52  {
53  label tgtI = tgtFaces.remove();
54  visitedFaces.append(tgtI);
55 
56  scalar dTest = magSqr(tgtCf[tgtI] - srcP);
57  if (dTest < d)
58  {
59  tgtFacei = tgtI;
60  d = dTest;
61 
62  this->appendNbrFaces
63  (
64  tgtFacei,
65  tgtPatch,
66  visitedFaces,
67  tgtFaces
68  );
69  }
70 
71  } while (tgtFaces.size() > 0);
72 }
73 
74 
75 template<class SourcePatch, class TargetPatch>
77 (
78  boolList& mapFlag,
79  label& startSeedI,
80  label& srcFacei,
81  label& tgtFacei
82 ) const
83 {
84  const labelList& srcNbr = this->srcPatch_.faceFaces()[srcFacei];
85 
86  srcFacei = -1;
87 
88  forAll(srcNbr, i)
89  {
90  label facei = srcNbr[i];
91  if (mapFlag[facei])
92  {
93  srcFacei = facei;
94  startSeedI = facei + 1;
95 
96  return;
97  }
98  }
99 
100  forAll(mapFlag, facei)
101  {
102  if (mapFlag[facei])
103  {
104  srcFacei = facei;
105  tgtFacei = this->findTargetFace(facei);
106 
107  if (tgtFacei == -1)
108  {
109  const vectorField& srcCf = this->srcPatch_.faceCentres();
110 
112  << "Unable to find target face for source face "
113  << srcFacei << " with face centre " << srcCf[srcFacei]
114  << abort(FatalError);
115  }
116 
117  break;
118  }
119  }
120 }
121 
122 
123 template<class SourcePatch, class TargetPatch>
125 (
126  const label tgtFacei,
127  const List<DynamicList<label>>& tgtToSrc
128 ) const
129 {
130  DynamicList<label> testFaces(10);
131  DynamicList<label> visitedFaces(10);
132 
133  testFaces.append(tgtFacei);
134 
135  do
136  {
137  // search target tgtFacei neighbours for match with source face
138  label tgtI = testFaces.remove();
139 
140  if (findIndex(visitedFaces, tgtI) == -1)
141  {
142  visitedFaces.append(tgtI);
143 
144  if (tgtToSrc[tgtI].size())
145  {
146  return tgtToSrc[tgtI][0];
147  }
148  else
149  {
150  const labelList& nbrFaces = this->tgtPatch_.faceFaces()[tgtI];
151 
152  forAll(nbrFaces, i)
153  {
154  if (findIndex(visitedFaces, nbrFaces[i]) == -1)
155  {
156  testFaces.append(nbrFaces[i]);
157  }
158  }
159  }
160  }
161  } while (testFaces.size());
162 
163  // did not find any match - should not be possible to get here!
164  return -1;
165 }
166 
167 
168 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
169 
170 template<class SourcePatch, class TargetPatch>
172 (
173  const SourcePatch& srcPatch,
174  const TargetPatch& tgtPatch,
175  const scalarField& srcMagSf,
176  const scalarField& tgtMagSf,
178  const bool reverseTarget,
179  const bool requireMatch
180 )
181 :
183  (
184  srcPatch,
185  tgtPatch,
186  srcMagSf,
187  tgtMagSf,
188  triMode,
189  reverseTarget,
190  requireMatch
191  )
192 {}
193 
194 
195 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
196 
197 template<class SourcePatch, class TargetPatch>
199 {}
200 
201 
202 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
203 
204 template<class SourcePatch, class TargetPatch>
206 (
207  labelListList& srcAddress,
208  scalarListList& srcWeights,
209  labelListList& tgtAddress,
210  scalarListList& tgtWeights,
211  label srcFacei,
212  label tgtFacei
213 )
214 {
215  bool ok =
216  this->initialise
217  (
218  srcAddress,
219  srcWeights,
220  tgtAddress,
221  tgtWeights,
222  srcFacei,
223  tgtFacei
224  );
225 
226  if (!ok)
227  {
228  return;
229  }
230 
231 
232  // temporary storage for addressing and weights
233  List<DynamicList<label>> srcAddr(this->srcPatch_.size());
234  List<DynamicList<label>> tgtAddr(this->tgtPatch_.size());
235 
236 
237  // construct weights and addressing
238  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
239 
240  // list to keep track of whether src face can be mapped
241  boolList mapFlag(srcAddr.size(), true);
242 
243  // reset starting seed
244  label startSeedI = 0;
245 
246  DynamicList<label> nonOverlapFaces;
247  do
248  {
249  findNearestFace(this->srcPatch_, this->tgtPatch_, srcFacei, tgtFacei);
250 
251  srcAddr[srcFacei].append(tgtFacei);
252  tgtAddr[tgtFacei].append(srcFacei);
253 
254  mapFlag[srcFacei] = false;
255 
256  // Do advancing front starting from srcFacei, tgtFacei
257  setNextNearestFaces
258  (
259  mapFlag,
260  startSeedI,
261  srcFacei,
262  tgtFacei
263  );
264  } while (srcFacei >= 0);
265 
266 
267  // for the case of multiple source faces per target face, select the
268  // nearest source face only and discard the others
269  const vectorField& srcCf = this->srcPatch_.faceCentres();
270  const vectorField& tgtCf = this->tgtPatch_.faceCentres();
271 
272  forAll(tgtAddr, targetFacei)
273  {
274  if (tgtAddr[targetFacei].size() > 1)
275  {
276  const vector& tgtC = tgtCf[tgtFacei];
277 
278  DynamicList<label>& srcFaces = tgtAddr[targetFacei];
279 
280  label srcFacei = srcFaces[0];
281  scalar d = magSqr(tgtC - srcCf[srcFacei]);
282 
283  for (label i = 1; i < srcFaces.size(); i++)
284  {
285  label srcI = srcFaces[i];
286  scalar dNew = magSqr(tgtC - srcCf[srcI]);
287  if (dNew < d)
288  {
289  d = dNew;
290  srcFacei = srcI;
291  }
292  }
293 
294  srcFaces.clear();
295  srcFaces.append(srcFacei);
296  }
297  }
298 
299  // If there are more target faces than source faces, some target faces
300  // might not yet be mapped
301  forAll(tgtAddr, tgtFacei)
302  {
303  if (tgtAddr[tgtFacei].empty())
304  {
305  label srcFacei = findMappedSrcFace(tgtFacei, tgtAddr);
306 
307  if (srcFacei >= 0)
308  {
309  // note - reversed search from src->tgt to tgt->src
310  findNearestFace
311  (
312  this->tgtPatch_,
313  this->srcPatch_,
314  tgtFacei,
315  srcFacei
316  );
317 
318  tgtAddr[tgtFacei].append(srcFacei);
319  }
320  }
321  }
322 
323 
324  // transfer data to persistent storage
325  forAll(srcAddr, i)
326  {
327  scalar magSf = this->srcMagSf_[i];
328  srcAddress[i].transfer(srcAddr[i]);
329  srcWeights[i] = scalarList(1, magSf);
330  }
331  forAll(tgtAddr, i)
332  {
333  scalar magSf = this->tgtMagSf_[i];
334  tgtAddress[i].transfer(tgtAddr[i]);
335  tgtWeights[i] = scalarList(1, magSf);
336  }
337 }
338 
339 
340 // ************************************************************************* //
Nearest-mapping Arbitrary Mesh Interface (AMI) method.
Definition: mapNearestAMI.H:49
#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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
volVectorField vectorField(fieldObject, mesh)
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
virtual ~mapNearestAMI()
Destructor.
virtual void calculate(labelListList &srcAddress, scalarListList &srcWeights, labelListList &tgtAddress, scalarListList &tgtWeights, label srcFacei=-1, label tgtFacei=-1)
Update addressing and weights.
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:50
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:292
List< label > labelList
A List of labels.
Definition: labelList.H:56
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Base class for Arbitrary Mesh Interface (AMI) methods.
Definition: AMIMethod.H:55
dimensioned< scalar > magSqr(const dimensioned< Type > &)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
void clear()
Clear the addressed list, i.e. set the size to zero.
Definition: DynamicListI.H:224
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:342