GenericTracerModel.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef OPM_GENERIC_TRACER_MODEL_HPP
29#define OPM_GENERIC_TRACER_MODEL_HPP
30
31#include <dune/istl/bcrsmatrix.hh>
32
33#include <opm/grid/common/CartesianIndexMapper.hpp>
34
35#include <opm/input/eclipse/EclipseState/Phase.hpp>
36
38
41
42#include <array>
43#include <cstddef>
44#include <functional>
45#include <memory>
46#include <string>
47#include <unordered_map>
48#include <vector>
49
50namespace Opm {
51
52class EclipseState;
53class Well;
54
55template<class Grid, class GridView, class DofMapper, class Stencil, class FluidSystem, class Scalar>
57public:
58 using TracerVectorSingle = Dune::BlockVector<Dune::FieldVector<Scalar, 1>>;
59 using TracerMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<Scalar, 2, 2>>;
60 using TracerVector = Dune::BlockVector<Dune::FieldVector<Scalar, 2>>;
62 static constexpr int dimWorld = Grid::dimensionworld;
66 int numTracers() const;
67
71 const std::string& name(int tracerIdx) const;
72 std::string fname(int tracerIdx) const;
73 std::string sname(int tracerIdx) const;
74 std::string wellfname(int tracerIdx) const;
75 std::string wellsname(int tracerIdx) const;
76
77 Phase phase(int tracerIdx) const;
78 const std::vector<bool>& enableSolTracers() const;
79
83 Scalar freeTracerConcentration(int tracerIdx, int globalDofIdx) const;
84 Scalar solTracerConcentration(int tracerIdx, int globalDofIdx) const;
85 void setFreeTracerConcentration(int tracerIdx, int globalDofIdx, Scalar value);
86 void setSolTracerConcentration(int tracerIdx, int globalDofIdx, Scalar value);
87 void setEnableSolTracers(int tracerIdx, bool enableSolTracer);
88
92 const std::unordered_map<int, std::vector<WellTracerRate<Scalar>>>&
94 { return wellTracerRate_; }
95
96 const std::unordered_map<int, std::vector<WellTracerRate<Scalar>>>&
98 { return wellFreeTracerRate_; }
99
100 const std::unordered_map<int, std::vector<WellTracerRate<Scalar>>>&
102 { return wellSolTracerRate_; }
103
104 const std::unordered_map<int, std::vector<MSWellTracerRate<Scalar>>>&
106
107 template<class Serializer>
108 void serializeOp(Serializer& serializer)
109 {
110 serializer(tracerConcentration_);
111 serializer(freeTracerConcentration_);
112 serializer(solTracerConcentration_);
113 serializer(wellTracerRate_);
114 serializer(wellFreeTracerRate_);
115 serializer(wellSolTracerRate_);
116 serializer(mSwTracerRate_);
117 }
118
119protected:
120 GenericTracerModel(const GridView& gridView,
121 const EclipseState& eclState,
122 const CartesianIndexMapper& cartMapper,
123 const DofMapper& dofMapper,
124 const std::function<std::array<double,dimWorld>(int)> centroids);
125
129 void doInit(bool rst,
130 std::size_t numGridDof,
131 std::size_t gasPhaseIdx,
132 std::size_t oilPhaseIdx,
133 std::size_t waterPhaseIdx);
134
135 bool linearSolve_(const TracerMatrix& M, TracerVector& x, TracerVector& b);
136
138 std::vector<TracerVector>& x,
139 std::vector<TracerVector>& b);
140
141 Scalar currentConcentration_(const Well& eclWell, const std::string& name) const;
142
145 Free = 0,
147 };
148
149 const GridView& gridView_;
150 const EclipseState& eclState_;
152 const DofMapper& dofMapper_;
153
154 std::vector<int> tracerPhaseIdx_;
155 std::vector<bool> enableSolTracers_;
156 std::vector<TracerVector> tracerConcentration_;
157 std::unique_ptr<TracerMatrix> tracerMatrix_;
158 std::vector<TracerVectorSingle> freeTracerConcentration_;
159 std::vector<TracerVectorSingle> solTracerConcentration_;
160
161 // well_index -> tracer rates
162 std::unordered_map<int, std::vector<WellTracerRate<Scalar>>> wellTracerRate_;
163 std::unordered_map<int, std::vector<WellTracerRate<Scalar>>> wellFreeTracerRate_;
164 std::unordered_map<int, std::vector<WellTracerRate<Scalar>>> wellSolTracerRate_;
165
166 std::unordered_map<int, std::vector<MSWellTracerRate<Scalar>>> mSwTracerRate_;
167
169 std::function<std::array<double,dimWorld>(int)> centroids_;
170};
171
172} // namespace Opm
173
174#endif // OPM_GENERIC_TRACER_MODEL_HPP
Definition: CollectDataOnIORank.hpp:49
Definition: GenericTracerModel.hpp:56
void setEnableSolTracers(int tracerIdx, bool enableSolTracer)
Definition: GenericTracerModel_impl.hpp:155
int numTracers() const
Return the number of tracers considered by the tracerModel.
Definition: GenericTracerModel_impl.hpp:162
Scalar currentConcentration_(const Well &eclWell, const std::string &name) const
Definition: GenericTracerModel_impl.hpp:211
const std::unordered_map< int, std::vector< WellTracerRate< Scalar > > > & getWellTracerRates() const
Return well tracer rates.
Definition: GenericTracerModel.hpp:93
void doInit(bool rst, std::size_t numGridDof, std::size_t gasPhaseIdx, std::size_t oilPhaseIdx, std::size_t waterPhaseIdx)
Initialize all internal data structures needed by the tracer module.
Definition: GenericTracerModel_impl.hpp:225
std::string sname(int tracerIdx) const
Definition: GenericTracerModel_impl.hpp:176
std::vector< TracerVector > tracerConcentration_
Definition: GenericTracerModel.hpp:156
std::vector< int > tracerPhaseIdx_
Definition: GenericTracerModel.hpp:154
std::unordered_map< int, std::vector< WellTracerRate< Scalar > > > wellFreeTracerRate_
Definition: GenericTracerModel.hpp:163
const DofMapper & dofMapper_
Definition: GenericTracerModel.hpp:152
std::unordered_map< int, std::vector< WellTracerRate< Scalar > > > wellTracerRate_
Definition: GenericTracerModel.hpp:162
const GridView & gridView_
Definition: GenericTracerModel.hpp:149
std::unordered_map< int, std::vector< MSWellTracerRate< Scalar > > > mSwTracerRate_
Definition: GenericTracerModel.hpp:166
const std::unordered_map< int, std::vector< WellTracerRate< Scalar > > > & getWellFreeTracerRates() const
Definition: GenericTracerModel.hpp:97
const EclipseState & eclState_
Definition: GenericTracerModel.hpp:150
void setSolTracerConcentration(int tracerIdx, int globalDofIdx, Scalar value)
Definition: GenericTracerModel_impl.hpp:147
std::string wellfname(int tracerIdx) const
Definition: GenericTracerModel_impl.hpp:183
GenericTracerModel(const GridView &gridView, const EclipseState &eclState, const CartesianIndexMapper &cartMapper, const DofMapper &dofMapper, const std::function< std::array< double, dimWorld >(int)> centroids)
Definition: GenericTracerModel_impl.hpp:102
const std::string & name(int tracerIdx) const
Return the tracer name.
Definition: GenericTracerModel_impl.hpp:218
const std::vector< bool > & enableSolTracers() const
Definition: GenericTracerModel_impl.hpp:204
Phase phase(int tracerIdx) const
Definition: GenericTracerModel_impl.hpp:197
const std::unordered_map< int, std::vector< WellTracerRate< Scalar > > > & getWellSolTracerRates() const
Definition: GenericTracerModel.hpp:101
std::string wellsname(int tracerIdx) const
Definition: GenericTracerModel_impl.hpp:190
const std::unordered_map< int, std::vector< MSWellTracerRate< Scalar > > > & getMswTracerRates() const
Definition: GenericTracerModel.hpp:105
std::string fname(int tracerIdx) const
Definition: GenericTracerModel_impl.hpp:169
std::vector< TracerVectorSingle > freeTracerConcentration_
Definition: GenericTracerModel.hpp:158
std::function< std::array< double, dimWorld >(int)> centroids_
Function returning the cell centers.
Definition: GenericTracerModel.hpp:169
bool linearSolve_(const TracerMatrix &M, TracerVector &x, TracerVector &b)
Definition: GenericTracerModel_impl.hpp:387
const CartesianIndexMapper & cartMapper_
Definition: GenericTracerModel.hpp:151
Dune::BlockVector< Dune::FieldVector< Scalar, 2 > > TracerVector
Definition: GenericTracerModel.hpp:60
Dune::CartesianIndexMapper< Grid > CartesianIndexMapper
Definition: GenericTracerModel.hpp:61
void setFreeTracerConcentration(int tracerIdx, int globalDofIdx, Scalar value)
Definition: GenericTracerModel_impl.hpp:139
std::vector< bool > enableSolTracers_
Definition: GenericTracerModel.hpp:155
Scalar solTracerConcentration(int tracerIdx, int globalDofIdx) const
Definition: GenericTracerModel_impl.hpp:128
TracerTypeIdx
Tracer type index.
Definition: GenericTracerModel.hpp:144
@ Free
Definition: GenericTracerModel.hpp:145
@ Solution
Definition: GenericTracerModel.hpp:146
static constexpr int dimWorld
Definition: GenericTracerModel.hpp:62
bool linearSolveBatchwise_(const TracerMatrix &M, std::vector< TracerVector > &x, std::vector< TracerVector > &b)
Definition: GenericTracerModel_impl.hpp:442
Dune::BCRSMatrix< Opm::MatrixBlock< Scalar, 2, 2 > > TracerMatrix
Definition: GenericTracerModel.hpp:59
void serializeOp(Serializer &serializer)
Definition: GenericTracerModel.hpp:108
std::vector< TracerVectorSingle > solTracerConcentration_
Definition: GenericTracerModel.hpp:159
Scalar freeTracerConcentration(int tracerIdx, int globalDofIdx) const
Return the tracer concentration for tracer index and global DofIdx.
Definition: GenericTracerModel_impl.hpp:117
Dune::BlockVector< Dune::FieldVector< Scalar, 1 > > TracerVectorSingle
Definition: GenericTracerModel.hpp:58
std::unique_ptr< TracerMatrix > tracerMatrix_
Definition: GenericTracerModel.hpp:157
std::unordered_map< int, std::vector< WellTracerRate< Scalar > > > wellSolTracerRate_
Definition: GenericTracerModel.hpp:164
Definition: blackoilboundaryratevector.hh:39