36#ifndef OPENRS_GRIDINTERFACEEULER_HEADER
37#define OPENRS_GRIDINTERFACEEULER_HEADER
39#include <dune/common/fvector.hh>
40#include <opm/core/utility/SparseTable.hpp>
41#include <opm/core/utility/StopWatch.hpp>
42#include <dune/grid/common/mcmgmapper.hh>
43#include <dune/grid/CpGrid.hpp>
45#include <boost/iterator/iterator_facade.hpp>
46#include <boost/scoped_ptr.hpp>
59 bool contains(Dune::GeometryType gt) {
return gt.dim() == dim; }
63 template <
class DuneGr
id>
66 typedef Dune::LeafMultipleCodimMultipleGeomTypeMapper<DuneGrid, AllCellsLayout>
Type;
75 template<
class EntityType>
76 int map (
const EntityType& e)
const
94 template <
class Gr
idInterface,
class EntityPo
interType>
101 typedef typename GI::GridType::template Codim<0>::EntityPointer
CellPtr;
102 typedef typename GI::GridType::ctype
Scalar;
103 typedef Dune::FieldVector<Scalar, GI::GridType::dimension>
Vector;
128 return iter_->geometry().volume();
137 return iter_->geometry().center();
146 return iter_->centerUnitOuterNormal();
154 return iter_->boundary();
162 return iter_->boundaryId();
194 if (
iter_->boundary()) {
222 return iter_->outside()->geometry().volume();
244 template <
class Gr
idInterface>
246 public boost::iterator_facade<FaceIterator<GridInterface>,
247 const Face<GridInterface>,
248 boost::forward_traversal_tag>,
249 public Face<GridInterface>
279 const int local_index)
317 template <
class Gr
idInterface,
class EntityPo
interType>
325 Cell(
const GridInterface& grid,
326 const EntityPointerType& it)
348 return iter_->geometry().volume();
357 return iter_->geometry().center();
370 template <
class Gr
idInterface>
372 :
public boost::iterator_facade<CellIterator<GridInterface>,
373 const Cell<GridInterface,
374 typename GridInterface::GridType::template Codim<0>::LeafIterator>,
375 boost::forward_traversal_tag>,
376 public Cell<GridInterface, typename GridInterface::GridType::template Codim<0>::LeafIterator>
379 typedef typename GridInterface::GridType::template Codim<0>::LeafIterator DuneCellIter;
410 template <
class DuneGr
id>
427 : pgrid_(0), num_faces_(0), max_faces_per_cell_(0)
431 : pgrid_(&
grid), pmapper_(new
Mapper(
grid)), num_faces_(0), max_faces_per_cell_(0)
437 void init(
const DuneGrid&
grid,
bool build_facemap =
true)
455 return grid().size(0);
459 assert(num_faces_ != 0);
464 assert(max_faces_per_cell_ != 0);
465 return max_faces_per_cell_;
482 assert(num_faces_ != 0);
483 return face_indices_[cell_index][local_face_index];
485 typedef Opm::SparseTable<int>::row_type
Indices;
488 assert(num_faces_ != 0);
489 return face_indices_[cell_index];
492 const DuneGrid* pgrid_;
493 boost::scoped_ptr<Mapper> pmapper_;
495 int max_faces_per_cell_;
496 Opm::SparseTable<int> face_indices_;
498 void buildFaceIndices()
501 std::cout <<
"Building unique face indices... " << std::flush;
502 Opm::time::StopWatch clock;
506 typedef typename CI::FaceIterator FI;
516 std::vector<int> cell(nc, -1);
517 std::vector<int> num_faces(nc);
518 std::vector<int> fpos; fpos .reserve(nc + 1);
519 std::vector<int> num_cf; num_cf.reserve(nc);
520 std::vector<int> faces ;
523 int cellno = 0; fpos.push_back(0);
524 int tot_ncf = 0, tot_ncf2 = 0, max_ncf = 0;
526 const int c0 = c->index();
527 assert((0 <= c0) && (c0 < nc) && (cell[c0] == -1));
530 int& ncf = num_cf.back();
531 for (FI f = c->facebegin(); f != c-> faceend(); ++f) {
532 if (!f->boundary()) {
533 const int c1 = f->neighbourCellIndex();
534 assert((0 <= c1) && (c1 < nc) && (c1 != c0));
536 if (cell[c1] == -1) {
544 fpos.push_back(
int(faces.size()));
545 max_ncf = std::max(max_ncf, ncf);
547 tot_ncf2 += ncf * ncf;
549 assert(cellno == nc);
554 cumul_num_faces[0] = 0;
555 std::partial_sum(num_faces.begin(), num_faces.end(), cumul_num_faces.begin() + 1);
558 std::vector<int> l2g;
559 l2g.reserve(max_ncf);
560 std::vector<double> cfdata(tot_ncf);
561 int total_num_faces = int(faces.size());
564 typedef std::vector<int>::iterator VII;
566 const int c0 = c->index();
567 assert ((0 <= c0 ) && ( c0 < nc) &&
568 (0 <= cell[c0]) && (cell[c0] < nc));
569 const int ncf = num_cf[cell[c0]];
571 for (FI f = c->facebegin(); f != c->faceend(); ++f) {
574 l2g[f->localIndex()] = total_num_faces++;
585 const int c1 = f->neighbourCellIndex();
586 assert ((0 <= c1 ) && ( c1 < nc) &&
587 (0 <= cell[c1]) && (cell[c1] < nc));
589 int t = c0, seek = c1;
590 if (cell[seek] < cell[t])
592 int s = fpos[cell[t]], e = fpos[cell[t] + 1];
593 VII p = std::find(faces.begin() + s, faces.begin() + e, seek);
594 assert(p != faces.begin() + e);
595 l2g[f->localIndex()] = p - faces.begin();
598 assert(
int(l2g.size()) == num_faces[c0]);
599 std::copy(l2g.begin(), l2g.end(), cfdata.begin() + cumul_num_faces[c0]);
601 num_faces_ = total_num_faces;
602 max_faces_per_cell_ = max_ncf;
603 face_indices_.assign(cfdata.begin(), cfdata.end(),
604 num_faces.begin(), num_faces.end());
608 double elapsed = clock.secsSinceStart();
609 std::cout <<
"done. Time elapsed: " << elapsed << std::endl;
Definition: GridInterfaceEuler.hpp:377
const CellIterator & dereference() const
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:391
CellType::Vector Vector
Definition: GridInterfaceEuler.hpp:382
CellIterator(const GridInterface &grid, DuneCellIter it)
Definition: GridInterfaceEuler.hpp:386
CellType::Scalar Scalar
Definition: GridInterfaceEuler.hpp:383
bool equal(const CellIterator &other) const
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:396
void increment()
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:401
CellType::Index Index
Definition: GridInterfaceEuler.hpp:384
Definition: GridInterfaceEuler.hpp:319
Scalar volume() const
Definition: GridInterfaceEuler.hpp:346
Vector centroid() const
Definition: GridInterfaceEuler.hpp:351
FaceIterator faceend() const
Definition: GridInterfaceEuler.hpp:340
Cell()
Definition: GridInterfaceEuler.hpp:321
FaceIterator::Scalar Scalar
Definition: GridInterfaceEuler.hpp:332
Cell(const GridInterface &grid, const EntityPointerType &it)
Definition: GridInterfaceEuler.hpp:325
FaceIterator::Index Index
Definition: GridInterfaceEuler.hpp:333
Index index() const
Definition: GridInterfaceEuler.hpp:360
FaceIterator facebegin() const
Definition: GridInterfaceEuler.hpp:335
FaceIterator::Vector Vector
Definition: GridInterfaceEuler.hpp:331
EntityPointerType iter_
Definition: GridInterfaceEuler.hpp:366
const GridInterface * pgrid_
Definition: GridInterfaceEuler.hpp:365
GIE::FaceIterator< GridInterface > FaceIterator
Definition: GridInterfaceEuler.hpp:330
Intersection (face) iterator for solver-near grid interface.
Definition: GridInterfaceEuler.hpp:250
bool equal(const FaceIterator &other) const
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:291
void increment()
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:299
Face< GridInterface >::DuneIntersectionIter DuneIntersectionIter
Type of low-level intersection iterator. Copied from the Dune grid.
Definition: GridInterfaceEuler.hpp:258
bool operator<(const FaceIterator &other) const
Gives an ordering of intersectionIterators.
Definition: GridInterfaceEuler.hpp:306
const FaceIterator & dereference() const
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:285
FaceIterator(const GridInterface &grid, const DuneIntersectionIter &it, const int local_index)
Constructor.
Definition: GridInterfaceEuler.hpp:277
FaceIterator()
Default constructor.
Definition: GridInterfaceEuler.hpp:261
Definition: GridInterfaceEuler.hpp:98
Dune::FieldVector< Scalar, GI::GridType::dimension-1 > LocalVector
Definition: GridInterfaceEuler.hpp:104
Index local_index_
Definition: GridInterfaceEuler.hpp:228
const GI * pgrid_
Definition: GridInterfaceEuler.hpp:226
Vector normal() const
Definition: GridInterfaceEuler.hpp:143
Vector centroid() const
Definition: GridInterfaceEuler.hpp:134
GI::DuneIntersectionIterator DuneIntersectionIter
Definition: GridInterfaceEuler.hpp:100
Face(const GI &grid, const DuneIntersectionIter &it, const Index loc_ind)
Definition: GridInterfaceEuler.hpp:116
GI::GridType::template Codim< 0 >::EntityPointer CellPtr
Definition: GridInterfaceEuler.hpp:101
Index neighbourCellIndex() const
Definition: GridInterfaceEuler.hpp:192
DuneIntersectionIter iter_
Definition: GridInterfaceEuler.hpp:227
bool boundary() const
Definition: GridInterfaceEuler.hpp:152
int Index
Definition: GridInterfaceEuler.hpp:105
@ LocalEndIndex
Definition: GridInterfaceEuler.hpp:109
@ BoundaryMarkerIndex
Definition: GridInterfaceEuler.hpp:108
Cell cell() const
Definition: GridInterfaceEuler.hpp:168
GIE::Cell< GI, CellPtr > Cell
Definition: GridInterfaceEuler.hpp:106
Index cellIndex() const
Definition: GridInterfaceEuler.hpp:176
GI::GridType::ctype Scalar
Definition: GridInterfaceEuler.hpp:102
Cell neighbourCell() const
Definition: GridInterfaceEuler.hpp:184
int boundaryId() const
Definition: GridInterfaceEuler.hpp:160
Index localIndex() const
Definition: GridInterfaceEuler.hpp:212
Scalar neighbourCellVolume() const
Definition: GridInterfaceEuler.hpp:220
Face()
Definition: GridInterfaceEuler.hpp:111
Dune::FieldVector< Scalar, GI::GridType::dimension > Vector
Definition: GridInterfaceEuler.hpp:103
Scalar area() const
Definition: GridInterfaceEuler.hpp:126
Index index() const
Definition: GridInterfaceEuler.hpp:204
Definition: GridInterfaceEuler.hpp:412
int numberOfFaces() const
Definition: GridInterfaceEuler.hpp:457
Opm::SparseTable< int >::row_type Indices
Definition: GridInterfaceEuler.hpp:485
void init(const DuneGrid &grid, bool build_facemap=true)
Definition: GridInterfaceEuler.hpp:437
Indices faceIndices(int cell_index) const
Definition: GridInterfaceEuler.hpp:486
const DuneGrid & grid() const
Definition: GridInterfaceEuler.hpp:467
@ Dimension
Definition: GridInterfaceEuler.hpp:424
const Mapper & mapper() const
Definition: GridInterfaceEuler.hpp:475
GridInterfaceEuler< DuneGrid > InterfaceType
Definition: GridInterfaceEuler.hpp:416
CellIterator::Index Index
Definition: GridInterfaceEuler.hpp:420
CellIterator cellend() const
Definition: GridInterfaceEuler.hpp:449
GIE::CellIterator< InterfaceType > CellIterator
Definition: GridInterfaceEuler.hpp:417
GridInterfaceEuler()
Definition: GridInterfaceEuler.hpp:426
CellIterator cellbegin() const
Definition: GridInterfaceEuler.hpp:445
GICellMapper< DuneGrid >::Type Mapper
Definition: GridInterfaceEuler.hpp:422
CellIterator::Scalar Scalar
Definition: GridInterfaceEuler.hpp:419
Index faceIndex(int cell_index, int local_face_index) const
Definition: GridInterfaceEuler.hpp:480
CellIterator::Vector Vector
Definition: GridInterfaceEuler.hpp:418
int numberOfCells() const
Definition: GridInterfaceEuler.hpp:453
GridInterfaceEuler(const DuneGrid &grid, bool build_facemap=true)
Definition: GridInterfaceEuler.hpp:430
DuneGrid GridType
Definition: GridInterfaceEuler.hpp:415
DuneGrid::LeafIntersectionIterator DuneIntersectionIterator
Definition: GridInterfaceEuler.hpp:414
int maxFacesPerCell() const
Definition: GridInterfaceEuler.hpp:462
Definition: BlackoilFluid.hpp:32
General cell layout.
Definition: GridInterfaceEuler.hpp:58
bool contains(Dune::GeometryType gt)
Definition: GridInterfaceEuler.hpp:59
A mapper for Dune::CpGrid cells only.
Definition: GridInterfaceEuler.hpp:71
int map(const EntityType &e) const
Definition: GridInterfaceEuler.hpp:76
CpGridCellMapper(const Dune::CpGrid &)
Definition: GridInterfaceEuler.hpp:72
CpGridCellMapper Type
Definition: GridInterfaceEuler.hpp:88
Mapper for general grids.
Definition: GridInterfaceEuler.hpp:65
Dune::LeafMultipleCodimMultipleGeomTypeMapper< DuneGrid, AllCellsLayout > Type
Definition: GridInterfaceEuler.hpp:66