From b167ed2ef6ae44ec3d718721c61cac179f8d0c3c Mon Sep 17 00:00:00 2001
From: kronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Thu, 23 Apr 2009 07:18:49 +0000
Subject: [PATCH] I had introduced some subtle errors in the SSOR
 preconditioner implementation. This should be fixed now.

git-svn-id: https://svn.dealii.org/trunk@18695 0785d39b-7218-0410-832d-ea1e28bc413d
---
 .../lac/include/lac/sparse_matrix.templates.h | 28 +++++++++++--------
 .../include/lac/sparse_matrix_ez.templates.h  | 14 ++++++----
 2 files changed, 24 insertions(+), 18 deletions(-)

diff --git a/deal.II/lac/include/lac/sparse_matrix.templates.h b/deal.II/lac/include/lac/sparse_matrix.templates.h
index 3fe59f41ff..e5874f41e0 100644
--- a/deal.II/lac/include/lac/sparse_matrix.templates.h
+++ b/deal.II/lac/include/lac/sparse_matrix.templates.h
@@ -1198,12 +1198,13 @@ SparseMatrix<number>::precondition_SSOR (Vector<somenumber>              &dst,
 	    pos_right_of_diagonal[row];
 	  Assert (first_right_of_diagonal_index <= *(rowstart_ptr+1),
 		  ExcInternalError());
-	  number s = *dst_ptr;
+	  number s = 0;
 	  for (unsigned int j=(*rowstart_ptr)+1; j<first_right_of_diagonal_index; ++j)
-	    s -= val[j] * dst(cols->colnums[j]);
+	    s += val[j] * dst(cols->colnums[j]);
 
 				   // divide by diagonal element
-	  *dst_ptr = s * om / val[*rowstart_ptr];
+	  *dst_ptr -= s * om;
+	  *dst_ptr /= val[*rowstart_ptr];
 	};
   
       rowstart_ptr = &cols->rowstart[0];
@@ -1219,11 +1220,12 @@ SparseMatrix<number>::precondition_SSOR (Vector<somenumber>              &dst,
 	  const unsigned int end_row = *(rowstart_ptr+1);
 	  const unsigned int first_right_of_diagonal_index
 	    = pos_right_of_diagonal[row];
-	  number s = *dst_ptr;
+	  number s = 0;
 	  for (unsigned int j=first_right_of_diagonal_index; j<end_row; ++j)
-	    s -= val[j] * dst(cols->colnums[j]);
+	    s += val[j] * dst(cols->colnums[j]);
       
-	  *dst_ptr = s * om / val[*rowstart_ptr];
+	  *dst_ptr -= s * om;
+	  *dst_ptr /= val[*rowstart_ptr];
 	};
       return;
     }
@@ -1249,12 +1251,13 @@ SparseMatrix<number>::precondition_SSOR (Vector<somenumber>              &dst,
 	   -
 	   &cols->colnums[0]);
 
-      number s = *dst_ptr;
+      number s = 0;
       for (unsigned int j=(*rowstart_ptr)+1; j<first_right_of_diagonal_index; ++j)
-	s -= val[j] * dst(cols->colnums[j]);
+	s += val[j] * dst(cols->colnums[j]);
 
 				       // divide by diagonal element
-      *dst_ptr = s * om / val[*rowstart_ptr];
+      *dst_ptr -= s * om;
+      *dst_ptr /= val[*rowstart_ptr];
     };
   
   rowstart_ptr = &cols->rowstart[0];
@@ -1273,10 +1276,11 @@ SparseMatrix<number>::precondition_SSOR (Vector<somenumber>              &dst,
 								   &cols->colnums[end_row],
 								   static_cast<unsigned int>(row)) -
 	   &cols->colnums[0]);
-      number s = *dst_ptr;
+      number s = 0;
       for (unsigned int j=first_right_of_diagonal_index; j<end_row; ++j)
-	s -= val[j] * dst(cols->colnums[j]);
-      *dst_ptr = s * om / val[*rowstart_ptr];
+	s += val[j] * dst(cols->colnums[j]);
+      *dst_ptr -= s * om;
+      *dst_ptr /= val[*rowstart_ptr];
     };
 }
 
diff --git a/deal.II/lac/include/lac/sparse_matrix_ez.templates.h b/deal.II/lac/include/lac/sparse_matrix_ez.templates.h
index 9ff6264c6a..64f553d672 100644
--- a/deal.II/lac/include/lac/sparse_matrix_ez.templates.h
+++ b/deal.II/lac/include/lac/sparse_matrix_ez.templates.h
@@ -342,12 +342,13 @@ SparseMatrixEZ<number>::precondition_SSOR (Vector<somenumber>       &dst,
   for (ri = row_info.begin(); ri != end; ++dst_ptr, ++src_ptr, ++ri)
     {
       Assert (ri->diagonal != RowInfo::invalid_diagonal, ExcNoDiagonal());
-      number s = *src_ptr;
+      number s = 0;
       const unsigned int end_row = ri->start + ri->diagonal;
       for (unsigned int i=ri->start;i<end_row;++i)
-	s -= data[i].value * dst(data[i].column);
+	s += data[i].value * dst(data[i].column);
       
-      *dst_ptr = s * om / data[ri->start + ri->diagonal].value;
+      *dst_ptr = *src_ptr - s * om;
+      *dst_ptr /= data[ri->start + ri->diagonal].value;
     }
 				   // Diagonal
   dst_ptr = dst.begin();
@@ -362,11 +363,12 @@ SparseMatrixEZ<number>::precondition_SSOR (Vector<somenumber>       &dst,
    for (rri = row_info.rbegin(); rri != rend; --dst_ptr, ++rri)
      {
       const unsigned int end_row = rri->start + rri->length;
-      number s = *dst_ptr;
+      number s = 0;
       for (unsigned int i=rri->start+rri->diagonal+1;i<end_row;++i)
-	s -= om * data[i].value * dst(data[i].column);
+	s += data[i].value * dst(data[i].column);
 
-      *dst_ptr = s * om / data[rri->start + rri->diagonal].value;
+      *dst_ptr -= s * om;
+      *dst_ptr /= data[rri->start + rri->diagonal].value;
      }
 }
 
-- 
2.39.5