added: StandardWellEquations::extract(SparseMatrixAdapter)

this adds the well matrices to a sparse matrix adapter.
this is the core of StandardWellEval::addWellContributions.
use the new method in the implementation.
This commit is contained in:
Arne Morten Kvarving 2022-11-11 21:41:24 +01:00
parent 65efe8e1fc
commit f2d8a073f9
3 changed files with 36 additions and 21 deletions

View File

@ -25,7 +25,9 @@
#include <opm/common/Exceptions.hpp>
#include <opm/simulators/linalg/bda/WellContributions.hpp>
#include <opm/simulators/linalg/istlsparsematrixadapter.hh>
#include <opm/simulators/linalg/matrixblock.hh>
#include <opm/simulators/linalg/SmallDenseMatrixUtils.hpp>
namespace Opm
{
@ -233,8 +235,36 @@ extract(const int numStaticWellEq,
colIndices.data(), nnzValues.data(), duneB_.nonzeroes());
}
template<class Scalar, int numEq>
template<class SparseMatrixAdapter>
void StandardWellEquations<Scalar,numEq>::
extract(SparseMatrixAdapter& jacobian) const
{
// We need to change matrx A as follows
// A -= C^T D^-1 B
// D is diagonal
// B and C have 1 row, nc colums and nonzero
// at (0,j) only if this well has a perforation at cell j.
typename SparseMatrixAdapter::MatrixBlock tmpMat;
Dune::DynamicMatrix<Scalar> tmp;
for (auto colC = duneC_[0].begin(),
endC = duneC_[0].end(); colC != endC; ++colC)
{
const auto row_index = colC.index();
for (auto colB = duneB_[0].begin(),
endB = duneB_[0].end(); colB != endB; ++colB)
{
detail::multMatrix(invDuneD_[0][0], (*colB), tmp);
detail::negativeMultMatrixTransposed((*colC), tmp, tmpMat);
jacobian.addToBlock(row_index, colB.index(), tmpMat);
}
}
}
#define INSTANCE(N) \
template class StandardWellEquations<double,N>;
template class StandardWellEquations<double,N>; \
template void StandardWellEquations<double,N>::extract(Linear::IstlSparseMatrixAdapter<MatrixBlock<double,N,N>>&) const;
INSTANCE(1)
INSTANCE(2)

View File

@ -94,6 +94,10 @@ public:
void extract(const int numStaticWellEq,
WellContributions& wellContribs) const;
//! \brief Add the matrices of this well to the sparse matrix adapter.
template<class SparseMatrixAdapter>
void extract(SparseMatrixAdapter& jacobian) const;
// two off-diagonal matrices
OffDiagMatWell duneB_;
OffDiagMatWell duneC_;

View File

@ -2129,26 +2129,7 @@ namespace Opm
void
StandardWell<TypeTag>::addWellContributions(SparseMatrixAdapter& jacobian) const
{
// We need to change matrx A as follows
// A -= C^T D^-1 B
// D is diagonal
// B and C have 1 row, nc colums and nonzero
// at (0,j) only if this well has a perforation at cell j.
typename SparseMatrixAdapter::MatrixBlock tmpMat;
Dune::DynamicMatrix<Scalar> tmp;
for (auto colC = this->linSys_.duneC_[0].begin(),
endC = this->linSys_.duneC_[0].end(); colC != endC; ++colC)
{
const auto row_index = colC.index();
for (auto colB = this->linSys_.duneB_[0].begin(),
endB = this->linSys_.duneB_[0].end(); colB != endB; ++colB)
{
detail::multMatrix(this->linSys_.invDuneD_[0][0], (*colB), tmp);
detail::negativeMultMatrixTransposed((*colC), tmp, tmpMat);
jacobian.addToBlock( row_index, colB.index(), tmpMat );
}
}
this->linSys_.extract(jacobian);
}