From d74d2a1ee0cf7e5b5ca646d68c71fb5f3ab53d4c Mon Sep 17 00:00:00 2001
From: guido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Fri, 6 Oct 2000 10:31:06 +0000
Subject: [PATCH] transposed SOR

git-svn-id: https://svn.dealii.org/trunk@3385 0785d39b-7218-0410-832d-ea1e28bc413d
---
 .../lac/include/lac/sparse_matrix.2.templates |  5 +++
 deal.II/lac/include/lac/sparse_matrix.h       | 19 ++++++++++
 .../lac/include/lac/sparse_matrix.templates.h | 38 +++++++++++++++++++
 3 files changed, 62 insertions(+)

diff --git a/deal.II/lac/include/lac/sparse_matrix.2.templates b/deal.II/lac/include/lac/sparse_matrix.2.templates
index 3e1d6fab19..60134fec6d 100644
--- a/deal.II/lac/include/lac/sparse_matrix.2.templates
+++ b/deal.II/lac/include/lac/sparse_matrix.2.templates
@@ -51,11 +51,16 @@ template void SparseMatrix<TYPEMAT>::precondition_SOR (Vector<TYPE2> &,
 						     const Vector<TYPE2> &,
 						     TYPEMAT) const;
 
+template void SparseMatrix<TYPEMAT>::precondition_TSOR (Vector<TYPE2> &,
+						     const Vector<TYPE2> &,
+						     TYPEMAT) const;
+
 template void SparseMatrix<TYPEMAT>::precondition_Jacobi (Vector<TYPE2> &,
 							const Vector<TYPE2> &,
 							TYPEMAT) const;
 
 template void SparseMatrix<TYPEMAT>::SOR (Vector<TYPE2> &, TYPEMAT) const;
+template void SparseMatrix<TYPEMAT>::TSOR (Vector<TYPE2> &, TYPEMAT) const;
 template void SparseMatrix<TYPEMAT>::SSOR (Vector<TYPE2> &, TYPEMAT) const;
 template void SparseMatrix<TYPEMAT>::SOR_step (Vector<TYPE2> &, const Vector<TYPE2> &, TYPEMAT) const;
 template void SparseMatrix<TYPEMAT>::TSOR_step (Vector<TYPE2> &, const Vector<TYPE2> &, TYPEMAT) const;
diff --git a/deal.II/lac/include/lac/sparse_matrix.h b/deal.II/lac/include/lac/sparse_matrix.h
index adb40a1c8b..9fe43678a1 100644
--- a/deal.II/lac/include/lac/sparse_matrix.h
+++ b/deal.II/lac/include/lac/sparse_matrix.h
@@ -526,6 +526,16 @@ class SparseMatrix : public Subscriptor
 			   const Vector<somenumber> &src,
  			   const number              om = 1.) const;
     
+				     /**
+				      * Apply transpose SOR preconditioning matrix to @p{src}.
+				      * The result of this method is
+				      * $dst = (om D - U)^{-1} src$.
+				      */
+    template <typename somenumber>
+    void precondition_TSOR (Vector<somenumber>       &dst,
+			    const Vector<somenumber> &src,
+			    const number              om = 1.) const;
+    
 				     /**
 				      * Perform SSOR preconditioning
 				      * in-place.  Apply the
@@ -545,6 +555,15 @@ class SparseMatrix : public Subscriptor
 				      */
     template <typename somenumber>
     void SOR (Vector<somenumber> &v,
+	      const number        om = 1.) const;
+
+				     /**
+				      * Perform a transpose SOR preconditioning in-place.
+				      * The result is $v = (\omega D - L)^{-1} v$.
+				      * @p{omega} is the damping parameter.
+				      */
+    template <typename somenumber>
+    void TSOR (Vector<somenumber> &v,
 	      const number        om = 1.) const;
 
 				     /**
diff --git a/deal.II/lac/include/lac/sparse_matrix.templates.h b/deal.II/lac/include/lac/sparse_matrix.templates.h
index 36a98a23ff..40fbc14c59 100644
--- a/deal.II/lac/include/lac/sparse_matrix.templates.h
+++ b/deal.II/lac/include/lac/sparse_matrix.templates.h
@@ -766,6 +766,21 @@ SparseMatrix<number>::precondition_SOR (Vector<somenumber>& dst, const Vector<so
 };
 
 
+template <typename number>
+template <typename somenumber>
+void
+SparseMatrix<number>::precondition_TSOR (Vector<somenumber>& dst, const Vector<somenumber>& src,
+			    const number om) const
+{
+  Assert (cols != 0, ExcMatrixNotInitialized());
+  Assert (val != 0, ExcMatrixNotInitialized());
+  Assert (m() == n(), ExcMatrixNotSquare());
+
+  dst = src;
+  TSOR(dst,om);
+};
+
+
 template <typename number>
 template <typename somenumber>
 void
@@ -788,6 +803,29 @@ SparseMatrix<number>::SOR (Vector<somenumber>& dst, const number om) const
 }
 
 
+template <typename number>
+template <typename somenumber>
+void
+SparseMatrix<number>::TSOR (Vector<somenumber>& dst, const number om) const
+{
+  Assert (cols != 0, ExcMatrixNotInitialized());
+  Assert (val != 0, ExcMatrixNotInitialized());
+  Assert (m() == n(), ExcMatrixNotSquare());
+  Assert (m() == dst.size(), ExcDimensionsDontMatch(m(),dst.size()));
+
+  for (unsigned int row=m(); row!=0;)
+    {
+      --row;
+      somenumber s = dst(row);
+      for (unsigned int j=cols->rowstart[row]; j<cols->rowstart[row+1]; ++j)
+	if (cols->colnums[j] > row)
+	  s -= val[j] * dst(cols->colnums[j]);
+
+      dst(row) = s * om / val[cols->rowstart[row]];
+    }
+}
+
+
 template <typename number>
 template <typename somenumber>
 void
-- 
2.39.5