My Project
setupGridAndProps.hpp
1//===========================================================================
2//
3// File: setupGridAndProps.hpp
4//
5// Created: Tue Aug 11 14:47:35 2009
6//
7// Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
8// Bård Skaflestad <bard.skaflestad@sintef.no>
9//
10// $Date$
11//
12// $Revision$
13//
14//===========================================================================
15
16/*
17 Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
18 Copyright 2009, 2010 Statoil ASA.
19
20 This file is part of The Open Reservoir Simulator Project (OpenRS).
21
22 OpenRS is free software: you can redistribute it and/or modify
23 it under the terms of the GNU General Public License as published by
24 the Free Software Foundation, either version 3 of the License, or
25 (at your option) any later version.
26
27 OpenRS is distributed in the hope that it will be useful,
28 but WITHOUT ANY WARRANTY; without even the implied warranty of
29 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 GNU General Public License for more details.
31
32 You should have received a copy of the GNU General Public License
33 along with OpenRS. If not, see <http://www.gnu.org/licenses/>.
34*/
35
36#ifndef OPM_SETUPGRIDANDPROPS_HEADER
37#define OPM_SETUPGRIDANDPROPS_HEADER
38
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>
45
46#include <opm/common/utility/platform_dependent/disable_warnings.h>
47
48#include <opm/common/utility/FileSystem.hpp>
49
50#include <opm/common/utility/platform_dependent/reenable_warnings.h>
51
52namespace Opm
53{
54
56 template <class RP>
57 bool useJ()
58 {
59 return false;
60 }
61
62 template<>
63 bool useJ< ReservoirPropertyCapillary<3> >();
64
68 template <template <int> class ResProp>
69 inline void setupGridAndProps(const Opm::ParameterGroup& param,
70 Dune::CpGrid& grid,
71 ResProp<3>& res_prop)
72 {
73 // Initialize grid and reservoir properties.
74 // Parts copied from Dune::CpGrid::init().
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");
83
84 Opm::Parser parser;
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";
88 }
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);
92 {
93 Opm::EclipseGrid inputGrid(deck);
94 grid.processEclipseFormat(&inputGrid, nullptr, periodic_extension, turn_normals, clip_z, true);
95 }
96 // Save EGRID file in case we are writing ECL output.
97 if (param.getDefault("output_ecl", false)) {
98 OPM_THROW(std::runtime_error, "Saving to EGRID files is not yet implemented");
99 /*
100 Opm::filesystem::path ecl_path(ecl_file);
101 const std::vector<int>& globalCell = grid.globalCell();
102 ecl_path.replace_extension(".EGRID");
103 parser.saveEGRID(ecl_path.string() , (int) globalCell.size() , &globalCell[0]);
104 */
105 }
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> >());
111 double sigma = 1.0;
112 double theta = 0.0;
113 if (use_j) {
114 sigma = param.getDefault("sigma", sigma);
115 theta = param.getDefault("theta", theta);
116 }
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);
121 }
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);
137 } else {
138 OPM_THROW(std::runtime_error, "Unknown file format string: " << fileformat);
139 }
140 if (param.getDefault("use_unique_boundary_ids", false)) {
141 grid.setUniqueBoundaryIds(true);
142 }
143 }
144
148 template <template <int> class ResProp>
149 inline void setupGridAndPropsEclipse(const Opm::Deck& deck,
150 bool periodic_extension,
151 bool turn_normals,
152 bool clip_z,
153 bool unique_bids,
154 double perm_threshold,
155 const std::string& rock_list,
156 bool use_jfunction_scaling,
157 double sigma,
158 double theta,
159 Dune::CpGrid& grid,
160 ResProp<3>& res_prop)
161 {
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);
166 if (unique_bids) {
167 grid.setUniqueBoundaryIds(true);
168 }
169 }
170
171} // namespace Opm
172
173
174#endif // OPENRS_SETUPGRIDANDPROPS_HEADER
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 &param, Dune::CpGrid &grid, ResProp< 3 > &res_prop)
Definition: setupGridAndProps.hpp:69