20#ifndef OPM_INITSTATE_IMPL_HEADER_INCLUDED
21#define OPM_INITSTATE_IMPL_HEADER_INCLUDED
23#include <opm/common/data/SimulationDataContainer.hpp>
24#include <opm/common/ErrorMacros.hpp>
26#include <opm/core/utility/parameters/ParameterGroup.hpp>
27#include <opm/core/grid.h>
28#include <opm/core/grid/GridHelpers.hpp>
29#include <opm/core/utility/MonotCubicInterpolator.hpp>
30#include <opm/parser/eclipse/Units/Units.hpp>
36#include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
37#include <opm/parser/eclipse/EclipseState/InitConfig/Equil.hpp>
46 template <
class Props>
47 static void initSaturation(
const std::vector<int>& cells ,
const Props& props , SimulationDataContainer& state ,
ExtremalSat satType) {
48 const int num_phases = state.numPhases();
49 std::vector<double> min_sat(num_phases * cells.size());
50 std::vector<double> max_sat(num_phases * cells.size());
51 props.satRange(cells.size() ,cells.data() , min_sat.data() , max_sat.data());
54 std::vector<double> second_sat(cells.size());
55 std::vector<double> first_sat(cells.size());
57 for (
size_t index = 0; index < cells.size(); index++) {
59 first_sat[index] = min_sat[ num_phases * index];
60 second_sat[index] = 1 - min_sat[ num_phases * index];
62 first_sat[index] = max_sat[ num_phases * index];
63 second_sat[index] = 1 - max_sat[ num_phases * index];
67 state.setCellDataComponent(
"SATURATION" , 0 , cells , first_sat );
68 state.setCellDataComponent(
"SATURATION" , 1 , cells , second_sat );
80#pragma clang diagnostic push
81#pragma clang diagnostic ignored "-Wunneeded-internal-declaration"
87 void cellsBelowAbove(
int number_of_cells,
88 T begin_cell_centroids,
91 std::vector<int>& below,
92 std::vector<int>& above)
94 below.reserve(number_of_cells);
95 above.reserve(number_of_cells);
96 for (
int c = 0; c < number_of_cells; ++c) {
97 const double z = UgGridHelpers
98 ::getCoordinate(UgGridHelpers::increment(begin_cell_centroids, c,
109#pragma clang diagnostic pop
113 enum WaterInit { WaterBelow, WaterAbove };
122 template <
class Props,
class State>
124 std::vector<double> min_sat(state.numPhases() * cells.size());
125 std::vector<double> max_sat(state.numPhases() * cells.size());
126 props.satRange(cells.size() ,cells.data() , min_sat.data() , max_sat.data());
129 std::vector<double> first_sat(cells.size());
130 std::vector<double> second_sat(cells.size());
132 for (
size_t index=0; index < cells.size(); index++) {
134 first_sat[index] = min_sat[index * state.numPhases()];
135 second_sat[index] = 1 - min_sat[index * state.numPhases()];
137 first_sat[index] = max_sat[index * state.numPhases()];
138 second_sat[index] = 1 - max_sat[index * state.numPhases()];
141 state.setCellDataComponent(
"SATURATION" , 0 , cells , first_sat );
142 state.setCellDataComponent(
"SATURATION" , 1 , cells , second_sat );
150 template <
class T,
class Props,
class State>
151 void initWaterOilContact(
int number_of_cells,
152 T begin_cell_centroids,
156 const WaterInit waterinit,
160 std::vector<int>
water;
161 std::vector<int>
oil;
169 cellsBelowAbove(number_of_cells, begin_cell_centroids, dimensions,
173 cellsBelowAbove(number_of_cells, begin_cell_centroids, dimensions,
193 template <
class T,
class State>
194 void initHydrostaticPressure(
int number_of_cells,
195 T begin_cell_centroids,
197 const double* densities,
199 const double gravity,
200 const double datum_z,
201 const double datum_p,
204 std::vector<double>& p = state.pressure();
206 const double rho_datum = datum_z > woc ? densities[0] : densities[1];
207 const double woc_p = datum_p + (woc - datum_z)*gravity*rho_datum;
208 for (
int c = 0; c < number_of_cells; ++c) {
210 const double z = UgGridHelpers::
211 getCoordinate(UgGridHelpers::increment(begin_cell_centroids, c, dimensions),
213 const double rho = z > woc ? densities[0] : densities[1];
214 p[c] = woc_p + (z - woc)*gravity*rho;
221 template <
class T,
class State>
222 void initHydrostaticPressure(
int number_of_cells,
223 T begin_cell_centroids,
225 const IncompPropertiesInterface& props,
227 const double gravity,
228 const double datum_z,
229 const double datum_p,
232 const double* densities = props.density();
233 initHydrostaticPressure(number_of_cells, begin_cell_centroids, dimensions,
234 densities, woc, gravity, datum_z, datum_p, state);
242 const BlackoilPropertiesInterface&
props_;
243 Density(
const BlackoilPropertiesInterface& props) :
props_(props) {}
244 double operator()(
const double pressure,
const double temperature,
const int phase)
246 assert(
props_.numPhases() == 2);
247 const double surfvol[2][2] = { { 1.0, 0.0 },
251 double A[4] = { 0.0 };
253 double rho[2] = { 0.0 };
254 props_.density(1, A, &cell, rho);
266 template <
class T,
class State>
267 void initHydrostaticPressure(
int number_of_cells,
268 T begin_cell_centroids,
270 const BlackoilPropertiesInterface& props,
272 const double gravity,
273 const double datum_z,
274 const double datum_p,
277 assert(props.numPhases() == 2);
280 double zlim[2] = { 1e100, -1e100 };
281 for (
int c = 0; c < number_of_cells; ++c) {
282 const double z = UgGridHelpers::
283 getCoordinate(UgGridHelpers::increment(begin_cell_centroids, c, dimensions),
285 zlim[0] = std::min(zlim[0], z);
286 zlim[1] = std::max(zlim[1], z);
290 const double hmin = (zlim[1] - zlim[0])/100.0;
293 std::map<double, double> press_by_z;
294 press_by_z[datum_z] = datum_p;
300 int phase = (datum_z > woc) ? 0 : 1;
301 int num_steps = int(std::ceil(std::fabs(woc - datum_z)/hmin));
302 double pval = datum_p;
303 double Tval = 273.15 + 20;
304 double zval = datum_z;
305 double h = (woc - datum_z)/
double(num_steps);
306 for (
int i = 0; i < num_steps; ++i) {
308 const double dp = rho(pval, Tval, phase)*gravity;
310 press_by_z[zval] = pval;
315 double z_end = (datum_z > woc) ? zlim[1] : zlim[0];
316 num_steps = int(std::ceil(std::fabs(z_end - datum_z)/hmin));
319 h = (z_end - datum_z)/
double(num_steps);
320 for (
int i = 0; i < num_steps; ++i) {
322 const double dp = rho(pval, Tval, phase)*gravity;
324 press_by_z[zval] = pval;
329 phase = (phase + 1) % 2;
330 z_end = (datum_z > woc) ? zlim[0] : zlim[1];
333 h = (z_end - datum_z)/
double(num_steps);
334 for (
int i = 0; i < num_steps; ++i) {
336 const double dp = rho(pval, Tval, phase)*gravity;
338 press_by_z[zval] = pval;
342 std::vector<double> zv, pv;
343 zv.reserve(press_by_z.size());
344 pv.reserve(press_by_z.size());
345 for (std::map<double, double>::const_iterator it = press_by_z.begin();
346 it != press_by_z.end(); ++it) {
347 zv.push_back(it->first);
348 pv.push_back(it->second);
350 MonotCubicInterpolator press(zv, pv);
353 std::vector<double>& p = state.pressure();
354 for (
int c = 0; c < number_of_cells; ++c) {
355 const double z = UgGridHelpers::
356 getCoordinate(UgGridHelpers::increment(begin_cell_centroids, c, dimensions),
363 template <
class State,
class FaceCells,
class FCI,
class CCI>
364 void initFacePressure(
int dimensions,
366 FaceCells face_cells,
367 FCI begin_face_centroids,
368 CCI begin_cell_centroids,
371 const std::vector<double>& cp = state.pressure();
372 std::vector<double>& fp = state.facepressure();
373 for (
int f = 0; f < number_of_faces; ++f) {
374 double dist[2] = { 0.0, 0.0 };
375 double press[2] = { 0.0, 0.0 };
377 for (
int j = 0; j < 2; ++j) {
378 const int c = face_cells(f, j);
381 for (
int dd = 0; dd < dimensions; ++dd) {
382 double diff = UgGridHelpers::
383 getCoordinate(UgGridHelpers::increment(begin_face_centroids, f,
387 getCoordinate(UgGridHelpers::increment(begin_cell_centroids, c,
390 dist[j] += diff*diff;
392 dist[j] = std::sqrt(dist[j]);
399 fp[f] = press[0]*(dist[1]/(dist[0] + dist[1])) + press[1]*(dist[0]/(dist[0] + dist[1]));
401 fp[f] = press[(bdy_idx + 1) % 2];
411 template <
class State>
414 const ParameterGroup& param,
415 const double gravity,
418 initStateBasic(grid.number_of_cells, grid.global_cell, grid.cartdims,
419 grid.number_of_faces, UgGridHelpers::faceCells(grid),
420 grid.face_centroids, grid.cell_centroids,
421 grid.dimensions, props, param,
425 template <
class FaceCells,
class CCI,
class FCI,
class State>
427 const int* global_cell,
430 FaceCells face_cells,
431 FCI begin_face_centroids,
432 CCI begin_cell_centroids,
435 const ParameterGroup& param,
436 const double gravity,
439 const int num_phases = props.
numPhases();
440 if (num_phases != 2) {
441 OPM_THROW(std::runtime_error,
"initStateTwophaseBasic(): currently handling only two-phase scenarios.");
443 const int num_cells = props.
numCells();
445 std::vector<int> all_cells(num_cells);
446 for (
int i = 0; i < num_cells; ++i) {
453 const bool convection_testcase = param.getDefault(
"convection_testcase",
false);
454 const bool segregation_testcase = param.getDefault(
"segregation_testcase",
false);
455 if (convection_testcase) {
457 std::vector<int> left_cells;
458 left_cells.reserve(num_cells/2);
459 for (
int cell = 0; cell < num_cells; ++cell) {
460 const int gc = global_cell == 0 ? cell : global_cell[cell];
461 bool left = (gc % cartdims[0]) < cartdims[0]/2;
463 left_cells.push_back(cell);
468 const double init_p = param.getDefault(
"ref_pressure", 100.0)*unit::barsa;
469 std::fill(state.pressure().begin(), state.pressure().end(), init_p);
470 }
else if (segregation_testcase) {
472 if (gravity == 0.0) {
473 std::cout <<
"**** Warning: running gravity segregation scenario, but gravity is zero." << std::endl;
475 if (cartdims[2] <= 1) {
476 std::cout <<
"**** Warning: running gravity segregation scenario, which expects nz > 1." << std::endl;
479 const double woc = param.get<
double>(
"water_oil_contact");
480 initWaterOilContact(number_of_cells, begin_cell_centroids, dimensions,
481 props, woc, WaterAbove, state);
483 const double ref_p = param.getDefault(
"ref_pressure", 100.0)*unit::barsa;
485 initHydrostaticPressure(number_of_cells, begin_cell_centroids, dimensions,
486 dens, woc, gravity, woc, ref_p, state);
487 }
else if (param.has(
"water_oil_contact")) {
489 if (gravity == 0.0) {
490 std::cout <<
"**** Warning: running gravity convection scenario, but gravity is zero." << std::endl;
492 if (cartdims[2] <= 1) {
493 std::cout <<
"**** Warning: running gravity convection scenario, which expects nz > 1." << std::endl;
496 const double woc = param.get<
double>(
"water_oil_contact");
497 initWaterOilContact(number_of_cells, begin_cell_centroids, dimensions,
498 props, woc, WaterBelow, state);
500 const double ref_p = param.getDefault(
"ref_pressure", 100.0)*unit::barsa;
501 initHydrostaticPressure(number_of_cells, begin_cell_centroids, dimensions,
502 props.
density(), woc, gravity, woc, ref_p, state);
503 }
else if (param.has(
"init_saturation")) {
505 const double init_saturation = param.get<
double>(
"init_saturation");
506 for (
int cell = 0; cell < num_cells; ++cell) {
507 state.saturation()[2*cell] = init_saturation;
508 state.saturation()[2*cell + 1] = 1.0 - init_saturation;
511 const double ref_p = param.getDefault(
"ref_pressure", 100.0)*unit::barsa;
512 const double rho = props.
density()[0]*init_saturation + props.
density()[1]*(1.0 - init_saturation);
513 const double dens[2] = { rho, rho };
514 const double ref_z = UgGridHelpers::getCoordinate(begin_cell_centroids, dimensions - 1);
515 initHydrostaticPressure(number_of_cells, begin_cell_centroids, dimensions,
516 dens, ref_z, gravity, ref_z, ref_p, state);
520 const double ref_p = param.getDefault(
"ref_pressure", 100.0)*unit::barsa;
521 const double rho = props.
density()[1];
522 const double dens[2] = { rho, rho };
523 const double ref_z = UgGridHelpers::getCoordinate(begin_cell_centroids, dimensions - 1);
524 initHydrostaticPressure(number_of_cells, begin_cell_centroids, dimensions,
525 dens, ref_z, gravity, ref_z, ref_p, state);
529 initFacePressure(dimensions, number_of_faces, face_cells, begin_face_centroids,
530 begin_cell_centroids, state);
535 template <
class State>
538 const ParameterGroup& param,
539 const double gravity,
542 initStateBasic(grid.number_of_cells, grid.global_cell, grid.cartdims,
543 grid.number_of_faces, UgGridHelpers::faceCells(grid),
544 grid.face_centroids, grid.cell_centroids, grid.dimensions,
545 props, param, gravity, state);
548 template <
class FaceCells,
class FCI,
class CCI,
class State>
550 const int* global_cell,
553 FaceCells face_cells,
554 FCI begin_face_centroids,
555 CCI begin_cell_centroids,
558 const ParameterGroup& param,
559 const double gravity,
563 const int num_phases = props.
numPhases();
564 if (num_phases != 2) {
565 OPM_THROW(std::runtime_error,
"initStateTwophaseBasic(): currently handling only two-phase scenarios.");
567 const int num_cells = props.
numCells();
569 std::vector<int> all_cells(num_cells);
570 for (
int i = 0; i < num_cells; ++i) {
574 const bool convection_testcase = param.getDefault(
"convection_testcase",
false);
575 if (convection_testcase) {
577 std::vector<int> left_cells;
578 left_cells.reserve(num_cells/2);
579 for (
int cell = 0; cell < num_cells; ++cell) {
580 const int gc = global_cell == 0 ? cell : global_cell[cell];
581 bool left = (gc % cartdims[0]) < cartdims[0]/2;
583 left_cells.push_back(cell);
587 const double init_p = param.getDefault(
"ref_pressure", 100.0)*unit::barsa;
588 std::fill(state.pressure().begin(), state.pressure().end(), init_p);
589 }
else if (param.has(
"water_oil_contact")) {
591 if (gravity == 0.0) {
592 std::cout <<
"**** Warning: running gravity convection scenario, but gravity is zero." << std::endl;
594 if (cartdims[2] <= 1) {
595 std::cout <<
"**** Warning: running gravity convection scenario, which expects nz > 1." << std::endl;
598 const double woc = param.get<
double>(
"water_oil_contact");
599 initWaterOilContact(number_of_cells, begin_cell_centroids, dimensions,
600 props, woc, WaterBelow, state);
602 const double ref_p = param.getDefault(
"ref_pressure", 100.0)*unit::barsa;
603 initHydrostaticPressure(number_of_cells, begin_cell_centroids, dimensions,
604 props, woc, gravity, woc, ref_p, state);
608 const double ref_p = param.getDefault(
"ref_pressure", 100.0)*unit::barsa;
609 const double ref_z = UgGridHelpers::getCoordinate(begin_cell_centroids, dimensions - 1);
610 const double woc = -1e100;
611 initHydrostaticPressure(number_of_cells, begin_cell_centroids, dimensions,
612 props, woc, gravity, ref_z, ref_p, state);
616 initFacePressure(dimensions, number_of_faces, face_cells, begin_face_centroids,
617 begin_cell_centroids, state);
622 template <
class Props,
class State>
625 const Opm::EclipseState& es,
626 const double gravity,
630 grid.number_of_faces, UgGridHelpers::faceCells(grid),
631 grid.face_centroids, grid.cell_centroids, grid.dimensions,
632 props, es, gravity, state);
635 template <
class FaceCells,
class FCI,
class CCI,
class Props,
class State>
637 const int* global_cell,
639 FaceCells face_cells,
640 FCI begin_face_centroids,
641 CCI begin_cell_centroids,
644 const Opm::EclipseState& es,
645 const double gravity,
648 const int num_phases = props.numPhases();
651 const auto& init_config = es.cfg().init();
652 const auto& grid_props = es.get3DProperties();
653 const auto has_equil = init_config.hasEquil();
654 const auto has_pressure = grid_props.hasDeckDoubleGridProperty(
"PRESSURE");
655 const auto has_sgas = grid_props.hasDeckDoubleGridProperty(
"SGAS");
656 const auto has_swat = grid_props.hasDeckDoubleGridProperty(
"SWAT");
659 OPM_THROW(std::runtime_error,
"initStateFromDeck(): user specified property object with " << num_phases <<
" phases, "
660 "found " << pu.
num_phases <<
" phases in deck.");
662 if (has_equil && has_pressure) {
663 OPM_THROW(std::runtime_error,
"initStateFromDeck(): The deck must either specify the initial "
664 "condition using the PRESSURE _or_ the EQUIL keyword (currently it has both)");
668 if (num_phases != 2) {
669 OPM_THROW(std::runtime_error,
"initStateFromDeck(): EQUIL-based init currently handling only two-phase scenarios.");
672 OPM_THROW(std::runtime_error,
"initStateFromDeck(): EQUIL-based init currently handling only oil-water scenario (no gas).");
675 const auto& equil = init_config.getEquil();
676 if (equil.size() != 1) {
677 OPM_THROW(std::runtime_error,
"initStateFromDeck(): No region support yet.");
679 const double woc = equil.getRecord( 0 ).waterOilContactDepth();
680 initWaterOilContact(number_of_cells, begin_cell_centroids, dimensions,
681 props, woc, WaterBelow, state);
683 const double datum_z = equil.getRecord( 0 ).datumDepth();
684 const double datum_p = equil.getRecord( 0 ).datumDepthPressure();
685 initHydrostaticPressure(number_of_cells, begin_cell_centroids, dimensions,
686 props, woc, gravity, datum_z, datum_p, state);
687 }
else if (has_pressure) {
689 std::vector<double>& s = state.saturation();
690 std::vector<double>& p = state.pressure();
691 const auto& p_deck = grid_props.getDoubleGridProperty(
"PRESSURE").getData();
692 const int num_cells = number_of_cells;
694 if (num_phases == 2) {
698 OPM_THROW(std::runtime_error,
"initStateFromDeck(): missing SGAS keyword in 2-phase init");
700 const auto& sg_deck = grid_props.getDoubleGridProperty(
"SGAS").getData();
703 for (
int c = 0; c < num_cells; ++c) {
704 int c_deck = (global_cell == NULL) ? c : global_cell[c];
705 s[2*c + gpos] = sg_deck[c_deck];
706 s[2*c + opos] = 1.0 - sg_deck[c_deck];
707 p[c] = p_deck[c_deck];
712 OPM_THROW(std::runtime_error,
"initStateFromDeck(): missing SWAT keyword in 2-phase init");
714 const auto& sw_deck = grid_props.getDoubleGridProperty(
"SWAT").getData();
716 const int nwpos = (wpos + 1) % 2;
717 for (
int c = 0; c < num_cells; ++c) {
718 int c_deck = (global_cell == NULL) ? c : global_cell[c];
719 s[2*c + wpos] = sw_deck[c_deck];
720 s[2*c + nwpos] = 1.0 - sw_deck[c_deck];
721 p[c] = p_deck[c_deck];
724 }
else if (num_phases == 3) {
725 const bool has_swat_sgas = has_swat && has_sgas;
726 if (!has_swat_sgas) {
727 OPM_THROW(std::runtime_error,
"initStateFromDeck(): missing SGAS or SWAT keyword in 3-phase init.");
732 const auto& sw_deck = grid_props.getDoubleGridProperty(
"SWAT").getData();
733 const auto& sg_deck = grid_props.getDoubleGridProperty(
"SGAS").getData();
734 for (
int c = 0; c < num_cells; ++c) {
735 int c_deck = (global_cell == NULL) ? c : global_cell[c];
736 s[3*c + wpos] = sw_deck[c_deck];
737 s[3*c + opos] = 1.0 - (sw_deck[c_deck] + sg_deck[c_deck]);
738 s[3*c + gpos] = sg_deck[c_deck];
739 p[c] = p_deck[c_deck];
742 OPM_THROW(std::runtime_error,
"initStateFromDeck(): init with SWAT etc. only available with 2 or 3 phases.");
745 OPM_THROW(std::runtime_error,
"initStateFromDeck(): we must either have EQUIL, or PRESSURE and SWAT/SOIL/SGAS.");
749 initFacePressure(dimensions, number_of_faces, face_cells, begin_face_centroids,
750 begin_cell_centroids, state);
756 template <
class Props,
class State>
764 template <
class Props,
class State>
769 state.surfacevol() = state.saturation();
770 const int np = props.numPhases();
771 const int nc = number_of_cells;
772 std::vector<double> allA(nc*np*np);
773 std::vector<int> allcells(nc);
774 for (
int c = 0; c < nc; ++c) {
780 props.matrix(nc, &state.pressure()[0], &state.temperature()[0], &state.surfacevol()[0], &allcells[0], &allA[0], 0);
781 for (
int c = 0; c < nc; ++c) {
783 double* z = &state.surfacevol()[c*np];
784 const double* s = &state.saturation()[c*np];
785 const double* A = &allA[c*np*np];
787 for (
int row = 0; row < np; ++row) { z[row] = 0.0; }
789 for (
int col = 0; col < np; ++col) {
790 for (
int row = 0; row < np; ++row) {
792 z[row] += A[row + np*col]*s[col];
800 template <
class Props,
class State>
808 template <
class Props,
class State>
813 const std::vector<double>& rs = state.gasoilratio();
814 const std::vector<double>& rv = state.rv();
817 state.surfacevol() = state.saturation();
820 const int np = props.numPhases();
822 std::vector<double> allA_a(number_of_cells*np*np);
823 std::vector<double> allA_l(number_of_cells*np*np);
824 std::vector<double> allA_v(number_of_cells*np*np);
826 std::vector<int> allcells(number_of_cells);
827 std::vector<double> z_init(number_of_cells*np);
828 std::fill(z_init.begin(),z_init.end(),0.0);
830 for (
int c = 0; c < number_of_cells; ++c) {
834 std::vector<double> capPressures(number_of_cells*np);
835 props.capPress(number_of_cells,&state.saturation()[0],&allcells[0],&capPressures[0],NULL);
841 for (
int c = 0; c < number_of_cells ; ++c){
842 for (
int p = 0; p < np ; ++p){
848 z_init[c*np + p] = z_tmp;
851 props.matrix(number_of_cells, &state.pressure()[0], &state.temperature()[0], &z_init[0], &allcells[0], &allA_a[0], 0);
855 for (
int c = 0; c < number_of_cells ; ++c){
856 for (
int p = 0; p < np ; ++p){
858 if(state.saturation()[np*c + p] > 0)
868 z_init[c*np + p] = z_tmp;
873 props.matrix(number_of_cells, &state.pressure()[0], &state.temperature()[0], &z_init[0], &allcells[0], &allA_l[0], 0);
876 for (
int c = 0; c < number_of_cells ; ++c){
877 for (
int p = 0; p < np ; ++p){
879 if(state.saturation()[np*c + p] > 0)
889 z_init[c*np + p] = z_tmp;
894 props.matrix(number_of_cells, &state.pressure()[0], &state.temperature()[0], &z_init[0], &allcells[0], &allA_v[0], 0);
896 for (
int c = 0; c < number_of_cells; ++c) {
898 double* z = &state.surfacevol()[c*np];
899 const double* s = &state.saturation()[c*np];
900 const double* A_a = &allA_a[c*np*np];
901 const double* A_l = &allA_l[c*np*np];
902 const double* A_v = &allA_v[c*np*np];
904 for (
int row = 0; row < np; ++row) { z[row] = 0.0; }
906 for (
int col = 0; col < np; ++col) {
907 z[0] += A_a[0 + np*col]*s[col];
908 z[1] += A_l[1 + np*col]*s[col];
910 z[2] += A_v[2 + np*col]*s[col];
923 template <
class Props,
class State>
926 const Opm::EclipseState& es,
927 const double gravity,
931 grid.number_of_faces, UgGridHelpers::faceCells(grid),
932 grid.face_centroids, grid.cell_centroids,grid.dimensions,
933 props, es, gravity, state);
937 template <
class FaceCells,
class FCI,
class CCI,
class Props,
class State>
939 const int* global_cell,
941 FaceCells face_cells,
942 FCI begin_face_centroids,
943 CCI begin_cell_centroids,
946 const Opm::EclipseState& es,
947 const double gravity,
951 face_cells, begin_face_centroids, begin_cell_centroids,
952 dimensions, props, es, gravity, state);
954 const auto& grid_props = es.get3DProperties();
956 if (grid_props.hasDeckDoubleGridProperty(
"RS")) {
957 const auto& rs_deck = grid_props.getDoubleGridProperty(
"RS").getData();
958 const int num_cells = number_of_cells;
959 for (
int c = 0; c < num_cells; ++c) {
960 int c_deck = (global_cell == NULL) ? c : global_cell[c];
961 state.gasoilratio()[c] = rs_deck[c_deck];
965 }
else if (grid_props.hasDeckDoubleGridProperty(
"RV")) {
966 const auto& rv_deck = grid_props.getDoubleGridProperty(
"RV").getData();
967 const int num_cells = number_of_cells;
968 for (
int c = 0; c < num_cells; ++c) {
969 int c_deck = (global_cell == NULL) ? c : global_cell[c];
970 state.rv()[c] = rv_deck[c_deck];
976 state.gasoilratio() = std::vector<double>(number_of_cells, 0.0);
977 state.rv() = std::vector<double>(number_of_cells, 0.0);
@ Liquid
Definition: BlackoilPhases.hpp:40
@ Aqua
Definition: BlackoilPhases.hpp:40
@ Vapour
Definition: BlackoilPhases.hpp:40
Definition: BlackoilPropertiesInterface.hpp:38
virtual int numCells() const =0
virtual int numPhases() const =0
Definition: IncompPropertiesInterface.hpp:36
virtual int numPhases() const =0
virtual const double * density() const =0
virtual int numCells() const =0
const BlackoilPropertiesInterface & props_
Definition: initState_impl.hpp:242
void water(const Grid &G, const Region ®, const std::array< double, 2 > &span, const double grav, double &po_woc, const CellRange &cells, std::vector< double > &press)
Definition: initStateEquil_impl.hpp:270
void oil(const Grid &G, const Region ®, const std::array< double, 2 > &span, const double grav, const CellRange &cells, std::vector< double > &press, double &po_woc, double &po_goc)
Definition: initStateEquil_impl.hpp:319
std::vector< double > temperature(const Grid &, const Region &, const CellRange &cells)
Definition: initStateEquil_impl.hpp:583
Definition: AnisotropicEikonal.hpp:44
void initBlackoilSurfvolUsingRSorRV(const UnstructuredGrid &grid, const Props &props, State &state)
Definition: initState_impl.hpp:801
void initBlackoilStateFromDeck(const UnstructuredGrid &grid, const Props &props, const Opm::EclipseState &es, const double gravity, State &state)
Initialize a blackoil state from input deck.
Definition: initState_impl.hpp:924
static void initSaturation(const std::vector< int > &cells, const Props &props, SimulationDataContainer &state, ExtremalSat satType)
void initStateBasic(const UnstructuredGrid &grid, const IncompPropertiesInterface &props, const ParameterGroup ¶m, const double gravity, State &state)
Initialize a twophase state from parameters.
Definition: initState_impl.hpp:412
ExtremalSat
Definition: initState.hpp:46
@ MaxSat
Definition: initState.hpp:46
@ MinSat
Definition: initState.hpp:46
void computeSaturation(const BlackoilPropertiesInterface &props, BlackoilState &state)
Computes saturation from surface volume densities.
void initStateFromDeck(const UnstructuredGrid &grid, const Props &props, const EclipseState &es, const double gravity, State &state)
void initBlackoilSurfvol(const UnstructuredGrid &grid, const Props &props, State &state)
Definition: initState_impl.hpp:757
PhaseUsage phaseUsageFromDeck(const Opm::EclipseState &eclipseState)
Definition: phaseUsageFromDeck.hpp:37
Definition: BlackoilPhases.hpp:44
int phase_pos[MaxNumPhases+NumCryptoPhases]
Definition: BlackoilPhases.hpp:47
int phase_used[MaxNumPhases+NumCryptoPhases]
Definition: BlackoilPhases.hpp:46
int num_phases
Definition: BlackoilPhases.hpp:45