]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add real/imag to PETSc vector. Adjust description of operator* for complex case
authorDenis Davydov <davydden@gmail.com>
Tue, 16 Jun 2015 21:09:05 +0000 (23:09 +0200)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 17 Jun 2015 19:33:33 +0000 (14:33 -0500)
include/deal.II/lac/petsc_vector_base.h
source/lac/petsc_vector_base.cc

index 38412c37abbfb033fe1f3567b1364152029db873..6765f92117841920bbc632bc9910189e70e3110f 100644 (file)
@@ -136,6 +136,19 @@ namespace PETScWrappers
        */
       const VectorReference &operator /= (const PetscScalar &s) const;
 
+      /**
+       * Return the real part of the value of the referenced element.
+       */
+      PetscReal real () const;
+
+      /**
+       * Return the imaginary part of the value of the referenced element.
+       *
+       * @note This operation is not defined for real numbers and an
+       * exception is thrown.
+       */
+      PetscReal imag () const;
+
       /**
        * Convert the reference to an actual value, i.e. return the value of
        * the referenced element of the vector.
@@ -443,6 +456,8 @@ namespace PETScWrappers
     /**
      * Return the scalar product of two vectors. The vectors must have the
      * same size.
+     *
+     * For complex valued vector, this gives$\left(v^\ast,vec\right)$.
      */
     PetscScalar operator * (const VectorBase &vec) const;
 
@@ -1022,10 +1037,35 @@ namespace PETScWrappers
 
       return *this;
     }
-  }
 
 
 
+    inline
+    PetscReal
+    VectorReference::real () const
+    {
+#ifndef PETSC_USE_COMPLEX
+      return static_cast<PetscScalar>(*this);
+#else
+      return PetscRealPart (static_cast<PetscScalar>(*this));
+#endif
+    }
+
+
+
+    inline
+    PetscReal
+    VectorReference::imag () const
+    {
+#ifndef PETSC_USE_COMPLEX
+      return PetscReal (0);
+#else
+      return PetscImaginaryPart (static_cast<PetscScalar>(*this));
+#endif
+    }
+
+  } // namespace internal
+
   inline
   bool
   VectorBase::in_local_range (const size_type index) const
index 0feeaef91a33185c0530bf845dd2ced36341a423..e3191cd202d708f06d81267a0e0a90672d503dbb 100644 (file)
@@ -389,7 +389,11 @@ namespace PETScWrappers
 
     PetscScalar result;
 
-    const int ierr = VecDot (vector, vec.vector, &result);
+    //For complex vectors, VecDot() computes
+    //    val = (x,y) = y^H x,
+    //where y^H denotes the conjugate transpose of y.
+    //Note that this corresponds to the usual "mathematicians" complex inner product where the SECOND argument gets the complex conjugate.
+    const int ierr = VecDot (vec.vector, vector, &result);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
     return result;

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.