From: Wolfgang Bangerth Date: Tue, 21 Jul 1998 08:53:37 +0000 (+0000) Subject: Implement a matrix scalar product for the dFMatrix class. X-Git-Tag: v8.0.0~22804 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d23ba50998419f852e7d81c4a4a685eab53eeb9f;p=dealii.git Implement a matrix scalar product for the dFMatrix class. git-svn-id: https://svn.dealii.org/trunk@450 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/dfmatrix.h b/deal.II/lac/include/lac/dfmatrix.h index c0254a2ed6..104513ee17 100644 --- a/deal.II/lac/include/lac/dfmatrix.h +++ b/deal.II/lac/include/lac/dfmatrix.h @@ -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$. 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 diff --git a/deal.II/lac/source/dfmatrix.cc b/deal.II/lac/source/dfmatrix.cc index 2f50134288..2b95838794 100644 --- a/deal.II/lac/source/dfmatrix.cc +++ b/deal.II/lac/source/dfmatrix.cc @@ -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