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

[bugfix] circumvent bug in dune-istl when matrix was assembled in

Browse Source 
row_wise mode. Also cleanup of output.
This commit is contained in:
Robert Kloefkorn
2017-05-18 18:23:55 +02:00
parent e4f1bf9955
commit 1fa04096c3
2 changed files with 15 additions and 7 deletions
Download Patch File Download Diff File Expand all files Collapse all files

1
opm/autodiff/ISTLSolver.hpp
View File

@@ -289,7 +289,6 @@ namespace Opm
} }
} }
//typedef Dune::SeqILU0<Matrix, Vector, Vector> SeqPreconditioner;
typedef ParallelOverlappingILU0<Matrix,Vector,Vector> SeqPreconditioner; typedef ParallelOverlappingILU0<Matrix,Vector,Vector> SeqPreconditioner;
template <class Operator> template <class Operator>

21
opm/autodiff/ParallelOverlappingILU0.hpp
View File

@@ -81,15 +81,25 @@ namespace Opm
template<class M, class CRS, class InvVector> template<class M, class CRS, class InvVector>
void convertToCRS(const M& A, CRS& lower, CRS& upper, InvVector& inv ) void convertToCRS(const M& A, CRS& lower, CRS& upper, InvVector& inv )
{ {
if( A.nonzeroes() == 0 ) { typedef typename M :: size_type size_type;
DUNE_THROW(Dune::ISTLError,"convertToCRS: ILU for matrix with zero entries");
// count non-zeros in case 0 is returned
// this is due to a bug in dune-istl < 2.5
size_type nonZeros = A.nonzeroes();
if( nonZeros == 0 )
{
// count non-zeros
const auto end = A.end();
for( auto row = A.begin(); row != end; ++row )
{
nonZeros += (*row).size();
}
} }
lower.resize( A.N(), A.nonzeroes() ); lower.resize( A.N(), nonZeros );
upper.resize( A.N(), A.nonzeroes() ); upper.resize( A.N(), nonZeros );
inv.resize( A.N() ); inv.resize( A.N() );
typedef typename M :: size_type size_type;
// implement left looking variant with stored inverse // implement left looking variant with stored inverse
const auto endi = A.end(); const auto endi = A.end();
size_type row = 0; size_type row = 0;
@@ -296,7 +306,6 @@ public:
typedef typename Domain::block_type vblock; typedef typename Domain::block_type vblock;
const size_type iEnd = lower_.rows(); const size_type iEnd = lower_.rows();
std::cout << "Rows = " << iEnd << std::endl;
const size_type lastRow = iEnd - 1; const size_type lastRow = iEnd - 1;
if( iEnd != upper_.rows() ) if( iEnd != upper_.rows() )
{ {