diff --git a/dune/solvers/iterationsteps/truncatedblockgsstep.hh b/dune/solvers/iterationsteps/truncatedblockgsstep.hh index 6069091bb5ea45c71260d51970fe5114cbb292e3..6d93b0570aafb8bca3579921ed3bdf32d0cc5386 100644 --- a/dune/solvers/iterationsteps/truncatedblockgsstep.hh +++ b/dune/solvers/iterationsteps/truncatedblockgsstep.hh @@ -293,19 +293,30 @@ public: { for(typename MBlock::size_type i=0; i<A.N(); ++i) { + assert(A.exists(i,i)); // No implementation for matrices w/o explicit diagonal elements is available here. + + bool rowIsZero = true; if (ignore[i]) { A[i] = 0.0; } + else + { + typename MBlock::row_type::Iterator inner_it = A[i].begin(); + typename MBlock::row_type::Iterator inner_end = A[i].end(); + for(; inner_it!=inner_end; ++inner_it) + if (*inner_it != 0.0) + { + rowIsZero = false; + break; + } + } - typename MBlock::row_type::Iterator inner_it = A[i].begin(); - typename MBlock::row_type::Iterator inner_end = A[i].end(); - for(; inner_it!=inner_end; ++inner_it) - if (inner_it.index()==i and *inner_it==0.0) - { - *inner_it = 1.0; - b[i] = x[i]; - } + if (rowIsZero) + { + A[i][i] = 1.0; + b[i] = x[i]; + } } A.solve(x,b);