From 5a9c473933063f090af65c393c4061172be514ac Mon Sep 17 00:00:00 2001
From: Eric Bechet <eric.bechet@ulg.ac.be>
Date: Thu, 5 Nov 2015 16:26:48 +0000
Subject: [PATCH]
---
Solver/linearSystemPETSc.cpp | 24 ++++++++++++++++++++++++
Solver/linearSystemPETSc.hpp | 27 ++++++++++-----------------
2 files changed, 34 insertions(+), 17 deletions(-)
diff --git a/Solver/linearSystemPETSc.cpp b/Solver/linearSystemPETSc.cpp
index 840e119f7e..e39af9541e 100644
--- a/Solver/linearSystemPETSc.cpp
+++ b/Solver/linearSystemPETSc.cpp
@@ -15,8 +15,32 @@
#include "linearSystemPETSc.hpp"
template class linearSystemPETSc<double>;
+
+
#ifdef PETSC_USE_COMPLEX
template class linearSystemPETSc<std::complex<double> >;
+
+// this specialization will cast to a double (take the real part) if "val" is a double wheras Petsc is built in complex mode.
+template <>
+void linearSystemPETSc<double>::getFromRightHandSide(int row, double &val) const
+{
+ PetscScalar *tmp;
+ _try(VecGetArray(_b, &tmp));
+ PetscScalar s = tmp[row];
+ _try(VecRestoreArray(_b, &tmp));
+ val = s.real();
+}
+
+// this specialization will cast to a double (take the real part) if "val" is a double wheras Petsc is built in complex mode.
+template <>
+void linearSystemPETSc<double>::getFromSolution(int row, double &val) const
+{
+ PetscScalar *tmp;
+ _try(VecGetArray(_x, &tmp));
+ PetscScalar s = tmp[row];
+ _try(VecRestoreArray(_x, &tmp));
+ val = s.real();
+}
#endif
template<>
diff --git a/Solver/linearSystemPETSc.hpp b/Solver/linearSystemPETSc.hpp
index d835b73ecb..63f0c524ba 100644
--- a/Solver/linearSystemPETSc.hpp
+++ b/Solver/linearSystemPETSc.hpp
@@ -266,20 +266,16 @@ void linearSystemPETSc<scalar>::addToRightHandSide(int row, const scalar &val)
PetscScalar s = val;
_try(VecSetValues(_b, 1, &i, &s, ADD_VALUES));
}
+#if defined(PETSC_USE_COMPLEX)
+// this specialization will cast to a double (take the real part) if "val" is a double wheras Petsc is built in complex mode.
+template <>
+void linearSystemPETSc<double>::getFromRightHandSide(int row, double &val) const;
+#endif
template <class scalar>
void linearSystemPETSc<scalar>::getFromRightHandSide(int row, scalar &val) const
{
-#if defined(PETSC_USE_COMPLEX)
- PetscScalar *tmp;
- _try(VecGetArray(_b, &tmp));
- PetscScalar s = tmp[row];
- _try(VecRestoreArray(_b, &tmp));
- // FIXME specialize this routine
- val = s.real();
-#else
_try(VecGetValues(_b, 1, &row, &val));
-#endif
}
template <class scalar>
@@ -310,19 +306,16 @@ void linearSystemPETSc<scalar>::addToSolution(int row, const scalar &val)
PetscScalar s = val;
_try(VecSetValues(_x, 1, &i, &s, ADD_VALUES));
}
+#if defined(PETSC_USE_COMPLEX)
+// this specialization will cast to a double (take the real part) if "val" is a double wheras Petsc is built in complex mode.
+template <>
+void linearSystemPETSc<double>::getFromSolution(int row, double &val) const;
+#endif
template <class scalar>
void linearSystemPETSc<scalar>::getFromSolution(int row, scalar &val) const
{
-#if defined(PETSC_USE_COMPLEX)
- PetscScalar *tmp;
- _try(VecGetArray(_x, &tmp));
- PetscScalar s = tmp[row];
- _try(VecRestoreArray(_x, &tmp));
- val = s.real();
-#else
_try(VecGetValues(_x, 1, &row, &val));
-#endif
}
template <class scalar>
--
GitLab