]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a function computing the norm with respect to a matrix.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 15 Apr 1998 15:55:57 +0000 (15:55 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 15 Apr 1998 15:55:57 +0000 (15:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@178 0785d39b-7218-0410-832d-ea1e28bc413d

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

index d8eef99fd68bafdb92edf7f74a22dbb71a081cd4..b7b5f8979f62d91aba8fd7f8b6bb66c7d7537e31 100644 (file)
@@ -54,7 +54,7 @@ public:
                                     //////////
     void compress ();
                                     //////////
-    int operator() (const unsigned int i, const unsigned int j);
+    int operator() (const unsigned int i, const unsigned int j) const;
                                     //////////
     void add (const unsigned int i, const unsigned int j);
                                     //////////
@@ -83,6 +83,12 @@ public:
                                      * access to the rowstart array, but
                                      * readonly.
                                      *
+                                     * Though the return value is declared
+                                     * #const#, you shoudl be aware that it
+                                     * may change if you call any nonconstant
+                                     * function of objects which operate on
+                                     * it.
+                                     *
                                      * You should use this interface very
                                      * carefully and only if you are absolutely
                                      * sure to know what you do. You should
@@ -100,6 +106,12 @@ public:
                                      * access to the colnums array, but
                                      * readonly.
                                      *
+                                     * Though the return value is declared
+                                     * #const#, you shoudl be aware that it
+                                     * may change if you call any nonconstant
+                                     * function of objects which operate on
+                                     * it.
+                                     *
                                      * You should use this interface very
                                      * carefully and only if you are absolutely
                                      * sure to know what you do. You should
@@ -146,7 +158,7 @@ CLASS
    */
 class dSMatrix
 {
-    dSMatrixStruct * cols;
+    const dSMatrixStruct * cols;
     double* val;
     unsigned int max_len;
   public:
@@ -287,6 +299,30 @@ class dSMatrix
                                     //
     void Tvmult (dVector& dst,const dVector& src) const;
   
+
+                                    /**
+                                     * Return the norm of the vector #v# with
+                                     * respect to the norm induced by this
+                                     * matrix, i.e. $\left<v,Mv\right>$. This
+                                     * is useful, e.g. in the finite element
+                                     * context, where the $L_2$ norm of a
+                                     * function equals the matrix norm with
+                                     * respect to the mass matrix of the vector
+                                     * representing the nodal values of the
+                                     * finite element function.
+                                     *
+                                     * Note the order in which the matrix
+                                     * appears. For non-symmetric matrices
+                                     * there is a difference whether the
+                                     * matrix operates on the first
+                                     * or on the second operand of the
+                                     * scalar product.
+                                     *
+                                     * Obviously, the matrix needs to be square
+                                     * for this operation.
+                                     */
+    double matrix_norm (const dVector &v) const;
+    
                                     //
     double residual (dVector& dst, const dVector& x, const dVector& b);
                                     //
@@ -309,6 +345,12 @@ class dSMatrix
                                      * Return a (constant) reference to the
                                      * underlying sparsity pattern of this
                                      * matrix.
+                                     *
+                                     * Though the return value is declared
+                                     * #const#, you shoudl be aware that it
+                                     * may change if you call any nonconstant
+                                     * function of objects which operate on
+                                     * it.
                                      */
     const dSMatrixStruct & get_sparsity_pattern () const;
 
index be393ca0a1c799bd2147bfbdaf3fb9af8e517909..b387aa2aa30680b75594e2062aced9574396828d 100644 (file)
@@ -232,8 +232,11 @@ dSMatrixStruct::compress ()
 
 
 int
-dSMatrixStruct::operator () (const unsigned int i, const unsigned int j)
+dSMatrixStruct::operator () (const unsigned int i, const unsigned int j) const
 {
+  Assert (i<rows, ExcInvalidIndex(i,rows));
+  Assert (j<cols, ExcInvalidIndex(j,cols));
+
   for (unsigned int k=rowstart[i] ; k<rowstart[i+1] ; k++)
     if (colnums[k] == (signed int)j) return k;
   return -1;
@@ -471,6 +474,29 @@ dSMatrix::Tvmult (dVector& dst, const dVector& src) const
     }
 }
 
+
+
+double
+dSMatrix::matrix_norm (const dVector& v) const
+{
+  Assert (cols != 0, ExcMatrixNotInitialized());
+  Assert(m() == v.n(), ExcDimensionsDontMatch(m(),v.n()));
+  Assert(n() == v.n(), ExcDimensionsDontMatch(n(),v.n()));
+
+  double sum = 0.;
+  for (unsigned int i=0;i<m();i++)
+    {
+      double s = 0.;
+      for (unsigned int j=cols->rowstart[i]; j<cols->rowstart[i+1]; ++j) 
+       s += val[j] * v(cols->colnums[j]);
+      sum += s*v(i);
+    };
+
+  return sum;
+}
+
+
+
 double
 dSMatrix::residual (dVector& dst, const dVector& u, const dVector& b)
 {

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.