Chombo + EB + MF  3.2
BlockBaseRegisterImplem.H
Go to the documentation of this file.
1 #ifdef CH_LANG_CC
2 /*
3  * _______ __
4  * / ___/ / ___ __ _ / / ___
5  * / /__/ _ \/ _ \/ V \/ _ \/ _ \
6  * \___/_//_/\___/_/_/_/_.__/\___/
7  * Please refer to Copyright.txt, in Chombo's root directory.
8  */
9 #endif
10 
11 #ifndef _BLOCKBASEREGISTERIMPLEM_H_
12 #define _BLOCKBASEREGISTERIMPLEM_H_
13 
14 #include "BlockBaseRegister.H"
15 #include "BoxIterator.H"
16 #include "CornerCopier.H"
17 #include "CH_Timer.H"
18 
19 #include "NamespaceHeader.H"
20 
21 inline int getBlockBdryIndex( const int& a_dir, const Side::LoHiSide& a_side );
22 
23 template <class T>
25  /// CELL-centered grids
26  const DisjointBoxLayout& a_grids,
27  /// ghost width for state variables
28  const int& a_ghost )
29  : m_verbosity(0),
30  m_hasData(false),
31  m_hasDefinedDataLayouts(false)
32 {
34  define( a_mblock, a_grids, a_ghost );
35 }
36 
37 
38 template <class T>
40 {
41  //do nothing for now
42 }
43 
44 
45 template <class T>
47  /// CELL-centered grids
48  const DisjointBoxLayout& a_grids,
49  /// ghost width for state variables
50  const int& a_ghost )
51 {
52  CH_TIME("BlockBaseRegister<T>::define");
53 
54  m_coordsys = (MultiBlockCoordSys*) a_mblock;
55  m_grids = a_grids;
56  m_ghost = a_ghost;
57  for (int dir(0); dir<CH_SPACEDIM; dir++)
58  {
59  m_ghostVect.setVal( dir, m_ghost );
60  }
61 
62  const ProblemDomain& inDomain( m_grids.physDomain() );
63  Box newBoundingBox( inDomain.domainBox() );
64  newBoundingBox.grow( m_ghostVect + IntVect::Unit );
65  m_domain.define( newBoundingBox );
66 
67  m_blockBdryIndexMap.resize( m_coordsys->numBlocks() );
68 
69  buildRegisterLayout();
70  buildSrcRegisterLayout();
71  buildExchangeRegisterLayout();
72 
73  m_isDefined = true;
74 }
75 
76 
77 template <class T>
78 void BlockBaseRegister<T>::defineDataLayouts( const int a_ncomp )
79 {
80  for (int j(0); j<2*CH_SPACEDIM; j++)
81  {
82  IntVect noDirGhost( m_ghost * IntVect::Unit );
83  int dir(j%CH_SPACEDIM);
84  noDirGhost[dir] = 0;
85 
86  m_register[j].define( m_registerLayout[j], a_ncomp, noDirGhost );
87 // m_excRegister[j].define( m_excRegisterLayout[j], a_ncomp );
88  m_excRegister[j].define( m_excRegisterLayout[j], a_ncomp, m_ghostVect );
89 // m_srcRegister[j].define( m_srcRegisterLayout[j], a_ncomp );
90  m_srcRegister[j].define( m_srcRegisterLayout[j], a_ncomp, noDirGhost );
91  }
92  m_hasDefinedDataLayouts = true;
93 }
94 
95 
96 template <class T>
98 {
99  for (int j(0); j<2*CH_SPACEDIM; j++)
100  {
101  m_register[j].clear();
102  m_excRegister[j].clear();
103  m_srcRegister[j].clear();
104  }
105  m_hasDefinedDataLayouts = false;
106 }
107 
108 
109 template <class T>
110 void BlockBaseRegister<T>::store( const T& a_data,
111  const DataIndex& a_dataIndex,
112  int a_dir,
113  Side::LoHiSide a_side )
114 {
115  CH_TIME("BlockBaseRegister<T>::store");
116 
117  if ( isDefined() )
118  {
119 
120  if (!hasStorage())
121  {
122  defineDataLayouts( a_data.nComp() );
123  }
124 
125  const int index( getBlockBdryIndex( a_dir, a_side ) );
126 
127  LevelData<T >& dataRegister( m_register[index] );
128 
129  LevelData<T >& dataExcRegister( m_excRegister[index] );
130  if (hasInterface( a_dataIndex, a_dir, a_side ))
131  {
132 
133  const DataIndex& regIndex( m_grdToDstMap[index][a_dataIndex] );
134  T& regData( dataRegister[regIndex] );
135 
136  // We cast away the const-ness in a_flux for the scope of this
137  // function so that we can shift without copying.
138  T& inData( const_cast<T&>(a_data) ); // Muhahaha.
139  inData.shiftHalf( a_dir, -sign(a_side) );
140 
141  Box srcBox( adjCellBox( m_grids[a_dataIndex], a_dir, a_side, -1 ) );
142  if (m_ghost>0)
143  {
144  for (int gdir(0); gdir<CH_SPACEDIM; ++gdir)
145  {
146  if (gdir!=a_dir)
147  {
148  srcBox.grow( gdir, m_ghost );
149  }
150  else
151  {
152  srcBox.growDir( gdir, flip(a_side), m_ghost );
153  }
154  }
155  }
156 
157  const Box window( srcBox & inData.box() );
158  regData.copy( inData, window & regData.box() );
159 
160  // Done being evil...for now ;-)
161  inData.shiftHalf( a_dir, sign(a_side) );
162 
163  const DataIndex& excIndex( m_dstToExcMap[index][regIndex] );
164  T& excData( dataExcRegister[excIndex] );
165  const Box excWindow( regData.box() & excData.box() );
166  excData.copy( regData, excWindow );
167 
168  m_hasData = true;
169  }
170  }
171  else
172  {
173  MayDay::Error(" BlockBaseRegister<T>::store: BlockBaseRegister<T> is not defined!");
174  }
175 }
176 
177 
178 template <class T>
179 void BlockBaseRegister<T>::increment( const T& a_data,
180  const DataIndex& a_dataIndex,
181  int a_dir,
182  Side::LoHiSide a_side )
183 {
184  CH_TIME("BlockBaseRegister<T>::increment");
185 
186  if ( isDefined() )
187  {
188 
189  if (!hasStorage())
190  {
191  defineDataLayouts( a_data.nComp() );
192  }
193 
194  const int index( getBlockBdryIndex( a_dir, a_side ) );
195 
196  BoxLayoutData<T >& dataRegister( m_register[index] );
197  LevelData<T >& dataExcRegister( m_excRegister[index] );
198  if (hasInterface( a_dataIndex, a_dir, a_side ))
199  {
200 
201  const DataIndex& regIndex( m_grdToDstMap[index][a_dataIndex] );
202  T& regData( dataRegister[regIndex] );
203 
204  // We cast away the const-ness in a_flux for the scope of this
205  // function so that we can shift without copying.
206  T& inData( const_cast<T&>(a_data) ); // Muhahaha.
207  inData.shiftHalf( a_dir, -sign(a_side) );
208 
209  Box srcBox( adjCellBox( m_grids[a_dataIndex], a_dir, a_side, -1 ) );
210  if (m_ghost>0)
211  {
212  for (int gdir(0); gdir<CH_SPACEDIM; ++gdir)
213  {
214  if (gdir!=a_dir)
215  {
216  srcBox.grow( gdir, m_ghost );
217  }
218  else
219  {
220  srcBox.growDir( gdir, flip(a_side), m_ghost );
221  }
222  }
223  }
224 
225  const Box window( srcBox & inData.box() );
226  const int n_comp( a_data.nComp() );
227  regData.plus( inData, window & regData.box(), 0, 0, n_comp );
228 
229  // Done being evil...for now ;-)
230  inData.shiftHalf( a_dir, sign(a_side) );
231 
232  const DataIndex& excIndex( m_dstToExcMap[index][regIndex] );
233  T& excData( dataExcRegister[excIndex] );
234  excData.copy( regData, regData.box() & excData.box() );
235 
236  m_hasData = true;
237  }
238  }
239  else
240  {
241  MayDay::Error(" BlockBaseRegister<T>::increment: BlockBaseRegister<T> is not defined!");
242  }
243 }
244 
245 
246 template <class T>
248 {
249  CH_TIME("BlockBaseRegister<T>::exchange");
250  if ( isDefined() )
251  {
252  if ( hasData() )
253  {
254  copyToSrcRegister();
255  copyToDstRegister();
256  if (m_verbosity>0) printRegisters();
257  }
258  }
259  else
260  {
261  MayDay::Error(" BlockBaseRegister<T>::exchange: BlockBaseRegister<T> is not defined!");
262  }
263 }
264 
265 
266 template <class T>
267 inline
269 {
270  for (DataIterator dit(a_layoutData.dataIterator()); dit.ok(); ++dit)
271  {
272  a_layoutData[dit].setVal( 0.0 );
273  }
274 }
275 
276 
277 template <class T>
278 void BlockBaseRegister<T>::zeroRegister( const int a_n_comp )
279 {
280  CH_TIME("BlockBaseRegister<T>::zeroRegister");
281  if ( isDefined() )
282  {
283  if (!hasStorage())
284  {
285  if (a_n_comp>0)
286  {
287  defineDataLayouts( a_n_comp );
288  }
289  else
290  {
291  MayDay::Error(" BlockRegister::zeroRegister: Attempt to zero register of unknown component length!");
292  }
293  }
294  else if ( a_n_comp>0 && a_n_comp != m_register[0].nComp() )
295  {
296  clearDataLayouts();
297  defineDataLayouts( a_n_comp );
298  }
299  for (int index(0); index<2*CH_SPACEDIM; index++)
300  {
301  zeroBoxLayoutData( m_register[index] );
302  zeroBoxLayoutData( m_excRegister[index] );
303  zeroBoxLayoutData( m_srcRegister[index] );
304  }
305  m_hasData = true;
306  }
307  else
308  {
309  MayDay::Error(" BlockRegister::zeroRegister: Attempt to zero undefined register!");
310  }
311 }
312 
313 
314 template <class T>
316  int a_dir,
317  Side::LoHiSide a_side ) const
318 {
319  CH_TIME("BlockBaseRegister<T>::hasInterface");
320  // JAFH: We assume the a_dataIndex is valid; should we check?
321 
322  // construct ghost box grown in a_dir direction on a_side
323  const Box& box( m_grids[a_dataIndex] );
324  const Box clipBox( adjCellBox( box, a_dir, a_side, 1 ) );
325 
326  // construct block box grown in a_dir direction on a_side
327  const int block( getBlockID( a_dataIndex ) );
328  const Box& blockBox( getBlockBox( block ) );
329  const Box blockBdryBox( adjCellBox( blockBox, a_dir, a_side, 1 ) );
330 
331  // Has interface if touches interfacial block boundary
332  if (blockBdryBox.intersectsNotEmpty( clipBox ))
333  {
334  if (isInterface( block, a_dir, a_side ))
335  {
336  return true;
337  }
338  }
339  return false;
340 }
341 
342 
343 template <class T>
345  const DataIndex& a_dataIndex,
346  int a_dir,
347  Side::LoHiSide a_side,
348  Side::LoHiSide a_sideData ) const
349 {
350  CH_TIME("BlockBaseRegister<T>::fill");
351 
352  if ( isDefined() )
353  {
354  if ( hasData() )
355  {
356  if ( hasInterface( a_dataIndex, a_dir, a_side ) )
357  {
358  const int index( getBlockBdryIndex( a_dir, a_side ) );
359  const DataIndex& regIndex( m_grdToDstMap[index][a_dataIndex] );
360  const T& regData( m_register[index][regIndex] );
361 
362  const int srcSide( sign( a_side ) * sign( a_sideData ) );
363  const Box& box( m_grids[a_dataIndex] );
364  Box srcBox( adjCellBox( box, a_dir, a_side, srcSide ) );
365  const int shift( sign( a_sideData ) );
366 
367  if (m_ghost>0)
368  {
369  for (int gdir(0); gdir<CH_SPACEDIM; ++gdir)
370  {
371  if (gdir!=a_dir)
372  {
373  srcBox.grow( m_ghost );
374  }
375  else
376  {
377  Side::LoHiSide side( (srcSide>0) ? Side::Lo : Side::Hi );
378  srcBox.growDir( gdir, side, m_ghost );
379  }
380  }
381  }
382 
383  a_data.shiftHalf( a_dir, shift );
384  const Box window( srcBox & a_data.box() & regData.box() );
385  a_data.copy( regData, window );
386  a_data.shiftHalf( a_dir, -shift );
387  }
388  }
389  }
390  else
391  {
392  MayDay::Error(" BlockBaseRegister<T>::fill: BlockBaseRegister<T> is not defined!");
393  }
394 }
395 
396 
397 //////////////////////////////////////////////////////////////////////////////
398 // Begin Private Methods
399 //////////////////////////////////////////////////////////////////////////////
400 
401 
402 template <class T>
403 inline
404 int BlockBaseRegister<T>::getBlockID( const DataIndex& a_dataIndex ) const
405 {
406  const Box& box( m_grids[a_dataIndex] );
407  return m_coordsys->whichBlock( box );
408 }
409 
410 
411 template <class T>
412 inline
413 const Box& BlockBaseRegister<T>::getBlockBox( const int& a_block ) const
414 {
415  const Vector<Box>& blocks( m_coordsys->mappingBlocks() );
416  return blocks[a_block];
417 }
418 
419 
420 template <class T>
421 inline
423  const int& a_blockBdryIndex ) const
424 {
426  boundary( m_coordsys->boundaries()[a_block] );
427  const int neighborBlock( boundary[a_blockBdryIndex].neighbor() );
428  return getBlockBox( neighborBlock );
429 }
430 
431 
432 inline
433 int getBlockBdryIndex( const int& a_dir, const Side::LoHiSide& a_side )
434 {
435  return ( a_side * CH_SPACEDIM + a_dir );
436 }
437 
438 
439 template <class T>
440 inline
442  const int& a_block,
443  const int& a_blockBdryIndex ) const
444 {
446  boundary( m_coordsys->boundaries()[a_block] );
447  return boundary[a_blockBdryIndex].getTransformation();
448 }
449 
450 
451 template <class T>
452 inline
454  const int& a_block,
455  const int& a_dir,
456  const Side::LoHiSide& a_side )
457 {
458  // assumes one neighbor per block boundary
459  // assumes connection only at block boundaries
460  const int index( getBlockBdryIndex( a_dir, a_side ) );
461  BlockBdryIndexMap::iterator it( m_blockBdryIndexMap[a_block].find( index ) );
462  if (it == m_blockBdryIndexMap[a_block].end())
463  {
464 
465  const IndicesTransformation& itrans( getTransformation( a_block, index ) );
466  Box mySrcBox( itrans.transformFwd( getBlockBox( a_block ) ) );
467  const int neighborDir( itrans.getPermutation()[a_dir] );
468  mySrcBox.growLo( neighborDir );
469  const Box& neighborBox( getNeighborBlockBox( a_block, index ) );
470  const Box intersectBox( mySrcBox & neighborBox );
471 
472  Side::LoHiSide neighborSide( Side::Hi );
473  if ( intersectBox.isEmpty() )
474  {
475  neighborSide = flip( neighborSide );
476  }
477  int src_index( getBlockBdryIndex( neighborDir, neighborSide ) );
478  m_blockBdryIndexMap[a_block].insert(
479  BlockBdryIndexMap::value_type( index, src_index ) );
480  }
481  return m_blockBdryIndexMap[a_block][index];
482 }
483 
484 
485 template <class T>
486 inline
487 bool BlockBaseRegister<T>::isInterface( const int& a_block,
488  const int& a_dir,
489  const Side::LoHiSide& a_side ) const
490 
491 {
493  boundary( m_coordsys->boundaries()[a_block] );
494  const int index( getBlockBdryIndex( a_dir, a_side ) );
495  return boundary[index].isInterface();
496 }
497 
498 
499 inline
501  const Vector<Box>& a_boxes,
502  const Vector<int>& a_procs )
503 {
504  int myRank(0);
505 #ifdef CH_MPI
506  MPI_Comm_rank( Chombo_MPI::comm, &myRank );
507 #endif
508  for (int i(0); i<a_boxes.size(); i++)
509  {
510  if (a_procs[i]==myRank)
511  {
512  a_list.append(i);
513  }
514  }
515 }
516 
517 
518 template <class T>
519 inline
521  LayoutData<DataIndex>& a_map,
522  const Vector<Box>& a_boxes,
523  const Vector<int>& a_procs,
524  const Vector<DataIndex>& a_indices )
525 {
526  a_layout.define( a_boxes, a_procs, m_domain );
527  a_layout.close();
528 
529  List<int> vi;
530  buildLocalIndexList( vi, a_boxes, a_procs );
531 
532  for (DataIterator dit(a_layout.dataIterator()); dit.ok(); ++dit)
533  {
534  const Box& box( a_layout[dit] );
535  for (ListIterator<int> li( vi.first() ); li.ok(); ++li)
536  {
537  int i( li() );
538  if ( box.eq( a_boxes[i] ) )
539  {
540  DataIndex foo( dit() );
541  a_map[a_indices[i]] = foo;
542 // a_map[a_indices[i]] = dit();
543 // a_map[a_indices[i]] = dit();
544  vi.remove( li() );
545  break;
546  }
547  }
548  }
549  if (vi.isNotEmpty())
550  {
551  MayDay::Error(" BlockBaseRegister<T>::createDstBoxLayout: One or more orphaned boxes!");
552  }
553 }
554 
555 
556 template <class T>
557 inline
559  DisjointBoxLayout& a_srcLayout,
561  const Vector<Box>& a_srcBoxes,
562  const Vector<int>& a_srcProcs,
563  const Vector<DataIndex>& a_srcDataIndices,
564  const Vector<int>& a_dstBlockBdryIndices,
565  const LayoutData<DataIndex> a_dstMap[2*CH_SPACEDIM] )
566 {
567  a_srcLayout.define( a_srcBoxes, a_srcProcs, m_domain );
568  a_srcLayout.close();
569 
570  List<int> vi;
571  buildLocalIndexList( vi, a_srcBoxes, a_srcProcs );
572 
573  for (DataIterator dit(a_srcLayout.dataIterator()); dit.ok(); ++dit)
574  {
575  const Box& box( a_srcLayout[dit] );
576  for (ListIterator<int> li( vi.first() ); li.ok(); ++li)
577  {
578  int i( li() );
579  if ( box.eq( a_srcBoxes[i] ) )
580  {
581  const int dstBlockBdryIndex( a_dstBlockBdryIndices[i] );
582  const DataIndex& gridDataIndex( a_srcDataIndices[i] );
583  const DataIndex& dstDataIndex( a_dstMap[dstBlockBdryIndex][gridDataIndex] );
584  a_srcMap[dstBlockBdryIndex][dstDataIndex] = dit();
585  vi.remove( li() );
586  break;
587  }
588  }
589  }
590  if (vi.isNotEmpty())
591  {
592  MayDay::Error(" BlockBaseRegister<T>::createSrcBoxLayout: One or more orphaned boxes!");
593  }
594 }
595 
596 
597 template <class T>
598 inline
600  DisjointBoxLayout& a_exchangeLayout,
601  LayoutData<DataIndex>& a_dstToExcMap,
602  const Vector<Box>& a_boxes,
603  const Vector<int>& a_procs,
604  const Vector<DataIndex>& a_dataIndices,
605  const LayoutData<DataIndex>& a_grdToDstMap )
606 {
607  a_exchangeLayout.define( a_boxes, a_procs, m_domain );
608  a_exchangeLayout.close();
609 
610  List<int> vi;
611  buildLocalIndexList( vi, a_boxes, a_procs );
612 
613  for (DataIterator dit(a_exchangeLayout.dataIterator()); dit.ok(); ++dit)
614  {
615  const Box& box( a_exchangeLayout[dit] );
616  for (ListIterator<int> li( vi.first() ); li.ok(); ++li)
617  {
618  int i( li() );
619  if ( box.eq( a_boxes[i] ) )
620  {
621  const DataIndex& gridDataIndex( a_dataIndices[i] );
622  const DataIndex& dstDataIndex( a_grdToDstMap[gridDataIndex] );
623  a_dstToExcMap[dstDataIndex] = dit();
624  vi.remove( li() );
625  break;
626  }
627  }
628  }
629  if (vi.isNotEmpty())
630  {
631  MayDay::Error(" BlockBaseRegister<T>::createSrcBoxLayout: One or more orphaned boxes!");
632  }
633 }
634 
635 
636 template <class T>
638 {
639  Vector<Box> boxes[2*CH_SPACEDIM];
640  Vector<int> boxProcMap[2*CH_SPACEDIM];
641  Vector<DataIndex> dataIndex[2*CH_SPACEDIM];
642  buildBoxVectors( boxes, boxProcMap, dataIndex, DST );
643 
644  // For each block boundary, create a BoxLayout for the destination registers
645  // and define a LayoutData that will be used to associate each destination
646  // register DataIndex with a source register DataIndex
647  for (int dir(0); dir<CH_SPACEDIM; ++dir)
648  {
649  for (SideIterator side; side.ok(); ++side)
650  {
651  const int index( getBlockBdryIndex( dir, side() ) );
652  m_grdToDstMap[index].define( m_grids );
653  createDstBoxLayout( m_registerLayout[index],
654  m_grdToDstMap[index],
655  boxes[index],
656  boxProcMap[index],
657  dataIndex[index] );
658  m_dstToSrcMap[index].define( m_registerLayout[index] );
659  m_dstToExcMap[index].define( m_registerLayout[index] );
660  }
661  }
662 }
663 
664 
665 template <class T>
667 {
668  Vector<Box> boxes[2*CH_SPACEDIM];
669  Vector<int> boxProcMap[2*CH_SPACEDIM];
670  Vector<DataIndex> dataIndex[2*CH_SPACEDIM];
671  buildBoxVectors( boxes, boxProcMap, dataIndex, SRC );
672 
673  Vector<int> blockBdryIndexMap[2*CH_SPACEDIM];
674  buildInverseBlockBdryIndexMap( blockBdryIndexMap );
675 
676  // For each block boundary, create a BoxLayout for the source registers.
677  // At the same time, we associate the source DataIndexes with destination
678  // DataIndexes.
679  for (int dir(0); dir<CH_SPACEDIM; ++dir)
680  {
681  for (SideIterator side; side.ok(); ++side)
682  {
683  const int index( getBlockBdryIndex( dir, side() ) );
684  createSrcBoxLayout( m_srcRegisterLayout[index],
685  m_dstToSrcMap,
686  boxes[index],
687  boxProcMap[index],
688  dataIndex[index],
689  blockBdryIndexMap[index],
690  m_grdToDstMap );
691  }
692  }
693 }
694 
695 
696 template <class T>
698 {
699  Vector<Box> boxes[2*CH_SPACEDIM];
700  Vector<int> boxProcMap[2*CH_SPACEDIM];
701  Vector<DataIndex> dataIndex[2*CH_SPACEDIM];
702  buildBoxVectors( boxes, boxProcMap, dataIndex, EXC );
703 
704  // For each block boundary, create a BoxLayout for the source registers.
705  // At the same time, we associate the source DataIndexes with destination
706  // DataIndexes.
707  for (int dir(0); dir<CH_SPACEDIM; ++dir)
708  {
709  for (SideIterator side; side.ok(); ++side)
710  {
711  const int index( getBlockBdryIndex( dir, side() ) );
712  createExchangeBoxLayout( m_excRegisterLayout[index],
713  m_dstToExcMap[index],
714  boxes[index],
715  boxProcMap[index],
716  dataIndex[index],
717  m_grdToDstMap[index] );
718  }
719  }
720 }
721 
722 
723 template <class T>
725  Vector<int> a_boxProcMap[2*CH_SPACEDIM],
726  Vector<DataIndex> a_dataIndex[2*CH_SPACEDIM],
727  const UseType a_type )
728 {
729  for (LayoutIterator lit(m_grids.layoutIterator()); lit.ok(); ++lit)
730  {
731  const Box& box( m_grids[lit] );
732  const int block( m_coordsys->whichBlock( box ) );
733  for (int dir(0); dir<CH_SPACEDIM; ++dir)
734  {
735  for (SideIterator si; si.ok(); ++si)
736  {
737  Side::LoHiSide side( si() );
738  if (hasInterface( DataIndex( lit() ), dir, side ))
739  {
740  int index( getBlockBdryIndex( dir, side ) );
741 
742  Box registerBox;
743 
744  if (a_type==DST)
745  {
746  registerBox = adjCellBox( box, dir, side, -1 );
747  registerBox.growDir( dir, side, 1 );
748 // registerBox.grow( m_ghost );
749  }
750 
751  else if (a_type==SRC)
752  {
753  registerBox = adjCellBox( box, dir, side, 1 );
754  if (m_ghost>0)
755  {
756 #if 0
757  for (int gdir(0); gdir<CH_SPACEDIM; gdir++)
758  {
759  if (gdir!=dir)
760  {
761 // registerBox.grow( gdir, m_ghost );
762  }
763  else
764  {
765 // registerBox.growDir( gdir, side, m_ghost );
766  }
767  }
768 #endif
769  }
770  const IndicesTransformation& itrans( getTransformation( block, index ) );
771  registerBox = itrans.transformFwd( registerBox );
772  index = getSrcBlockBdryIndex( block, dir, side );
773  }
774 
775  else if (a_type==EXC)
776  {
777  registerBox = adjCellBox( box, dir, side, -(1+m_ghost) );
778  if (m_ghost>0)
779  {
780  for (int gdir(0); gdir<CH_SPACEDIM; ++gdir)
781  {
782  if (gdir!=dir)
783  {
784 // registerBox.grow( gdir, m_ghost );
785 #if 0
786  if (hasInterface( DataIndex( lit() ), gdir, Side::Lo ))
787  {
788  registerBox.growLo( gdir, m_ghost );
789  }
790  if (hasInterface( DataIndex( lit() ), gdir, Side::Hi ))
791  {
792  registerBox.growHi( gdir, m_ghost );
793  }
794 #endif
795  }
796  }
797  }
798  }
799 
800  else
801  {
802  MayDay::Error(" BlockBaseRegister<T>::buildBoxVectors: Unknown use type!");
803  }
804 
805  a_boxes[index].push_back( registerBox );
806  a_boxProcMap[index].push_back( m_grids.procID( lit() ) );
807  a_dataIndex[index].push_back( DataIndex( lit() ) );
808  }
809  }
810  }
811  }
812 }
813 
814 
815 template <class T>
817  Vector<int> a_map[2*CH_SPACEDIM] )
818 {
819  for (LayoutIterator lit(m_grids.layoutIterator()); lit.ok(); ++lit)
820  {
821  const Box& box( m_grids[lit] );
822  const int block( m_coordsys->whichBlock( box ) );
823  for (int dir(0); dir<CH_SPACEDIM; ++dir)
824  {
825  for (SideIterator si; si.ok(); ++si)
826  {
827  Side::LoHiSide side( si() );
828  if (hasInterface( DataIndex( lit() ), dir, side ))
829  {
830  const int index( getBlockBdryIndex( dir, side ) );
831  const int transformedIndex( getSrcBlockBdryIndex( block, dir, side ) );
832  a_map[transformedIndex].push_back( index );
833  }
834  }
835  }
836  }
837 }
838 
839 
840 template <class T>
841 inline
843 {
844  if (m_hasData)
845  {
846  for (int dir(0); dir<CH_SPACEDIM; dir++)
847  {
848  for (SideIterator side; side.ok(); ++side)
849  {
850  const int index( getBlockBdryIndex( dir, side() ) );
851  // petermc, 23 May 2013: removed constness in order to do exchange
852  LevelData<T >& srcLevelData( m_excRegister[index] );
853 // BoxLayoutData<T >& dstBoxLayoutData( m_srcRegister[index] );
854  LevelData<T >& dstLevelData( m_srcRegister[index] );
855  if (m_ghost==0)
856  {
857  srcLevelData.copyTo( dstLevelData );
858  }
859  else
860  {
861  Copier copier;
862  copier.ghostDefine( srcLevelData.disjointBoxLayout(),
863  dstLevelData.disjointBoxLayout(),
864  m_domain,
865  srcLevelData.ghostVect(),
866  dstLevelData.ghostVect() );
867  srcLevelData.exchange(); // added by petermc, 23 May 2013
868  srcLevelData.copyTo( dstLevelData, copier );
869  }
870  }
871  }
872  }
873 }
874 
875 
876 inline
877 IntVect permute( const IntVect& a_vec, const IntVect& a_permutation )
878 {
879  IntVect tmp_vec;
880  for (int dir(0); dir<SpaceDim; ++dir)
881  {
882  tmp_vec[dir] = a_vec[a_permutation[dir]];
883  }
884  return tmp_vec;
885 }
886 
887 template <class T>
888 inline
890 {
891  if (m_hasData)
892  {
893  for (int dir(0); dir<CH_SPACEDIM; dir++)
894  {
895  for (SideIterator side; side.ok(); ++side)
896  {
897 
898  const int blockBdryIndex( getBlockBdryIndex( dir, side() ) );
899 
900  DataIterator dit( m_grids.dataIterator() );
901  for (dit.begin(); dit.ok(); ++dit)
902  {
903 
904  if (hasInterface( dit(), dir, side() ))
905  {
906  const int block( getBlockID( dit() ) );
907 
908  const int srcBlockBdryIndex( getSrcBlockBdryIndex( block, dir, side() ) );
909  const IndicesTransformation& itrans( getTransformation( block, blockBdryIndex ) );
910 
911  const LayoutData<DataIndex>& grdToDst( m_grdToDstMap[blockBdryIndex] );
912  const DataIndex& regDataIndex( grdToDst[dit()] );
913 
914  const LayoutData<DataIndex>& dstToSrc( m_dstToSrcMap[blockBdryIndex] );
915  const DataIndex& srcDataIndex( dstToSrc[regDataIndex] );
916 
917  LevelData<T>& dstRegister( m_register[blockBdryIndex] );
918  const LevelData<T>& srcRegister( m_srcRegister[srcBlockBdryIndex] );
919 
920  T& dst( dstRegister[regDataIndex] );
921  const T& src( srcRegister[srcDataIndex] );
922 
923  const Box& srcWindow( src.box() );
924  const Box& regWindow( itrans.transformBack( srcWindow ) );
925  const IntVect perm( itrans.getPermutation() );
926  const IntVect sign( itrans.getSign() );
927  const IntVect trans( itrans.getTranslation() );
928 
929  for (int n(0); n<dst.nComp(); n++)
930  {
931  BoxIterator bit(regWindow);
932  for (bit.begin(); bit.ok(); ++bit)
933  {
934  const IntVect dst_loc( bit() );
935  IntVect src_loc( permute( dst_loc, perm ) );
936  src_loc *= sign;
937  src_loc += trans;
938  dst(dst_loc,n) = src(src_loc,n);
939  }
940  }
941  }
942  }
943  m_register[blockBdryIndex].exchange();
944  }
945  }
946  }
947 }
948 
949 
950 #if 0
951 template <class T>
952 inline
953 void BlockBaseRegister<T>::print( const BoxLayoutData<T > a_register[2*CH_SPACEDIM] ) const
954 {
955  for (int dir(0); dir<CH_SPACEDIM; dir++)
956  {
957  for (SideIterator side; side.ok(); ++side)
958  {
959  pout() << "========= dir: " << dir << " ";
960  pout() << "side: " << side() << " ";
961  const int index( getBlockBdryIndex( dir, side() ) );
962  pout() << "block boundary index: " << index << std::endl;
963  const BoxLayoutData<T >& thisRegister( a_register[index] );
964  if (thisRegister.isDefined())
965  {
966  DataIterator dit( thisRegister.dataIterator() );
967  for (dit.begin(); dit.ok(); ++dit)
968  {
969  const T& data( thisRegister[dit()] );
970  const Box& box( data.box() );
971  pout() << "--------------------------------------" << std::endl;
972  pout() << "Box: " << box << endl;
973  for (int ncomp(0); ncomp<data.nComp(); ++ncomp)
974  {
975  pout() << "Component: " << ncomp << endl;
976  BoxIterator bit( box );
977  for (bit.begin(); bit.ok(); ++bit)
978  {
979  pout() << bit() << "\t" << data.get( bit(), ncomp ) << std::endl;;
980  }
981  }
982  pout() << "--------------------------------------" << std::endl;
983  }
984  }
985  }
986  }
987 }
988 #endif
989 
990 
991 template <class T>
992 inline
993 void BlockBaseRegister<T>::print( const LevelData<T > a_register[2*CH_SPACEDIM] ) const
994 {
995  for (int dir(0); dir<CH_SPACEDIM; dir++)
996  {
997  for (SideIterator side; side.ok(); ++side)
998  {
999  pout() << "========= dir: " << dir << " ";
1000  pout() << "side: " << side() << " ";
1001  const int index( getBlockBdryIndex( dir, side() ) );
1002  pout() << "block boundary index: " << index << std::endl;
1003  const LevelData<T >& thisRegister( a_register[index] );
1004  if (thisRegister.isDefined())
1005  {
1006  DataIterator dit( thisRegister.dataIterator() );
1007  for (dit.begin(); dit.ok(); ++dit)
1008  {
1009  const T& data( thisRegister[dit()] );
1010  const Box& box( data.box() );
1011  pout() << "--------------------------------------" << std::endl;
1012  pout() << "Box: " << box << endl;
1013  for (int ncomp(0); ncomp<data.nComp(); ++ncomp)
1014  {
1015  pout() << "Component: " << ncomp << endl;
1016  BoxIterator bit( box );
1017  for (bit.begin(); bit.ok(); ++bit)
1018  {
1019  pout() << bit() << "\t" << data( bit(), ncomp ) << std::endl;;
1020  }
1021  }
1022  pout() << "--------------------------------------" << std::endl;
1023  }
1024  }
1025  }
1026  }
1027 }
1028 
1029 
1030 template <class T>
1032 {
1033  pout() << ":::::::::::::::::::::::::::::::::::::::::::::::" << std::endl;
1034  pout() << ":: Source " << m_name << " Registers" << std::endl;
1035  pout() << ":::::::::::::::::::::::::::::::::::::::::::::::" << std::endl;
1036  print( m_srcRegister );
1037 
1038  pout() << ":::::::::::::::::::::::::::::::::::::::::::::::" << std::endl;
1039  pout() << ":: Destination " << m_name << " Registers" << std::endl;
1040  pout() << ":::::::::::::::::::::::::::::::::::::::::::::::" << std::endl;
1041  print( m_register );
1042 }
1043 
1044 
1045 #include "NamespaceFooter.H"
1046 
1047 #endif
std::ostream & pout()
Use this in place of std::cout for program output.
bool ok()
Definition: BoxIterator.H:281
void ghostDefine(const DisjointBoxLayout &a_src, const DisjointBoxLayout &a_dest, const ProblemDomain &a_domain, const IntVect &a_srcGhost)
performs the computation necessary for moving from ghost+valid cells to all valid cells in a_dest ...
IntVect transformBack(const IntVect &a_ivNew) const
void buildSrcRegisterLayout()
Build all box layouts for the source boxes.
Definition: BlockBaseRegisterImplem.H:666
int getBlockBdryIndex(const int &a_dir, const Side::LoHiSide &a_side)
Definition: BlockBaseRegisterImplem.H:433
IntVect getTranslation() const
Definition: IndicesTransformation.H:425
const Box & getBlockBox(const int &a_block) const
Returns the box corresponding to the block specified by the block integer index.
Definition: BlockBaseRegisterImplem.H:413
#define CH_SPACEDIM
Definition: SPACE.H:51
A class to facilitate interaction with physical boundary conditions.
Definition: ProblemDomain.H:141
const Box & getNeighborBlockBox(const int &a_block, const int &a_blockBdryIndex) const
Returns the box corresponding to the block opposite to the block specified by the block integer index...
Definition: BlockBaseRegisterImplem.H:422
IntVect getPermutation() const
Definition: IndicesTransformation.H:417
bool isNotEmpty() const
Returns true if the List<T> is not empty.
Definition: List.H:627
void copyToDstRegister()
Local copy of data from exterior source registers to exterior full register.
Definition: BlockBaseRegisterImplem.H:889
void remove(const T &value)
Removes all objects in the List<T> equal to value.
Definition: ListImplem.H:356
Box copy() const
Definition: Box.H:571
void zeroBoxLayoutData(BoxLayoutData< T > &a_layoutData)
Definition: BlockBaseRegisterImplem.H:268
IntVect transformFwd(const IntVect &a_ivOld) const
Class to implement a single data register a mapped multiblock domain block boundary.
Definition: BlockBaseRegister.H:27
void createExchangeBoxLayout(DisjointBoxLayout &a_layout, LayoutData< DataIndex > &a_dstToExcMap, const Vector< Box > &a_boxes, const Vector< int > &a_procs, const Vector< DataIndex > &a_indices, const LayoutData< DataIndex > &a_grdToDstMap)
Build the box layout and associated box mappings for the data exchange.
Definition: BlockBaseRegisterImplem.H:599
A strange but true thing to make copying from one boxlayoutdata to another fast.
Definition: Copier.H:152
void append(const T &value)
Adds a copy of the value to the end of the List<T>.
Definition: List.H:582
IntVect getSign() const
Definition: IndicesTransformation.H:421
void defineDataLayouts(const int a_ncomp)
Define each of the three internally used data layouts.
Definition: BlockBaseRegisterImplem.H:78
BlockBaseRegister()
Empty constructor.
Definition: BlockBaseRegister.H:39
bool hasInterface(const DataIndex &a_dataIndex, int a_dir, Side::LoHiSide a_side) const
Returns whether or not a box has a block-boundary interface on one of its faces.
Definition: BlockBaseRegisterImplem.H:315
Definition: DataIterator.H:190
Box & shiftHalf(int dir, int num_halfs)
iterates through the IntVects of a Box
Definition: BoxIterator.H:37
void clearDataLayouts()
Undefine each of the three internally used data layouts.
Definition: BlockBaseRegisterImplem.H:97
void define(const IntVect &small, const IntVect &big)
void buildExchangeRegisterLayout()
Build all disjoint box layouts for the data exchange.
Definition: BlockBaseRegisterImplem.H:697
void printRegisters() const
Print all registers to pout().
Definition: BlockBaseRegisterImplem.H:1031
An Iterator based on a BoxLayout object.
Definition: LayoutIterator.H:35
Definition: LoHiSide.H:31
void store(const T &a_data, const DataIndex &a_dataIndex, int a_dir, Side::LoHiSide a_side)
Store data on block-boundary faces.
Definition: BlockBaseRegisterImplem.H:110
int getBlockID(const DataIndex &a_dataIndex) const
Returns the integer index of the block in which the box corresponding to the DataIndex resides...
Definition: BlockBaseRegisterImplem.H:404
void zeroRegister(const int a_n_comp=-1)
Definition: BlockBaseRegisterImplem.H:278
bool intersectsNotEmpty(const Box &b) const
const int SpaceDim
Definition: SPACE.H:38
void buildBoxVectors(Vector< Box > a_boxes[2 *CH_SPACEDIM], Vector< int > a_boxProcMap[2 *CH_SPACEDIM], Vector< DataIndex > a_dataIndex[2 *CH_SPACEDIM], const UseType a_type)
Definition: BlockBaseRegisterImplem.H:724
virtual void exchange(void)
Simplest case – do all components.
Definition: LevelDataI.H:451
virtual ~BlockBaseRegister()
virtual destructor
Definition: BlockBaseRegisterImplem.H:39
Definition: EBInterface.H:45
virtual void close()
virtual void copyTo(const Interval &srcComps, BoxLayoutData< T > &dest, const Interval &destComps) const
Definition: LevelDataI.H:202
#define CH_TIME(name)
Definition: CH_Timer.H:82
A Doubly-Linked List Class.
Definition: List.H:21
static const IntVect Unit
Definition: IntVect.H:663
new code
Definition: BoxLayoutData.H:170
void createSrcBoxLayout(DisjointBoxLayout &a_layout, LayoutData< DataIndex > a_dstToSrcMap[2 *CH_SPACEDIM], const Vector< Box > &a_boxes, const Vector< int > &a_procs, const Vector< DataIndex > &a_indices, const Vector< int > &a_dstBlockBdryIndices, const LayoutData< DataIndex > a_grdToDstMap[2 *CH_SPACEDIM])
Build the box layout and associated box mappings for the source boxes.
Definition: BlockBaseRegisterImplem.H:558
virtual bool isDefined() const
Definition: BoxLayoutDataI.H:39
Iterator over a List.
Definition: List.H:20
int sign(const Side::LoHiSide &a_side)
Data on a BoxLayout.
Definition: BoxLayoutData.H:97
Box & growDir(int a_idir, const Side::LoHiSide &a_sd, int a_cell)
Definition: Box.H:2363
Ordered Tuples for Types T.
Definition: Tuple.H:30
const DisjointBoxLayout & disjointBoxLayout() const
Definition: LevelData.H:209
A BoxLayout that has a concept of disjointedness.
Definition: DisjointBoxLayout.H:30
bool isInterface(const int &a_block, const int &a_dir, const Side::LoHiSide &a_side) const
Returns true if the block-boundary in this direction on this side is an interface between blocks...
Definition: BlockBaseRegisterImplem.H:487
LoHiSide
Definition: LoHiSide.H:27
void fill(T &a_data, const DataIndex &a_dataIndex, int a_dir, Side::LoHiSide a_side, Side::LoHiSide a_sideData) const
Return data on a block-boundary interface of a box, if hasInterface(a_dataIndex, a_dir, a_side) is true.
Definition: BlockBaseRegisterImplem.H:344
void buildLocalIndexList(List< int > &a_list, const Vector< Box > &a_boxes, const Vector< int > &a_procs)
Definition: BlockBaseRegisterImplem.H:500
virtual interface class encapsulating multi-block mapping API
Definition: MultiBlockCoordSys.H:34
static void Error(const char *const a_msg=m_nullString, int m_exitCode=CH_DEFAULT_ERROR_CODE)
Print out message to cerr and exit with the specified exit code.
int getSrcBlockBdryIndex(const int &a_block, const int &a_dir, const Side::LoHiSide &a_side)
Returns the index of the neighboring block opposite the block in the direction and side spceified...
Definition: BlockBaseRegisterImplem.H:453
bool isEmpty() const
{ Comparison Functions}
Definition: Box.H:1862
DataIterator dataIterator() const
Parallel iterator.
void print(const LevelData< T > a_register[2 *CH_SPACEDIM]) const
Print the specified register.
Definition: BlockBaseRegisterImplem.H:993
A Rectangular Domain on an Integer Lattice.
Definition: Box.H:469
Iterator for low and high side.
Definition: LoHiSide.H:74
Definition: DataIndex.H:114
void begin()
Definition: BoxIterator.H:150
void buildRegisterLayout()
Build all box layouts for the data register.
Definition: BlockBaseRegisterImplem.H:637
void copyToSrcRegister()
Copy, including parallel exchange, of data from interior exchange registers to exterior source regist...
Definition: BlockBaseRegisterImplem.H:842
ListIterator< T > first() const
Returns a ListIterator<T> to the first object in this List<T>.
Definition: List.H:697
void increment(const T &a_flux, const DataIndex &a_dataIndex, int a_dir, Side::LoHiSide a_side)
Increment data on block-boundary faces.
Definition: BlockBaseRegisterImplem.H:179
void createDstBoxLayout(DisjointBoxLayout &a_layout, LayoutData< DataIndex > &a_map, const Vector< Box > &a_boxes, const Vector< int > &a_procs, const Vector< DataIndex > &a_indices)
Build a box layout for the data register given a vector of boxes with processor mappings and related ...
Definition: BlockBaseRegisterImplem.H:520
IndicesTransformation getTransformation(const int &a_block, const int &a_blockBdryIndex) const
Returns the index transformation box corresponding to the Blockboundary specified by the block intege...
Definition: BlockBaseRegisterImplem.H:441
An integer Vector in SpaceDim-dimensional space.
Definition: CHArray.H:42
bool ok() const
Return true if the iterator is not past the end of the list.
Definition: List.H:465
Definition: LoHiSide.H:30
void exchange()
Exchange data to fill exterior register regions.
Definition: BlockBaseRegisterImplem.H:247
void define(const MultiBlockCoordSys *a_mblock, const DisjointBoxLayout &a_grids, const int &a_ghost=0)
Define function that matches constructor.
Definition: BlockBaseRegisterImplem.H:46
DataIterator dataIterator() const
Definition: LayoutDataI.H:78
Class to describe transformation of indices from one block to another.
Definition: IndicesTransformation.H:25
size_t size() const
Definition: Vector.H:192
Box & grow(int i)
grow functions
Definition: Box.H:2263
virtual bool ok() const
return true if this iterator is still in its Layout
Definition: LayoutIterator.H:117
UseType
Builds an array of vectors of boxes, one for each direction-side combination, that corresponds to eit...
Definition: BlockBaseRegister.H:223
virtual void define(const Vector< Box > &a_boxes, const Vector< int > &a_procIDs)
Box & growLo(int idir, int n_cell=1)
Definition: Box.H:2349
bool ok() const
bool eq(const Box &b) const
void buildInverseBlockBdryIndexMap(Vector< int > a_map[2 *CH_SPACEDIM])
Builds a mapping from the block-boundary index of neighboring block to the local block-boundary index...
Definition: BlockBaseRegisterImplem.H:816
Box & growHi(int idir, int n_cell=1)
Definition: Box.H:2380
IntVect permute(const IntVect &a_vec, const IntVect &a_permutation)
Definition: BlockBaseRegisterImplem.H:877
const IntVect & ghostVect() const
Definition: LevelData.H:170