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-2018 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_.neighbPatch().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).neighbPatch();
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 tensorField& rotTensor,
461  const label nFaces,
462  List<Type>& faceInfo
463 )
464 {
465  // Transform. Implementation referred to Type
466 
467  if (rotTensor.size() == 1)
468  {
469  const tensor& T = rotTensor[0];
470 
471  for (label facei = 0; facei < nFaces; facei++)
472  {
473  faceInfo[facei].transform(mesh_, T, td_);
474  }
475  }
476  else
477  {
478  for (label facei = 0; facei < nFaces; facei++)
479  {
480  faceInfo[facei].transform(mesh_, rotTensor[facei], td_);
481  }
482  }
483 }
484 
485 
486 template<class Type, class TrackingData>
488 (
489  const vectorTensorTransform& trans,
490  const label nFaces,
491  List<Type>& faceInfo
492 )
493 {
494  // Transform. Implementation referred to Type
495 
496  if (trans.hasR())
497  {
498  for (label facei = 0; facei < nFaces; facei++)
499  {
500  faceInfo[facei].transform(mesh_, trans.R(), td_);
501  }
502  }
503 }
504 
505 
506 template<class Type, class TrackingData>
508 (
509  const polyPatch&,
510  const label cycOffset,
511  const label nFaces,
512  labelList& faces
513 )
514 {
515  // Offset mesh face. Used for transferring from one cyclic half to the
516  // other.
517 
518  for (label facei = 0; facei < nFaces; facei++)
519  {
520  faces[facei] += cycOffset;
521  }
522 }
523 
524 
525 template<class Type, class TrackingData>
527 {
528  // Transfer all the information to/from neighbouring processors
529 
530  const globalMeshData& pData = mesh_.globalData();
531 
532  // Which patches are processor patches
533  const labelList& procPatches = pData.processorPatches();
534 
535  // Send all
536 
538 
539  forAll(procPatches, i)
540  {
541  label patchi = procPatches[i];
542 
543  const processorPolyPatch& procPatch =
544  refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchi]);
545 
546  // Allocate buffers
547  label nSendFaces;
548  labelList sendFaces(procPatch.size());
549  List<Type> sendFacesInfo(procPatch.size());
550 
551  // Determine which faces changed on current patch
552  nSendFaces = getChangedPatchFaces
553  (
554  procPatch,
555  0,
556  procPatch.size(),
557  sendFaces,
558  sendFacesInfo
559  );
560 
561  // Adapt wallInfo for leaving domain
562  leaveDomain
563  (
564  procPatch,
565  nSendFaces,
566  sendFaces,
567  sendFacesInfo
568  );
569 
570  if (debug & 2)
571  {
572  Pout<< " Processor patch " << patchi << ' ' << procPatch.name()
573  << " communicating with " << procPatch.neighbProcNo()
574  << " Sending:" << nSendFaces
575  << endl;
576  }
577 
578  UOPstream toNeighbour(procPatch.neighbProcNo(), pBufs);
579  // writeFaces(nSendFaces, sendFaces, sendFacesInfo, toNeighbour);
580  toNeighbour
581  << SubList<label>(sendFaces, nSendFaces)
582  << SubList<Type>(sendFacesInfo, nSendFaces);
583  }
584 
585  pBufs.finishedSends();
586 
587  // Receive all
588 
589  forAll(procPatches, i)
590  {
591  label patchi = procPatches[i];
592 
593  const processorPolyPatch& procPatch =
594  refCast<const processorPolyPatch>(mesh_.boundaryMesh()[patchi]);
595 
596  // Allocate buffers
597  labelList receiveFaces;
598  List<Type> receiveFacesInfo;
599 
600  {
601  UIPstream fromNeighbour(procPatch.neighbProcNo(), pBufs);
602  fromNeighbour >> receiveFaces >> receiveFacesInfo;
603  }
604 
605  if (debug & 2)
606  {
607  Pout<< " Processor patch " << patchi << ' ' << procPatch.name()
608  << " communicating with " << procPatch.neighbProcNo()
609  << " Receiving:" << receiveFaces.size()
610  << endl;
611  }
612 
613  // Apply transform to received data for non-parallel planes
614  if (!procPatch.parallel())
615  {
616  transform
617  (
618  procPatch.forwardT(),
619  receiveFaces.size(),
620  receiveFacesInfo
621  );
622  }
623 
624  // Adapt wallInfo for entering domain
625  enterDomain
626  (
627  procPatch,
628  receiveFaces.size(),
629  receiveFaces,
630  receiveFacesInfo
631  );
632 
633  // Merge received info
634  mergeFaceInfo
635  (
636  procPatch,
637  receiveFaces.size(),
638  receiveFaces,
639  receiveFacesInfo
640  );
641  }
642 }
643 
644 
645 template<class Type, class TrackingData>
647 {
648  // Transfer information across cyclic halves.
649 
650  forAll(mesh_.boundaryMesh(), patchi)
651  {
652  const polyPatch& patch = mesh_.boundaryMesh()[patchi];
653 
654  if (isA<cyclicPolyPatch>(patch))
655  {
656  const cyclicPolyPatch& nbrPatch =
657  refCast<const cyclicPolyPatch>(patch).neighbPatch();
658 
659  // Allocate buffers
660  label nReceiveFaces;
661  labelList receiveFaces(patch.size());
662  List<Type> receiveFacesInfo(patch.size());
663 
664  // Determine which faces changed
665  nReceiveFaces = getChangedPatchFaces
666  (
667  nbrPatch,
668  0,
669  nbrPatch.size(),
670  receiveFaces,
671  receiveFacesInfo
672  );
673 
674  // Adapt wallInfo for leaving domain
675  leaveDomain
676  (
677  nbrPatch,
678  nReceiveFaces,
679  receiveFaces,
680  receiveFacesInfo
681  );
682 
683  const cyclicPolyPatch& cycPatch =
684  refCast<const cyclicPolyPatch>(patch);
685 
686  if (!cycPatch.parallel())
687  {
688  // received data from other half
689  transform
690  (
691  cycPatch.forwardT(),
692  nReceiveFaces,
693  receiveFacesInfo
694  );
695  }
696 
697  if (debug & 2)
698  {
699  Pout<< " Cyclic patch " << patchi << ' ' << cycPatch.name()
700  << " Changed : " << nReceiveFaces
701  << endl;
702  }
703 
704  // Half2: Adapt wallInfo for entering domain
705  enterDomain
706  (
707  cycPatch,
708  nReceiveFaces,
709  receiveFaces,
710  receiveFacesInfo
711  );
712 
713  // Merge into global storage
714  mergeFaceInfo
715  (
716  cycPatch,
717  nReceiveFaces,
718  receiveFaces,
719  receiveFacesInfo
720  );
721 
722  if (debug)
723  {
724  checkCyclic(cycPatch);
725  }
726  }
727  }
728 }
729 
730 
731 template<class Type, class TrackingData>
733 {
734  // Transfer information across cyclicAMI halves.
735 
736  forAll(mesh_.boundaryMesh(), patchi)
737  {
738  const polyPatch& patch = mesh_.boundaryMesh()[patchi];
739 
740  if (isA<cyclicAMIPolyPatch>(patch))
741  {
742  const cyclicAMIPolyPatch& cycPatch =
743  refCast<const cyclicAMIPolyPatch>(patch);
744 
745  List<Type> receiveInfo;
746 
747  {
748  const cyclicAMIPolyPatch& nbrPatch =
749  refCast<const cyclicAMIPolyPatch>(patch).neighbPatch();
750 
751  // Get nbrPatch data (so not just changed faces)
752  typename List<Type>::subList sendInfo
753  (
754  nbrPatch.patchSlice
755  (
756  allFaceInfo_
757  )
758  );
759 
760  if (!nbrPatch.parallel() || nbrPatch.separated())
761  {
762  // Adapt sendInfo for leaving domain
763  const vectorField::subField fc = nbrPatch.faceCentres();
764  forAll(sendInfo, i)
765  {
766  sendInfo[i].leaveDomain(mesh_, nbrPatch, i, fc[i], td_);
767  }
768  }
769 
770  // Transfer sendInfo to cycPatch
771  combine<Type, TrackingData> cmb(*this, cycPatch);
772 
773  List<Type> defVals;
774  if (cycPatch.applyLowWeightCorrection())
775  {
776  defVals = cycPatch.patchInternalList(allCellInfo_);
777  }
778 
779  if (cycPatch.owner())
780  {
781  forAll(cycPatch.AMIs(), i)
782  {
783  List<Type> sendInfoT(sendInfo);
784  transform
785  (
786  cycPatch.AMITransforms()[i],
787  sendInfoT.size(),
788  sendInfoT
789  );
790  cycPatch.AMIs()[i].interpolateToSource
791  (
792  sendInfoT,
793  cmb,
794  receiveInfo,
795  defVals
796  );
797  }
798  }
799  else
800  {
801  forAll(cycPatch.neighbPatch().AMIs(), i)
802  {
803  List<Type> sendInfoT(sendInfo);
804  transform
805  (
806  cycPatch.neighbPatch().AMITransforms()[i],
807  sendInfoT.size(),
808  sendInfoT
809  );
810  cycPatch.neighbPatch().AMIs()[i].interpolateToTarget
811  (
812  sendInfoT,
813  cmb,
814  receiveInfo,
815  defVals
816  );
817  }
818  }
819  }
820 
821  // Apply transform to received data for non-parallel planes
822  if (!cycPatch.parallel())
823  {
824  transform
825  (
826  cycPatch.forwardT(),
827  receiveInfo.size(),
828  receiveInfo
829  );
830  }
831 
832  if (!cycPatch.parallel() || cycPatch.separated())
833  {
834  // Adapt receiveInfo for entering domain
835  const vectorField::subField fc = cycPatch.faceCentres();
836  forAll(receiveInfo, i)
837  {
838  receiveInfo[i].enterDomain(mesh_, cycPatch, i, fc[i], td_);
839  }
840  }
841 
842  // Merge into global storage
843  forAll(receiveInfo, i)
844  {
845  label meshFacei = cycPatch.start()+i;
846 
847  Type& currentWallInfo = allFaceInfo_[meshFacei];
848 
849  if
850  (
851  receiveInfo[i].valid(td_)
852  && !currentWallInfo.equal(receiveInfo[i], td_)
853  )
854  {
855  updateFace
856  (
857  meshFacei,
858  receiveInfo[i],
859  propagationTol_,
860  currentWallInfo
861  );
862  }
863  }
864  }
865  }
866 }
867 
868 
869 template<class Type, class TrackingData>
871 {
872  // Collect changed information
873 
874  DynamicList<label> f0Baffle(explicitConnections_.size());
875  DynamicList<Type> f0Info(explicitConnections_.size());
876 
877  DynamicList<label> f1Baffle(explicitConnections_.size());
878  DynamicList<Type> f1Info(explicitConnections_.size());
879 
880  forAll(explicitConnections_, connI)
881  {
882  const labelPair& baffle = explicitConnections_[connI];
883 
884  label f0 = baffle[0];
885  if (changedFace_[f0])
886  {
887  f0Baffle.append(connI);
888  f0Info.append(allFaceInfo_[f0]);
889  }
890 
891  label f1 = baffle[1];
892  if (changedFace_[f1])
893  {
894  f1Baffle.append(connI);
895  f1Info.append(allFaceInfo_[f1]);
896  }
897  }
898 
899 
900  // Update other side with changed information
901 
902  forAll(f1Info, index)
903  {
904  const labelPair& baffle = explicitConnections_[f1Baffle[index]];
905 
906  label f0 = baffle[0];
907  Type& currentWallInfo = allFaceInfo_[f0];
908 
909  if (!currentWallInfo.equal(f1Info[index], td_))
910  {
911  updateFace
912  (
913  f0,
914  f1Info[index],
915  propagationTol_,
916  currentWallInfo
917  );
918  }
919  }
920 
921  forAll(f0Info, index)
922  {
923  const labelPair& baffle = explicitConnections_[f0Baffle[index]];
924 
925  label f1 = baffle[1];
926  Type& currentWallInfo = allFaceInfo_[f1];
927 
928  if (!currentWallInfo.equal(f0Info[index], td_))
929  {
930  updateFace
931  (
932  f1,
933  f0Info[index],
934  propagationTol_,
935  currentWallInfo
936  );
937  }
938  }
939 }
940 
941 
942 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
943 
944 template<class Type, class TrackingData>
946 (
947  const polyMesh& mesh,
948  UList<Type>& allFaceInfo,
949  UList<Type>& allCellInfo,
950  TrackingData& td
951 )
952 :
953  mesh_(mesh),
954  explicitConnections_(0),
955  allFaceInfo_(allFaceInfo),
956  allCellInfo_(allCellInfo),
957  td_(td),
958  changedFace_(mesh_.nFaces(), false),
959  changedFaces_(mesh_.nFaces()),
960  changedCell_(mesh_.nCells(), false),
961  changedCells_(mesh_.nCells()),
962  hasCyclicPatches_(hasPatch<cyclicPolyPatch>()),
963  hasCyclicAMIPatches_
964  (
965  returnReduce(hasPatch<cyclicAMIPolyPatch>(), orOp<bool>())
966  ),
967  nEvals_(0),
968  nUnvisitedCells_(mesh_.nCells()),
969  nUnvisitedFaces_(mesh_.nFaces())
970 {
971  if
972  (
973  allFaceInfo.size() != mesh_.nFaces()
974  || allCellInfo.size() != mesh_.nCells()
975  )
976  {
978  << "face and cell storage not the size of mesh faces, cells:"
979  << endl
980  << " allFaceInfo :" << allFaceInfo.size() << endl
981  << " mesh_.nFaces():" << mesh_.nFaces() << endl
982  << " allCellInfo :" << allCellInfo.size() << endl
983  << " mesh_.nCells():" << mesh_.nCells()
984  << exit(FatalError);
985  }
986 }
987 
988 
989 template<class Type, class TrackingData>
991 (
992  const polyMesh& mesh,
993  const labelList& changedFaces,
994  const List<Type>& changedFacesInfo,
995  UList<Type>& allFaceInfo,
996  UList<Type>& allCellInfo,
997  const label maxIter,
998  TrackingData& td
999 )
1000 :
1001  mesh_(mesh),
1002  explicitConnections_(0),
1003  allFaceInfo_(allFaceInfo),
1004  allCellInfo_(allCellInfo),
1005  td_(td),
1006  changedFace_(mesh_.nFaces(), false),
1007  changedFaces_(mesh_.nFaces()),
1008  changedCell_(mesh_.nCells(), false),
1009  changedCells_(mesh_.nCells()),
1010  hasCyclicPatches_(hasPatch<cyclicPolyPatch>()),
1011  hasCyclicAMIPatches_
1012  (
1013  returnReduce(hasPatch<cyclicAMIPolyPatch>(), orOp<bool>())
1014  ),
1015  nEvals_(0),
1016  nUnvisitedCells_(mesh_.nCells()),
1017  nUnvisitedFaces_(mesh_.nFaces())
1018 {
1019  if
1020  (
1021  allFaceInfo.size() != mesh_.nFaces()
1022  || allCellInfo.size() != mesh_.nCells()
1023  )
1024  {
1026  << "face and cell storage not the size of mesh faces, cells:"
1027  << endl
1028  << " allFaceInfo :" << allFaceInfo.size() << endl
1029  << " mesh_.nFaces():" << mesh_.nFaces() << endl
1030  << " allCellInfo :" << allCellInfo.size() << endl
1031  << " mesh_.nCells():" << mesh_.nCells()
1032  << exit(FatalError);
1033  }
1034 
1035  // Copy initial changed faces data
1036  setFaceInfo(changedFaces, changedFacesInfo);
1037 
1038  // Iterate until nothing changes
1039  label iter = iterate(maxIter);
1040 
1041  if ((maxIter > 0) && (iter >= maxIter))
1042  {
1044  << "Maximum number of iterations reached. Increase maxIter." << endl
1045  << " maxIter:" << maxIter << endl
1046  << " nChangedCells:" << changedCells_.size() << endl
1047  << " nChangedFaces:" << changedFaces_.size() << endl
1048  << exit(FatalError);
1049  }
1050 }
1051 
1052 
1053 template<class Type, class TrackingData>
1056  const polyMesh& mesh,
1057  const labelPairList& explicitConnections,
1058  const bool handleCyclicAMI,
1059  const labelList& changedFaces,
1060  const List<Type>& changedFacesInfo,
1061  UList<Type>& allFaceInfo,
1062  UList<Type>& allCellInfo,
1063  const label maxIter,
1064  TrackingData& td
1065 )
1066 :
1067  mesh_(mesh),
1068  explicitConnections_(explicitConnections),
1069  allFaceInfo_(allFaceInfo),
1070  allCellInfo_(allCellInfo),
1071  td_(td),
1072  changedFace_(mesh_.nFaces(), false),
1073  changedFaces_(mesh_.nFaces()),
1074  changedCell_(mesh_.nCells(), false),
1075  changedCells_(mesh_.nCells()),
1076  hasCyclicPatches_(hasPatch<cyclicPolyPatch>()),
1077  hasCyclicAMIPatches_
1078  (
1079  handleCyclicAMI
1080  && returnReduce(hasPatch<cyclicAMIPolyPatch>(), orOp<bool>())
1081  ),
1082  nEvals_(0),
1083  nUnvisitedCells_(mesh_.nCells()),
1084  nUnvisitedFaces_(mesh_.nFaces())
1085 {
1086  if
1087  (
1088  allFaceInfo.size() != mesh_.nFaces()
1089  || allCellInfo.size() != mesh_.nCells()
1090  )
1091  {
1093  << "face and cell storage not the size of mesh faces, cells:"
1094  << endl
1095  << " allFaceInfo :" << allFaceInfo.size() << endl
1096  << " mesh_.nFaces():" << mesh_.nFaces() << endl
1097  << " allCellInfo :" << allCellInfo.size() << endl
1098  << " mesh_.nCells():" << mesh_.nCells()
1099  << exit(FatalError);
1100  }
1101 
1102  // Copy initial changed faces data
1103  setFaceInfo(changedFaces, changedFacesInfo);
1104 
1105  // Iterate until nothing changes
1106  label iter = iterate(maxIter);
1107 
1108  if ((maxIter > 0) && (iter >= maxIter))
1109  {
1111  << "Maximum number of iterations reached. Increase maxIter." << endl
1112  << " maxIter:" << maxIter << endl
1113  << " nChangedCells:" << changedCells_.size() << endl
1114  << " nChangedFaces:" << changedFaces_.size() << endl
1115  << exit(FatalError);
1116  }
1117 }
1118 
1119 
1120 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
1121 
1122 template<class Type, class TrackingData>
1124 {
1125  return nUnvisitedCells_;
1126 }
1127 
1128 
1129 template<class Type, class TrackingData>
1131 {
1132  return nUnvisitedFaces_;
1133 }
1134 
1135 
1136 template<class Type, class TrackingData>
1138 {
1139  // Propagate face to cell
1140 
1141  const labelList& owner = mesh_.faceOwner();
1142  const labelList& neighbour = mesh_.faceNeighbour();
1143  label nInternalFaces = mesh_.nInternalFaces();
1144 
1145  forAll(changedFaces_, changedFacei)
1146  {
1147  label facei = changedFaces_[changedFacei];
1148  if (!changedFace_[facei])
1149  {
1151  << "Face " << facei
1152  << " not marked as having been changed"
1153  << abort(FatalError);
1154  }
1155 
1156 
1157  const Type& neighbourWallInfo = allFaceInfo_[facei];
1158 
1159  // Evaluate all connected cells
1160 
1161  // Owner
1162  label celli = owner[facei];
1163  Type& currentWallInfo = allCellInfo_[celli];
1164 
1165  if (!currentWallInfo.equal(neighbourWallInfo, td_))
1166  {
1167  updateCell
1168  (
1169  celli,
1170  facei,
1171  neighbourWallInfo,
1172  propagationTol_,
1173  currentWallInfo
1174  );
1175  }
1176 
1177  // Neighbour.
1178  if (facei < nInternalFaces)
1179  {
1180  celli = neighbour[facei];
1181  Type& currentWallInfo2 = allCellInfo_[celli];
1182 
1183  if (!currentWallInfo2.equal(neighbourWallInfo, td_))
1184  {
1185  updateCell
1186  (
1187  celli,
1188  facei,
1189  neighbourWallInfo,
1190  propagationTol_,
1191  currentWallInfo2
1192  );
1193  }
1194  }
1195 
1196  // Reset status of face
1197  changedFace_[facei] = false;
1198  }
1199 
1200  // Handled all changed faces by now
1201  changedFaces_.clear();
1202 
1203  if (debug & 2)
1204  {
1205  Pout<< " Changed cells : " << changedCells_.size() << endl;
1206  }
1207 
1208  // Sum changedCells over all procs
1209  label totNChanged = changedCells_.size();
1210 
1211  reduce(totNChanged, sumOp<label>());
1212 
1213  return totNChanged;
1214 }
1215 
1216 
1217 template<class Type, class TrackingData>
1219 {
1220  // Propagate cell to face
1221 
1222  const cellList& cells = mesh_.cells();
1223 
1224  forAll(changedCells_, changedCelli)
1225  {
1226  label celli = changedCells_[changedCelli];
1227  if (!changedCell_[celli])
1228  {
1230  << "Cell " << celli << " not marked as having been changed"
1231  << abort(FatalError);
1232  }
1233 
1234  const Type& neighbourWallInfo = allCellInfo_[celli];
1235 
1236  // Evaluate all connected faces
1237 
1238  const labelList& faceLabels = cells[celli];
1239  forAll(faceLabels, faceLabelI)
1240  {
1241  label facei = faceLabels[faceLabelI];
1242  Type& currentWallInfo = allFaceInfo_[facei];
1243 
1244  if (!currentWallInfo.equal(neighbourWallInfo, td_))
1245  {
1246  updateFace
1247  (
1248  facei,
1249  celli,
1250  neighbourWallInfo,
1251  propagationTol_,
1252  currentWallInfo
1253  );
1254  }
1255  }
1256 
1257  // Reset status of cell
1258  changedCell_[celli] = false;
1259  }
1260 
1261  // Handled all changed cells by now
1262  changedCells_.clear();
1263 
1264 
1265  // Transfer across any explicitly provided internal connections
1266  handleExplicitConnections();
1267 
1268  if (hasCyclicPatches_)
1269  {
1270  // Transfer changed faces across cyclic halves
1271  handleCyclicPatches();
1272  }
1273 
1274  if (hasCyclicAMIPatches_)
1275  {
1276  handleAMICyclicPatches();
1277  }
1278 
1279  if (Pstream::parRun())
1280  {
1281  // Transfer changed faces from neighbouring processors.
1282  handleProcPatches();
1283  }
1284 
1285  if (debug & 2)
1286  {
1287  Pout<< " Changed faces : " << changedFaces_.size() << endl;
1288  }
1289 
1290  // Sum nChangedFaces over all procs
1291  label totNChanged = changedFaces_.size();
1292 
1293  reduce(totNChanged, sumOp<label>());
1294 
1295  return totNChanged;
1296 }
1297 
1298 
1299 // Iterate
1300 template<class Type, class TrackingData>
1302 {
1303  if (hasCyclicPatches_)
1304  {
1305  // Transfer changed faces across cyclic halves
1306  handleCyclicPatches();
1307  }
1308 
1309  if (hasCyclicAMIPatches_)
1310  {
1311  handleAMICyclicPatches();
1312  }
1313 
1314  if (Pstream::parRun())
1315  {
1316  // Transfer changed faces from neighbouring processors.
1317  handleProcPatches();
1318  }
1319 
1320  label iter = 0;
1321 
1322  while (iter < maxIter)
1323  {
1324  if (debug)
1325  {
1326  Info<< " Iteration " << iter << endl;
1327  }
1328 
1329  nEvals_ = 0;
1330 
1331  label nCells = faceToCell();
1332 
1333  if (debug)
1334  {
1335  Info<< " Total changed cells : " << nCells << endl;
1336  }
1337 
1338  if (nCells == 0)
1339  {
1340  break;
1341  }
1342 
1343  label nFaces = cellToFace();
1344 
1345  if (debug)
1346  {
1347  Info<< " Total changed faces : " << nFaces << nl
1348  << " Total evaluations : " << nEvals_ << nl
1349  << " Remaining unvisited cells: " << nUnvisitedCells_ << nl
1350  << " Remaining unvisited faces: " << nUnvisitedFaces_ << endl;
1351  }
1352 
1353  if (nFaces == 0)
1354  {
1355  break;
1356  }
1357 
1358  ++iter;
1359  }
1360 
1361  return iter;
1362 }
1363 
1364 
1365 // ************************************************************************* //
void handleAMICyclicPatches()
Merge data from across AMI cyclics.
Definition: FaceCellWave.C:732
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:870
#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.
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.
virtual bool separated() const
Are the planes separated.
const UIndirectList< T > patchInternalList(const UList< T > &internalValues) const
Extract face cell data.
Definition: polyPatch.H:327
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:319
Vector-tensor class used to perform translations and rotations in 3D space.
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:278
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
static scalar propagationTol()
Access to tolerance.
Definition: FaceCellWave.H:259
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
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:946
bool applyLowWeightCorrection() const
Return true if applying the low weight correction.
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:526
scalar f1
Definition: createFields.H:28
virtual const tensorField & forwardT() const
Return face transformation tensor.
Neighbour processor patch.
const List< vectorTensorTransform > & AMITransforms() const
Return a reference to the AMI transforms.
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
virtual bool parallel() const
Are the cyclic planes parallel.
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
volScalarField & e
Elementary charge.
Definition: createFields.H:11
bool hasPatch() const
Has cyclic patch?
Definition: FaceCellWave.C:302
const PtrList< AMIInterpolation > & AMIs() const
Return a reference to the AMI interpolators.
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:177
void transform(const tensorField &rotTensor, const label nFaces, List< Type > &faceInfo)
Apply transformation to Type.
Definition: FaceCellWave.C:459
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
static const char nl
Definition: Ostream.H:265
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 T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
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
virtual const cyclicAMIPolyPatch & neighbPatch() const
Return a reference to the neighbour patch.
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:303
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:646
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:299
static void offset(const polyPatch &patch, const label off, const label nFaces, labelList &faces)
Offset face labels by constant value.
Definition: FaceCellWave.C:508
const List< T >::subList patchSlice(const UList< T > &l) const
Slice list to patch.
Definition: polyPatch.H:336
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:284
label getUnsetCells() const
Get number of unvisited cells, i.e. cells that were not (yet)