Skip to content
Snippets Groups Projects
Commit 6e9d5b16 authored by Max Kahnt's avatar Max Kahnt
Browse files

Fix TruncatedBlockGSStep localExactSolve.

Implementation mistakenly changed the solved system for indefinite
matrices with diagonal entries of value 0.

The idea is to also be able to tackle systems that are semi-definite
due to trivial rows by simply not altering the initial x in these components.
This was accomplished by setting the corresponding diagonal entry to 1.

Updated the check for trivial rows such that it is not triggered
if only the diagonal element is zero but only the entire row.
parent 5750cf1b
No related branches found
No related tags found
No related merge requests found
...@@ -293,19 +293,30 @@ public: ...@@ -293,19 +293,30 @@ public:
{ {
for(typename MBlock::size_type i=0; i<A.N(); ++i) 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]) if (ignore[i])
{ {
A[i] = 0.0; 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(); if (rowIsZero)
typename MBlock::row_type::Iterator inner_end = A[i].end(); {
for(; inner_it!=inner_end; ++inner_it) A[i][i] = 1.0;
if (inner_it.index()==i and *inner_it==0.0) b[i] = x[i];
{ }
*inner_it = 1.0;
b[i] = x[i];
}
} }
A.solve(x,b); A.solve(x,b);
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment