22#ifndef OPM_MULTISEGMENTWELL_EQUATIONS_HEADER_INCLUDED
23#define OPM_MULTISEGMENTWELL_EQUATIONS_HEADER_INCLUDED
25#include <dune/common/fmatrix.hh>
26#include <dune/common/fvector.hh>
27#include <dune/istl/bcrsmatrix.hh>
28#include <dune/istl/bvector.hh>
33template<
class M>
class UMFPack;
39template<
class Scalar,
int numWellEq,
int numEq>
class MultisegmentWellEquationAccess;
40template<
class Scalar>
class MultisegmentWellGeneric;
41class WellContributions;
42class WellInterfaceGeneric;
45template<
class Scalar,
int numWellEq,
int numEq>
54 using VectorBlockWellType = Dune::FieldVector<Scalar,numWellEq>;
55 using BVectorWell = Dune::BlockVector<VectorBlockWellType>;
57 using VectorBlockType = Dune::FieldVector<Scalar,numEq>;
58 using BVector = Dune::BlockVector<VectorBlockType>;
61 using DiagMatrixBlockWellType = Dune::FieldMatrix<Scalar,numWellEq,numWellEq>;
62 using DiagMatWell = Dune::BCRSMatrix<DiagMatrixBlockWellType>;
65 using OffDiagMatrixBlockWellType = Dune::FieldMatrix<Scalar,numWellEq,numEq>;
66 using OffDiagMatWell = Dune::BCRSMatrix<OffDiagMatrixBlockWellType>;
76 void init(
const int num_cells,
78 const std::vector<int>& cells,
79 const std::vector<std::vector<int>>& segment_inlets,
80 const std::vector<std::vector<int>>& segment_perforations);
86 void apply(
const BVector& x, BVector& Ax)
const;
89 void apply(BVector& r)
const;
95 BVectorWell
solve()
const;
105 template<
class SparseMatrixAdapter>
106 void extract(SparseMatrixAdapter& jacobian)
const;
109 template<
class PressureMatrix>
111 const BVector& weights,
112 const int pressureVarIndex,
115 const int seg_pressure_var_ind,
127 OffDiagMatWell duneB_;
128 OffDiagMatWell duneC_;
135 mutable std::shared_ptr<Dune::UMFPack<DiagMatWell>> duneDSolver_;
138 BVectorWell resWell_;
Class administering assembler access to equation system.
Definition: MultisegmentWellAssemble.cpp:45
Definition: MultisegmentWellEquations.hpp:47
void extract(WellContributions &wellContribs) const
Add the matrices of this well to the WellContributions object.
Definition: MultisegmentWellEquations.cpp:194
void clear()
Set all coefficients to 0.
Definition: MultisegmentWellEquations.cpp:125
void extractCPRPressureMatrix(PressureMatrix &jacobian, const BVector &weights, const int pressureVarIndex, const bool, const WellInterfaceGeneric &well, const int seg_pressure_var_ind, const WellState &well_state) const
Extract CPR pressure matrix.
Definition: MultisegmentWellEquations.cpp:299
void apply(const BVector &x, BVector &Ax) const
Apply linear operator to vector.
Definition: MultisegmentWellEquations.cpp:136
void recoverSolutionWell(const BVector &x, BVectorWell &xw) const
Recover well solution.
Definition: MultisegmentWellEquations.cpp:183
void createSolver()
Compute the LU-decomposition of D matrix.
Definition: MultisegmentWellEquations.cpp:160
BVectorWell solve() const
Apply inverted D matrix to residual and return result.
Definition: MultisegmentWellEquations.cpp:176
const BVectorWell & residual() const
Returns a const reference to the residual.
Definition: MultisegmentWellEquations.hpp:119
void init(const int num_cells, const int numPerfs, const std::vector< int > &cells, const std::vector< std::vector< int > > &segment_inlets, const std::vector< std::vector< int > > &segment_perforations)
Setup sparsity pattern for the matrices.
Definition: MultisegmentWellEquations.cpp:53
Definition: MultisegmentWellGeneric.hpp:42
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
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27