Go to the documentation of this file.
36#ifndef OPM_SETUPGRIDANDPROPS_HEADER
37#define OPM_SETUPGRIDANDPROPS_HEADER
39#include <opm/core/utility/parameters/ParameterGroup.hpp>
40#include <opm/core/utility/Units.hpp>
41#include <dune/grid/CpGrid.hpp>
42#include <dune/grid/sgrid.hh>
44#include <opm/parser/eclipse/Parser/Parser.hpp>
45#include <opm/parser/eclipse/Parser/ParseContext.hpp>
46#include <opm/parser/eclipse/Deck/Deck.hpp>
47#include <boost/filesystem.hpp>
60 bool useJ< ReservoirPropertyCapillary<3> >();
65 template < template < int> class ResProp>
72 std::string fileformat = param.getDefault<std::string>( "fileformat", "cartesian");
73 if (fileformat == "sintef_legacy") {
74 std::string grid_prefix = param.get<std::string>( "grid_prefix");
75 grid.readSintefLegacyFormat(grid_prefix);
76 OPM_MESSAGE( "Warning: We do not yet read legacy reservoir properties. Using defaults.");
77 res_prop.init(grid.size(0));
78 } else if (fileformat == "eclipse") {
79 std::string ecl_file = param.get<std::string>( "filename");
81 Opm::ParseContext parseContext;
82 Opm::ParserPtr parser( new Opm::Parser());
83 Opm::DeckConstPtr deck(parser->parseFile(ecl_file , parseContext));
84 if (param.has( "z_tolerance")) {
85 std::cerr << "****** Warning: z_tolerance parameter is obsolete, use PINCH in deck input instead\n";
87 bool periodic_extension = param.getDefault< bool>( "periodic_extension", false);
88 bool turn_normals = param.getDefault< bool>( "turn_normals", false);
89 grid.processEclipseFormat(deck, periodic_extension, turn_normals);
91 if (param.getDefault( "output_ecl", false)) {
92 OPM_THROW(std::runtime_error, "Saving to EGRID files is not yet implemented");
100 double perm_threshold_md = param.getDefault( "perm_threshold_md", 0.0);
101 double perm_threshold = Opm::unit::convert::from(perm_threshold_md, Opm::prefix::milli*Opm::unit::darcy);
102 std::string rock_list = param.getDefault<std::string>( "rock_list", "no_list");
103 std::string* rl_ptr = (rock_list == "no_list") ? 0 : &rock_list;
104 bool use_j = param.getDefault( "use_jfunction_scaling", useJ<ResProp<3> >());
108 sigma = param.getDefault( "sigma", sigma);
109 theta = param.getDefault( "theta", theta);
111 if (param.has( "viscosity1") || param.has( "viscosity2")) {
112 double v1 = param.getDefault( "viscosity1", 0.001);
113 double v2 = param.getDefault( "viscosity2", 0.003);
114 res_prop.setViscosities(v1, v2);
116 res_prop.init(deck, grid.globalCell(), perm_threshold, rl_ptr,
117 use_j, sigma, theta);
118 } else if (fileformat == "cartesian") {
119 std::array<int, 3> dims = {{ param.getDefault< int>( "nx", 1),
120 param.getDefault< int>( "ny", 1),
121 param.getDefault< int>( "nz", 1) }};
122 std::array<double, 3> cellsz = {{ param.getDefault< double>( "dx", 1.0),
123 param.getDefault< double>( "dy", 1.0),
124 param.getDefault< double>( "dz", 1.0) }};
125 grid.createCartesian(dims, cellsz);
126 double default_poro = param.getDefault( "default_poro", 0.2);
127 double default_perm_md = param.getDefault( "default_perm_md", 100.0);
128 double default_perm = Opm::unit::convert::from(default_perm_md, Opm::prefix::milli*Opm::unit::darcy);
129 OPM_MESSAGE( "Warning: For generated cartesian grids, we use uniform reservoir properties.");
130 res_prop.init(grid.size(0), default_poro, default_perm);
132 OPM_THROW(std::runtime_error, "Unknown file format string: " << fileformat);
134 if (param.getDefault( "use_unique_boundary_ids", false)) {
135 grid.setUniqueBoundaryIds( true);
142 template < template < int> class ResProp>
144 bool periodic_extension,
148 double perm_threshold,
149 const std::string& rock_list,
150 bool use_jfunction_scaling,
154 ResProp<3>& res_prop)
156 grid.processEclipseFormat(deck, periodic_extension, turn_normals, clip_z);
157 const std::string* rl_ptr = (rock_list == "no_list") ? 0 : &rock_list;
158 res_prop.init(deck, grid.globalCell(), perm_threshold, rl_ptr, use_jfunction_scaling, sigma, theta);
160 grid.setUniqueBoundaryIds( true);
167 template < template < int> class ResProp>
169 Dune::SGrid<3, 3>& grid,
170 ResProp<3>& res_prop)
174 std::string fileformat = param.getDefault<std::string>( "fileformat", "cartesian");
175 if (fileformat == "cartesian") {
176 std::array<int, 3> dims = {{ param.getDefault< int>( "nx", 1),
177 param.getDefault< int>( "ny", 1),
178 param.getDefault< int>( "nz", 1) }};
179 std::array<double, 3> cellsz = {{ param.getDefault< double>( "dx", 1.0),
180 param.getDefault< double>( "dy", 1.0),
181 param.getDefault< double>( "dz", 1.0) }};
183 new (&grid) Dune::SGrid<3, 3>(&dims[0], &cellsz[0]);
184 double default_poro = param.getDefault( "default_poro", 0.2);
185 double default_perm = param.getDefault( "default_perm", 100.0*Opm::prefix::milli*Opm::unit::darcy);
186 OPM_MESSAGE( "Warning: For generated cartesian grids, we use uniform reservoir properties.");
187 res_prop.init(grid.size(0), default_poro, default_perm);
189 OPM_THROW(std::runtime_error, "Dune::SGrid can only handle cartesian grids, unsupported file format string: " << fileformat);
Definition: BlackoilFluid.hpp:32
void setupGridAndProps(const Opm::parameter::ParameterGroup ¶m, Dune::CpGrid &grid, ResProp< 3 > &res_prop) Definition: setupGridAndProps.hpp:66
bool useJ() Helper for determining whether we should. Definition: setupGridAndProps.hpp:54
void setupGridAndPropsEclipse(Opm::DeckConstPtr 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:143
|