24 #ifndef OPM_BLACKOILDETAILS_HEADER_INCLUDED
25 #define OPM_BLACKOILDETAILS_HEADER_INCLUDED
27 #include <opm/core/linalg/ParallelIstlInformation.hpp>
35 buildAllCells(
const int nc)
37 std::vector<int> all_cells(nc);
39 for (
int c = 0; c < nc; ++c) { all_cells[c] = c; }
48 activePhases(
const PU& pu)
50 const int maxnp = Opm::BlackoilPhases::MaxNumPhases;
51 std::vector<bool> active(maxnp,
false);
53 for (
int p = 0; p < pu.MaxNumPhases; ++p) {
54 active[ p ] = pu.phase_used[ p ] != 0;
64 active2Canonical(
const PU& pu)
66 const int maxnp = Opm::BlackoilPhases::MaxNumPhases;
67 std::vector<int> act2can(maxnp, -1);
69 for (
int phase = 0; phase < maxnp; ++phase) {
70 if (pu.phase_used[ phase ]) {
71 act2can[ pu.phase_pos[ phase ] ] = phase;
81 double getGravity(
const double* g,
const int dim) {
85 for (
int dd = 0; dd < dim - 1; ++dd) {
103 template <
class Iterator>
105 double euclidianNormSquared( Iterator it,
const Iterator end,
int num_components,
const boost::any& pinfo = boost::any() )
107 static_cast<void>(num_components);
108 static_cast<void>(pinfo);
110 if ( pinfo.type() ==
typeid(ParallelISTLInformation) )
112 const ParallelISTLInformation& info =
113 boost::any_cast<
const ParallelISTLInformation&>(pinfo);
114 typedef typename Iterator::value_type Scalar;
115 Scalar product = 0.0;
116 int size_per_component = (end - it);
117 size_per_component /= num_components;
118 assert((end - it) == num_components * size_per_component);
120 if( num_components == 1 )
122 auto component_container =
123 boost::make_iterator_range(it, end);
124 info.computeReduction(component_container,
125 Opm::Reduction::makeInnerProductFunctor<double>(),
130 auto& maskContainer = info.getOwnerMask();
131 auto mask = maskContainer.begin();
132 assert(static_cast<int>(maskContainer.size()) == size_per_component);
134 for(
int cell = 0; cell < size_per_component; ++cell, ++mask)
136 Scalar cell_product = (*it) * (*it);
138 for(
int component=1; component < num_components;
141 cell_product += (*it) * (*it);
143 product += cell_product * (*mask);
146 return info.communicator().sum(product);
151 double product = 0.0 ;
152 for( ; it != end; ++it ) {
153 product += ( *it * *it );
164 std::size_t countLocalInteriorCells(
const Grid& grid)
166 if ( grid.comm().size() == 1)
170 std::size_t count = 0;
171 const auto& gridView = grid.leafGridView();
172 for(
auto cell = gridView.template begin<0, Dune::Interior_Partition>(),
173 endCell = gridView.template end<0, Dune::Interior_Partition>();
174 cell != endCell; ++cell)
189 std::size_t countGlobalCells(
const Grid& grid)
191 if ( grid.comm().size() == 1)
195 std::size_t count = countLocalInteriorCells(grid);
196 return grid.comm().sum(count);
203 #endif // OPM_BLACKOILDETAILS_HEADER_INCLUDED