OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

added: StandardWellEquations::apply(x,Ax)

Browse Source 
this applies the equation system to a vector.
use the new method in the well implementation.
This commit is contained in:
Arne Morten Kvarving 2022-11-11 21:41:24 +01:00
parent 97d158da38
commit 21bb4dc955
3 changed files with 23 additions and 13 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

19
opm/simulators/wells/StandardWellEquations.cpp
View File

@ -111,6 +111,25 @@ void StandardWellEquations<Scalar,numEq>::clear()
resWell_ = 0.0;
}
template<class Scalar, int numEq>
void StandardWellEquations<Scalar,numEq>::apply(const BVector& x, BVector& Ax) const
{
assert(Bx_.size() == duneB_.N());
assert(invDrw_.size() == invDuneD_.N());
// Bx_ = duneB_ * x
parallelB_.mv(x, Bx_);
// invDBx = invDuneD_ * Bx_
// TODO: with this, we modified the content of the invDrw_.
// Is it necessary to do this to save some memory?
auto& invDBx = invDrw_;
invDuneD_.mv(Bx_, invDBx);
// Ax = Ax - duneC_^T * invDBx
duneC_.mmtv(invDBx, Ax);
}
#define INSTANCE(N) \
template class StandardWellEquations<double,N>;

3
opm/simulators/wells/StandardWellEquations.hpp
View File

@ -73,6 +73,9 @@ public:
//! \brief Set all coefficients to 0.
void clear();
//! \brief Apply linear operator to vector.
void apply(const BVector& x, BVector& Ax) const;
// two off-diagonal matrices
OffDiagMatWell duneB_;
OffDiagMatWell duneC_;

14
opm/simulators/wells/StandardWell_impl.hpp
View File

@ -1722,20 +1722,8 @@ namespace Opm
// Contributions are already in the matrix itself
return;
}
assert(this->linSys_.Bx_.size() == this->linSys_.duneB_.N());
assert(this->linSys_.invDrw_.size() == this->linSys_.invDuneD_.N());
// Bx_ = duneB_ * x
this->linSys_.parallelB_.mv(x, this->linSys_.Bx_);
// invDBx = invDuneD_ * Bx_
// TODO: with this, we modified the content of the invDrw_.
// Is it necessary to do this to save some memory?
BVectorWell& invDBx = this->linSys_.invDrw_;
this->linSys_.invDuneD_.mv(this->linSys_.Bx_, invDBx);
// Ax = Ax - duneC_^T * invDBx
this->linSys_.duneC_.mmtv(invDBx,Ax);
this->linSys_.apply(x, Ax);
}