FaceCellWave.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-2020 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 "FaceCellWave.H"
27 #include "polyMesh.H"
28 #include "processorPolyPatch.H"
29 #include "cyclicPolyPatch.H"
30 #include "cyclicAMIPolyPatch.H"
31 #include "OPstream.H"
32 #include "IPstream.H"
33 #include "PstreamReduceOps.H"
34 #include "debug.H"
35 #include "typeInfo.H"
36 #include "SubField.H"
37 #include "globalMeshData.H"
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 template<class Type, class TrackingData>
43 
44 template<class Type, class TrackingData>
46 
47 template<class Type, class TrackingData>
49 
50 namespace Foam
51 {
52  template<class Type, class TrackingData>
53  class combine
54  {
55  //- Combine operator for AMIInterpolation
56 
58 
59  const cyclicAMIPolyPatch& patch_;
60 
61  public:
62 
63  combine
64  (
66  const cyclicAMIPolyPatch& patch
67  )
68  :
69  solver_(solver),
70  patch_(patch)
71  {}
72 
73 
74  void operator()
75  (
76  Type& x,
77  const label facei,
78  const Type& y,
79  const scalar weight
80  ) const
81  {
82  if (y.valid(solver_.data()))
83  {
84  label meshFacei = -1;
85  if (patch_.owner())
86  {
87  meshFacei = patch_.start() + facei;
88  }
89  else
90  {
91  meshFacei = patch_.nbrPatch().start() + facei;
92  }
93  x.updateFace
94  (
95  solver_.mesh(),
96  meshFacei,
97  y,
98  solver_.propagationTol(),
99  solver_.data()
100  );
101  }
102  }
103  };
104 }
105 
106 
107 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
108 
109 template<class Type, class TrackingData>
111 (
112  const label celli,
113  const label neighbourFacei,
114  const Type& neighbourInfo,
115  const scalar tol,
116  Type& cellInfo
117 )
118 {
119  // Update info for celli, at position pt, with information from
120  // neighbouring face/cell.
121  // Updates:
122  // - changedCell_, changedCells_
123  // - statistics: nEvals_, nUnvisitedCells_
124 
125  nEvals_++;
126 
127  bool wasValid = cellInfo.valid(td_);
128 
129  bool propagate =
130  cellInfo.updateCell
131  (
132  mesh_,
133  celli,
134  neighbourFacei,
135  neighbourInfo,
136  tol,
137  td_
138  );
139 
140  if (propagate)
141  {
142  if (!changedCell_[celli])
143  {
144  changedCell_[celli] = true;
145  changedCells_.append(celli);
146  }
147  }
148 
149  if (!wasValid && cellInfo.valid(td_))
150  {
151  --nUnvisitedCells_;
152  }
153 
154  return propagate;
155 }
156 
157 
158 template<class Type, class TrackingData>
160 (
161  const label facei,
162  const label neighbourCelli,
163  const Type& neighbourInfo,
164  const scalar tol,
165  Type& faceInfo
166 )
167 {
168  // Update info for facei, at position pt, with information from
169  // neighbouring face/cell.
170  // Updates:
171  // - changedFace_, changedFaces_,
172  // - statistics: nEvals_, nUnvisitedFaces_
173 
174  nEvals_++;
175 
176  bool wasValid = faceInfo.valid(td_);
177 
178  bool propagate =
179  faceInfo.updateFace
180  (
181  mesh_,
182  facei,
183  neighbourCelli,
184  neighbourInfo,
185  tol,
186  td_
187  );
188 
189  if (propagate)
190  {
191  if (!changedFace_[facei])
192  {
193  changedFace_[facei] = true;
194  changedFaces_.append(facei);
195  }
196  }
197 
198  if (!wasValid && faceInfo.valid(td_))
199  {
200  --nUnvisitedFaces_;
201  }
202 
203  return propagate;
204 }
205 
206 
207 template<class Type, class TrackingData>
209 (
210  const label facei,
211  const Type& neighbourInfo,
212  const scalar tol,
213  Type& faceInfo
214 )
215 {
216  // Update info for facei, at position pt, with information from
217  // same face.
218  // Updates:
219  // - changedFace_, changedFaces_,
220  // - statistics: nEvals_, nUnvisitedFaces_
221 
222  nEvals_++;
223 
224  bool wasValid = faceInfo.valid(td_);
225 
226  bool propagate =
227  faceInfo.updateFace
228  (
229  mesh_,
230  facei,
231  neighbourInfo,
232  tol,
233  td_
234  );
235 
236  if (propagate)
237  {
238  if (!changedFace_[facei])
239  {
240  changedFace_[facei] = true;
241  changedFaces_.append(facei);
242  }
243  }
244 
245  if (!wasValid && faceInfo.valid(td_))
246  {
247  --nUnvisitedFaces_;
248  }
249 
250  return propagate;
251 }
252 
253 
254 template<class Type, class TrackingData>
256 (
257  const polyPatch& patch
258 ) const
259 {
260  // For debugging: check status on both sides of cyclic
261 
262  const cyclicPolyPatch& nbrPatch =
263  refCast<const cyclicPolyPatch>(patch).nbrPatch();
264 
265  forAll(patch, patchFacei)
266  {
267  label i1 = patch.start() + patchFacei;
268  label i2 = nbrPatch.start() + patchFacei;
269 
270  if
271  (
272  !allFaceInfo_[i1].sameGeometry
273  (
274  mesh_,
275  allFaceInfo_[i2],
276  geomTol_,
277  td_
278  )
279  )
280  {
282  << " faceInfo:" << allFaceInfo_[i1]
283  << " otherfaceInfo:" << allFaceInfo_[i2]
284  << abort(FatalError);
285  }
286 
287  if (changedFace_[i1] != changedFace_[i2])
288  {
290  << " faceInfo:" << allFaceInfo_[i1]
291  << " otherfaceInfo:" << allFaceInfo_[i2]
292  << " changedFace:" << changedFace_[i1]
293  << " otherchangedFace:" << changedFace_[i2]
294  << abort(FatalError);
295  }
296  }
297 }
298 
299 
300 template<class Type, class TrackingData>
301 template<class PatchType>
303 {
304  forAll(mesh_.boundaryMesh(), patchi)
305  {
306  if (isA<PatchType>(mesh_.boundaryMesh()[patchi]))
307  {
308  return true;
309  }
310  }
311  return false;
312 }
313 
314 
315 template<class Type, class TrackingData>
317 (
318  const labelList& changedFaces,
319  const List<Type>& changedFacesInfo
320 )
321 {
322  forAll(changedFaces, changedFacei)
323  {
324  label facei = changedFaces[changedFacei];
325 
326  bool wasValid = allFaceInfo_[facei].valid(td_);
327 
328  // Copy info for facei
329  allFaceInfo_[facei] = changedFacesInfo[changedFacei];
330 
331  // Maintain count of unset faces
332  if (!wasValid && allFaceInfo_[facei].valid(td_))
333  {
334  --nUnvisitedFaces_;
335  }
336 
337  // Mark facei as changed, both on list and on face itself.
338 
339  changedFace_[facei] = true;
340  changedFaces_.append(facei);
341  }
342 }
343 
344 
345 template<class Type, class TrackingData>
347 (
348  const polyPatch& patch,
349  const label nFaces,
350  const labelList& changedFaces,
351  const List<Type>& changedFacesInfo
352 )
353 {
354  // Merge face information into member data
355 
356  for (label changedFacei = 0; changedFacei < nFaces; changedFacei++)
357  {
358  const Type& neighbourWallInfo = changedFacesInfo[changedFacei];
359  label patchFacei = changedFaces[changedFacei];
360 
361  label meshFacei = patch.start() + patchFacei;
362 
363  Type& currentWallInfo = allFaceInfo_[meshFacei];
364 
365  if (!currentWallInfo.equal(neighbourWallInfo, td_))
366  {
367  updateFace
368  (
369  meshFacei,
370  neighbourWallInfo,
371  propagationTol_,
372  currentWallInfo
373  );
374  }
375  }
376 }
377 
378 
379 template<class Type, class TrackingData>
381 (
382  const polyPatch& patch,
383  const label startFacei,
384  const label nFaces,
385  labelList& changedPatchFaces,
386  List<Type>& changedPatchFacesInfo
387 ) const
388 {
389  // Construct compact patchFace change arrays for a (slice of a) single
390  // patch. changedPatchFaces in local patch numbering.
391  // Return length of arrays.
392  label nChangedPatchFaces = 0;
393 
394  for (label i = 0; i < nFaces; i++)
395  {
396  label patchFacei = i + startFacei;
397 
398  label meshFacei = patch.start() + patchFacei;
399 
400  if (changedFace_[meshFacei])
401  {
402  changedPatchFaces[nChangedPatchFaces] = patchFacei;
403  changedPatchFacesInfo[nChangedPatchFaces] = allFaceInfo_[meshFacei];
404  nChangedPatchFaces++;
405  }
406  }
407  return nChangedPatchFaces;
408 }
409 
410 
411 template<class Type, class TrackingData>
413 (
414  const polyPatch& patch,
415  const label nFaces,
416  const labelList& faceLabels,
417  List<Type>& faceInfo
418 ) const
419 {
420  // Handle leaving domain. Implementation referred to Type
421 
422  const vectorField& fc = mesh_.faceCentres();
423 
424  for (label i = 0; i < nFaces; i++)
425  {
426  label patchFacei = faceLabels[i];
427 
428  label meshFacei = patch.start() + patchFacei;
429  faceInfo[i].leaveDomain(mesh_, patch, patchFacei, fc[meshFacei], td_);
430  }
431 }
432 
433 
434 template<class Type, class TrackingData>
436 (
437  const polyPatch& patch,
438  const label nFaces,
439  const labelList& faceLabels,
440  List<Type>& faceInfo
441 ) const
442 {
443  // Handle entering domain. Implementation referred to Type
444 
445  const vectorField& fc = mesh_.faceCentres();
446 
447  for (label i = 0; i < nFaces; i++)
448  {
449  label patchFacei = faceLabels[i];
450 
451  label meshFacei = patch.start() + patchFacei;
452  faceInfo[i].enterDomain(mesh_, patch, patchFacei, fc[meshFacei], td_);
453  }
454 }
455 
456 
457 template<class Type, class TrackingData>
459 (
460  const tensor& rotTensor,
461  const label nFaces,
462  List<Type>& faceInfo
463 )
464 {
465  for (label facei = 0; facei < nFaces; facei++)
466  {
467  faceInfo[facei].transform(mesh_, rotTensor, td_);
468  }
469 }
470 
471 
472 template<class Type, class TrackingData>
474 (
475  const transformer& trans,
476  const label nFaces,
477  List<Type>& faceInfo
478 )
479 {
480  // Transform. Implementation referred to Type
481 
482  if (trans.transforms())
483  {
484  for (label facei = 0; facei < nFaces; facei++)
485  {
486  faceInfo[facei].transform(mesh_, trans.T(), td_);
487  }
488  }
489 }
490 
491 
492 template<class Type, class TrackingData>
494 (
495  const polyPatch&,
496  const label cycOffset,
497  const label nFaces,
498  labelList& faces
499 )
500 {
501  // Offset mesh face. Used for transferring from one cyclic half to the
502  // other.
503 
504  for (label facei = 0; facei < nFaces; facei++)
505  {
506  faces[facei] += cycOffset;
507  }
508 }
509 
510 
511 template<class Type, class TrackingData>
513 {
514  // Transfer all the information to/from neighbouring processors
515 
516  const globalMeshData& pData = mesh_.globalData();
517 
518  // Which patches are processor patches
519  const labelList& procPatches = pData.processorPatches();
520 
521  // Send all
522 
524 
525  forAll(procPatches, i)
526  {
527  label patchi = procPatches[i];
528 
529  const processorPolyPatch& procPatch =
530  refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchi]);
531 
532  // Allocate buffers
533  label nSendFaces;
534  labelList sendFaces(procPatch.size());
535  List<Type> sendFacesInfo(procPatch.size());
536 
537  // Determine which faces changed on current patch
538  nSendFaces = getChangedPatchFaces
539  (
540  procPatch,
541  0,
542  procPatch.size(),
543  sendFaces,
544  sendFacesInfo
545  );
546 
547  // Adapt wallInfo for leaving domain
548  leaveDomain
549  (
550  procPatch,
551  nSendFaces,
552  sendFaces,
553  sendFacesInfo
554  );
555 
556  if (debug & 2)
557  {
558  Pout<< " Processor patch " << patchi << ' ' << procPatch.name()
559  << " communicating with " << procPatch.neighbProcNo()
560  << " Sending:" << nSendFaces
561  << endl;
562  }
563 
564  UOPstream toNeighbour(procPatch.neighbProcNo(), pBufs);
565  // writeFaces(nSendFaces, sendFaces, sendFacesInfo, toNeighbour);
566  toNeighbour
567  << SubList<label>(sendFaces, nSendFaces)
568  << SubList<Type>(sendFacesInfo, nSendFaces);
569  }
570 
571  pBufs.finishedSends();
572 
573  // Receive all
574 
575  forAll(procPatches, i)
576  {
577  label patchi = procPatches[i];
578 
579  const processorPolyPatch& procPatch =
580  refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchi]);
581 
582  // Allocate buffers
583  labelList receiveFaces;
584  List<Type> receiveFacesInfo;
585 
586  {
587  UIPstream fromNeighbour(procPatch.neighbProcNo(), pBufs);
588  fromNeighbour >> receiveFaces >> receiveFacesInfo;
589  }
590 
591  if (debug & 2)
592  {
593  Pout<< " Processor patch " << patchi << ' ' << procPatch.name()
594  << " communicating with " << procPatch.neighbProcNo()
595  << " Receiving:" << receiveFaces.size()
596  << endl;
597  }
598 
599  // Apply transform to received data for non-parallel planes
600  if (procPatch.transform().transforms())
601  {
602  transform
603  (
604  procPatch.transform().T(),
605  receiveFaces.size(),
606  receiveFacesInfo
607  );
608  }
609 
610  // Adapt wallInfo for entering domain
611  enterDomain
612  (
613  procPatch,
614  receiveFaces.size(),
615  receiveFaces,
616  receiveFacesInfo
617  );
618 
619  // Merge received info
620  mergeFaceInfo
621  (
622  procPatch,
623  receiveFaces.size(),
624  receiveFaces,
625  receiveFacesInfo
626  );
627  }
628 }
629 
630 
631 template<class Type, class TrackingData>
633 {
634  // Transfer information across cyclic halves.
635 
636  forAll(mesh_.boundaryMesh(), patchi)
637  {
638  const polyPatch& patch = mesh_.boundaryMesh()[patchi];
639 
640  if (isA<cyclicPolyPatch>(patch))
641  {
642  const cyclicPolyPatch& nbrPatch =
643  refCast<const cyclicPolyPatch>(patch).nbrPatch();
644 
645  // Allocate buffers
646  label nReceiveFaces;
647  labelList receiveFaces(patch.size());
648  List<Type> receiveFacesInfo(patch.size());
649 
650  // Determine which faces changed
651  nReceiveFaces = getChangedPatchFaces
652  (
653  nbrPatch,
654  0,
655  nbrPatch.size(),
656  receiveFaces,
657  receiveFacesInfo
658  );
659 
660  // Adapt wallInfo for leaving domain
661  leaveDomain
662  (
663  nbrPatch,
664  nReceiveFaces,
665  receiveFaces,
666  receiveFacesInfo
667  );
668 
669  const cyclicPolyPatch& cycPatch =
670  refCast<const cyclicPolyPatch>(patch);
671 
672  if (cycPatch.transform().transforms())
673  {
674  // received data from other half
675  transform
676  (
677  cycPatch.transform().T(),
678  nReceiveFaces,
679  receiveFacesInfo
680  );
681  }
682 
683  if (debug & 2)
684  {
685  Pout<< " Cyclic patch " << patchi << ' ' << cycPatch.name()
686  << " Changed : " << nReceiveFaces
687  << endl;
688  }
689 
690  // Half2: Adapt wallInfo for entering domain
691  enterDomain
692  (
693  cycPatch,
694  nReceiveFaces,
695  receiveFaces,
696  receiveFacesInfo
697  );
698 
699  // Merge into global storage
700  mergeFaceInfo
701  (
702  cycPatch,
703  nReceiveFaces,
704  receiveFaces,
705  receiveFacesInfo
706  );
707 
708  if (debug)
709  {
710  checkCyclic(cycPatch);
711  }
712  }
713  }
714 }
715 
716 
717 template<class Type, class TrackingData>
719 {
720  // Transfer information across cyclicAMI halves.
721 
722  forAll(mesh_.boundaryMesh(), patchi)
723  {
724  const polyPatch& patch = mesh_.boundaryMesh()[patchi];
725 
726  if (isA<cyclicAMIPolyPatch>(patch))
727  {
728  const cyclicAMIPolyPatch& cycPatch =
729  refCast<const cyclicAMIPolyPatch>(patch);
730 
731  List<Type> receiveInfo;
732 
733  {
734  const cyclicAMIPolyPatch& nbrPatch =
735  refCast<const cyclicAMIPolyPatch>(patch).nbrPatch();
736 
737  // Get nbrPatch data (so not just changed faces)
738  typename List<Type>::subList sendInfo
739  (
740  nbrPatch.patchSlice
741  (
742  allFaceInfo_
743  )
744  );
745 
746  if (nbrPatch.transform().transformsPosition())
747  {
748  // Adapt sendInfo for leaving domain
749  const vectorField::subField fc = nbrPatch.faceCentres();
750  forAll(sendInfo, i)
751  {
752  sendInfo[i].leaveDomain(mesh_, nbrPatch, i, fc[i], td_);
753  }
754  }
755 
756  // Transfer sendInfo to cycPatch
757  combine<Type, TrackingData> cmb(*this, cycPatch);
758 
759  List<Type> defVals;
760  if (cycPatch.applyLowWeightCorrection())
761  {
762  defVals = cycPatch.patchInternalList(allCellInfo_);
763  }
764 
765  if (cycPatch.owner())
766  {
767  forAll(cycPatch.AMIs(), i)
768  {
769  List<Type> sendInfoT(sendInfo);
770  transform
771  (
772  cycPatch.AMITransforms()[i],
773  sendInfoT.size(),
774  sendInfoT
775  );
776  cycPatch.AMIs()[i].interpolateToSource
777  (
778  sendInfoT,
779  cmb,
780  receiveInfo,
781  defVals
782  );
783  }
784  }
785  else
786  {
787  forAll(cycPatch.nbrPatch().AMIs(), i)
788  {
789  List<Type> sendInfoT(sendInfo);
790  transform
791  (
792  cycPatch.nbrPatch().AMITransforms()[i],
793  sendInfoT.size(),
794  sendInfoT
795  );
796  cycPatch.nbrPatch().AMIs()[i].interpolateToTarget
797  (
798  sendInfoT,
799  cmb,
800  receiveInfo,
801  defVals
802  );
803  }
804  }
805  }
806 
807  // Apply transform to received data for non-parallel planes
808  if (cycPatch.transform().transforms())
809  {
810  transform
811  (
812  cycPatch.transform().T(),
813  receiveInfo.size(),
814  receiveInfo
815  );
816  }
817 
818  if (cycPatch.transform().transformsPosition())
819  {
820  // Adapt receiveInfo for entering domain
821  const vectorField::subField fc = cycPatch.faceCentres();
822  forAll(receiveInfo, i)
823  {
824  receiveInfo[i].enterDomain(mesh_, cycPatch, i, fc[i], td_);
825  }
826  }
827 
828  // Merge into global storage
829  forAll(receiveInfo, i)
830  {
831  label meshFacei = cycPatch.start()+i;
832 
833  Type& currentWallInfo = allFaceInfo_[meshFacei];
834 
835  if
836  (
837  receiveInfo[i].valid(td_)
838  && !currentWallInfo.equal(receiveInfo[i], td_)
839  )
840  {
841  updateFace
842  (
843  meshFacei,
844  receiveInfo[i],
845  propagationTol_,
846  currentWallInfo
847  );
848  }
849  }
850  }
851  }
852 }
853 
854 
855 template<class Type, class TrackingData>
857 {
858  // Collect changed information
859 
860  DynamicList<label> f0Baffle(explicitConnections_.size());
861  DynamicList<Type> f0Info(explicitConnections_.size());
862 
863  DynamicList<label> f1Baffle(explicitConnections_.size());
864  DynamicList<Type> f1Info(explicitConnections_.size());
865 
866  forAll(explicitConnections_, connI)
867  {
868  const labelPair& baffle = explicitConnections_[connI];
869 
870  label f0 = baffle[0];
871  if (changedFace_[f0])
872  {
873  f0Baffle.append(connI);
874  f0Info.append(allFaceInfo_[f0]);
875  }
876 
877  label f1 = baffle[1];
878  if (changedFace_[f1])
879  {
880  f1Baffle.append(connI);
881  f1Info.append(allFaceInfo_[f1]);
882  }
883  }
884 
885 
886  // Update other side with changed information
887 
888  forAll(f1Info, index)
889  {
890  const labelPair& baffle = explicitConnections_[f1Baffle[index]];
891 
892  label f0 = baffle[0];
893  Type& currentWallInfo = allFaceInfo_[f0];
894 
895  if (!currentWallInfo.equal(f1Info[index], td_))
896  {
897  updateFace
898  (
899  f0,
900  f1Info[index],
901  propagationTol_,
902  currentWallInfo
903  );
904  }
905  }
906 
907  forAll(f0Info, index)
908  {
909  const labelPair& baffle = explicitConnections_[f0Baffle[index]];
910 
911  label f1 = baffle[1];
912  Type& currentWallInfo = allFaceInfo_[f1];
913 
914  if (!currentWallInfo.equal(f0Info[index], td_))
915  {
916  updateFace
917  (
918  f1,
919  f0Info[index],
920  propagationTol_,
921  currentWallInfo
922  );
923  }
924  }
925 }
926 
927 
928 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
929 
930 template<class Type, class TrackingData>
932 (
933  const polyMesh& mesh,
934  UList<Type>& allFaceInfo,
935  UList<Type>& allCellInfo,
936  TrackingData& td
937 )
938 :
939  mesh_(mesh),
940  explicitConnections_(0),
941  allFaceInfo_(allFaceInfo),
942  allCellInfo_(allCellInfo),
943  td_(td),
944  changedFace_(mesh_.nFaces(), false),
945  changedFaces_(mesh_.nFaces()),
946  changedCell_(mesh_.nCells(), false),
947  changedCells_(mesh_.nCells()),
948  hasCyclicPatches_(hasPatch<cyclicPolyPatch>()),
949  hasCyclicAMIPatches_
950  (
951  returnReduce(hasPatch<cyclicAMIPolyPatch>(), orOp<bool>())
952  ),
953  nEvals_(0),
954  nUnvisitedCells_(mesh_.nCells()),
955  nUnvisitedFaces_(mesh_.nFaces())
956 {
957  if
958  (
959  allFaceInfo.size() != mesh_.nFaces()
960  || allCellInfo.size() != mesh_.nCells()
961  )
962  {
964  << "face and cell storage not the size of mesh faces, cells:"
965  << endl
966  << " allFaceInfo :" << allFaceInfo.size() << endl
967  << " mesh_.nFaces():" << mesh_.nFaces() << endl
968  << " allCellInfo :" << allCellInfo.size() << endl
969  << " mesh_.nCells():" << mesh_.nCells()
970  << exit(FatalError);
971  }
972 }
973 
974 
975 template<class Type, class TrackingData>
977 (
978  const polyMesh& mesh,
979  const labelList& changedFaces,
980  const List<Type>& changedFacesInfo,
981  UList<Type>& allFaceInfo,
982  UList<Type>& allCellInfo,
983  const label maxIter,
984  TrackingData& td
985 )
986 :
987  mesh_(mesh),
988  explicitConnections_(0),
989  allFaceInfo_(allFaceInfo),
990  allCellInfo_(allCellInfo),
991  td_(td),
992  changedFace_(mesh_.nFaces(), false),
993  changedFaces_(mesh_.nFaces()),
994  changedCell_(mesh_.nCells(), false),
995  changedCells_(mesh_.nCells()),
996  hasCyclicPatches_(hasPatch<cyclicPolyPatch>()),
997  hasCyclicAMIPatches_
998  (
999  returnReduce(hasPatch<cyclicAMIPolyPatch>(), orOp<bool>())
1000  ),
1001  nEvals_(0),
1002  nUnvisitedCells_(mesh_.nCells()),
1003  nUnvisitedFaces_(mesh_.nFaces())
1004 {
1005  if
1006  (
1007  allFaceInfo.size() != mesh_.nFaces()
1008  || allCellInfo.size() != mesh_.nCells()
1009  )
1010  {
1012  << "face and cell storage not the size of mesh faces, cells:"
1013  << endl
1014  << " allFaceInfo :" << allFaceInfo.size() << endl
1015  << " mesh_.nFaces():" << mesh_.nFaces() << endl
1016  << " allCellInfo :" << allCellInfo.size() << endl
1017  << " mesh_.nCells():" << mesh_.nCells()
1018  << exit(FatalError);
1019  }
1020 
1021  // Copy initial changed faces data
1022  setFaceInfo(changedFaces, changedFacesInfo);
1023 
1024  // Iterate until nothing changes
1025  label iter = iterate(maxIter);
1026 
1027  if ((maxIter > 0) && (iter >= maxIter))
1028  {
1030  << "Maximum number of iterations reached. Increase maxIter." << endl
1031  << " maxIter:" << maxIter << endl
1032  << " nChangedCells:" << changedCells_.size() << endl
1033  << " nChangedFaces:" << changedFaces_.size() << endl
1034  << exit(FatalError);
1035  }
1036 }
1037 
1038 
1039 template<class Type, class TrackingData>
1042  const polyMesh& mesh,
1043  const labelPairList& explicitConnections,
1044  const bool handleCyclicAMI,
1045  const labelList& changedFaces,
1046  const List<Type>& changedFacesInfo,
1047  UList<Type>& allFaceInfo,
1048  UList<Type>& allCellInfo,
1049  const label maxIter,
1050  TrackingData& td
1051 )
1052 :
1053  mesh_(mesh),
1054  explicitConnections_(explicitConnections),
1055  allFaceInfo_(allFaceInfo),
1056  allCellInfo_(allCellInfo),
1057  td_(td),
1058  changedFace_(mesh_.nFaces(), false),
1059  changedFaces_(mesh_.nFaces()),
1060  changedCell_(mesh_.nCells(), false),
1061  changedCells_(mesh_.nCells()),
1062  hasCyclicPatches_(hasPatch<cyclicPolyPatch>()),
1063  hasCyclicAMIPatches_
1064  (
1065  handleCyclicAMI
1066  && returnReduce(hasPatch<cyclicAMIPolyPatch>(), orOp<bool>())
1067  ),
1068  nEvals_(0),
1069  nUnvisitedCells_(mesh_.nCells()),
1070  nUnvisitedFaces_(mesh_.nFaces())
1071 {
1072  if
1073  (
1074  allFaceInfo.size() != mesh_.nFaces()
1075  || allCellInfo.size() != mesh_.nCells()
1076  )
1077  {
1079  << "face and cell storage not the size of mesh faces, cells:"
1080  << endl
1081  << " allFaceInfo :" << allFaceInfo.size() << endl
1082  << " mesh_.nFaces():" << mesh_.nFaces() << endl
1083  << " allCellInfo :" << allCellInfo.size() << endl
1084  << " mesh_.nCells():" << mesh_.nCells()
1085  << exit(FatalError);
1086  }
1087 
1088  // Copy initial changed faces data
1089  setFaceInfo(changedFaces, changedFacesInfo);
1090 
1091  // Iterate until nothing changes
1092  label iter = iterate(maxIter);
1093 
1094  if ((maxIter > 0) && (iter >= maxIter))
1095  {
1097  << "Maximum number of iterations reached. Increase maxIter." << endl
1098  << " maxIter:" << maxIter << endl
1099  << " nChangedCells:" << changedCells_.size() << endl
1100  << " nChangedFaces:" << changedFaces_.size() << endl
1101  << exit(FatalError);
1102  }
1103 }
1104 
1105 
1106 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
1107 
1108 template<class Type, class TrackingData>
1110 {
1111  return nUnvisitedCells_;
1112 }
1113 
1114 
1115 template<class Type, class TrackingData>
1117 {
1118  return nUnvisitedFaces_;
1119 }
1120 
1121 
1122 template<class Type, class TrackingData>
1124 {
1125  // Propagate face to cell
1126 
1127  const labelList& owner = mesh_.faceOwner();
1128  const labelList& neighbour = mesh_.faceNeighbour();
1129  label nInternalFaces = mesh_.nInternalFaces();
1130 
1131  forAll(changedFaces_, changedFacei)
1132  {
1133  label facei = changedFaces_[changedFacei];
1134  if (!changedFace_[facei])
1135  {
1137  << "Face " << facei
1138  << " not marked as having been changed"
1139  << abort(FatalError);
1140  }
1141 
1142 
1143  const Type& neighbourWallInfo = allFaceInfo_[facei];
1144 
1145  // Evaluate all connected cells
1146 
1147  // Owner
1148  label celli = owner[facei];
1149  Type& currentWallInfo = allCellInfo_[celli];
1150 
1151  if (!currentWallInfo.equal(neighbourWallInfo, td_))
1152  {
1153  updateCell
1154  (
1155  celli,
1156  facei,
1157  neighbourWallInfo,
1158  propagationTol_,
1159  currentWallInfo
1160  );
1161  }
1162 
1163  // Neighbour.
1164  if (facei < nInternalFaces)
1165  {
1166  celli = neighbour[facei];
1167  Type& currentWallInfo2 = allCellInfo_[celli];
1168 
1169  if (!currentWallInfo2.equal(neighbourWallInfo, td_))
1170  {
1171  updateCell
1172  (
1173  celli,
1174  facei,
1175  neighbourWallInfo,
1176  propagationTol_,
1177  currentWallInfo2
1178  );
1179  }
1180  }
1181 
1182  // Reset status of face
1183  changedFace_[facei] = false;
1184  }
1185 
1186  // Handled all changed faces by now
1187  changedFaces_.clear();
1188 
1189  if (debug & 2)
1190  {
1191  Pout<< " Changed cells : " << changedCells_.size() << endl;
1192  }
1193 
1194  // Sum changedCells over all procs
1195  label totNChanged = changedCells_.size();
1196 
1197  reduce(totNChanged, sumOp<label>());
1198 
1199  return totNChanged;
1200 }
1201 
1202 
1203 template<class Type, class TrackingData>
1205 {
1206  // Propagate cell to face
1207 
1208  const cellList& cells = mesh_.cells();
1209 
1210  forAll(changedCells_, changedCelli)
1211  {
1212  label celli = changedCells_[changedCelli];
1213  if (!changedCell_[celli])
1214  {
1216  << "Cell " << celli << " not marked as having been changed"
1217  << abort(FatalError);
1218  }
1219 
1220  const Type& neighbourWallInfo = allCellInfo_[celli];
1221 
1222  // Evaluate all connected faces
1223 
1224  const labelList& faceLabels = cells[celli];
1225  forAll(faceLabels, faceLabelI)
1226  {
1227  label facei = faceLabels[faceLabelI];
1228  Type& currentWallInfo = allFaceInfo_[facei];
1229 
1230  if (!currentWallInfo.equal(neighbourWallInfo, td_))
1231  {
1232  updateFace
1233  (
1234  facei,
1235  celli,
1236  neighbourWallInfo,
1237  propagationTol_,
1238  currentWallInfo
1239  );
1240  }
1241  }
1242 
1243  // Reset status of cell
1244  changedCell_[celli] = false;
1245  }
1246 
1247  // Handled all changed cells by now
1248  changedCells_.clear();
1249 
1250 
1251  // Transfer across any explicitly provided internal connections
1252  handleExplicitConnections();
1253 
1254  if (hasCyclicPatches_)
1255  {
1256  // Transfer changed faces across cyclic halves
1257  handleCyclicPatches();
1258  }
1259 
1260  if (hasCyclicAMIPatches_)
1261  {
1262  handleAMICyclicPatches();
1263  }
1264 
1265  if (Pstream::parRun())
1266  {
1267  // Transfer changed faces from neighbouring processors.
1268  handleProcPatches();
1269  }
1270 
1271  if (debug & 2)
1272  {
1273  Pout<< " Changed faces : " << changedFaces_.size() << endl;
1274  }
1275 
1276  // Sum nChangedFaces over all procs
1277  label totNChanged = changedFaces_.size();
1278 
1279  reduce(totNChanged, sumOp<label>());
1280 
1281  return totNChanged;
1282 }
1283 
1284 
1285 // Iterate
1286 template<class Type, class TrackingData>
1288 {
1289  if (hasCyclicPatches_)
1290  {
1291  // Transfer changed faces across cyclic halves
1292  handleCyclicPatches();
1293  }
1294 
1295  if (hasCyclicAMIPatches_)
1296  {
1297  handleAMICyclicPatches();
1298  }
1299 
1300  if (Pstream::parRun())
1301  {
1302  // Transfer changed faces from neighbouring processors.
1303  handleProcPatches();
1304  }
1305 
1306  label iter = 0;
1307 
1308  while (iter < maxIter)
1309  {
1310  if (debug)
1311  {
1312  Info<< " Iteration " << iter << endl;
1313  }
1314 
1315  nEvals_ = 0;
1316 
1317  label nCells = faceToCell();
1318 
1319  if (debug)
1320  {
1321  Info<< " Total changed cells : " << nCells << endl;
1322  }
1323 
1324  if (nCells == 0)
1325  {
1326  break;
1327  }
1328 
1329  label nFaces = cellToFace();
1330 
1331  if (debug)
1332  {
1333  Info<< " Total changed faces : " << nFaces << nl
1334  << " Total evaluations : " << nEvals_ << nl
1335  << " Remaining unvisited cells: " << nUnvisitedCells_ << nl
1336  << " Remaining unvisited faces: " << nUnvisitedFaces_ << endl;
1337  }
1338 
1339  if (nFaces == 0)
1340  {
1341  break;
1342  }
1343 
1344  ++iter;
1345  }
1346 
1347  return iter;
1348 }
1349 
1350 
1351 // ************************************************************************* //
void handleAMICyclicPatches()
Merge data from across AMI cyclics.
Definition: FaceCellWave.C:718
void mergeFaceInfo(const polyPatch &patch, const label nFaces, const labelList &, const List< Type > &)
Merge received patch data into global data.
Definition: FaceCellWave.C:347
void handleExplicitConnections()
Merge data across explicitly provided local connections (usually.
Definition: FaceCellWave.C:856
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const labelList & processorPatches() const
Return list of processor patch labels.
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
const word & name() const
Return name.
Inter-processor communication reduction functions.
Vector-tensor class used to perform translations, rotations and scaling operations in 3D space...
Definition: transformer.H:83
void finishedSends(const bool block=true)
Mark all sends as having been done. This will start receives.
virtual label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
const UIndirectList< T > patchInternalList(const UList< T > &internalValues) const
Extract face cell data.
Definition: polyPatch.H:333
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void setFaceInfo(const labelList &changedFaces, const List< Type > &changedFacesInfo)
Set initial changed faces.
Definition: FaceCellWave.C:317
virtual label cellToFace()
Propagate from cell to face. Returns total number of faces.
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
bool transforms() const
Return true if the transformer transforms a type.
Definition: transformerI.H:133
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:291
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static scalar propagationTol()
Access to tolerance.
Definition: FaceCellWave.H:259
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Wave propagation of information through grid. Every iteration information goes through one layer of c...
Definition: FaceCellWave.H:76
FaceCellWave(const polyMesh &, UList< Type > &allFaceInfo, UList< Type > &allCellInfo, TrackingData &td=dummyTrackData_)
Definition: FaceCellWave.C:932
bool applyLowWeightCorrection() const
Return true if applying the low weight correction.
virtual const transformer & transform() const
Return transformation between the coupled patches.
const polyMesh & mesh() const
Access mesh.
Definition: FaceCellWave.H:346
Pre-declare related SubField type.
Definition: Field.H:60
void handleProcPatches()
Merge data from across processor boundaries.
Definition: FaceCellWave.C:512
scalar f1
Definition: createFields.H:15
Neighbour processor patch.
const List< transformer > & AMITransforms() const
Return a reference to the AMI transforms.
virtual const cyclicAMIPolyPatch & nbrPatch() const
Return a reference to the neighbour patch.
A topoSetSource to select a faceSet from cells.
Definition: cellToFace.H:53
Input inter-processor communications stream operating on external buffer.
Definition: UIPstream.H:53
scalar y
label getUnsetFaces() const
Get number of unvisited faces.
A List obtained as a section of another List.
Definition: SubList.H:53
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
void enterDomain(const polyPatch &patch, const label nFaces, const labelList &faceLabels, List< Type > &faceInfo) const
Handle leaving domain. Implementation referred to Type.
Definition: FaceCellWave.C:436
const cellShapeList & cells
bool hasPatch() const
Has cyclic patch?
Definition: FaceCellWave.C:302
const PtrList< AMIInterpolation > & AMIs() const
Return a reference to the AMI interpolators.
void transform(const tensor &rotTensor, const label nFaces, List< Type > &faceInfo)
Apply transformation to Type.
Definition: FaceCellWave.C:459
Cyclic plane patch.
A topoSetSource to select cells based on usage in faces.
Definition: faceToCell.H:49
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:178
virtual const transformer & transform() const
Return transformation between the coupled patches.
Cyclic patch for Arbitrary Mesh Interface (AMI)
int neighbProcNo() const
Return neighbour processor number.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Output inter-processor communications stream operating on external buffer.
Definition: UOPstream.H:54
const tensor & T() const
Return the transformation tensor.
Definition: transformerI.H:94
static direction size()
Return the number of elements in the VectorSpace = Ncmpts.
Definition: VectorSpaceI.H:83
static const char nl
Definition: Ostream.H:260
bool updateFace(const label facei, const label neighbourCelli, const Type &neighbourInfo, const scalar tol, Type &faceInfo)
Updates faceInfo with information from neighbour. Updates all.
Definition: FaceCellWave.C:160
bool updateCell(const label celli, const label neighbourFacei, const Type &neighbourInfo, const scalar tol, Type &cellInfo)
Updates cellInfo with information from neighbour. Updates all.
Definition: FaceCellWave.C:111
Buffers for inter-processor communications streams (UOPstream, UIPstream).
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
bool transformsPosition() const
Return true if the transformer transforms a point.
Definition: transformerI.H:146
virtual const transformer & transform() const
Return null transform between processor patches.
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
label patchi
label getChangedPatchFaces(const polyPatch &patch, const label startFacei, const label nFaces, labelList &changedPatchFaces, List< Type > &changedPatchFacesInfo) const
Extract info for single patch only.
Definition: FaceCellWave.C:381
combine(FaceCellWave< Type, TrackingData > &solver, const cyclicAMIPolyPatch &patch)
Definition: FaceCellWave.C:64
const TrackingData & data() const
Additional data to be passed into container.
Definition: FaceCellWave.H:340
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
virtual bool owner() const
Does this side own the patch?
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:309
messageStream Info
void checkCyclic(const polyPatch &pPatch) const
Debugging: check info on both sides of cyclic.
Definition: FaceCellWave.C:256
void handleCyclicPatches()
Merge data from across cyclics.
Definition: FaceCellWave.C:632
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
void leaveDomain(const polyPatch &patch, const label nFaces, const labelList &faceLabels, List< Type > &faceInfo) const
Handle leaving domain. Implementation referred to Type.
Definition: FaceCellWave.C:413
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
label size() const
Return the number of elements in the UList.
Definition: UListI.H:311
static void offset(const polyPatch &patch, const label off, const label nFaces, labelList &faces)
Offset face labels by constant value.
Definition: FaceCellWave.C:494
const dimensionedScalar e
Elementary charge.
Definition: doubleScalar.H:105
const List< T >::subList patchSlice(const UList< T > &l) const
Slice list to patch.
Definition: polyPatch.H:342
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477
virtual label faceToCell()
Propagate from face to cell. Returns total number of cells.
const vectorField::subField faceCentres() const
Return face centres.
Definition: polyPatch.C:297
label getUnsetCells() const
Get number of unvisited cells, i.e. cells that were not (yet)