My Project
Loading...
Searching...
No Matches
StandardWellEquations.hpp
1/*
2 Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2017 Statoil ASA.
4 Copyright 2016 - 2017 IRIS AS.
5
6 This file is part of the Open Porous Media project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22
23#ifndef OPM_STANDARDWELL_EQUATIONS_HEADER_INCLUDED
24#define OPM_STANDARDWELL_EQUATIONS_HEADER_INCLUDED
25
26#include <opm/simulators/utils/ParallelCommunication.hpp>
27#include <opm/simulators/wells/WellHelpers.hpp>
28#include <opm/common/TimingMacros.hpp>
29#include <dune/common/dynmatrix.hh>
30#include <dune/common/dynvector.hh>
31#include <dune/istl/bcrsmatrix.hh>
32#include <dune/istl/bvector.hh>
33
34namespace Opm
35{
36
37template<class Scalar> class ParallelWellInfo;
38template<class Scalar, int numEq> class StandardWellEquationAccess;
39#if COMPILE_BDA_BRIDGE
40template<class Scalar> class WellContributions;
41#endif
42template<class Scalar> class WellInterfaceGeneric;
43template<class Scalar> class WellState;
44
45template<class Scalar, int numEq>
47{
48public:
49 // sparsity pattern for the matrices
50 //[A C^T [x = [ res
51 // B D ] x_well] res_well]
52
53 // the vector type for the res_well and x_well
54 using VectorBlockWellType = Dune::DynamicVector<Scalar>;
55 using BVectorWell = Dune::BlockVector<VectorBlockWellType>;
56
57 // the matrix type for the diagonal matrix D
58 using DiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
59 using DiagMatWell = Dune::BCRSMatrix<DiagMatrixBlockWellType>;
60
61 // the matrix type for the non-diagonal matrix B and C^T
62 using OffDiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
63 using OffDiagMatWell = Dune::BCRSMatrix<OffDiagMatrixBlockWellType>;
64
65 // block vector type
66 using BVector = Dune::BlockVector<Dune::FieldVector<Scalar,numEq>>;
67
69
74 void init(const int numWellEq,
75 const int numPerfs,
76 const std::vector<int>& cells);
77
79 void clear();
80
82 void apply(const BVector& x, BVector& Ax) const;
83
85 void apply(BVector& r) const;
86
88 void solve(BVectorWell& dx_well) const;
89
91 void solve(const BVectorWell& rhs_well, BVectorWell& x_well) const;
92
94 void invert();
95
98 void recoverSolutionWell(const BVector& x, BVectorWell& xw) const;
99
100#if COMPILE_BDA_BRIDGE
102 void extract(const int numStaticWellEq,
104#endif
105
107 template<class SparseMatrixAdapter>
108 void extract(SparseMatrixAdapter& jacobian) const;
109
111 template<class PressureMatrix>
112 void extractCPRPressureMatrix(PressureMatrix& jacobian,
113 const BVector& weights,
114 const int pressureVarIndex,
115 const bool use_well_weights,
117 const int bhp_var_index,
118 const WellState<Scalar>& well_state) const;
119
121 unsigned int getNumBlocks() const;
122
124 void sumDistributed(Parallel::Communication comm);
125
127 const BVectorWell& residual() const
128 {
129 return resWell_;
130 }
131
132private:
133 friend class StandardWellEquationAccess<Scalar,numEq>;
134
135 // two off-diagonal matrices
136 OffDiagMatWell duneB_;
137 OffDiagMatWell duneC_;
138 // diagonal matrix for the well
139 DiagMatWell invDuneD_;
140 DiagMatWell duneD_;
141
142 // Wrapper for the parallel application of B for distributed wells
144
145 // residuals of the well equations
146 BVectorWell resWell_;
147
148 // several vector used in the matrix calculation
149 mutable BVectorWell Bx_;
150 mutable BVectorWell invDrw_;
151
152 // Store the global index of well perforated cells
153 std::vector<int> cells_;
154};
155
156}
157
158#endif // OPM_STANDARDWELL_EQUATIONS_HEADER_INCLUDED
Class encapsulating some information about parallel wells.
Definition ParallelWellInfo.hpp:186
Class administering assembler access to equation system.
Definition StandardWellAssemble.cpp:45
Definition StandardWellEquations.hpp:47
const BVectorWell & residual() const
Returns a const reference to the residual.
Definition StandardWellEquations.hpp:127
void extract(SparseMatrixAdapter &jacobian) const
Add the matrices of this well to the sparse matrix adapter.
Definition StandardWellEquations.cpp:254
void invert()
Invert D matrix.
Definition StandardWellEquations.cpp:158
void extractCPRPressureMatrix(PressureMatrix &jacobian, const BVector &weights, const int pressureVarIndex, const bool use_well_weights, const WellInterfaceGeneric< Scalar > &well, const int bhp_var_index, const WellState< Scalar > &well_state) const
Extract CPR pressure matrix.
Definition StandardWellEquations.cpp:291
void apply(const BVector &x, BVector &Ax) const
Apply linear operator to vector.
Definition StandardWellEquations.cpp:128
unsigned int getNumBlocks() const
Get the number of blocks of the C and B matrices.
Definition StandardWellEquations.cpp:283
void sumDistributed(Parallel::Communication comm)
Sum with off-process contribution.
Definition StandardWellEquations.cpp:404
void solve(BVectorWell &dx_well) const
Apply inverted D matrix to residual and store in vector.
Definition StandardWellEquations.cpp:173
void recoverSolutionWell(const BVector &x, BVectorWell &xw) const
Recover well solution.
Definition StandardWellEquations.cpp:186
void init(const int numWellEq, const int numPerfs, const std::vector< int > &cells)
Setup sparsity pattern for the matrices.
Definition StandardWellEquations.cpp:55
void clear()
Set all coefficients to 0.
Definition StandardWellEquations.cpp:119
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition WellContributions.hpp:53
Definition WellInterfaceGeneric.hpp:52
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:62
A wrapper around the B matrix for distributed wells.
Definition WellHelpers.hpp:51
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242