My Project
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
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
37class ParallelWellInfo;
38template<class Scalar, int numEq> class StandardWellEquationAccess;
39class WellContributions;
40class WellInterfaceGeneric;
41class WellState;
42
43template<class Scalar, int numEq>
45{
46public:
47 // sparsity pattern for the matrices
48 //[A C^T [x = [ res
49 // B D ] x_well] res_well]
50
51 // the vector type for the res_well and x_well
52 using VectorBlockWellType = Dune::DynamicVector<Scalar>;
53 using BVectorWell = Dune::BlockVector<VectorBlockWellType>;
54
55 // the matrix type for the diagonal matrix D
56 using DiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
57 using DiagMatWell = Dune::BCRSMatrix<DiagMatrixBlockWellType>;
58
59 // the matrix type for the non-diagonal matrix B and C^T
60 using OffDiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
61 using OffDiagMatWell = Dune::BCRSMatrix<OffDiagMatrixBlockWellType>;
62
63 // block vector type
64 using BVector = Dune::BlockVector<Dune::FieldVector<Scalar,numEq>>;
65
66 StandardWellEquations(const ParallelWellInfo& parallel_well_info);
67
73 void init(const int num_cells,
74 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 invert();
92
95 void recoverSolutionWell(const BVector& x, BVectorWell& xw) const;
96
98 void extract(const int numStaticWellEq,
99 WellContributions& wellContribs) const;
100
102 template<class SparseMatrixAdapter>
103 void extract(SparseMatrixAdapter& jacobian) const;
104
106 template<class PressureMatrix>
107 void extractCPRPressureMatrix(PressureMatrix& jacobian,
108 const BVector& weights,
109 const int pressureVarIndex,
110 const bool use_well_weights,
111 const WellInterfaceGeneric& well,
112 const int bhp_var_index,
113 const WellState& well_state) const;
114
116 unsigned int getNumBlocks() const;
117
119 void sumDistributed(Parallel::Communication comm);
120
122 const BVectorWell& residual() const
123 {
124 return resWell_;
125 }
126
127private:
128 friend class StandardWellEquationAccess<Scalar,numEq>;
129
130 // two off-diagonal matrices
131 OffDiagMatWell duneB_;
132 OffDiagMatWell duneC_;
133 // diagonal matrix for the well
134 DiagMatWell invDuneD_;
135 DiagMatWell duneD_;
136
137 // Wrapper for the parallel application of B for distributed wells
139
140 // residuals of the well equations
141 BVectorWell resWell_;
142
143 // several vector used in the matrix calculation
144 mutable BVectorWell Bx_;
145 mutable BVectorWell invDrw_;
146};
147
148}
149
150#endif // OPM_STANDARDWELL_EQUATIONS_HEADER_INCLUDED
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:184
Class administering assembler access to equation system.
Definition: StandardWellAssemble.cpp:45
Definition: StandardWellEquations.hpp:45
const BVectorWell & residual() const
Returns a const reference to the residual.
Definition: StandardWellEquations.hpp:122
void invert()
Invert D matrix.
Definition: StandardWellEquations.cpp:156
void apply(const BVector &x, BVector &Ax) const
Apply linear operator to vector.
Definition: StandardWellEquations.cpp:126
void init(const int num_cells, const int numWellEq, const int numPerfs, const std::vector< int > &cells)
Setup sparsity pattern for the matrices.
Definition: StandardWellEquations.cpp:51
void extract(const int numStaticWellEq, WellContributions &wellContribs) const
Add the matrices of this well to the WellContributions object.
Definition: StandardWellEquations.cpp:190
unsigned int getNumBlocks() const
Get the number of blocks of the C and B matrices.
Definition: StandardWellEquations.cpp:271
void sumDistributed(Parallel::Communication comm)
Sum with off-process contribution.
Definition: StandardWellEquations.cpp:390
void solve(BVectorWell &dx_well) const
Apply inverted D matrix to residual and store in vector.
Definition: StandardWellEquations.cpp:171
void recoverSolutionWell(const BVector &x, BVectorWell &xw) const
Recover well solution.
Definition: StandardWellEquations.cpp:179
void clear()
Set all coefficients to 0.
Definition: StandardWellEquations.cpp:117
void extractCPRPressureMatrix(PressureMatrix &jacobian, const BVector &weights, const int pressureVarIndex, const bool use_well_weights, const WellInterfaceGeneric &well, const int bhp_var_index, const WellState &well_state) const
Extract CPR pressure matrix.
Definition: StandardWellEquations.cpp:279
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:52
Definition: WellInterfaceGeneric.hpp:51
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellState.hpp:60
A wrapper around the B matrix for distributed wells.
Definition: WellHelpers.hpp:53
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27