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

Merge pull request #9 from andlaus/improve_cpr_mod

Browse Source 
make the branch self contained
This commit is contained in:
Atgeirr Flø Rasmussen 2019-03-20 10:21:34 +01:00 committed by GitHub
commit dc1135324f
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
1 changed files with 13 additions and 2 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

15
opm/autodiff/ISTLSolverEbos.hpp
View File

@ -61,6 +61,17 @@ END_PROPERTIES
namespace Opm namespace Opm
{ {
template <class DenseMatrix>
DenseMatrix transposeDenseMatrix(const DenseMatrix& M)
{
DenseMatrix tmp;
for (int i = 0; i < M.rows; ++i)
for (int j = 0; j < M.cols; ++j)
tmp[j][i] = M[i][j];
return tmp;
}
//===================================================================== //=====================================================================
// Implementation for ISTL-matrix based operator // Implementation for ISTL-matrix based operator
//===================================================================== //=====================================================================
@ -662,7 +673,7 @@ protected:
} }
} }
BlockVector bweights; BlockVector bweights;
MatrixBlockType block_transpose = block.transpose(); MatrixBlockType block_transpose = Opm::transposeDenseMatrix(block);
block_transpose.solve(bweights, rhs); block_transpose.solve(bweights, rhs);
bweights /= 1000.0; // given normal densities this scales weights to about 1. bweights /= 1000.0; // given normal densities this scales weights to about 1.
weights[index] = bweights; weights[index] = bweights;
@ -731,7 +742,7 @@ protected:
} }
} }
BlockVector bweights; BlockVector bweights;
auto diag_block_transpose = diag_block.transpose(); auto diag_block_transpose = Opm::transposeDenseMatrix(diag_block);
diag_block_transpose.solve(bweights, rhs); diag_block_transpose.solve(bweights, rhs);
double abs_max = double abs_max =
*std::max_element(bweights.begin(), bweights.end(), [](double a, double b){ return std::abs(a) < std::abs(b); } ); *std::max_element(bweights.begin(), bweights.end(), [](double a, double b){ return std::abs(a) < std::abs(b); } );