]> https://gitweb.dealii.org/ - dealii.git/commitdiff
modify methods description in PETSc sparse matrix for complex numbers
authorDenis Davydov <davydden@gmail.com>
Tue, 16 Jun 2015 21:06:07 +0000 (23:06 +0200)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 17 Jun 2015 19:33:33 +0000 (14:33 -0500)
include/deal.II/lac/petsc_parallel_sparse_matrix.h
source/lac/petsc_parallel_sparse_matrix.cc

index 13499cbf3a51ad37b37dde23436728edfe56b38d..b10c92f6e9165aa4d8e6ca354cc9646a328082ce 100644 (file)
@@ -353,7 +353,7 @@ namespace PETScWrappers
 
       /**
        * Return the square of the norm of the vector $v$ with respect to the
-       * norm induced by this matrix, i.e. $\left(v,Mv\right)$. This is
+       * norm induced by this matrix, i.e. $\left(v^\ast,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
@@ -369,7 +369,7 @@ namespace PETScWrappers
       PetscScalar matrix_norm_square (const Vector &v) const;
 
       /**
-       * Compute the matrix scalar product $\left(u,Mv\right)$.
+       * Compute the matrix scalar product $\left(u^\ast,Mv\right)$.
        *
        * The implementation of this function is not as efficient as the one in
        * the @p MatrixBase class used in deal.II (i.e. the original one, not
index 4e1f2fa63acb1d6a9c011312b9642f364d5abed0..ba9d281a72d440b44733d2f9d66b00926ceddefc 100644 (file)
@@ -858,7 +858,8 @@ namespace PETScWrappers
     {
       Vector tmp (v);
       vmult (tmp, v);
-      return tmp*v;
+      // note, that v*tmp returns  sum_i conjugate(v)_i * tmp_i
+      return v*tmp;
     }
 
     PetscScalar
@@ -867,6 +868,7 @@ namespace PETScWrappers
     {
       Vector tmp (v);
       vmult (tmp, v);
+      // note, that v*tmp returns  sum_i conjugate(v)_i * tmp_i
       return u*tmp;
     }
 

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.