]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement a matrix scalar product for the dFMatrix class.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 21 Jul 1998 08:53:37 +0000 (08:53 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 21 Jul 1998 08:53:37 +0000 (08:53 +0000)
git-svn-id: https://svn.dealii.org/trunk@450 0785d39b-7218-0410-832d-ea1e28bc413d

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

index c0254a2ed6c6c341fd6e72767433332a5e6f18d7..104513ee17a37feb2256d0e768069e42341c67e6 100644 (file)
@@ -237,6 +237,38 @@ class dFMatrix
                                      */
     void Tvmult (dVector& w, const dVector& v, const bool adding=false) 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;
+
+                                    /**
+                                     * Build the matrix scalar product
+                                     * #u^T M v#. This function is mostly
+                                     * useful when building the cellwise
+                                     * scalar product of two functions in
+                                     * the finite element context.
+                                     */
+    double matrix_scalar_product (const dVector &u, const dVector &v) const;
+    
                                     /**
                                      * A=Inverse(A). Inversion of this by
                                      * Gauss-Jordan-algorithm
index 2f501342881ddbec8562ad421d687fb2edc40c8d..2b958387947416e47a15697b63362b72827cd63a 100644 (file)
@@ -571,6 +571,59 @@ void dFMatrix::Tmmult (dFMatrix& dst, const dFMatrix& src) const
       }
 }*/
 
+
+
+double dFMatrix::matrix_norm (const dVector &v) const {
+  Assert(m() == v.size(), ExcDimensionMismatch(m(),v.size()));
+  Assert(n() == v.size(), ExcDimensionMismatch(n(),v.size()));
+
+  double sum = 0.;
+  const unsigned int n_rows = m();
+  const double *val_ptr = &val[0];
+  const double *v_ptr;
+  
+  for (unsigned int row=0; row<n_rows; ++row)
+    {
+      double s = 0.;
+      const double * const val_end_of_row = val_ptr+n_rows;
+      v_ptr = v.begin();
+      while (val_ptr != val_end_of_row)
+       s += *val_ptr++ * *v_ptr++;
+
+      sum += s* v(row);
+    };
+
+  return sum;
+};
+
+
+
+double dFMatrix::matrix_scalar_product (const dVector &u, const dVector &v) const {
+  Assert(m() == u.size(), ExcDimensionMismatch(m(),v.size()));
+  Assert(n() == v.size(), ExcDimensionMismatch(n(),v.size()));
+
+  double sum = 0.;
+  const unsigned int n_rows = m();
+  const unsigned int n_cols = n();
+  const double *val_ptr = &val[0];
+  const double *v_ptr;
+  
+  for (unsigned int row=0; row<n_rows; ++row)
+    {
+      double s = 0.;
+      const double * const val_end_of_row = val_ptr+n_cols;
+      v_ptr = v.begin();
+      while (val_ptr != val_end_of_row)
+       s += *val_ptr++ * *v_ptr++;
+
+      sum += s* u(row);
+    };
+
+  return sum;
+};
+
+
+
 void dFMatrix::print (FILE* f, const char* format) const
 {
   if (!format) format = " %5.2f";

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.