00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024 #ifndef OPM_BLACKOILDETAILS_HEADER_INCLUDED
00025 #define OPM_BLACKOILDETAILS_HEADER_INCLUDED
00026
00027 #include <opm/core/linalg/ParallelIstlInformation.hpp>
00028
00029 namespace Opm {
00030 namespace detail {
00031
00032
00033 inline
00034 std::vector<int>
00035 buildAllCells(const int nc)
00036 {
00037 std::vector<int> all_cells(nc);
00038
00039 for (int c = 0; c < nc; ++c) { all_cells[c] = c; }
00040
00041 return all_cells;
00042 }
00043
00044
00045
00046 template <class PU>
00047 std::vector<bool>
00048 activePhases(const PU& pu)
00049 {
00050 const int maxnp = Opm::BlackoilPhases::MaxNumPhases;
00051 std::vector<bool> active(maxnp, false);
00052
00053 for (int p = 0; p < pu.MaxNumPhases; ++p) {
00054 active[ p ] = pu.phase_used[ p ] != 0;
00055 }
00056
00057 return active;
00058 }
00059
00060
00061
00062 template <class PU>
00063 std::vector<int>
00064 active2Canonical(const PU& pu)
00065 {
00066 const int maxnp = Opm::BlackoilPhases::MaxNumPhases;
00067 std::vector<int> act2can(maxnp, -1);
00068
00069 for (int phase = 0; phase < maxnp; ++phase) {
00070 if (pu.phase_used[ phase ]) {
00071 act2can[ pu.phase_pos[ phase ] ] = phase;
00072 }
00073 }
00074
00075 return act2can;
00076 }
00077
00078
00079
00080 inline
00081 double getGravity(const double* g, const int dim) {
00082 double grav = 0.0;
00083 if (g) {
00084
00085 for (int dd = 0; dd < dim - 1; ++dd) {
00086 assert(g[dd] == 0.0);
00087 }
00088 grav = g[dim - 1];
00089 }
00090 return grav;
00091 }
00092
00093
00103 template <class Iterator>
00104 inline
00105 double euclidianNormSquared( Iterator it, const Iterator end, int num_components, const boost::any& pinfo = boost::any() )
00106 {
00107 static_cast<void>(num_components);
00108 static_cast<void>(pinfo);
00109 #if HAVE_MPI
00110 if ( pinfo.type() == typeid(ParallelISTLInformation) )
00111 {
00112 const ParallelISTLInformation& info =
00113 boost::any_cast<const ParallelISTLInformation&>(pinfo);
00114 typedef typename Iterator::value_type Scalar;
00115 Scalar product = 0.0;
00116 int size_per_component = (end - it);
00117 size_per_component /= num_components;
00118 assert((end - it) == num_components * size_per_component);
00119
00120 if( num_components == 1 )
00121 {
00122 auto component_container =
00123 boost::make_iterator_range(it, end);
00124 info.computeReduction(component_container,
00125 Opm::Reduction::makeInnerProductFunctor<double>(),
00126 product);
00127 }
00128 else
00129 {
00130 auto& maskContainer = info.getOwnerMask();
00131 auto mask = maskContainer.begin();
00132 assert(static_cast<int>(maskContainer.size()) == size_per_component);
00133
00134 for(int cell = 0; cell < size_per_component; ++cell, ++mask)
00135 {
00136 Scalar cell_product = (*it) * (*it);
00137 ++it;
00138 for(int component=1; component < num_components;
00139 ++component, ++it)
00140 {
00141 cell_product += (*it) * (*it);
00142 }
00143 product += cell_product * (*mask);
00144 }
00145 }
00146 return info.communicator().sum(product);
00147 }
00148 else
00149 #endif
00150 {
00151 double product = 0.0 ;
00152 for( ; it != end; ++it ) {
00153 product += ( *it * *it );
00154 }
00155 return product;
00156 }
00157 }
00163 template<class Grid>
00164 std::size_t countLocalInteriorCells(const Grid& grid)
00165 {
00166 if ( grid.comm().size() == 1)
00167 {
00168 return grid.size(0);
00169 }
00170 std::size_t count = 0;
00171 const auto& gridView = grid.leafGridView();
00172 for(auto cell = gridView.template begin<0, Dune::Interior_Partition>(),
00173 endCell = gridView.template end<0, Dune::Interior_Partition>();
00174 cell != endCell; ++cell)
00175 {
00176 ++count;
00177 }
00178 return count;
00179 }
00180
00188 template<class Grid>
00189 std::size_t countGlobalCells(const Grid& grid)
00190 {
00191 if ( grid.comm().size() == 1)
00192 {
00193 return grid.size(0);
00194 }
00195 std::size_t count = countLocalInteriorCells(grid);
00196 return grid.comm().sum(count);
00197 }
00198
00199
00200 }
00201 }
00202
00203 #endif // OPM_BLACKOILDETAILS_HEADER_INCLUDED