]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add some 'const' declarations and change the naming of some functions a bit.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 3 Jul 1998 06:50:14 +0000 (06:50 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 3 Jul 1998 06:50:14 +0000 (06:50 +0000)
git-svn-id: https://svn.dealii.org/trunk@433 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/dsmatrix.h
deal.II/lac/source/dsmatrix.cc

index 9291ff083399b7ffb9906629d4168a7aa59884e5..e778499b96f249c8fe3aceacd2ff7b8ccb06d5db 100644 (file)
@@ -593,22 +593,23 @@ class dSMatrix
     double matrix_norm (const dVector &v) const;
     
                                     //
-    double residual (dVector& dst, const dVector& x, const dVector& b);
+    double residual (dVector& dst, const dVector& x,
+                    const dVector& b) const;
                                     //
-    void Jacobi_precond (dVector& dst, const dVector& src,
-                        const double om = 1.);
+    void precondition_Jacobi (dVector& dst, const dVector& src,
+                             const double om = 1.) const;
                                     //
-    void SSOR_precond (dVector& dst, const dVector& src,
-                      const double om = 1.);
+    void precondition_SSOR (dVector& dst, const dVector& src,
+                           const double om = 1.) const;
                                     //
-    void SOR_precond (dVector& dst, const dVector& src,
-                     const double om = 1.);
+    void precondition_SOR (dVector& dst, const dVector& src,
+                          const double om = 1.) const;
                                     //
-    void SSOR (dVector& dst, const double om = 1.);
+    void SSOR (dVector& dst, const double om = 1.) const;
                                     //
-    void SOR (dVector& dst, const double om = 1.);
+    void SOR (dVector& dst, const double om = 1.) const;
                                     //
-    void precondition (dVector& dst, const dVector& src) { dst=src; }
+    void precondition (dVector& dst, const dVector& src) const;
 
                                     /**
                                      * Return a (constant) reference to the
index 77b5b7bf82bc57fa8548ef27e2af14df951d5880..32ce7ee08b3e009cd34f8cfe7b5ac4a9ba129b37 100644 (file)
@@ -539,7 +539,7 @@ dSMatrix::matrix_norm (const dVector& v) const
 
 
 double
-dSMatrix::residual (dVector& dst, const dVector& u, const dVector& b)
+dSMatrix::residual (dVector& dst, const dVector& u, const dVector& b) const
 {
   Assert (cols != 0, ExcMatrixNotInitialized());
   Assert (val != 0, ExcMatrixNotInitialized());
@@ -564,21 +564,27 @@ dSMatrix::residual (dVector& dst, const dVector& u, const dVector& b)
 }
 
 void
-dSMatrix::Jacobi_precond (dVector& dst, const dVector& src, const double om)
+dSMatrix::precondition_Jacobi (dVector& dst, const dVector& src,
+                              const double om) const
 {
   Assert (cols != 0, ExcMatrixNotInitialized());
   Assert (val != 0, ExcMatrixNotInitialized());
+  Assert (m() == n(), ExcMatrixNotSquare());
 
   const unsigned int n = src.size();
 
-  for (unsigned int i=0;i<n;++i)
-    {
-      dst(i) = om * src(i) * val[cols->rowstart[i]];
-    }
-}
+  for (unsigned int i=0; i<n; ++i)
+                                    // note that for square matrices,
+                                    // the diagonal entry is the first
+                                    // in each row, i.e. at index
+                                    // rowstart[i]
+    dst(i) = om * src(i) * val[cols->rowstart[i]];
+};
+
 
 void
-dSMatrix::SSOR_precond (dVector& dst, const dVector& src, const double om)
+dSMatrix::precondition_SSOR (dVector& dst, const dVector& src,
+                            const double om) const
 {
   Assert (cols != 0, ExcMatrixNotInitialized());
   Assert (val != 0, ExcMatrixNotInitialized());
@@ -613,16 +619,23 @@ dSMatrix::SSOR_precond (dVector& dst, const dVector& src, const double om)
 }
 
 void
-dSMatrix::SOR_precond (dVector& dst, const dVector& src, const double om)
+dSMatrix::precondition_SOR (dVector& dst, const dVector& src,
+                           const double om) const
 {
   Assert (cols != 0, ExcMatrixNotInitialized());
   Assert (val != 0, ExcMatrixNotInitialized());
   dst = src;
   SOR(dst,om);
-}
+};
+
+
+void dSMatrix::precondition (dVector &dst, const dVector &src) const {
+  dst=src;
+};
+
 
 void
-dSMatrix::SOR (dVector& dst, const double om)
+dSMatrix::SOR (dVector& dst, const double om) const
 {
   Assert (cols != 0, ExcMatrixNotInitialized());
   Assert (val != 0, ExcMatrixNotInitialized());
@@ -643,7 +656,7 @@ dSMatrix::SOR (dVector& dst, const double om)
 }
 
 void
-dSMatrix::SSOR (dVector& dst, const double om)
+dSMatrix::SSOR (dVector& dst, const double om) const
 {
   Assert (cols != 0, ExcMatrixNotInitialized());
   Assert (val != 0, ExcMatrixNotInitialized());

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.