dune-grid  2.2.1
gridptr.hh
Go to the documentation of this file.
1 #ifndef DUNE_DGF_GRIDPTR_HH
2 #define DUNE_DGF_GRIDPTR_HH
3 
4 #include <iostream>
5 #include <string>
6 #include <vector>
7 #include <memory>
8 #include <map>
9 #include <assert.h>
10 
11 //- Dune includes
12 #include <dune/common/mpihelper.hh>
15 
19 
21 
22 namespace Dune
23 {
24  // forward declarations
25  // --------------------
26  template < class G >
27  struct DGFGridFactory;
28 
29  template< class GridImp, template < class > class IntersectionImp >
30  class Intersection;
31 
44  template< class GridType >
45  struct GridPtr
46  {
47  typedef MPIHelper::MPICommunicator MPICommunicatorType;
48  static const int dimension = GridType::dimension;
49 
51  explicit GridPtr ( const std::string &filename,
52  MPICommunicatorType comm = MPIHelper::getCommunicator() )
53  : gridPtr_( 0 ),
54  elParam_(),
55  vtxParam_(),
56  bndParam_(),
57  bndId_(),
58  emptyParam_(),
59  nofElParam_( 0 ),
60  nofVtxParam_( 0 ),
61  haveBndParam_( false )
62  {
63  DGFGridFactory< GridType > dgfFactory( filename, comm );
64  initialize( dgfFactory );
65  }
66 
68  explicit GridPtr ( std::istream &input,
69  MPICommunicatorType comm = MPIHelper::getCommunicator() )
70  : gridPtr_( 0 ),
71  elParam_(),
72  vtxParam_(),
73  bndParam_(),
74  bndId_(),
75  emptyParam_(),
76  nofElParam_( 0 ),
77  nofVtxParam_( 0 ),
78  haveBndParam_( false )
79  {
80  DGFGridFactory< GridType > dgfFactory( input, comm );
81  initialize( dgfFactory );
82  }
83 
86  : gridPtr_(0),
87  elParam_(),
88  vtxParam_(),
89  bndParam_(),
90  bndId_(),
91  emptyParam_(),
92  nofElParam_(0),
93  nofVtxParam_(0),
94  haveBndParam_( false )
95  {}
96 
98  explicit GridPtr( GridType *grd )
99  : gridPtr_(grd),
100  elParam_(),
101  vtxParam_(),
102  bndParam_(),
103  bndId_(),
104  emptyParam_(),
105  nofElParam_(0),
106  nofVtxParam_(0),
107  haveBndParam_( false )
108  {}
109 
111  GridPtr( const GridPtr &org )
112  : gridPtr_(org.gridPtr_),
113  elParam_(org.elParam_),
114  vtxParam_(org.vtxParam_),
115  bndParam_(org.bndParam_),
116  bndId_(org.bndId_),
117  emptyParam_( org.emptyParam_ ),
121  {}
122 
124  GridPtr &operator= ( const GridPtr &org )
125  {
126  gridPtr_ = org.gridPtr_;
127  elParam_ = org.elParam_;
128  vtxParam_ = org.vtxParam_;
129  bndParam_ = org.bndParam_;
130  bndId_ = org.bndId_;
131  emptyParam_ = org.emptyParam_;
132 
133  nofElParam_ = org.nofElParam_;
134  nofVtxParam_ = org.nofVtxParam_;
136  return *this;
137  }
138 
140  GridPtr & operator = (GridType * grd)
141  {
142  gridPtr_ = std::auto_ptr<GridType>(grd);
143  elParam_.resize(0);
144  vtxParam_.resize(0);
145  bndParam_.resize(0);
146  bndId_.resize(0);
147  emptyParam_.resize(0);
148 
149  nofVtxParam_ = 0;
150  nofElParam_ = 0;
151  haveBndParam_ = false;
152  return *this;
153  }
154 
156  GridType& operator*() {
157  return *gridPtr_;
158  }
159 
161  GridType* operator->() {
162  return gridPtr_.operator -> ();
163  }
164 
166  const GridType& operator*() const {
167  return *gridPtr_;
168  }
169 
171  const GridType* operator->() const {
172  return gridPtr_.operator -> ();
173  }
174 
176  GridType* release () {
177  return gridPtr_.release();
178  }
179 
181  int nofParameters(int cdim) const {
182  switch (cdim) {
183  case 0: return nofElParam_; break;
184  case GridType::dimension: return nofVtxParam_; break;
185  }
186  return 0;
187  }
188 
190  template <class Entity>
191  int nofParameters ( const Entity & ) const
192  {
193  return nofParameters( (int) Entity::codimension );
194  }
195 
197  template< class GridImp, template< class > class IntersectionImp >
198  int nofParameters ( const Intersection< GridImp, IntersectionImp > & intersection ) const
199  {
200  return parameters( intersection ).size();
201  }
202 
204  template <class Entity>
205  const std::vector< double > &parameters ( const Entity &entity ) const
206  {
207  typedef typename GridType::LevelGridView GridView;
208  GridView gridView = gridPtr_->levelView( 0 );
209  switch( (int)Entity::codimension )
210  {
211  case 0:
212  if( nofElParam_ > 0 )
213  {
214  assert( (unsigned int)gridView.indexSet().index( entity ) < elParam_.size() );
215  return elParam_[ gridView.indexSet().index( entity ) ];
216  }
217  break;
218  case GridType::dimension:
219  if( nofVtxParam_ > 0 )
220  {
221  assert( (unsigned int)gridView.indexSet().index( entity ) < vtxParam_.size() );
222  return vtxParam_[ gridView.indexSet().index( entity ) ];
223  }
224  break;
225  }
226  return emptyParam_;
227  }
228 
230  template< class GridImp, template< class > class IntersectionImp >
232  {
233  // if no parameters given return empty vecto
234  if ( !haveBndParam_ )
236 
237  return bndParam_[ intersection.boundarySegmentIndex() ];
238  }
239 
240  void loadBalance()
241  {
242  if ( gridPtr_->comm().size() == 1 )
243  return;
245  if ( haveBndParam_ )
246  params += 1;
247  if ( gridPtr_->comm().max( params ) > 0 )
248  {
249  DataHandle dh(*this);
250  gridPtr_->loadBalance( dh.interface() );
252  } else
253  {
254  gridPtr_->loadBalance();
255  }
256  }
257 
258  protected:
260  {
261  gridPtr_ = std::auto_ptr< GridType >( dgfFactory.grid() );
262 
263  typedef typename GridType::LevelGridView GridView;
264  GridView gridView = gridPtr_->levelView( 0 );
265  const typename GridView::IndexSet &indexSet = gridView.indexSet();
266 
267  nofElParam_ = dgfFactory.template numParameters< 0 >();
268  nofVtxParam_ = dgfFactory.template numParameters< dimension >();
269  haveBndParam_ = dgfFactory.haveBoundaryParameters();
270 
271  if ( nofElParam_ > 0 )
272  elParam_.resize( indexSet.size(0) );
273  if ( nofVtxParam_ > 0 )
274  vtxParam_.resize( indexSet.size(dimension) );
275  bndId_.resize( indexSet.size(1) );
276  if ( haveBndParam_ )
277  bndParam_.resize( gridPtr_->numBoundarySegments() );
278 
279  const PartitionIteratorType partType = Interior_Partition;
280  typedef typename GridView::template Codim< 0 >::template Partition< partType >::Iterator Iterator;
281  const Iterator enditer = gridView.template end< 0, partType >();
282  for( Iterator iter = gridView.template begin< 0, partType >(); iter != enditer; ++iter )
283  {
284  const typename Iterator::Entity &el = *iter;
285  if ( nofElParam_ > 0 ) {
286  std::swap( elParam_[ indexSet.index(el) ], dgfFactory.parameter(el) );
287  assert( elParam_[ indexSet.index(el) ].size() == (size_t)nofElParam_ );
288  }
289  if ( nofVtxParam_ > 0 )
290  {
291  for ( int v = 0; v < el.template count<dimension>(); ++v)
292  {
293  typename GridView::IndexSet::IndexType index = indexSet.subIndex(el,v,dimension);
294  if ( vtxParam_[ index ].empty() )
295  std::swap( vtxParam_[ index ], dgfFactory.parameter(*el.template subEntity<dimension>(v) ) );
296  assert( vtxParam_[ index ].size() == (size_t)nofVtxParam_ );
297  }
298  }
299  if ( el.hasBoundaryIntersections() )
300  {
302  const IntersectionIterator iend = gridView.iend(el);
303  for( IntersectionIterator iiter = gridView.ibegin(el); iiter != iend; ++iiter )
304  {
305  const typename IntersectionIterator::Intersection &inter = *iiter;
306  // dirty hack: check for "none" to make corner point grid work
307  if ( inter.boundary() && !inter.type().isNone() )
308  {
309  const int k = indexSet.subIndex(el,inter.indexInInside(),1);
310  bndId_[ k ] = dgfFactory.boundaryId( inter );
311  if ( haveBndParam_ )
312  bndParam_[ inter.boundarySegmentIndex() ] = dgfFactory.boundaryParameter( inter );
313  }
314  }
315  }
316  }
317  }
318 
319  template <class Entity>
320  std::vector< double > &params ( const Entity &entity )
321  {
322  typedef typename GridType::LevelGridView GridView;
323  GridView gridView = gridPtr_->levelView( 0 );
324  switch( (int)Entity::codimension )
325  {
326  case 0:
327  if( nofElParam_ > 0 ) {
328  if ( gridView.indexSet().index( entity ) >= elParam_.size() )
329  elParam_.resize( gridView.indexSet().index( entity ) );
330  return elParam_[ gridView.indexSet().index( entity ) ];
331  }
332  break;
333  case GridType::dimension:
334  if( nofVtxParam_ > 0 ) {
335  if ( gridView.indexSet().index( entity ) >= vtxParam_.size() )
336  vtxParam_.resize( gridView.indexSet().index( entity ) );
337  return vtxParam_[ gridView.indexSet().index( entity ) ];
338  }
339  break;
340  }
341  return emptyParam_;
342  }
343 
344  void setNofParams( int cdim, int nofP )
345  {
346  switch (cdim) {
347  case 0: nofElParam_ = nofP; break;
348  case GridType::dimension: nofVtxParam_ = nofP; break;
349  }
350  }
351 
352  struct DataHandle : public CommDataHandleIF<DataHandle,double>
353  {
354  DataHandle( GridPtr& gridPtr) :
355  gridPtr_(gridPtr),
356  idSet_(gridPtr->localIdSet())
357  {
358  typedef typename GridType::LevelGridView GridView;
359  GridView gridView = gridPtr_->levelView( 0 );
360  const typename GridView::IndexSet &indexSet = gridView.indexSet();
361 
362  const PartitionIteratorType partType = Interior_Partition;
363  typedef typename GridView::template Codim< 0 >::template Partition< partType >::Iterator Iterator;
364  const Iterator enditer = gridView.template end< 0, partType >();
365  for( Iterator iter = gridView.template begin< 0, partType >(); iter != enditer; ++iter )
366  {
367  const typename Iterator::Entity &el = *iter;
368  if ( gridPtr_.nofElParam_ > 0 )
369  std::swap( gridPtr_.elParam_[ indexSet.index(el) ], elData_[ idSet_.id(el) ] );
370  if ( gridPtr_.nofVtxParam_ > 0 )
371  {
372  for ( int v = 0; v < el.template count<dimension>(); ++v)
373  {
374  typename GridView::IndexSet::IndexType index = indexSet.subIndex(el,v,dimension);
375  if ( ! gridPtr_.vtxParam_[ index ].empty() )
376  std::swap( gridPtr_.vtxParam_[ index ], vtxData_[ idSet_.subId(el,v,dimension) ] );
377  }
378  }
379  }
380  }
381 
383  {
384  typedef typename GridType::LevelGridView GridView;
385  GridView gridView = gridPtr_->levelView( 0 );
386  const typename GridView::IndexSet &indexSet = gridView.indexSet();
387 
388  if ( gridPtr_.nofElParam_ > 0 )
389  gridPtr_.elParam_.resize( indexSet.size(0) );
390  if ( gridPtr_.nofVtxParam_ > 0 )
391  gridPtr_.vtxParam_.resize( indexSet.size(dimension) );
392 
393  const PartitionIteratorType partType = All_Partition;
394  typedef typename GridView::template Codim< 0 >::template Partition< partType >::Iterator Iterator;
395  const Iterator enditer = gridView.template end< 0, partType >();
396  for( Iterator iter = gridView.template begin< 0, partType >(); iter != enditer; ++iter )
397  {
398  const typename Iterator::Entity &el = *iter;
399  if ( gridPtr_.nofElParam_ > 0 )
400  {
401  std::swap( gridPtr_.elParam_[ indexSet.index(el) ], elData_[ idSet_.id(el) ] );
402  assert( gridPtr_.elParam_[ indexSet.index(el) ].size() == (unsigned int)gridPtr_.nofElParam_ );
403  }
404  if ( gridPtr_.nofVtxParam_ > 0 )
405  {
406  for ( int v = 0; v < el.template count<dimension>(); ++v)
407  {
408  typename GridView::IndexSet::IndexType index = indexSet.subIndex(el,v,dimension);
409  if ( gridPtr_.vtxParam_[ index ].empty() )
410  std::swap( gridPtr_.vtxParam_[ index ], vtxData_[ idSet_.subId(el,v,dimension) ] );
411  assert( gridPtr_.vtxParam_[ index ].size() == (unsigned int)gridPtr_.nofVtxParam_ );
412  }
413  }
414  }
415  }
416 
418  {
419  return *this;
420  }
421 
422  bool contains (int dim, int codim) const
423  {
424  return (codim==dim || codim==0);
425  }
426 
427  bool fixedsize (int dim, int codim) const
428  {
429  return false;
430  }
431 
432  template<class EntityType>
433  size_t size (const EntityType& e) const
434  {
435  return gridPtr_.nofParameters( (int) e.codimension);
436  }
437 
438  template<class MessageBufferImp, class EntityType>
439  void gather (MessageBufferImp& buff, const EntityType& e) const
440  {
441  const std::vector<double> &v = (e.codimension==0)? elData_[idSet_.id(e)]:vtxData_[idSet_.id(e)];
442  const size_t s = v.size();
443  for (size_t i=0;i<s;++i)
444  buff.write( v[i] );
445  assert( s == (size_t)gridPtr_.nofParameters(e.codimension) );
446  }
447 
448  template<class MessageBufferImp, class EntityType>
449  void scatter (MessageBufferImp& buff, const EntityType& e, size_t n)
450  {
451  std::vector<double> &v = (e.codimension==0)? elData_[idSet_.id(e)]:vtxData_[idSet_.id(e)];
452  v.resize( n );
453  gridPtr_.setNofParams( e.codimension, n );
454  for (size_t i=0;i<n;++i)
455  buff.read( v[i] );
456  }
457 
458  private:
459  typedef typename GridType::LocalIdSet IdSet;
460  GridPtr &gridPtr_;
461  const IdSet &idSet_;
462  mutable std::map< typename IdSet::IdType, std::vector<double> > elData_, vtxData_;
463  };
464 
465  // grid auto pointer
466  mutable std::auto_ptr<GridType> gridPtr_;
467  // element and vertex parameters
468  std::vector< std::vector< double > > elParam_;
469  std::vector< std::vector< double > > vtxParam_;
470  std::vector< DGFBoundaryParameter::type > bndParam_;
471  std::vector< int > bndId_;
472  std::vector < double > emptyParam_;
473 
477  }; // end of class GridPtr
478 
479 } // end namespace Dune
480 
481 #endif