FreeFOAM The Cross-Platform CFD Toolkit
refinementHistory.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) 1991-2010 OpenCFD Ltd.
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 <OpenFOAM/DynamicList.H>
27 #include "refinementHistory.H"
28 #include <OpenFOAM/ListOps.H>
29 #include <OpenFOAM/mapPolyMesh.H>
31 #include <OpenFOAM/polyMesh.H>
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 
38 defineTypeNameAndDebug(refinementHistory, 0);
39 
40 }
41 
42 
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
44 
45 void Foam::refinementHistory::writeEntry
46 (
47  const List<splitCell8>& splitCells,
48  const splitCell8& split
49 )
50 {
51  // Write me:
52  if (split.addedCellsPtr_.valid())
53  {
54  Pout<< "parent:" << split.parent_
55  << " subCells:" << split.addedCellsPtr_()
56  << endl;
57  }
58  else
59  {
60  Pout<< "parent:" << split.parent_
61  << " no subcells"
62  << endl;
63  }
64 
65  if (split.parent_ >= 0)
66  {
67  Pout<< "parent data:" << endl;
68  // Write my parent
69  string oldPrefix = Pout.prefix();
70  Pout.prefix() = " " + oldPrefix;
71  writeEntry(splitCells, splitCells[split.parent_]);
72  Pout.prefix() = oldPrefix;
73  }
74 }
75 
76 
78 (
79  const labelList& visibleCells,
80  const List<splitCell8>& splitCells
81 )
82 {
83  string oldPrefix = Pout.prefix();
84  Pout.prefix() = "";
85 
86  forAll(visibleCells, cellI)
87  {
88  label index = visibleCells[cellI];
89 
90  if (index >= 0)
91  {
92  Pout<< "Cell from refinement:" << cellI << " index:" << index
93  << endl;
94 
95  string oldPrefix = Pout.prefix();
96  Pout.prefix() = " " + oldPrefix;
97  writeEntry(splitCells, splitCells[index]);
98  Pout.prefix() = oldPrefix;
99  }
100  else
101  {
102  Pout<< "Unrefined cell:" << cellI << " index:" << index << endl;
103  }
104  }
105  Pout.prefix() = oldPrefix;
106 }
107 
108 
109 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
110 
111 //- Construct null
113 :
114  parent_(-1),
115  addedCellsPtr_(NULL)
116 {}
117 
118 
119 //- Construct as child element of parent
121 :
122  parent_(parent),
123  addedCellsPtr_(NULL)
124 {}
125 
126 
127 //- Construct from Istream
129 {
130  is >> *this;
131 }
132 
133 
134 //- Construct as (deep) copy.
136 :
137  parent_(sc.parent_),
138  addedCellsPtr_
139  (
140  sc.addedCellsPtr_.valid()
141  ? new FixedList<label, 8>(sc.addedCellsPtr_())
142  : NULL
143  )
144 {}
145 
146 
147 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
148 
150 {
151  labelList addedCells;
152 
153  is >> sc.parent_ >> addedCells;
154 
155  if (addedCells.size())
156  {
157  sc.addedCellsPtr_.reset(new FixedList<label, 8>(addedCells));
158  }
159  else
160  {
161  sc.addedCellsPtr_.reset(NULL);
162  }
163 
164  return is;
165 }
166 
167 
168 Foam::Ostream& Foam::operator<<
169 (
170  Ostream& os,
171  const refinementHistory::splitCell8& sc
172 )
173 {
174  // Output as labelList so we can have 0 sized lists. Alternative is to
175  // output as fixedlist with e.g. -1 elements and check for this upon
176  // reading. However would cause much more data to be transferred.
177 
178  if (sc.addedCellsPtr_.valid())
179  {
180  return os
181  << sc.parent_
182  << token::SPACE
183  << labelList(sc.addedCellsPtr_());
184  }
185  else
186  {
187  return os << sc.parent_ << token::SPACE << labelList(0);
188  }
189 }
190 
191 
192 void Foam::refinementHistory::checkIndices() const
193 {
194  // Check indices.
195  forAll(visibleCells_, i)
196  {
197  if (visibleCells_[i] < 0 && visibleCells_[i] >= splitCells_.size())
198  {
199  FatalErrorIn("refinementHistory::checkIndices() const")
200  << "Illegal entry " << visibleCells_[i]
201  << " in visibleCells at location" << i << nl
202  << "It points outside the range of splitCells : 0.."
203  << splitCells_.size()-1
204  << abort(FatalError);
205  }
206  }
207 }
208 
209 
210 Foam::label Foam::refinementHistory::allocateSplitCell
211 (
212  const label parent,
213  const label i
214 )
215 {
216  label index = -1;
217 
218  if (freeSplitCells_.size())
219  {
220  index = freeSplitCells_.remove();
221 
222  splitCells_[index] = splitCell8(parent);
223  }
224  else
225  {
226  index = splitCells_.size();
227 
228  splitCells_.append(splitCell8(parent));
229  }
230 
231 
232  // Update the parent field
233  if (parent >= 0)
234  {
235  splitCell8& parentSplit = splitCells_[parent];
236 
237  if (parentSplit.addedCellsPtr_.empty())
238  {
239  // Allocate storage on parent for the 8 subcells.
240  parentSplit.addedCellsPtr_.reset(new FixedList<label, 8>(-1));
241  }
242 
243 
244  // Store me on my parent
245  FixedList<label, 8>& parentSplits = parentSplit.addedCellsPtr_();
246 
247  parentSplits[i] = index;
248  }
249 
250  return index;
251 }
252 
253 
254 void Foam::refinementHistory::freeSplitCell(const label index)
255 {
256  splitCell8& split = splitCells_[index];
257 
258  // Make sure parent does not point to me anymore.
259  if (split.parent_ >= 0)
260  {
261  autoPtr<FixedList<label, 8> >& subCellsPtr =
262  splitCells_[split.parent_].addedCellsPtr_;
263 
264  if (subCellsPtr.valid())
265  {
266  FixedList<label, 8>& subCells = subCellsPtr();
267 
268  label myPos = findIndex(subCells, index);
269 
270  if (myPos == -1)
271  {
272  FatalErrorIn("refinementHistory::freeSplitCell")
273  << "Problem: cannot find myself in"
274  << " parents' children" << abort(FatalError);
275  }
276  else
277  {
278  subCells[myPos] = -1;
279  }
280  }
281  }
282 
283  // Mark splitCell as free
284  split.parent_ = -2;
285 
286  // Add to cache of free splitCells
287  freeSplitCells_.append(index);
288 }
289 
290 
291 // Mark entry in splitCells. Recursively mark its parent and subs.
292 void Foam::refinementHistory::markSplit
293 (
294  const label index,
295  labelList& oldToNew,
296  DynamicList<splitCell8>& newSplitCells
297 ) const
298 {
299  if (oldToNew[index] == -1)
300  {
301  // Not yet compacted.
302 
303  const splitCell8& split = splitCells_[index];
304 
305  oldToNew[index] = newSplitCells.size();
306  newSplitCells.append(split);
307 
308  if (split.parent_ >= 0)
309  {
310  markSplit(split.parent_, oldToNew, newSplitCells);
311  }
312  if (split.addedCellsPtr_.valid())
313  {
314  const FixedList<label, 8>& splits = split.addedCellsPtr_();
315 
316  forAll(splits, i)
317  {
318  if (splits[i] >= 0)
319  {
320  markSplit(splits[i], oldToNew, newSplitCells);
321  }
322  }
323  }
324  }
325 }
326 
327 
328 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
329 
331 :
332  regIOobject(io)
333 {
334  if
335  (
337  || (io.readOpt() == IOobject::READ_IF_PRESENT && headerOk())
338  )
339  {
340  readStream(typeName) >> *this;
341  close();
342  }
343 
344  if (debug)
345  {
346  Pout<< "refinementHistory::refinementHistory :"
347  << " constructed history from IOobject :"
348  << " splitCells:" << splitCells_.size()
349  << " visibleCells:" << visibleCells_.size()
350  << endl;
351  }
352 }
353 
354 
355 //- Read or construct
357 (
358  const IOobject& io,
359  const List<splitCell8>& splitCells,
360  const labelList& visibleCells
361 )
362 :
363  regIOobject(io),
364  splitCells_(splitCells),
365  freeSplitCells_(0),
366  visibleCells_(visibleCells)
367 {
368  if
369  (
371  || (io.readOpt() == IOobject::READ_IF_PRESENT && headerOk())
372  )
373  {
374  readStream(typeName) >> *this;
375  close();
376  }
377 
378  // Check indices.
379  checkIndices();
380 
381  if (debug)
382  {
383  Pout<< "refinementHistory::refinementHistory :"
384  << " constructed history from IOobject or components :"
385  << " splitCells:" << splitCells_.size()
386  << " visibleCells:" << visibleCells_.size()
387  << endl;
388  }
389 }
390 
391 
392 // Construct from initial number of cells (all visible)
394 (
395  const IOobject& io,
396  const label nCells
397 )
398 :
399  regIOobject(io),
400  freeSplitCells_(0)
401 {
402  if
403  (
405  || (io.readOpt() == IOobject::READ_IF_PRESENT && headerOk())
406  )
407  {
408  readStream(typeName) >> *this;
409  close();
410  }
411  else
412  {
413  visibleCells_.setSize(nCells);
414  splitCells_.setCapacity(nCells);
415 
416  for (label cellI = 0; cellI < nCells; cellI++)
417  {
418  visibleCells_[cellI] = cellI;
419  splitCells_.append(splitCell8());
420  }
421  }
422 
423  // Check indices.
424  checkIndices();
425 
426  if (debug)
427  {
428  Pout<< "refinementHistory::refinementHistory :"
429  << " constructed history from IOobject or initial size :"
430  << " splitCells:" << splitCells_.size()
431  << " visibleCells:" << visibleCells_.size()
432  << endl;
433  }
434 }
435 
436 
437 // Construct as copy
439 (
440  const IOobject& io,
441  const refinementHistory& rh
442 )
443 :
444  regIOobject(io),
445  splitCells_(rh.splitCells()),
446  freeSplitCells_(rh.freeSplitCells()),
447  visibleCells_(rh.visibleCells())
448 {
449  if (debug)
450  {
451  Pout<< "refinementHistory::refinementHistory : constructed initial"
452  << " history." << endl;
453  }
454 }
455 
456 
457 // Construct from Istream
459 :
460  regIOobject(io),
461  splitCells_(is),
462  freeSplitCells_(0),
463  visibleCells_(is)
464 {
465  // Check indices.
466  checkIndices();
467 
468  if (debug)
469  {
470  Pout<< "refinementHistory::refinementHistory :"
471  << " constructed history from Istream"
472  << " splitCells:" << splitCells_.size()
473  << " visibleCells:" << visibleCells_.size()
474  << endl;
475  }
476 }
477 
478 
479 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
480 
481 void Foam::refinementHistory::resize(const label size)
482 {
483  label oldSize = visibleCells_.size();
484 
485  if (debug)
486  {
487  Pout<< "refinementHistory::resize from " << oldSize << " to " << size
488  << " cells" << endl;
489  }
490 
491  visibleCells_.setSize(size);
492 
493  // Set additional elements to -1.
494  for (label i = oldSize; i < visibleCells_.size(); i++)
495  {
496  visibleCells_[i] = -1;
497  }
498 }
499 
500 
502 {
503  if (active())
504  {
505  const labelList& reverseCellMap = map.reverseCellMap();
506 
507  // Note that only the live cells need to be renumbered.
508 
509  labelList newVisibleCells(map.cellMap().size(), -1);
510 
511  forAll(visibleCells_, cellI)
512  {
513  if (visibleCells_[cellI] != -1)
514  {
515  label index = visibleCells_[cellI];
516 
517  // Check not already set
518  if (splitCells_[index].addedCellsPtr_.valid())
519  {
521  (
522  "refinementHistory::updateMesh(const mapPolyMesh&)"
523  ) << "Problem" << abort(FatalError);
524  }
525 
526  label newCellI = reverseCellMap[cellI];
527 
528  if (newCellI >= 0)
529  {
530  newVisibleCells[newCellI] = index;
531  }
532  }
533  }
534 
535  if (debug)
536  {
537  Pout<< "refinementHistory::updateMesh : from "
538  << visibleCells_.size()
539  << " to " << newVisibleCells.size()
540  << " cells" << endl;
541  }
542 
543  visibleCells_.transfer(newVisibleCells);
544  }
545 }
546 
547 
548 // Update numbering for subsetting
550 (
551  const labelList& pointMap,
552  const labelList& faceMap,
553  const labelList& cellMap
554 )
555 {
556  if (active())
557  {
558  labelList newVisibleCells(cellMap.size(), -1);
559 
560  forAll(newVisibleCells, cellI)
561  {
562  label oldCellI = cellMap[cellI];
563 
564  label index = visibleCells_[oldCellI];
565 
566  // Check that cell is live (so its parent has no refinement)
567  if (index >= 0 && splitCells_[index].addedCellsPtr_.valid())
568  {
570  (
571  "refinementHistory::subset"
572  "(const labelList&, const labelList&, const labelList&)"
573  ) << "Problem" << abort(FatalError);
574  }
575 
576  newVisibleCells[cellI] = index;
577  }
578 
579  if (debug)
580  {
581  Pout<< "refinementHistory::updateMesh : from "
582  << visibleCells_.size()
583  << " to " << newVisibleCells.size()
584  << " cells" << endl;
585  }
586 
587  visibleCells_.transfer(newVisibleCells);
588  }
589 }
590 
591 
592 void Foam::refinementHistory::countProc
593 (
594  const label index,
595  const label newProcNo,
596  labelList& splitCellProc,
597  labelList& splitCellNum
598 ) const
599 {
600  if (splitCellProc[index] != newProcNo)
601  {
602  // Different destination processor from other cells using this
603  // parent. Reset count.
604  splitCellProc[index] = newProcNo;
605  splitCellNum[index] = 1;
606  }
607  else
608  {
609  splitCellNum[index]++;
610 
611  // Increment parent if whole splitCell moves to same processor
612  if (splitCellNum[index] == 8)
613  {
614  Pout<< "Moving " << splitCellNum[index]
615  << " cells originating from cell " << index
616  << " from processor " << Pstream::myProcNo()
617  << " to processor " << splitCellProc[index]
618  << endl;
619 
620  label parent = splitCells_[index].parent_;
621 
622  if (parent >= 0)
623  {
624  string oldPrefix = Pout.prefix();
625  Pout.prefix() = " " + oldPrefix;
626 
627  countProc(parent, newProcNo, splitCellProc, splitCellNum);
628 
629  Pout.prefix() = oldPrefix;
630  }
631  }
632  }
633 }
634 
635 
637 {
638  if (!active())
639  {
641  (
642  "refinementHistory::distribute(const mapDistributePolyMesh&)"
643  ) << "Calling distribute on inactive history" << abort(FatalError);
644  }
645 
646 
647  if (!Pstream::parRun())
648  {
649  return;
650  }
651 
652  // Remove unreferenced history.
653  compact();
654 
655  Pout<< nl << "--BEFORE:" << endl;
656  writeDebug();
657  Pout<< "---------" << nl << endl;
658 
659 
660  // Distribution is only partially functional.
661  // If all 8 cells resulting from a single parent are sent across in one
662  // go it will also send across that part of the refinement history.
663  // If however e.g. first 1 and then the other 7 are sent across the
664  // history will not be reconstructed.
665 
666  // Determine clusters. This is per every entry in splitCells_ (that is
667  // a parent of some refinement) a label giving the processor it goes to
668  // if all its children are going to the same processor.
669 
670  // Per visible cell the processor it goes to.
671  labelList destination(visibleCells_.size());
672 
673  const labelListList& subCellMap = map.cellMap().subMap();
674 
675  forAll(subCellMap, procI)
676  {
677  const labelList& newToOld = subCellMap[procI];
678 
679  forAll(newToOld, i)
680  {
681  label oldCellI = newToOld[i];
682 
683  destination[oldCellI] = procI;
684  }
685  }
686 
687 //Pout<< "refinementHistory::distribute :"
688 // << " destination:" << destination << endl;
689 
690  // Per splitCell entry the processor it moves to
691  labelList splitCellProc(splitCells_.size(), -1);
692  // Per splitCell entry the number of live cells that move to that processor
693  labelList splitCellNum(splitCells_.size(), 0);
694 
695  forAll(visibleCells_, cellI)
696  {
697  label index = visibleCells_[cellI];
698 
699  if (index >= 0)
700  {
701  countProc
702  (
703  splitCells_[index].parent_,
704  destination[cellI],
705  splitCellProc,
706  splitCellNum
707  );
708  }
709  }
710 
711 Pout<< "refinementHistory::distribute :"
712  << " splitCellProc:" << splitCellProc << endl;
713 
714 Pout<< "refinementHistory::distribute :"
715  << " splitCellNum:" << splitCellNum << endl;
716 
717 
718  // Create subsetted refinement tree consisting of all parents that
719  // move in their whole to other processor.
720  for (label procI = 0; procI < Pstream::nProcs(); procI++)
721  {
722  Pout<< "-- Subetting for processor " << procI << endl;
723 
724  // From uncompacted to compacted splitCells.
725  labelList oldToNew(splitCells_.size(), -1);
726 
727  // Compacted splitCells. Similar to subset routine below.
728  DynamicList<splitCell8> newSplitCells(splitCells_.size());
729 
730  // Loop over all entries. Note: could recurse like countProc so only
731  // visit used entries but is probably not worth it.
732 
733  forAll(splitCells_, index)
734  {
735 // Pout<< "oldCell:" << index
736 // << " proc:" << splitCellProc[index]
737 // << " nCells:" << splitCellNum[index]
738 // << endl;
739 
740  if (splitCellProc[index] == procI && splitCellNum[index] == 8)
741  {
742  // Entry moves in its whole to procI
743  oldToNew[index] = newSplitCells.size();
744  newSplitCells.append(splitCells_[index]);
745 
746  Pout<< "Added oldCell " << index
747  << " info " << newSplitCells[newSplitCells.size()-1]
748  << " at position " << newSplitCells.size()-1
749  << endl;
750  }
751  }
752 
753  // Add live cells that are subsetted.
754  forAll(visibleCells_, cellI)
755  {
756  label index = visibleCells_[cellI];
757 
758  if (index >= 0 && destination[cellI] == procI)
759  {
760  label parent = splitCells_[index].parent_;
761 
762  Pout<< "Adding refined cell " << cellI
763  << " since moves to "
764  << procI << " old parent:" << parent << endl;
765 
766  // Create new splitCell with parent
767  oldToNew[index] = newSplitCells.size();
768  newSplitCells.append(splitCell8(parent));
769  }
770  }
771 
772  //forAll(oldToNew, index)
773  //{
774  // Pout<< "old:" << index << " new:" << oldToNew[index]
775  // << endl;
776  //}
777 
778  newSplitCells.shrink();
779 
780  // Renumber contents of newSplitCells
781  forAll(newSplitCells, index)
782  {
783  splitCell8& split = newSplitCells[index];
784 
785  if (split.parent_ >= 0)
786  {
787  split.parent_ = oldToNew[split.parent_];
788  }
789  if (split.addedCellsPtr_.valid())
790  {
791  FixedList<label, 8>& splits = split.addedCellsPtr_();
792 
793  forAll(splits, i)
794  {
795  if (splits[i] >= 0)
796  {
797  splits[i] = oldToNew[splits[i]];
798  }
799  }
800  }
801  }
802 
803 
804  const labelList& subMap = subCellMap[procI];
805 
806  // New visible cells.
807  labelList newVisibleCells(subMap.size(), -1);
808 
809  forAll(subMap, newCellI)
810  {
811  label oldCellI = subMap[newCellI];
812 
813  label oldIndex = visibleCells_[oldCellI];
814 
815  if (oldIndex >= 0)
816  {
817  newVisibleCells[newCellI] = oldToNew[oldIndex];
818  }
819  }
820 
821  //Pout<< nl << "--Subset for domain:" << procI << endl;
822  //writeDebug(newVisibleCells, newSplitCells);
823  //Pout<< "---------" << nl << endl;
824 
825 
826  // Send to neighbours
827  OPstream toNbr(Pstream::blocking, procI);
828  toNbr << newSplitCells << newVisibleCells;
829  }
830 
831 
832  // Receive from neighbours and merge
833  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
834 
835  // Remove all entries. Leave storage intact.
836  splitCells_.clear();
837 
838  visibleCells_.setSize(map.mesh().nCells());
839  visibleCells_ = -1;
840 
841  for (label procI = 0; procI < Pstream::nProcs(); procI++)
842  {
843  IPstream fromNbr(Pstream::blocking, procI);
844  List<splitCell8> newSplitCells(fromNbr);
845  labelList newVisibleCells(fromNbr);
846 
847  //Pout<< nl << "--Received from domain:" << procI << endl;
848  //writeDebug(newVisibleCells, newSplitCells);
849  //Pout<< "---------" << nl << endl;
850 
851 
852  // newSplitCells contain indices only into newSplitCells so
853  // renumbering can be done here.
854  label offset = splitCells_.size();
855 
856  Pout<< "**Renumbering data from proc " << procI << " with offset "
857  << offset << endl;
858 
859  forAll(newSplitCells, index)
860  {
861  splitCell8& split = newSplitCells[index];
862 
863  if (split.parent_ >= 0)
864  {
865  split.parent_ += offset;
866  }
867  if (split.addedCellsPtr_.valid())
868  {
869  FixedList<label, 8>& splits = split.addedCellsPtr_();
870 
871  forAll(splits, i)
872  {
873  if (splits[i] >= 0)
874  {
875  splits[i] += offset;
876  }
877  }
878  }
879 
880  splitCells_.append(split);
881  }
882 
883 
884  // Combine visibleCell.
885  const labelList& constructMap = map.cellMap().constructMap()[procI];
886 
887  forAll(newVisibleCells, i)
888  {
889  visibleCells_[constructMap[i]] = newVisibleCells[i] + offset;
890  }
891  }
892  splitCells_.shrink();
893 
894  Pout<< nl << "--AFTER:" << endl;
895  writeDebug();
896  Pout<< "---------" << nl << endl;
897 }
898 
899 
901 {
902  if (debug)
903  {
904  Pout<< "refinementHistory::compact() Entering with:"
905  << " freeSplitCells_:" << freeSplitCells_.size()
906  << " splitCells_:" << splitCells_.size()
907  << " visibleCells_:" << visibleCells_.size()
908  << endl;
909 
910  // Check all free splitCells are marked as such
911  forAll(freeSplitCells_, i)
912  {
913  label index = freeSplitCells_[i];
914 
915  if (splitCells_[index].parent_ != -2)
916  {
917  FatalErrorIn("refinementHistory::compact()")
918  << "Problem index:" << index
919  << abort(FatalError);
920  }
921  }
922 
923  // Check none of the visible cells are marked as free
924  forAll(visibleCells_, cellI)
925  {
926  if
927  (
928  visibleCells_[cellI] >= 0
929  && splitCells_[visibleCells_[cellI]].parent_ == -2
930  )
931  {
932  FatalErrorIn("refinementHistory::compact()")
933  << "Problem : visible cell:" << cellI
934  << " is marked as being free." << abort(FatalError);
935  }
936  }
937  }
938 
939  DynamicList<splitCell8> newSplitCells(splitCells_.size());
940 
941  // From uncompacted to compacted splitCells.
942  labelList oldToNew(splitCells_.size(), -1);
943 
944  // Mark all used splitCell entries. These are either indexed by visibleCells
945  // or indexed from other splitCell entries.
946 
947  // Mark from visibleCells
948  forAll(visibleCells_, cellI)
949  {
950  label index = visibleCells_[cellI];
951 
952  if (index >= 0)
953  {
954  // Make sure we only mark visible indices if they either have a
955  // parent or subsplits.
956  if
957  (
958  splitCells_[index].parent_ != -1
959  || splitCells_[index].addedCellsPtr_.valid()
960  )
961  {
962  markSplit(index, oldToNew, newSplitCells);
963  }
964  }
965  }
966 
967  // Mark from splitCells
968  forAll(splitCells_, index)
969  {
970  if (splitCells_[index].parent_ == -2)
971  {
972  // freed cell.
973  }
974  else if
975  (
976  splitCells_[index].parent_ == -1
977  && splitCells_[index].addedCellsPtr_.empty()
978  )
979  {
980  // recombined cell. No need to keep since no parent and no subsplits
981  // Note that gets marked if reachable from other index!
982  }
983  else
984  {
985  // Is used element.
986  markSplit(index, oldToNew, newSplitCells);
987  }
988  }
989 
990 
991  // Now oldToNew is fully complete and compacted elements are in
992  // newSplitCells.
993  // Renumber contents of newSplitCells and visibleCells.
994  forAll(newSplitCells, index)
995  {
996  splitCell8& split = newSplitCells[index];
997 
998  if (split.parent_ >= 0)
999  {
1000  split.parent_ = oldToNew[split.parent_];
1001  }
1002  if (split.addedCellsPtr_.valid())
1003  {
1004  FixedList<label, 8>& splits = split.addedCellsPtr_();
1005 
1006  forAll(splits, i)
1007  {
1008  if (splits[i] >= 0)
1009  {
1010  splits[i] = oldToNew[splits[i]];
1011  }
1012  }
1013  }
1014  }
1015 
1016 
1017  if (debug)
1018  {
1019  Pout<< "refinementHistory::compact : compacted splitCells from "
1020  << splitCells_.size() << " to " << newSplitCells.size() << endl;
1021  }
1022 
1023  splitCells_.transfer(newSplitCells);
1024  freeSplitCells_.clearStorage();
1025 
1026 
1027  if (debug)
1028  {
1029  Pout<< "refinementHistory::compact() NOW:"
1030  << " freeSplitCells_:" << freeSplitCells_.size()
1031  << " splitCells_:" << splitCells_.size()
1032  << " newSplitCells:" << newSplitCells.size()
1033  << " visibleCells_:" << visibleCells_.size()
1034  << endl;
1035  }
1036 
1037 
1038  // Adapt indices in visibleCells_
1039  forAll(visibleCells_, cellI)
1040  {
1041  label index = visibleCells_[cellI];
1042 
1043  if (index >= 0)
1044  {
1045  // Note that oldToNew can be -1 so it resets newVisibleCells.
1046  visibleCells_[cellI] = oldToNew[index];
1047  }
1048  else
1049  {
1050  // Keep -1 value.
1051  }
1052  }
1053 }
1054 
1055 
1057 {
1058  writeDebug(visibleCells_, splitCells_);
1059 }
1060 
1061 
1064  const label cellI,
1065  const labelList& addedCells
1066 )
1067 {
1068  label parentIndex = -1;
1069 
1070  if (visibleCells_[cellI] != -1)
1071  {
1072  // Was already live. The current live cell becomes the
1073  // parent of the cells split off from it.
1074 
1075  parentIndex = visibleCells_[cellI];
1076 
1077  // It is no longer live (note that actually cellI gets alive
1078  // again below since is addedCells[0])
1079  visibleCells_[cellI] = -1;
1080  }
1081  else
1082  {
1083  // Create 0th level. -1 parent to denote this.
1084  parentIndex = allocateSplitCell(-1, -1);
1085  }
1086 
1087  // Create live entries for added cells that point to the
1088  // cell they were created from (parentIndex)
1089  forAll(addedCells, i)
1090  {
1091  label addedCellI = addedCells[i];
1092 
1093  // Create entries for the split off cells. All of them
1094  // are visible.
1095  visibleCells_[addedCellI] = allocateSplitCell(parentIndex, i);
1096  }
1097 }
1098 
1099 
1102  const label masterCellI,
1103  const labelList& combinedCells
1104 )
1105 {
1106  // Save the parent structure
1107  label parentIndex = splitCells_[visibleCells_[masterCellI]].parent_;
1108 
1109  // Remove the information for the combined cells
1110  forAll(combinedCells, i)
1111  {
1112  label cellI = combinedCells[i];
1113 
1114  freeSplitCell(visibleCells_[cellI]);
1115  visibleCells_[cellI] = -1;
1116  }
1117 
1118  splitCell8& parentSplit = splitCells_[parentIndex];
1119  parentSplit.addedCellsPtr_.reset(NULL);
1120  visibleCells_[masterCellI] = parentIndex;
1121 }
1122 
1123 
1125 {
1126  is >> *this;
1127  return !is.bad();
1128 }
1129 
1130 
1132 {
1133  os << *this;
1134 
1135  return os.good();
1136 }
1137 
1138 
1139 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
1140 
1142 {
1143  rh.freeSplitCells_.clearStorage();
1144 
1145  is >> rh.splitCells_ >> rh.visibleCells_;
1146 
1147  // Check indices.
1148  rh.checkIndices();
1149 
1150  return is;
1151 }
1152 
1153 
1154 Foam::Ostream& Foam::operator<<(Ostream& os, const refinementHistory& rh)
1155 {
1156  const_cast<refinementHistory&>(rh).compact();
1157 
1158  return os << "// splitCells" << nl
1159  << rh.splitCells_ << nl
1160  << "// visibleCells" << nl
1161  << rh.visibleCells_;
1162 }
1163 
1164 
1165 // ************************ vim: set sw=4 sts=4 et: ************************ //