My Project
SimulatorBase.hpp
1//===========================================================================
2//
3// File: SimulatorBase.hpp
4//
5// Created: Tue Aug 11 15:01:48 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 OPENRS_SIMULATORBASE_HEADER
37#define OPENRS_SIMULATORBASE_HEADER
38
39
40#include <opm/common/utility/parameters/ParameterGroup.hpp>
41
42#include <opm/common/utility/numeric/SparseVector.hpp>
43#include <opm/grid/utility/SparseTable.hpp>
44#include <opm/input/eclipse/Units/Units.hpp>
45
46#include <opm/grid/common/Volumes.hpp>
47#include <opm/grid/CpGrid.hpp>
48
49#include <opm/porsol/common/GridInterfaceEuler.hpp>
50#include <opm/porsol/common/ReservoirPropertyCapillary.hpp>
51#include <opm/porsol/common/BoundaryConditions.hpp>
52#include <opm/porsol/common/setupGridAndProps.hpp>
53#include <opm/porsol/common/setupBoundaryConditions.hpp>
54#include <opm/porsol/common/SimulatorUtilities.hpp>
55
56#include <opm/porsol/euler/EulerUpstream.hpp>
57#include <opm/porsol/euler/ImplicitCapillarity.hpp>
58
59#include <opm/porsol/mimetic/MimeticIPEvaluator.hpp>
60#include <opm/porsol/mimetic/IncompFlowSolverHybrid.hpp>
61
62
63#include <opm/common/utility/platform_dependent/disable_warnings.h>
64
65#include <dune/grid/yaspgrid.hh>
66
67#include <opm/common/utility/platform_dependent/reenable_warnings.h>
68
69
70#include <fstream>
71#include <iterator>
72#include <iostream>
73
74
75namespace Opm
76{
77
78
79
80
84 template <class SimTraits>
86 {
87 public:
88
92 : simulation_steps_(1),
93 stepsize_(1.0), // init() expects units of days! Yes, but now the meaning of stepsize_ changes
94 // from days (here) to seconds (after init()). Solution to that?
95 residual_tolerance_(1e-8),
96 linsolver_verbosity_(1),
97 linsolver_type_(1)
98 {
99 }
100
103 void init(const Opm::ParameterGroup& param)
104 {
105 initControl(param);
106 initGridAndProps(param);
107 initInitialConditions(param);
108 initBoundaryConditions(param);
109 initSources(param);
110 initSolvers(param);
111
112 // Write any unused parameters.
113 std::cout << "==================== Unused parameters: ====================\n";
114 param.displayUsage();
115 std::cout << "================================================================\n";
116 }
117
118 protected:
119 typedef Dune::CpGrid GridType;
120 enum { Dimension = GridType::dimension };
121 typedef Dune::FieldVector<double, Dimension> Vector;
122 typedef typename SimTraits::template ResProp<Dimension>::Type ResProp;
123 typedef GridInterfaceEuler<GridType> GridInterface;
124 typedef GridInterface::CellIterator CellIter;
125 typedef CellIter::FaceIterator FaceIter;
126 typedef BasicBoundaryConditions<true, true> BCs;
127 typedef typename SimTraits::template FlowSolver<GridInterface, BCs>::Type FlowSolver;
128 typedef typename SimTraits::template TransportSolver<GridInterface, BCs>::Type TransportSolver;
129
130 int simulation_steps_;
131 double stepsize_;
132 std::vector<double> init_saturation_;
133 Vector gravity_;
134 double residual_tolerance_;
135 int linsolver_verbosity_;
136 int linsolver_type_;
137
138 GridType grid_;
139 GridInterface ginterf_;
140 ResProp res_prop_;
141 BCs bcond_;
142 Opm::SparseVector<double> injection_rates_;
143 std::vector<double> injection_rates_psolver_; // Should modify psolver to take SparseVector
144 FlowSolver flow_solver_;
145 TransportSolver transport_solver_;
146
147
148 virtual void initControl(const Opm::ParameterGroup& param)
149 {
150 simulation_steps_ = param.getDefault("simulation_steps", simulation_steps_);
151 stepsize_ = Opm::unit::convert::from(param.getDefault("stepsize", stepsize_),
152 Opm::unit::day);
153 }
154
155 virtual void initGridAndProps(const Opm::ParameterGroup& param)
156 {
157 setupGridAndProps(param, grid_, res_prop_);
158 ginterf_.init(grid_);
159
160 gravity_[0] = param.getDefault("gx", 0.0);
161 gravity_[1] = param.getDefault("gy", 0.0);
162 gravity_[2] = param.getDefault("gz", 0.0); //Dune::unit::gravity);
163 }
164
165 virtual void initInitialConditions(const Opm::ParameterGroup& param)
166 {
167 if (param.getDefault("init_saturation_from_file", false)) {
168 std::string filename = param.get<std::string>("init_saturation_filename");
169 std::ifstream satfile(filename.c_str());
170 if (!satfile) {
171 OPM_THROW(std::runtime_error, "Could not open initial saturation file: " << filename);
172 }
173 int num_sats;
174 satfile >> num_sats;
175 if (num_sats != ginterf_.numberOfCells()) {
176 OPM_THROW(std::runtime_error, "Number of saturation values claimed different from number of grid cells: "
177 << num_sats << " vs. " << ginterf_.numberOfCells());
178 }
179 std::istream_iterator<double> beg(satfile);
180 std::istream_iterator<double> end;
181 init_saturation_.assign(beg, end);
182 if (int(init_saturation_.size()) != num_sats) {
183 OPM_THROW(std::runtime_error, "Number of saturation values claimed different from actual file content: "
184 << num_sats << " vs. " << init_saturation_.size());
185 }
186 } else {
187 double init_s = param.getDefault("init_saturation", 0.0);
188 init_saturation_.clear();
189 init_saturation_.resize(ginterf_.numberOfCells(), init_s);
190 }
191 }
192
193 virtual void initBoundaryConditions(const Opm::ParameterGroup& param)
194 {
195 setupBoundaryConditions(param, ginterf_, bcond_);
196 }
197
198 virtual void initSources(const Opm::ParameterGroup& /* param */)
199 {
200 int nc = ginterf_.numberOfCells();
201 injection_rates_ = Opm::SparseVector<double>(nc);
202 injection_rates_psolver_.resize(nc, 0.0);
203// injection_rates_.addElement(1.0, 0);
204// injection_rates_.addElement(-1.0, nc - 1);
205// injection_rates_psolver_[0] = 1.0;
206// injection_rates_psolver_[nc - 1] = -1.0;
207 }
208
209 virtual void initSolvers(const Opm::ParameterGroup& param)
210 {
211 // Initialize flow solver.
212 flow_solver_.init(ginterf_, res_prop_, gravity_, bcond_);
213 residual_tolerance_ = param.getDefault("residual_tolerance", residual_tolerance_);
214 linsolver_verbosity_ = param.getDefault("linsolver_verbosity", linsolver_verbosity_);
215 linsolver_type_ = param.getDefault("linsolver_type", linsolver_type_);
216 //flow_solver_.assembleStatic(ginterf_, res_prop_);
217 // Initialize transport solver.
218 transport_solver_.init(param, ginterf_, res_prop_, bcond_);
219 }
220
221
222 };
223
224
225
226} // namespace Opm
227
228
229
230#endif // OPENRS_SIMULATORBASE_HEADER
Definition: SimulatorBase.hpp:86
void init(const Opm::ParameterGroup &param)
Initialization from parameters.
Definition: SimulatorBase.hpp:103
SimulatorBase()
Definition: SimulatorBase.hpp:91
Inverting small matrices.
Definition: ImplicitAssembly.hpp:43
void setupBoundaryConditions(const Opm::ParameterGroup &param, const GridInterface &g, BCs &bcs)
Setup boundary conditions for a simulation.
Definition: setupBoundaryConditions.hpp:51
void setupGridAndProps(const Opm::ParameterGroup &param, Dune::CpGrid &grid, ResProp< 3 > &res_prop)
Definition: setupGridAndProps.hpp:69