From 30f202279f595b85e71182e1ed6a4c18f6920b53 Mon Sep 17 00:00:00 2001
From: Christoph Lehmann <christoph.lehmann@ufz.de>
Date: Tue, 1 Mar 2016 11:23:26 +0100
Subject: [PATCH] [BL] added Euclidean norm
---
MathLib/LinAlg/BLAS.cpp | 37 +++++++++++++++++++++++++++++++++++++
MathLib/LinAlg/BLAS.h | 4 ++++
2 files changed, 41 insertions(+)
diff --git a/MathLib/LinAlg/BLAS.cpp b/MathLib/LinAlg/BLAS.cpp
index 93813d36138..bf2d8e3aba5 100644
--- a/MathLib/LinAlg/BLAS.cpp
+++ b/MathLib/LinAlg/BLAS.cpp
@@ -10,6 +10,27 @@
#include "BLAS.h"
+// Dense Eigen matrices and vectors ////////////////////////////////////////
+#ifdef OGS_USE_EIGEN
+
+#include <Eigen/Core>
+
+namespace MathLib { namespace BLAS
+{
+
+// Explicit specialization
+// Computes the Euclidean norm of x
+template<>
+double norm2(Eigen::VectorXd const& x)
+{
+ return x.norm();
+}
+
+} } // namespaces
+
+#endif
+
+
// Global PETScMatrix/PETScVector //////////////////////////////////////////
#ifdef USE_PETSC
@@ -52,6 +73,14 @@ void axpby(PETScVector& y, double const a, double const b, PETScVector const& x)
VecAXPBY(y.getRawVector(), a, b, x.getRawVector());
}
+// Explicit specialization
+// Computes the Euclidean norm of x
+template<>
+double norm2(PETScVector const& x)
+{
+ return x.getNorm(MathLib::VecNormType::NORM2);
+}
+
// Matrix
@@ -150,6 +179,14 @@ void axpby(EigenVector& y, double const a, double const b, EigenVector const& x)
y.getRawVector() = a*x.getRawVector() + b*y.getRawVector();
}
+// Explicit specialization
+// Euclidean norm
+template<>
+double norm2(EigenVector const& x)
+{
+ return x.getRawVector().norm();
+}
+
// Matrix
diff --git a/MathLib/LinAlg/BLAS.h b/MathLib/LinAlg/BLAS.h
index b8b96dd8d70..3c3623b2a8a 100644
--- a/MathLib/LinAlg/BLAS.h
+++ b/MathLib/LinAlg/BLAS.h
@@ -66,6 +66,10 @@ void axpby(MatrixOrVector& y, double const a, double const b, MatrixOrVector con
y = a*x + b*y;
}
+//! Computes the Euclidean norm of \c x.
+template<typename MatrixOrVector>
+double norm2(MatrixOrVector const& x);
+
// Matrix and Vector
--
GitLab