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/parser/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 
75 namespace 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