diff --git a/src/algebra/SmallVector.hpp b/src/algebra/SmallVector.hpp
index d20bb48f4cdb700f61a0881898fb5273217693aa..d4a0bf84268bdf1d164cf6216153117adafaab56 100644
--- a/src/algebra/SmallVector.hpp
+++ b/src/algebra/SmallVector.hpp
@@ -30,6 +30,18 @@ class SmallVector // LCOV_EXCL_LINE
return os << x.m_values;
}
+ [[nodiscard]] friend std::remove_const_t<DataType>
+ min(const SmallVector& x)
+ {
+ return min(x.m_values);
+ }
+
+ [[nodiscard]] friend std::remove_const_t<DataType>
+ max(const SmallVector& x)
+ {
+ return max(x.m_values);
+ }
+
friend SmallVector<std::remove_const_t<DataType>>
copy(const SmallVector& x)
{
@@ -40,7 +52,18 @@ class SmallVector // LCOV_EXCL_LINE
return x_copy;
}
- friend SmallVector operator*(const DataType& a, const SmallVector& x)
+ PUGS_INLINE
+ SmallVector<std::remove_const_t<DataType>>
+ operator-() const
+ {
+ SmallVector<std::remove_const_t<DataType>> opposite(this->size());
+ parallel_for(
+ opposite.size(), PUGS_LAMBDA(index_type i) { opposite.m_values[i] = -m_values[i]; });
+ return opposite;
+ }
+
+ friend SmallVector
+ operator*(const DataType& a, const SmallVector& x)
{
SmallVector<std::remove_const_t<DataType>> y = copy(x);
return y *= a;
@@ -51,14 +74,8 @@ class SmallVector // LCOV_EXCL_LINE
dot(const SmallVector& x, const SmallVector<DataType2>& y)
{
Assert(x.size() == y.size(), "cannot compute dot product: incompatible vector sizes");
- // Quite ugly, TODO: fix me in C++20
- auto promoted = [] {
- DataType a{0};
- DataType2 b{0};
- return a * b;
- }();
- decltype(promoted) sum = 0;
+ decltype(std::declval<DataType>() * std::declval<DataType2>()) sum = 0;
// Does not use parallel_for to preserve sum order
for (index_type i = 0; i < x.size(); ++i) {
@@ -136,7 +153,8 @@ class SmallVector // LCOV_EXCL_LINE
}
PUGS_INLINE
- DataType& operator[](index_type i) const noexcept(NO_ASSERT)
+ DataType&
+ operator[](index_type i) const noexcept(NO_ASSERT)
{
return m_values[i];
}
diff --git a/src/algebra/TinyVector.hpp b/src/algebra/TinyVector.hpp
index ae616a4849d0e681bf09e613ad4d19263d0f035a..e39dcec6e9834e5cb5eb80e1e845de1054beb1b2 100644
--- a/src/algebra/TinyVector.hpp
+++ b/src/algebra/TinyVector.hpp
@@ -42,14 +42,12 @@ class [[nodiscard]] TinyVector
return opposite;
}
- PUGS_INLINE
- constexpr size_t dimension() const
+ [[nodiscard]] PUGS_INLINE constexpr size_t dimension() const
{
return N;
}
- PUGS_INLINE
- constexpr bool operator==(const TinyVector& v) const
+ [[nodiscard]] PUGS_INLINE constexpr bool operator==(const TinyVector& v) const
{
for (size_t i = 0; i < N; ++i) {
if (m_values[i] != v.m_values[i])
@@ -58,14 +56,12 @@ class [[nodiscard]] TinyVector
return true;
}
- PUGS_INLINE
- constexpr bool operator!=(const TinyVector& v) const
+ [[nodiscard]] PUGS_INLINE constexpr bool operator!=(const TinyVector& v) const
{
return not this->operator==(v);
}
- PUGS_INLINE
- constexpr friend T dot(const TinyVector& u, const TinyVector& v)
+ [[nodiscard]] PUGS_INLINE constexpr friend T dot(const TinyVector& u, const TinyVector& v)
{
T t = u.m_values[0] * v.m_values[0];
for (size_t i = 1; i < N; ++i) {
@@ -242,7 +238,7 @@ class [[nodiscard]] TinyVector
};
template <size_t N, typename T>
-PUGS_INLINE constexpr T
+[[nodiscard]] PUGS_INLINE constexpr T
l2Norm(const TinyVector<N, T>& x)
{
static_assert(std::is_arithmetic<T>(), "Cannot compute L2 norm for non-arithmetic types");
@@ -250,13 +246,39 @@ l2Norm(const TinyVector<N, T>& x)
return std::sqrt(dot(x, x));
}
-// Cross product is only defined for K^3 vectors
+template <size_t N, typename T>
+[[nodiscard]] PUGS_INLINE constexpr T
+min(const TinyVector<N, T>& x)
+{
+ T m = x[0];
+ for (size_t i = 1; i < N; ++i) {
+ if (x[i] < m) {
+ m = x[i];
+ }
+ }
+ return m;
+}
+
+template <size_t N, typename T>
+[[nodiscard]] PUGS_INLINE constexpr T
+max(const TinyVector<N, T>& x)
+{
+ T m = x[0];
+ for (size_t i = 1; i < N; ++i) {
+ if (x[i] > m) {
+ m = x[i];
+ }
+ }
+ return m;
+}
+
+// Cross product is only defined for dimension 3 vectors
template <typename T>
-PUGS_INLINE constexpr TinyVector<3, T>
+[[nodiscard]] PUGS_INLINE constexpr TinyVector<3, T>
crossProduct(const TinyVector<3, T>& u, const TinyVector<3, T>& v)
{
TinyVector<3, T> cross_product(u[1] * v[2] - u[2] * v[1], u[2] * v[0] - u[0] * v[2], u[0] * v[1] - u[1] * v[0]);
return cross_product;
}
-#endif // TINYVECTOR_HPP
+#endif // TINY_VECTOR_HPP
diff --git a/src/algebra/Vector.hpp b/src/algebra/Vector.hpp
index c918ed409746a3c212720ad8ae6769cfb0882ba9..749a795a489f12a671f4ecfd452104702ba68147 100644
--- a/src/algebra/Vector.hpp
+++ b/src/algebra/Vector.hpp
@@ -31,6 +31,18 @@ class Vector // LCOV_EXCL_LINE
return os << x.m_values;
}
+ [[nodiscard]] friend std::remove_const_t<DataType>
+ min(const Vector& x)
+ {
+ return min(x.m_values);
+ }
+
+ [[nodiscard]] friend std::remove_const_t<DataType>
+ max(const Vector& x)
+ {
+ return max(x.m_values);
+ }
+
friend Vector<std::remove_const_t<DataType>>
copy(const Vector& x)
{
@@ -41,6 +53,16 @@ class Vector // LCOV_EXCL_LINE
return x_copy;
}
+ PUGS_INLINE
+ Vector<std::remove_const_t<DataType>>
+ operator-() const
+ {
+ Vector<std::remove_const_t<DataType>> opposite(this->size());
+ parallel_for(
+ opposite.size(), PUGS_LAMBDA(index_type i) { opposite.m_values[i] = -m_values[i]; });
+ return opposite;
+ }
+
friend Vector
operator*(const DataType& a, const Vector& x)
{
@@ -53,14 +75,7 @@ class Vector // LCOV_EXCL_LINE
dot(const Vector& x, const Vector<DataType2>& y)
{
Assert(x.size() == y.size(), "cannot compute dot product: incompatible vector sizes");
- // Quite ugly, TODO: fix me in C++20
- auto promoted = [] {
- DataType a{0};
- DataType2 b{0};
- return a * b;
- }();
-
- decltype(promoted) sum = 0;
+ decltype(std::declval<DataType>() * std::declval<DataType2>()) sum = 0;
// Does not use parallel_for to preserve sum order
for (index_type i = 0; i < x.size(); ++i) {
diff --git a/tests/test_SmallVector.cpp b/tests/test_SmallVector.cpp
index 3b499c3ab3e984247b7f20f96b626f8ee413763a..b5ffc6819736b2e55d33706aeb18f1e1d28dabca 100644
--- a/tests/test_SmallVector.cpp
+++ b/tests/test_SmallVector.cpp
@@ -360,6 +360,60 @@ TEST_CASE("SmallVector", "[algebra]")
REQUIRE(z[4] == 4);
}
+ SECTION("unary minus")
+ {
+ SmallVector<int> x{5};
+ x[0] = 2;
+ x[1] = 3;
+ x[2] = 5;
+ x[3] = 7;
+ x[4] = 8;
+
+ SmallVector<const int> y = x;
+
+ SmallVector<int> z = -x;
+
+ REQUIRE(z[0] == -2);
+ REQUIRE(z[1] == -3);
+ REQUIRE(z[2] == -5);
+ REQUIRE(z[3] == -7);
+ REQUIRE(z[4] == -8);
+
+ z = -y;
+ REQUIRE(z[0] == -2);
+ REQUIRE(z[1] == -3);
+ REQUIRE(z[2] == -5);
+ REQUIRE(z[3] == -7);
+ REQUIRE(z[4] == -8);
+ }
+
+ SECTION("min/max")
+ {
+ SmallVector<int> u(5);
+ u[0] = 1;
+ u[1] = 7;
+ u[2] = 6;
+ u[3] = 2;
+ u[4] = 9;
+
+ REQUIRE(min(u) == 1);
+ REQUIRE(max(u) == 9);
+ REQUIRE(min(-u) == -9);
+ REQUIRE(max(-u) == -1);
+
+ SmallVector<int> v(5);
+ v[0] = 1;
+ v[1] = 11;
+ v[2] = 6;
+ v[3] = -2;
+ v[4] = 9;
+
+ REQUIRE(min(v) == -2);
+ REQUIRE(max(v) == 11);
+ REQUIRE(min(-v) == -11);
+ REQUIRE(max(-v) == 2);
+ }
+
SECTION("output")
{
SmallVector<int> x{5};
diff --git a/tests/test_TinyVector.cpp b/tests/test_TinyVector.cpp
index 9c4dc2e5a2c4b185904fcf04a706fa70b572a197..6a3942e84b1c95680d8c63080b5c8b37330f139b 100644
--- a/tests/test_TinyVector.cpp
+++ b/tests/test_TinyVector.cpp
@@ -87,6 +87,23 @@ TEST_CASE("TinyVector", "[algebra]")
REQUIRE(crossProduct(a, b) == TinyVector<3, int>(-16, 6, 7));
}
+ SECTION("min/max")
+ {
+ TinyVector<5, int> u{1, 7, 6, 2, 9};
+
+ REQUIRE(min(u) == 1);
+ REQUIRE(max(u) == 9);
+ REQUIRE(min(-u) == -9);
+ REQUIRE(max(-u) == -1);
+
+ TinyVector<5, int> v{1, 11, 6, -2, 9};
+
+ REQUIRE(min(v) == -2);
+ REQUIRE(max(v) == 11);
+ REQUIRE(min(-v) == -11);
+ REQUIRE(max(-v) == 2);
+ }
+
#ifndef NDEBUG
SECTION("output with signaling NaN")
{
diff --git a/tests/test_Vector.cpp b/tests/test_Vector.cpp
index e31bf166713b70f69ca108d7dace6dc07179dc5b..0d8d0e898321070f74b3cf99f8fb0e04686c9b0d 100644
--- a/tests/test_Vector.cpp
+++ b/tests/test_Vector.cpp
@@ -359,6 +359,60 @@ TEST_CASE("Vector", "[algebra]")
REQUIRE(z[4] == 4);
}
+ SECTION("unary minus")
+ {
+ Vector<int> x{5};
+ x[0] = 2;
+ x[1] = 3;
+ x[2] = 5;
+ x[3] = 7;
+ x[4] = 8;
+
+ Vector<const int> y = x;
+
+ Vector<int> z = -x;
+
+ REQUIRE(z[0] == -2);
+ REQUIRE(z[1] == -3);
+ REQUIRE(z[2] == -5);
+ REQUIRE(z[3] == -7);
+ REQUIRE(z[4] == -8);
+
+ z = -y;
+ REQUIRE(z[0] == -2);
+ REQUIRE(z[1] == -3);
+ REQUIRE(z[2] == -5);
+ REQUIRE(z[3] == -7);
+ REQUIRE(z[4] == -8);
+ }
+
+ SECTION("min/max")
+ {
+ Vector<int> u(5);
+ u[0] = 1;
+ u[1] = 7;
+ u[2] = 6;
+ u[3] = 2;
+ u[4] = 9;
+
+ REQUIRE(min(u) == 1);
+ REQUIRE(max(u) == 9);
+ REQUIRE(min(-u) == -9);
+ REQUIRE(max(-u) == -1);
+
+ Vector<int> v(5);
+ v[0] = 1;
+ v[1] = 11;
+ v[2] = 6;
+ v[3] = -2;
+ v[4] = 9;
+
+ REQUIRE(min(v) == -2);
+ REQUIRE(max(v) == 11);
+ REQUIRE(min(-v) == -11);
+ REQUIRE(max(-v) == 2);
+ }
+
SECTION("output")
{
Vector<int> x{5};