36#ifndef OPENRS_GRIDINTERFACEEULER_HEADER
37#define OPENRS_GRIDINTERFACEEULER_HEADER
39#include <opm/grid/utility/SparseTable.hpp>
40#include <opm/grid/utility/StopWatch.hpp>
41#include <opm/grid/CpGrid.hpp>
43#include <opm/common/utility/platform_dependent/disable_warnings.h>
45#include <dune/common/fvector.hh>
46#include <dune/grid/common/mcmgmapper.hh>
47#include <dune/grid/common/defaultgridview.hh>
49#include <boost/iterator/iterator_facade.hpp>
51#include <opm/common/utility/platform_dependent/reenable_warnings.h>
63 template <
class DuneGr
id>
66 using Type = Dune::MultipleCodimMultipleGeomTypeMapper<Dune::GridView<Dune::DefaultLeafGridViewTraits<DuneGrid>>>;
75 template<
class EntityType>
76 int map (
const EntityType& e)
const
84 template <
class Gr
idInterface,
class EntityPo
interType>
91 typedef typename GI::GridType::Traits::template Codim<0>::Entity
CellPtr;
92 typedef typename GI::GridType::ctype
Scalar;
93 typedef Dune::FieldVector<Scalar, GI::GridType::dimension>
Vector;
118 return iter_->geometry().volume();
127 return iter_->geometry().center();
136 return iter_->centerUnitOuterNormal();
144 return iter_->boundary();
152 return iter_->boundaryId();
184 if (
iter_->boundary()) {
187 return pgrid_->mapper().index(
213 return iter_->outside()->geometry().volume();
235 template <
class Gr
idInterface>
237 public boost::iterator_facade<FaceIterator<GridInterface>,
238 const Face<GridInterface>,
239 boost::forward_traversal_tag>,
240 public Face<GridInterface>
270 const int local_index)
308 template <
class Gr
idInterface,
class EntityType>
316 Cell(
const GridInterface& grid,
317 const EntityType& cell)
339 return cell_.geometry().volume();
348 return cell_.geometry().center();
361 template <
class Gr
idInterface>
363 :
public boost::iterator_facade<CellIterator<GridInterface>,
364 const Cell<GridInterface,
365 typename GridInterface::GridType::template Codim<0>::LeafIterator>,
366 boost::forward_traversal_tag>,
367 public Cell<GridInterface, typename GridInterface::GridType::template Codim<0>::LeafIterator>
370 typedef typename GridInterface::GridType::template Codim<0>::LeafIterator DuneCellIter;
401 template <
class DuneGr
id>
418 : pgrid_(0), num_faces_(0), max_faces_per_cell_(0)
422 : pgrid_(&
grid), pmapper_(new
Mapper(
grid.leafGridView(), Dune::mcmgElementLayout())), num_faces_(0), max_faces_per_cell_(0)
428 void init(
const DuneGrid&
grid,
bool build_facemap =
true)
431 pmapper_.reset(
new Mapper(
grid.leafGridView(), Dune::mcmgElementLayout()));
446 return grid().size(0);
450 assert(num_faces_ != 0);
455 assert(max_faces_per_cell_ != 0);
456 return max_faces_per_cell_;
473 assert(num_faces_ != 0);
474 return face_indices_[cell_index][local_face_index];
476 typedef Opm::SparseTable<int>::row_type
Indices;
479 assert(num_faces_ != 0);
480 return face_indices_[cell_index];
483 const DuneGrid* pgrid_;
484 std::unique_ptr<Mapper> pmapper_;
486 int max_faces_per_cell_;
487 Opm::SparseTable<int> face_indices_;
489 void buildFaceIndices()
492 std::cout <<
"Building unique face indices... " << std::flush;
493 Opm::time::StopWatch clock;
497 typedef typename CI::FaceIterator FI;
507 std::vector<int> cell(nc, -1);
508 std::vector<int> num_faces(nc);
509 std::vector<int> fpos; fpos .reserve(nc + 1);
510 std::vector<int> num_cf; num_cf.reserve(nc);
511 std::vector<int> faces ;
514 int cellno = 0; fpos.push_back(0);
515 int tot_ncf = 0, max_ncf = 0;
517 const int c0 = c->index();
518 assert((0 <= c0) && (c0 < nc) && (cell[c0] == -1));
521 int& ncf = num_cf.back();
522 for (FI f = c->facebegin(); f != c-> faceend(); ++f) {
523 if (!f->boundary()) {
524 const int c1 = f->neighbourCellIndex();
525 assert((0 <= c1) && (c1 < nc) && (c1 != c0));
527 if (cell[c1] == -1) {
535 fpos.push_back(
int(faces.size()));
536 max_ncf = std::max(max_ncf, ncf);
539 assert(cellno == nc);
544 cumul_num_faces[0] = 0;
545 std::partial_sum(num_faces.begin(), num_faces.end(), cumul_num_faces.begin() + 1);
548 std::vector<int> l2g;
549 l2g.reserve(max_ncf);
550 std::vector<double> cfdata(tot_ncf);
551 int total_num_faces = int(faces.size());
554 typedef std::vector<int>::iterator VII;
556 const int c0 = c->index();
557 assert ((0 <= c0 ) && ( c0 < nc) &&
558 (0 <= cell[c0]) && (cell[c0] < nc));
559 const int ncf = num_cf[cell[c0]];
561 for (FI f = c->facebegin(); f != c->faceend(); ++f) {
564 l2g[f->localIndex()] = total_num_faces++;
575 const int c1 = f->neighbourCellIndex();
576 assert ((0 <= c1 ) && ( c1 < nc) &&
577 (0 <= cell[c1]) && (cell[c1] < nc));
579 int t = c0, seek = c1;
580 if (cell[seek] < cell[t])
582 int s = fpos[cell[t]], e = fpos[cell[t] + 1];
583 VII p = std::find(faces.begin() + s, faces.begin() + e, seek);
584 assert(p != faces.begin() + e);
585 l2g[f->localIndex()] = p - faces.begin();
588 assert(
int(l2g.size()) == num_faces[c0]);
589 std::copy(l2g.begin(), l2g.end(), cfdata.begin() + cumul_num_faces[c0]);
591 num_faces_ = total_num_faces;
592 max_faces_per_cell_ = max_ncf;
593 face_indices_.assign(cfdata.begin(), cfdata.end(),
594 num_faces.begin(), num_faces.end());
598 double elapsed = clock.secsSinceStart();
599 std::cout <<
"done. Time elapsed: " << elapsed << std::endl;
Definition: GridInterfaceEuler.hpp:368
const CellIterator & dereference() const
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:382
CellType::Vector Vector
Definition: GridInterfaceEuler.hpp:373
CellIterator(const GridInterface &grid, DuneCellIter it)
Definition: GridInterfaceEuler.hpp:377
CellType::Scalar Scalar
Definition: GridInterfaceEuler.hpp:374
bool equal(const CellIterator &other) const
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:387
void increment()
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:392
CellType::Index Index
Definition: GridInterfaceEuler.hpp:375
Definition: GridInterfaceEuler.hpp:310
Cell()
Definition: GridInterfaceEuler.hpp:312
FaceIterator::Scalar Scalar
Definition: GridInterfaceEuler.hpp:323
Scalar volume() const
Definition: GridInterfaceEuler.hpp:337
const GridInterface * pgrid_
Definition: GridInterfaceEuler.hpp:356
Vector centroid() const
Definition: GridInterfaceEuler.hpp:342
GIE::FaceIterator< GridInterface > FaceIterator
Definition: GridInterfaceEuler.hpp:321
Cell(const GridInterface &grid, const EntityType &cell)
Definition: GridInterfaceEuler.hpp:316
FaceIterator facebegin() const
Definition: GridInterfaceEuler.hpp:326
FaceIterator faceend() const
Definition: GridInterfaceEuler.hpp:331
EntityType cell_
Definition: GridInterfaceEuler.hpp:357
Index index() const
Definition: GridInterfaceEuler.hpp:351
FaceIterator::Index Index
Definition: GridInterfaceEuler.hpp:324
FaceIterator::Vector Vector
Definition: GridInterfaceEuler.hpp:322
Intersection (face) iterator for solver-near grid interface.
Definition: GridInterfaceEuler.hpp:241
bool equal(const FaceIterator &other) const
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:282
void increment()
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:290
Face< GridInterface >::DuneIntersectionIter DuneIntersectionIter
Type of low-level intersection iterator. Copied from the Dune grid.
Definition: GridInterfaceEuler.hpp:249
bool operator<(const FaceIterator &other) const
Gives an ordering of intersectionIterators.
Definition: GridInterfaceEuler.hpp:297
const FaceIterator & dereference() const
Used by iterator facade.
Definition: GridInterfaceEuler.hpp:276
FaceIterator(const GridInterface &grid, const DuneIntersectionIter &it, const int local_index)
Constructor.
Definition: GridInterfaceEuler.hpp:268
FaceIterator()
Default constructor.
Definition: GridInterfaceEuler.hpp:252
Definition: GridInterfaceEuler.hpp:88
Dune::FieldVector< Scalar, GI::GridType::dimension-1 > LocalVector
Definition: GridInterfaceEuler.hpp:94
Index local_index_
Definition: GridInterfaceEuler.hpp:219
const GI * pgrid_
Definition: GridInterfaceEuler.hpp:217
Vector normal() const
Definition: GridInterfaceEuler.hpp:133
Vector centroid() const
Definition: GridInterfaceEuler.hpp:124
GI::DuneIntersectionIterator DuneIntersectionIter
Definition: GridInterfaceEuler.hpp:90
Face(const GI &grid, const DuneIntersectionIter &it, const Index loc_ind)
Definition: GridInterfaceEuler.hpp:106
Index neighbourCellIndex() const
Definition: GridInterfaceEuler.hpp:182
@ LocalEndIndex
Definition: GridInterfaceEuler.hpp:99
@ BoundaryMarkerIndex
Definition: GridInterfaceEuler.hpp:98
DuneIntersectionIter iter_
Definition: GridInterfaceEuler.hpp:218
bool boundary() const
Definition: GridInterfaceEuler.hpp:142
int Index
Definition: GridInterfaceEuler.hpp:95
Cell cell() const
Definition: GridInterfaceEuler.hpp:158
GIE::Cell< GI, CellPtr > Cell
Definition: GridInterfaceEuler.hpp:96
Index cellIndex() const
Definition: GridInterfaceEuler.hpp:166
GI::GridType::ctype Scalar
Definition: GridInterfaceEuler.hpp:92
Cell neighbourCell() const
Definition: GridInterfaceEuler.hpp:174
int boundaryId() const
Definition: GridInterfaceEuler.hpp:150
Index localIndex() const
Definition: GridInterfaceEuler.hpp:203
Scalar neighbourCellVolume() const
Definition: GridInterfaceEuler.hpp:211
Face()
Definition: GridInterfaceEuler.hpp:101
GI::GridType::Traits::template Codim< 0 >::Entity CellPtr
Definition: GridInterfaceEuler.hpp:91
Dune::FieldVector< Scalar, GI::GridType::dimension > Vector
Definition: GridInterfaceEuler.hpp:93
Scalar area() const
Definition: GridInterfaceEuler.hpp:116
Index index() const
Definition: GridInterfaceEuler.hpp:195
Definition: GridInterfaceEuler.hpp:403
int numberOfFaces() const
Definition: GridInterfaceEuler.hpp:448
Opm::SparseTable< int >::row_type Indices
Definition: GridInterfaceEuler.hpp:476
void init(const DuneGrid &grid, bool build_facemap=true)
Definition: GridInterfaceEuler.hpp:428
Indices faceIndices(int cell_index) const
Definition: GridInterfaceEuler.hpp:477
const DuneGrid & grid() const
Definition: GridInterfaceEuler.hpp:458
const Mapper & mapper() const
Definition: GridInterfaceEuler.hpp:466
GridInterfaceEuler< DuneGrid > InterfaceType
Definition: GridInterfaceEuler.hpp:407
CellIterator::Index Index
Definition: GridInterfaceEuler.hpp:411
CellIterator cellend() const
Definition: GridInterfaceEuler.hpp:440
GIE::CellIterator< InterfaceType > CellIterator
Definition: GridInterfaceEuler.hpp:408
GridInterfaceEuler()
Definition: GridInterfaceEuler.hpp:417
CellIterator cellbegin() const
Definition: GridInterfaceEuler.hpp:436
GICellMapper< DuneGrid >::Type Mapper
Definition: GridInterfaceEuler.hpp:413
CellIterator::Scalar Scalar
Definition: GridInterfaceEuler.hpp:410
Index faceIndex(int cell_index, int local_face_index) const
Definition: GridInterfaceEuler.hpp:471
@ Dimension
Definition: GridInterfaceEuler.hpp:415
CellIterator::Vector Vector
Definition: GridInterfaceEuler.hpp:409
int numberOfCells() const
Definition: GridInterfaceEuler.hpp:444
GridInterfaceEuler(const DuneGrid &grid, bool build_facemap=true)
Definition: GridInterfaceEuler.hpp:421
DuneGrid GridType
Definition: GridInterfaceEuler.hpp:406
DuneGrid::LeafIntersectionIterator DuneIntersectionIterator
Definition: GridInterfaceEuler.hpp:405
int maxFacesPerCell() const
Definition: GridInterfaceEuler.hpp:453
Definition: ImplicitAssembly.hpp:43
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
Mapper for general grids.
Definition: GridInterfaceEuler.hpp:65
Dune::MultipleCodimMultipleGeomTypeMapper< Dune::GridView< Dune::DefaultLeafGridViewTraits< DuneGrid > > > Type
Definition: GridInterfaceEuler.hpp:66