From 93ff97c7a33166b4be2755ff16e98532dfa59a2c Mon Sep 17 00:00:00 2001
From: Christoph Lehmann <christoph.lehmann@ufz.de>
Date: Sun, 10 Apr 2016 15:33:04 +0200
Subject: [PATCH] [MaL] added Manhattan and max norm to BLAS
---
MathLib/LinAlg/BLAS.cpp | 48 +++++++++++++++++++++++++++++++++++++++++
MathLib/LinAlg/BLAS.h | 8 +++++++
2 files changed, 56 insertions(+)
diff --git a/MathLib/LinAlg/BLAS.cpp b/MathLib/LinAlg/BLAS.cpp
index caf6b3326c6..bffbbe9a722 100644
--- a/MathLib/LinAlg/BLAS.cpp
+++ b/MathLib/LinAlg/BLAS.cpp
@@ -19,6 +19,14 @@
namespace MathLib { namespace BLAS
{
+// Explicit specialization
+// Computes the Manhattan norm of x
+template<>
+double norm1(Eigen::VectorXd const& x)
+{
+ return x.lpNorm<1>();
+}
+
// Explicit specialization
// Computes the Euclidean norm of x
template<>
@@ -27,6 +35,14 @@ double norm2(Eigen::VectorXd const& x)
return x.norm();
}
+// Explicit specialization
+// Computes the Maximum norm of x
+template<>
+double normMax(Eigen::VectorXd const& x)
+{
+ return x.lpNorm<Eigen::Infinity>();
+}
+
} } // namespaces
#endif
@@ -74,6 +90,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 Manhattan norm of x
+template<>
+double norm1(PETScVector const& x)
+{
+ return x.getNorm(MathLib::VecNormType::NORM1);
+}
+
// Explicit specialization
// Computes the Euclidean norm of x
template<>
@@ -82,6 +106,14 @@ double norm2(PETScVector const& x)
return x.getNorm(MathLib::VecNormType::NORM2);
}
+// Explicit specialization
+// Computes the Maximum norm of x
+template<>
+double normMax(PETScVector const& x)
+{
+ return x.getNorm(MathLib::VecNormType::INFINITY_N);
+}
+
// Matrix
@@ -192,6 +224,14 @@ void axpby(EigenVector& y, double const a, double const b, EigenVector const& x)
y.getRawVector() = a*x.getRawVector() + b*y.getRawVector();
}
+// Explicit specialization
+// Computes the Manhattan norm of x
+template<>
+double norm1(EigenVector const& x)
+{
+ return x.getRawVector().lpNorm<1>();
+}
+
// Explicit specialization
// Euclidean norm
template<>
@@ -200,6 +240,14 @@ double norm2(EigenVector const& x)
return x.getRawVector().norm();
}
+// Explicit specialization
+// Computes the Maximum norm of x
+template<>
+double normMax(EigenVector const& x)
+{
+ return x.getRawVector().lpNorm<Eigen::Infinity>();
+}
+
// Matrix
diff --git a/MathLib/LinAlg/BLAS.h b/MathLib/LinAlg/BLAS.h
index 58755c9892c..355af434f79 100644
--- a/MathLib/LinAlg/BLAS.h
+++ b/MathLib/LinAlg/BLAS.h
@@ -66,10 +66,18 @@ void axpby(MatrixOrVector& y, double const a, double const b, MatrixOrVector con
y = a*x + b*y;
}
+//! Computes the Manhattan norm of \c x.
+template<typename MatrixOrVector>
+double norm1(MatrixOrVector const& x);
+
//! Computes the Euclidean norm of \c x.
template<typename MatrixOrVector>
double norm2(MatrixOrVector const& x);
+//! Computes the maximum norm of \c x.
+template<typename MatrixOrVector>
+double normMax(MatrixOrVector const& x);
+
template<typename Matrix>
void finalizeAssembly(Matrix& /*A*/)
{
--
GitLab