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-2014 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  (
113  "void Foam::mapNearestAMI<SourcePatch, TargetPatch>::"
114  "setNextNearestFaces"
115  "("
116  "boolList&, "
117  "label&, "
118  "label&, "
119  "label&"
120  ") const"
121  )
122  << "Unable to find target face for source face "
123  << srcFaceI << " with face centre " << srcCf[srcFaceI]
124  << abort(FatalError);
125  }
126 
127  break;
128  }
129  }
130 }
131 
132 
133 template<class SourcePatch, class TargetPatch>
135 (
136  const label tgtFaceI,
137  const List<DynamicList<label> >& tgtToSrc
138 ) const
139 {
140  DynamicList<label> testFaces(10);
141  DynamicList<label> visitedFaces(10);
142 
143  testFaces.append(tgtFaceI);
144 
145  do
146  {
147  // search target tgtFaceI neighbours for match with source face
148  label tgtI = testFaces.remove();
149 
150  if (findIndex(visitedFaces, tgtI) == -1)
151  {
152  visitedFaces.append(tgtI);
153 
154  if (tgtToSrc[tgtI].size())
155  {
156  return tgtToSrc[tgtI][0];
157  }
158  else
159  {
160  const labelList& nbrFaces = this->tgtPatch_.faceFaces()[tgtI];
161 
162  forAll(nbrFaces, i)
163  {
164  if (findIndex(visitedFaces, nbrFaces[i]) == -1)
165  {
166  testFaces.append(nbrFaces[i]);
167  }
168  }
169  }
170  }
171  } while (testFaces.size());
172 
173  // did not find any match - should not be possible to get here!
174  return -1;
175 }
176 
177 
178 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
179 
180 template<class SourcePatch, class TargetPatch>
182 (
183  const SourcePatch& srcPatch,
184  const TargetPatch& tgtPatch,
185  const scalarField& srcMagSf,
186  const scalarField& tgtMagSf,
188  const bool reverseTarget,
189  const bool requireMatch
190 )
191 :
193  (
194  srcPatch,
195  tgtPatch,
196  srcMagSf,
197  tgtMagSf,
198  triMode,
199  reverseTarget,
200  requireMatch
201  )
202 {}
203 
204 
205 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
206 
207 template<class SourcePatch, class TargetPatch>
209 {}
210 
211 
212 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
213 
214 template<class SourcePatch, class TargetPatch>
216 (
217  labelListList& srcAddress,
218  scalarListList& srcWeights,
219  labelListList& tgtAddress,
220  scalarListList& tgtWeights,
221  label srcFaceI,
222  label tgtFaceI
223 )
224 {
225  bool ok =
226  this->initialise
227  (
228  srcAddress,
229  srcWeights,
230  tgtAddress,
231  tgtWeights,
232  srcFaceI,
233  tgtFaceI
234  );
235 
236  if (!ok)
237  {
238  return;
239  }
240 
241 
242  // temporary storage for addressing and weights
243  List<DynamicList<label> > srcAddr(this->srcPatch_.size());
244  List<DynamicList<label> > tgtAddr(this->tgtPatch_.size());
245 
246 
247  // construct weights and addressing
248  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
249 
250  // list to keep track of whether src face can be mapped
251  boolList mapFlag(srcAddr.size(), true);
252 
253  // reset starting seed
254  label startSeedI = 0;
255 
256  DynamicList<label> nonOverlapFaces;
257  do
258  {
259  findNearestFace(this->srcPatch_, this->tgtPatch_, srcFaceI, tgtFaceI);
260 
261  srcAddr[srcFaceI].append(tgtFaceI);
262  tgtAddr[tgtFaceI].append(srcFaceI);
263 
264  mapFlag[srcFaceI] = false;
265 
266  // Do advancing front starting from srcFaceI, tgtFaceI
267  setNextNearestFaces
268  (
269  mapFlag,
270  startSeedI,
271  srcFaceI,
272  tgtFaceI
273  );
274  } while (srcFaceI >= 0);
275 
276 
277  // for the case of multiple source faces per target face, select the
278  // nearest source face only and discard the others
279  const vectorField& srcCf = this->srcPatch_.faceCentres();
280  const vectorField& tgtCf = this->tgtPatch_.faceCentres();
281 
282  forAll(tgtAddr, targetFaceI)
283  {
284  if (tgtAddr[targetFaceI].size() > 1)
285  {
286  const vector& tgtC = tgtCf[tgtFaceI];
287 
288  DynamicList<label>& srcFaces = tgtAddr[targetFaceI];
289 
290  label srcFaceI = srcFaces[0];
291  scalar d = magSqr(tgtC - srcCf[srcFaceI]);
292 
293  for (label i = 1; i < srcFaces.size(); i++)
294  {
295  label srcI = srcFaces[i];
296  scalar dNew = magSqr(tgtC - srcCf[srcI]);
297  if (dNew < d)
298  {
299  d = dNew;
300  srcFaceI = srcI;
301  }
302  }
303 
304  srcFaces.clear();
305  srcFaces.append(srcFaceI);
306  }
307  }
308 
309  // If there are more target faces than source faces, some target faces
310  // might not yet be mapped
311  forAll(tgtAddr, tgtFaceI)
312  {
313  if (tgtAddr[tgtFaceI].empty())
314  {
315  label srcFaceI = findMappedSrcFace(tgtFaceI, tgtAddr);
316 
317  if (srcFaceI >= 0)
318  {
319  // note - reversed search from src->tgt to tgt->src
320  findNearestFace
321  (
322  this->tgtPatch_,
323  this->srcPatch_,
324  tgtFaceI,
325  srcFaceI
326  );
327 
328  tgtAddr[tgtFaceI].append(srcFaceI);
329  }
330  }
331  }
332 
333 
334  // transfer data to persistent storage
335  forAll(srcAddr, i)
336  {
337  scalar magSf = this->srcMagSf_[i];
338  srcAddress[i].transfer(srcAddr[i]);
339  srcWeights[i] = scalarList(1, magSf);
340  }
341  forAll(tgtAddr, i)
342  {
343  scalar magSf = this->tgtMagSf_[i];
344  tgtAddress[i].transfer(tgtAddr[i]);
345  tgtWeights[i] = scalarList(1, magSf);
346  }
347 }
348 
349 
350 // ************************************************************************* //
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:390
dimensioned< scalar > magSqr(const dimensioned< Type > &)
Base class for Arbitrary Mesh Interface (AMI) methods.
Definition: AMIMethod.H:55
void clear()
Clear the addressed list, i.e. set the size to zero.
Definition: DynamicListI.H:242
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:310
volVectorField vectorField(fieldObject, mesh)
#define forAll(list, i)
Definition: UList.H:421
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:50
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Nearest-mapping Arbitrary Mesh Interface (AMI) method.
Definition: mapNearestAMI.H:49
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
virtual ~mapNearestAMI()
Destructor.
error FatalError
List< label > labelList
A List of labels.
Definition: labelList.H:56
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
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