OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2024-11-22 09:16:27 -06:00
Code Issues Packages Projects Releases Wiki Activity

addressing the first part of comments form PR#1680

Browse Source 
there should be no functional change.
This commit is contained in:
Kai Bao 2018-12-13 13:29:59 +01:00
parent 1ddcfcb679
commit c4254240e2
6 changed files with 39 additions and 34 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
opm/autodiff/BlackoilModelEbos.hpp
View File

@ -775,7 +775,7 @@ namespace Opm {
assert(has_polymer_); assert(has_polymer_);
B_avg[contiPolymerMWEqIdx] += 1.0 / fs.invB(FluidSystem::waterPhaseIdx).value(); B_avg[contiPolymerMWEqIdx] += 1.0 / fs.invB(FluidSystem::waterPhaseIdx).value();
// the residual of the polymer molecular equatinon is scaled down by a 100, since molecular weight // the residual of the polymer molecular equation is scaled down by a 100, since molecular weight
// can be much bigger than 1, and this equation shares the same tolerance with other mass balance equations // can be much bigger than 1, and this equation shares the same tolerance with other mass balance equations
// TODO: there should be a more general way to determine the scaling-down coefficient // TODO: there should be a more general way to determine the scaling-down coefficient
const auto R2 = ebosResid[cell_idx][contiPolymerMWEqIdx] / 100.; const auto R2 = ebosResid[cell_idx][contiPolymerMWEqIdx] / 100.;

47
opm/autodiff/MatrixBlock.hpp
View File

@ -32,7 +32,7 @@ namespace Dune
namespace FMatrixHelp { namespace FMatrixHelp {
//! invert 4x4 Matrix without changing the original matrix //! invert 4x4 Matrix without changing the original matrix
template <typename K> template <typename K>
static inline K invertMatrix (const FieldMatrix<K,4,4> &matrix, FieldMatrix<K,4,4> &inverse) static inline K invertMatrix(const FieldMatrix<K,4,4>& matrix, FieldMatrix<K,4,4>& inverse)
{ {
inverse[0][0] = matrix[1][1] * matrix[2][2] * matrix[3][3] - inverse[0][0] = matrix[1][1] * matrix[2][2] * matrix[3][3] -
matrix[1][1] * matrix[2][3] * matrix[3][2] - matrix[1][1] * matrix[2][3] * matrix[3][2] -
@ -146,7 +146,8 @@ static inline K invertMatrix (const FieldMatrix<K,4,4> &matrix, FieldMatrix<K,4,
matrix[2][0] * matrix[0][1] * matrix[1][2] - matrix[2][0] * matrix[0][1] * matrix[1][2] -
matrix[2][0] * matrix[0][2] * matrix[1][1]; matrix[2][0] * matrix[0][2] * matrix[1][1];
K det = matrix[0][0] * inverse[0][0] + matrix[0][1] * inverse[1][0] + matrix[0][2] * inverse[2][0] + matrix[0][3] * inverse[3][0]; K det = matrix[0][0] * inverse[0][0] + matrix[0][1] * inverse[1][0] +
matrix[0][2] * inverse[2][0] + matrix[0][3] * inverse[3][0];
// return identity for singular or nearly singular matrices. // return identity for singular or nearly singular matrices.
if (std::abs(det) < 1e-40) { if (std::abs(det) < 1e-40) {
@ -166,7 +167,7 @@ namespace ISTLUtility {
//! invert matrix by calling FMatrixHelp::invert //! invert matrix by calling FMatrixHelp::invert
template <typename K> template <typename K>
static inline void invertMatrix (FieldMatrix<K,1,1> &matrix) static inline void invertMatrix(FieldMatrix<K,1,1>& matrix)
{ {
FieldMatrix<K,1,1> A ( matrix ); FieldMatrix<K,1,1> A ( matrix );
FMatrixHelp::invertMatrix(A, matrix ); FMatrixHelp::invertMatrix(A, matrix );
@ -174,7 +175,7 @@ static inline void invertMatrix (FieldMatrix<K,1,1> &matrix)
//! invert matrix by calling FMatrixHelp::invert //! invert matrix by calling FMatrixHelp::invert
template <typename K> template <typename K>
static inline void invertMatrix (FieldMatrix<K,2,2> &matrix) static inline void invertMatrix(FieldMatrix<K,2,2>& matrix)
{ {
FieldMatrix<K,2,2> A ( matrix ); FieldMatrix<K,2,2> A ( matrix );
FMatrixHelp::invertMatrix(A, matrix ); FMatrixHelp::invertMatrix(A, matrix );
@ -182,7 +183,7 @@ static inline void invertMatrix (FieldMatrix<K,2,2> &matrix)
//! invert matrix by calling FMatrixHelp::invert //! invert matrix by calling FMatrixHelp::invert
template <typename K> template <typename K>
static inline void invertMatrix (FieldMatrix<K,3,3> &matrix) static inline void invertMatrix(FieldMatrix<K,3,3>& matrix)
{ {
FieldMatrix<K,3,3> A ( matrix ); FieldMatrix<K,3,3> A ( matrix );
FMatrixHelp::invertMatrix(A, matrix ); FMatrixHelp::invertMatrix(A, matrix );
@ -190,7 +191,7 @@ static inline void invertMatrix (FieldMatrix<K,3,3> &matrix)
//! invert matrix by calling FMatrixHelp::invert //! invert matrix by calling FMatrixHelp::invert
template <typename K> template <typename K>
static inline void invertMatrix (FieldMatrix<K,4,4> &matrix) static inline void invertMatrix(FieldMatrix<K,4,4>& matrix)
{ {
FieldMatrix<K,4,4> A ( matrix ); FieldMatrix<K,4,4> A ( matrix );
FMatrixHelp::invertMatrix(A, matrix ); FMatrixHelp::invertMatrix(A, matrix );
@ -198,7 +199,7 @@ static inline void invertMatrix (FieldMatrix<K,4,4> &matrix)
//! invert matrix by calling matrix.invert //! invert matrix by calling matrix.invert
template <typename K, int n> template <typename K, int n>
static inline void invertMatrix (FieldMatrix<K,n,n> &matrix) static inline void invertMatrix(FieldMatrix<K,n,n>& matrix)
{ {
#if ! DUNE_VERSION_NEWER( DUNE_COMMON, 2, 7 ) #if ! DUNE_VERSION_NEWER( DUNE_COMMON, 2, 7 )
Dune::FMatrixPrecision<K>::set_singular_limit(1.e-20); Dune::FMatrixPrecision<K>::set_singular_limit(1.e-20);
@ -208,7 +209,7 @@ static inline void invertMatrix (FieldMatrix<K,n,n> &matrix)
//! invert matrix by calling matrix.invert //! invert matrix by calling matrix.invert
template <typename K> template <typename K>
static inline void invertMatrix (Dune::DynamicMatrix<K> &matrix) static inline void invertMatrix(Dune::DynamicMatrix<K>& matrix)
{ {
#if ! DUNE_VERSION_NEWER( DUNE_COMMON, 2, 7 ) #if ! DUNE_VERSION_NEWER( DUNE_COMMON, 2, 7 )
Dune::FMatrixPrecision<K>::set_singular_limit(1.e-30); Dune::FMatrixPrecision<K>::set_singular_limit(1.e-30);
@ -238,17 +239,17 @@ public:
template<class K, int n, int m> template<class K, int n, int m>
void void
print_row (std::ostream& s, const MatrixBlock<K,n,m>& A, print_row(std::ostream& s, const MatrixBlock<K,n,m>& A,
typename FieldMatrix<K,n,m>::size_type I, typename FieldMatrix<K,n,m>::size_type I,
typename FieldMatrix<K,n,m>::size_type J, typename FieldMatrix<K,n,m>::size_type J,
typename FieldMatrix<K,n,m>::size_type therow, int width, typename FieldMatrix<K,n,m>::size_type therow, int width,
int precision) int precision)
{ {
print_row(s, A.asBase(), I, J, therow, width, precision); print_row(s, A.asBase(), I, J, therow, width, precision);
} }
template<class K, int n, int m> template<class K, int n, int m>
K& firstmatrixelement (MatrixBlock<K,n,m>& A) K& firstmatrixelement(MatrixBlock<K,n,m>& A)
{ {
return firstmatrixelement( A.asBase() ); return firstmatrixelement( A.asBase() );
} }
@ -313,9 +314,9 @@ namespace Detail
{ {
//! calculates ret = A^T * B //! calculates ret = A^T * B
template< class K, int m, int n, int p > template< class K, int m, int n, int p >
static inline void multMatrixTransposed ( const Dune::FieldMatrix< K, n, m > &A, static inline void multMatrixTransposed(const Dune::FieldMatrix< K, n, m >& A,
const Dune::FieldMatrix< K, n, p > &B, const Dune::FieldMatrix< K, n, p >& B,
Dune::FieldMatrix< K, m, p > &ret ) Dune::FieldMatrix< K, m, p >& ret)
{ {
typedef typename Dune::FieldMatrix< K, m, p > :: size_type size_type; typedef typename Dune::FieldMatrix< K, m, p > :: size_type size_type;
@ -332,9 +333,9 @@ namespace Detail
//! calculates ret = A * B //! calculates ret = A * B
template< class K> template< class K>
static inline void multMatrix( const Dune::DynamicMatrix<K> &A, static inline void multMatrix(const Dune::DynamicMatrix<K>& A,
const Dune::DynamicMatrix<K> &B, const Dune::DynamicMatrix<K>& B,
Dune::DynamicMatrix<K> &ret ) Dune::DynamicMatrix<K>& ret )
{ {
typedef typename Dune::DynamicMatrix<K> :: size_type size_type; typedef typename Dune::DynamicMatrix<K> :: size_type size_type;
@ -361,9 +362,9 @@ namespace Detail
//! calculates ret = A^T * B //! calculates ret = A^T * B
template< class K, int m, int p > template< class K, int m, int p >
static inline void multMatrixTransposed ( const Dune::DynamicMatrix<K> &A, static inline void multMatrixTransposed(const Dune::DynamicMatrix<K>& A,
const Dune::DynamicMatrix<K> &B, const Dune::DynamicMatrix<K>& B,
Dune::FieldMatrix< K, m, p> &ret ) Dune::FieldMatrix< K, m, p>& ret )
{ {
typedef typename Dune::DynamicMatrix<K> :: size_type size_type; typedef typename Dune::DynamicMatrix<K> :: size_type size_type;

2
opm/autodiff/MultisegmentWell_impl.hpp
View File

@ -1911,7 +1911,7 @@ namespace Opm
template<typename TypeTag> template<typename TypeTag>
void void
MultisegmentWell<TypeTag>:: MultisegmentWell<TypeTag>::
updateWaterThroughput(const double dt, WellState &well_state) const updateWaterThroughput(const double dt OPM_UNUSED, WellState& well_state OPM_UNUSED) const
{ {
} }

10
opm/autodiff/StandardWell.hpp
View File

@ -89,8 +89,8 @@ namespace Opm
// TODO: we should have indices for the well equations and well primary variables separately // TODO: we should have indices for the well equations and well primary variables separately
static const int Bhp = numStaticWellEq - numWellControlEq; static const int Bhp = numStaticWellEq - numWellControlEq;
// total number of the welll equations and primary variables // total number of the well equations and primary variables
// for StandardWell, there is no extra well equations will be used // for StandardWell, no extra well equations will be used.
static const int numWellEq = numStaticWellEq; static const int numWellEq = numStaticWellEq;
using typename Base::Scalar; using typename Base::Scalar;
@ -312,8 +312,10 @@ namespace Opm
const std::vector<EvalWell>& mob, const std::vector<EvalWell>& mob,
const EvalWell& bhp, const EvalWell& bhp,
const int perf, const int perf,
const bool allow_cf, std::vector<EvalWell>& cq_s, const bool allow_cf,
double& perf_dis_gas_rate, double& perf_vap_oil_rate) const; std::vector<EvalWell>& cq_s,
double& perf_dis_gas_rate,
double& perf_vap_oil_rate) const;
// TODO: maybe we should provide a light version of computePerfRate, which does not include the // TODO: maybe we should provide a light version of computePerfRate, which does not include the
// calculation of the derivatives // calculation of the derivatives

8
opm/autodiff/StandardWellV.hpp
View File

@ -94,7 +94,7 @@ namespace Opm
// TODO: we should have indices for the well equations and well primary variables separately // TODO: we should have indices for the well equations and well primary variables separately
static const int Bhp = numStaticWellEq - numWellControlEq; static const int Bhp = numStaticWellEq - numWellControlEq;
// total number of the welll equations and primary variables // total number of the well equations and primary variables
// there might be extra equations be used, numWellEq will be updated during the initialization // there might be extra equations be used, numWellEq will be updated during the initialization
int numWellEq = numStaticWellEq; int numWellEq = numStaticWellEq;
@ -317,8 +317,10 @@ namespace Opm
const std::vector<EvalWell>& mob, const std::vector<EvalWell>& mob,
const EvalWell& bhp, const EvalWell& bhp,
const int perf, const int perf,
const bool allow_cf, std::vector<EvalWell>& cq_s, const bool allow_cf,
double& perf_dis_gas_rate, double& perf_vap_oil_rate) const; std::vector<EvalWell>& cq_s,
double& perf_dis_gas_rate,
double& perf_vap_oil_rate) const;
// TODO: maybe we should provide a light version of computePerfRate, which does not include the // TODO: maybe we should provide a light version of computePerfRate, which does not include the
// calculation of the derivatives // calculation of the derivatives

4
opm/autodiff/StandardWell_impl.hpp
View File

@ -357,7 +357,7 @@ namespace Opm
const EvalWell cqt_i = - Tw * (total_mob_dense * drawdown); const EvalWell cqt_i = - Tw * (total_mob_dense * drawdown);
// surface volume fraction of fluids within wellbore // surface volume fraction of fluids within wellbore
std::vector<EvalWell> cmix_s(num_components_, EvalWell{numWellEq + numEq}); std::vector<EvalWell> cmix_s(num_components_, 0.);
for (int componentIdx = 0; componentIdx < num_components_; ++componentIdx) { for (int componentIdx = 0; componentIdx < num_components_; ++componentIdx) {
cmix_s[componentIdx] = wellSurfaceVolumeFraction(componentIdx); cmix_s[componentIdx] = wellSurfaceVolumeFraction(componentIdx);
} }
@ -2807,7 +2807,7 @@ namespace Opm
template<typename TypeTag> template<typename TypeTag>
void void
StandardWell<TypeTag>:: StandardWell<TypeTag>::
updateWaterThroughput(const double dt, WellState &well_state) const updateWaterThroughput(const double dt OPM_UNUSED, WellState& well_state OPM_UNUSED) const
{ {
} }
} }