36 #ifndef OPM_SETUPGRIDANDPROPS_HEADER
37 #define OPM_SETUPGRIDANDPROPS_HEADER
39 #include <opm/common/utility/parameters/ParameterGroup.hpp>
40 #include <opm/input/eclipse/Units/Units.hpp>
41 #include <opm/grid/CpGrid.hpp>
42 #include <opm/porsol/common/ReservoirPropertyCapillary.hpp>
43 #include <opm/input/eclipse/Parser/Parser.hpp>
44 #include <opm/input/eclipse/Deck/Deck.hpp>
46 #include <opm/common/utility/platform_dependent/disable_warnings.h>
48 #include <opm/common/utility/FileSystem.hpp>
50 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
63 bool useJ< ReservoirPropertyCapillary<3> >();
68 template <
template <
int>
class ResProp>
75 std::string fileformat = param.getDefault<std::string>(
"fileformat",
"cartesian");
76 if (fileformat ==
"sintef_legacy") {
77 std::string grid_prefix = param.get<std::string>(
"grid_prefix");
78 grid.readSintefLegacyFormat(grid_prefix);
79 OPM_MESSAGE(
"Warning: We do not yet read legacy reservoir properties. Using defaults.");
80 res_prop.init(grid.size(0));
81 }
else if (fileformat ==
"eclipse") {
82 std::string ecl_file = param.get<std::string>(
"filename");
85 auto deck = parser.parseFile(ecl_file);
86 if (param.has(
"z_tolerance")) {
87 std::cerr <<
"****** Warning: z_tolerance parameter is obsolete, use PINCH in deck input instead\n";
89 bool periodic_extension = param.getDefault<
bool>(
"periodic_extension",
false);
90 bool clip_z = param.getDefault<
bool>(
"clip_z",
false);
91 bool turn_normals = param.getDefault<
bool>(
"turn_normals",
false);
93 Opm::EclipseGrid inputGrid(deck);
94 grid.processEclipseFormat(&inputGrid,
nullptr, periodic_extension, turn_normals, clip_z,
true);
97 if (param.getDefault(
"output_ecl",
false)) {
98 OPM_THROW(std::runtime_error,
"Saving to EGRID files is not yet implemented");
106 double perm_threshold_md = param.getDefault(
"perm_threshold_md", 0.0);
107 double perm_threshold = Opm::unit::convert::from(perm_threshold_md, Opm::prefix::milli*Opm::unit::darcy);
108 std::string rock_list = param.getDefault<std::string>(
"rock_list",
"no_list");
109 std::string* rl_ptr = (rock_list ==
"no_list") ? 0 : &rock_list;
110 bool use_j = param.getDefault(
"use_jfunction_scaling",
useJ<ResProp<3> >());
114 sigma = param.getDefault(
"sigma", sigma);
115 theta = param.getDefault(
"theta", theta);
117 if (param.has(
"viscosity1") || param.has(
"viscosity2")) {
118 double v1 = param.getDefault(
"viscosity1", 0.001);
119 double v2 = param.getDefault(
"viscosity2", 0.003);
120 res_prop.setViscosities(v1, v2);
122 res_prop.init(deck, grid.globalCell(), perm_threshold, rl_ptr,
123 use_j, sigma, theta);
124 }
else if (fileformat ==
"cartesian") {
125 std::array<int, 3> dims = {{ param.getDefault<
int>(
"nx", 1),
126 param.getDefault<
int>(
"ny", 1),
127 param.getDefault<
int>(
"nz", 1) }};
128 std::array<double, 3> cellsz = {{ param.getDefault<
double>(
"dx", 1.0),
129 param.getDefault<
double>(
"dy", 1.0),
130 param.getDefault<
double>(
"dz", 1.0) }};
131 grid.createCartesian(dims, cellsz);
132 double default_poro = param.getDefault(
"default_poro", 0.2);
133 double default_perm_md = param.getDefault(
"default_perm_md", 100.0);
134 double default_perm = Opm::unit::convert::from(default_perm_md, Opm::prefix::milli*Opm::unit::darcy);
135 OPM_MESSAGE(
"Warning: For generated cartesian grids, we use uniform reservoir properties.");
136 res_prop.init(grid.size(0), default_poro, default_perm);
138 OPM_THROW(std::runtime_error,
"Unknown file format string: " << fileformat);
140 if (param.getDefault(
"use_unique_boundary_ids",
false)) {
141 grid.setUniqueBoundaryIds(
true);
148 template <
template <
int>
class ResProp>
150 bool periodic_extension,
154 double perm_threshold,
155 const std::string& rock_list,
156 bool use_jfunction_scaling,
160 ResProp<3>& res_prop)
162 Opm::EclipseGrid eg(deck);
163 const std::string* rl_ptr = (rock_list ==
"no_list") ? 0 : &rock_list;
164 grid.processEclipseFormat(&eg,
nullptr, periodic_extension, turn_normals, clip_z,
true);
165 res_prop.init(deck, grid.globalCell(), perm_threshold, rl_ptr, use_jfunction_scaling, sigma, theta);
167 grid.setUniqueBoundaryIds(
true);
Inverting small matrices.
Definition: ImplicitAssembly.hpp:43
bool useJ()
Helper for determining whether we should.
Definition: setupGridAndProps.hpp:57
void setupGridAndPropsEclipse(const Opm::Deck &deck, bool periodic_extension, bool turn_normals, bool clip_z, bool unique_bids, double perm_threshold, const std::string &rock_list, bool use_jfunction_scaling, double sigma, double theta, Dune::CpGrid &grid, ResProp< 3 > &res_prop)
Definition: setupGridAndProps.hpp:149
void setupGridAndProps(const Opm::ParameterGroup ¶m, Dune::CpGrid &grid, ResProp< 3 > &res_prop)
Definition: setupGridAndProps.hpp:69