From 14c4c5537f1be3258971ab9b60890cc78492978d Mon Sep 17 00:00:00 2001 From: Roy Stogner Date: Tue, 11 Aug 2026 16:45:15 -0500 Subject: [PATCH 1/9] libMesh::isfinite(DenseMatrix/SubMatrix/Vector) --- include/numerics/dense_matrix.h | 11 +++++++++++ include/numerics/dense_submatrix.h | 13 +++++++++++++ include/numerics/dense_vector.h | 14 ++++++++++++++ 3 files changed, 38 insertions(+) diff --git a/include/numerics/dense_matrix.h b/include/numerics/dense_matrix.h index 0cfa608e357..bd1a5549497 100644 --- a/include/numerics/dense_matrix.h +++ b/include/numerics/dense_matrix.h @@ -1218,6 +1218,17 @@ T DenseMatrix::transpose (const unsigned int i, // } +template +bool isfinite (const DenseMatrix & var) +{ + using std::isfinite; + const auto m = var.m(), n = var.n(); + for (auto i : make_range(m)) + for (auto j : make_range(n)) + if (!isfinite(var(i,j))) + return false; + return true; +} } // namespace libMesh diff --git a/include/numerics/dense_submatrix.h b/include/numerics/dense_submatrix.h index c4454ca3209..4b3a5dc428d 100644 --- a/include/numerics/dense_submatrix.h +++ b/include/numerics/dense_submatrix.h @@ -229,6 +229,19 @@ T & DenseSubMatrix::operator () (const unsigned int i, } +template +bool isfinite (const DenseSubMatrix & var) +{ + using std::isfinite; + const auto m = var.m(), n = var.n(); + for (auto i : make_range(m)) + for (auto j : make_range(n)) + if (!isfinite(var(i,j))) + return false; + return true; +} + + } // namespace libMesh diff --git a/include/numerics/dense_vector.h b/include/numerics/dense_vector.h index a5e3d39741a..35103b60d76 100644 --- a/include/numerics/dense_vector.h +++ b/include/numerics/dense_vector.h @@ -721,6 +721,20 @@ void DenseVector::get_principal_subvector (unsigned int sub_n, dest(i) = _val[i]; } + +template +bool isfinite (const DenseVector & var) +{ + using std::isfinite; + for (auto i : make_range(var.size())) + if (!isfinite(var(i))) + return false; + return true; +} + + + + } // namespace libMesh #ifdef LIBMESH_HAVE_METAPHYSICL From 462ac4f721f7d850132668c9b6171c7944ef10fd Mon Sep 17 00:00:00 2001 From: Roy Stogner Date: Tue, 11 Aug 2026 16:45:55 -0500 Subject: [PATCH 2/9] isfinite(Dense*) unit test coverage --- tests/numerics/dense_matrix_test.C | 30 ++++++++++++++++++++++++++++++ 1 file changed, 30 insertions(+) diff --git a/tests/numerics/dense_matrix_test.C b/tests/numerics/dense_matrix_test.C index ad89c27effe..7ee26a783d7 100644 --- a/tests/numerics/dense_matrix_test.C +++ b/tests/numerics/dense_matrix_test.C @@ -1,5 +1,6 @@ // libmesh includes #include +#include #include #ifdef LIBMESH_HAVE_PETSC @@ -8,6 +9,7 @@ #include "libmesh_cppunit.h" +#include using namespace libMesh; @@ -20,6 +22,7 @@ public: LIBMESH_CPPUNIT_TEST_SUITE(DenseMatrixTest); + CPPUNIT_TEST(testIsFinite); CPPUNIT_TEST(testOuterProduct); CPPUNIT_TEST(testSVD); CPPUNIT_TEST(testEVDreal); @@ -32,6 +35,33 @@ public: private: + void testIsFinite() + { + LOG_UNIT_TEST; + + DenseMatrix a(2, 2, {1, 2, + 3, 4}); + CPPUNIT_ASSERT(isfinite(a)); + + DenseVector diag = a.diagonal(); + CPPUNIT_ASSERT(isfinite(diag)); + + DenseSubMatrix suba1(a,0,0,1,1); + CPPUNIT_ASSERT(isfinite(suba1)); + + DenseSubMatrix suba2(a,1,1,1,1); + CPPUNIT_ASSERT(isfinite(suba2)); + + a(0, 0) = std::numeric_limits::infinity(); + CPPUNIT_ASSERT(!isfinite(a)); + CPPUNIT_ASSERT(!isfinite(suba1)); + CPPUNIT_ASSERT(isfinite(suba2)); + + CPPUNIT_ASSERT(isfinite(diag)); + diag = a.diagonal(); + CPPUNIT_ASSERT(!isfinite(diag)); + } + void testOuterProduct() { LOG_UNIT_TEST; From 5e7653d2c6575cff175cb6fc4b2c8eef6bd0847b Mon Sep 17 00:00:00 2001 From: Roy Stogner Date: Tue, 11 Aug 2026 17:37:38 -0500 Subject: [PATCH 3/9] Add libMesh::isfinite(TypeVector/TypeTensor) --- include/numerics/type_tensor.h | 13 +++++++++++++ include/numerics/type_vector.h | 12 ++++++++++++ 2 files changed, 25 insertions(+) diff --git a/include/numerics/type_tensor.h b/include/numerics/type_tensor.h index 470b745f120..4bc82666085 100644 --- a/include/numerics/type_tensor.h +++ b/include/numerics/type_tensor.h @@ -1435,6 +1435,19 @@ void TypeTensor::print(std::ostream & os) const #endif } + +template +bool isfinite (const TypeTensor & var) +{ + using std::isfinite; + for (unsigned int i=0; i inline TypeTensor::supertype> diff --git a/include/numerics/type_vector.h b/include/numerics/type_vector.h index 65bf6c3353f..b3dca47b3db 100644 --- a/include/numerics/type_vector.h +++ b/include/numerics/type_vector.h @@ -1186,6 +1186,18 @@ void TypeVector::print(std::ostream & os) const #endif } + +template +bool isfinite (const TypeVector & var) +{ + using std::isfinite; + for (unsigned int i=0; i struct CompareTypes, TypeVector> { From 0493dc13e788c1fb3c90395f50b9d51ccd9370bf Mon Sep 17 00:00:00 2001 From: Roy Stogner Date: Tue, 11 Aug 2026 17:38:02 -0500 Subject: [PATCH 4/9] isfinite(TypeVector/TypeTensor) unit tests --- tests/numerics/type_tensor_test.C | 16 ++++++++++++++++ tests/numerics/type_vector_test.h | 22 ++++++++++++++++++++++ 2 files changed, 38 insertions(+) diff --git a/tests/numerics/type_tensor_test.C b/tests/numerics/type_tensor_test.C index 5c2d88e94b7..fa9875cecbc 100644 --- a/tests/numerics/type_tensor_test.C +++ b/tests/numerics/type_tensor_test.C @@ -23,6 +23,7 @@ public: CPPUNIT_TEST(testRotation); #endif CPPUNIT_TEST(testRowCol); + CPPUNIT_TEST(testIsFinite); CPPUNIT_TEST(testIsZero); CPPUNIT_TEST(testIsHPD); #ifdef LIBMESH_HAVE_METAPHYSICL @@ -91,6 +92,21 @@ private: LIBMESH_ASSERT_FP_EQUAL(28, product(2, 2), tol); } + void testIsFinite() + { + LOG_UNIT_TEST; + + { + TensorValue tensor; + CPPUNIT_ASSERT(isfinite(tensor)); + } + { + TensorValue tensor(0,1,0,2,3); + tensor(0,0) = std::numeric_limits::infinity(); + CPPUNIT_ASSERT(!isfinite(tensor)); + } + } + void testIsZero() { LOG_UNIT_TEST; diff --git a/tests/numerics/type_vector_test.h b/tests/numerics/type_vector_test.h index 1a8e6afb8ad..9ea217a3f52 100644 --- a/tests/numerics/type_vector_test.h +++ b/tests/numerics/type_vector_test.h @@ -30,6 +30,7 @@ CPPUNIT_TEST( testVectorAddAssign ); \ CPPUNIT_TEST( testVectorSubAssign ); \ \ + CPPUNIT_TEST( testIsFinite ); \ CPPUNIT_TEST( testIsZero ); \ CPPUNIT_TEST( testSolidAngle ); \ CPPUNIT_TEST( testCircumcenter ); \ @@ -322,6 +323,27 @@ class TypeVectorTestBase : public CppUnit::TestCase { CPPUNIT_ASSERT_EQUAL( T(0), avector(i)); } + void testIsFinite() + { + LOG_UNIT_TEST; + + { +#if LIBMESH_DIM > 2 + DerivedClass avector(0,0,0); +#elif LIBMESH_DIM > 1 + DerivedClass avector(0,0); +#else + DerivedClass avector(0); +#endif + CPPUNIT_ASSERT(isfinite(avector)); + + avector(0) = + std::numeric_limits::infinity(); + CPPUNIT_ASSERT(!isfinite(avector)); + } + } + + void testIsZero() { LOG_UNIT_TEST; From 1308dca59c75361eb50ae4bbdee8f16da60e15d5 Mon Sep 17 00:00:00 2001 From: Roy Stogner Date: Tue, 11 Aug 2026 21:47:54 -0500 Subject: [PATCH 5/9] libMesh::isfinite(complex, Foo) You'd think std::isfinite would already support std::complex, right? But apparently the idea only got floated once back in 2022, the C++ people want to see it in C _Complex first, and I can't find a C equivalent to the C++ std-proposals mailing list so I'm not sure where it went from there. --- include/base/libmesh_common.h | 9 +++++++++ include/numerics/dense_matrix.h | 1 + include/numerics/dense_submatrix.h | 1 + include/numerics/dense_vector.h | 1 + include/numerics/type_tensor.h | 1 + include/numerics/type_vector.h | 1 + 6 files changed, 14 insertions(+) diff --git a/include/base/libmesh_common.h b/include/base/libmesh_common.h index 57c65bb38d9..95356b7110b 100644 --- a/include/base/libmesh_common.h +++ b/include/base/libmesh_common.h @@ -212,6 +212,15 @@ template inline bool libmesh_isinf(std::complex a) { return (std::isinf(std::real(a)) || std::isinf(std::imag(a))); } +// std::isfinite() doesn't support complex, but we really want an +// isfinite to support Number, at least in our own namespace +template +inline bool isfinite(std::complex a) +{ + using std::isfinite; + return (isfinite(std::real(a)) && isfinite(std::imag(a))); +} + // Define the value type for unknowns in simulations. // This is either Real or Complex, depending on how // the library was configures diff --git a/include/numerics/dense_matrix.h b/include/numerics/dense_matrix.h index bd1a5549497..65367e384bb 100644 --- a/include/numerics/dense_matrix.h +++ b/include/numerics/dense_matrix.h @@ -1222,6 +1222,7 @@ template bool isfinite (const DenseMatrix & var) { using std::isfinite; + using libMesh::isfinite; // for T==complex const auto m = var.m(), n = var.n(); for (auto i : make_range(m)) for (auto j : make_range(n)) diff --git a/include/numerics/dense_submatrix.h b/include/numerics/dense_submatrix.h index 4b3a5dc428d..39860ebc404 100644 --- a/include/numerics/dense_submatrix.h +++ b/include/numerics/dense_submatrix.h @@ -233,6 +233,7 @@ template bool isfinite (const DenseSubMatrix & var) { using std::isfinite; + using libMesh::isfinite; // for T==complex const auto m = var.m(), n = var.n(); for (auto i : make_range(m)) for (auto j : make_range(n)) diff --git a/include/numerics/dense_vector.h b/include/numerics/dense_vector.h index 35103b60d76..199c96c5e30 100644 --- a/include/numerics/dense_vector.h +++ b/include/numerics/dense_vector.h @@ -726,6 +726,7 @@ template bool isfinite (const DenseVector & var) { using std::isfinite; + using libMesh::isfinite; // for T==complex for (auto i : make_range(var.size())) if (!isfinite(var(i))) return false; diff --git a/include/numerics/type_tensor.h b/include/numerics/type_tensor.h index 4bc82666085..173dc451a24 100644 --- a/include/numerics/type_tensor.h +++ b/include/numerics/type_tensor.h @@ -1440,6 +1440,7 @@ template bool isfinite (const TypeTensor & var) { using std::isfinite; + using libMesh::isfinite; // for T==complex for (unsigned int i=0; i bool isfinite (const TypeVector & var) { using std::isfinite; + using libMesh::isfinite; // for T==complex for (unsigned int i=0; i Date: Tue, 11 Aug 2026 21:53:36 -0500 Subject: [PATCH 6/9] Avoid numeric_limits> Surely they can't blame *this* lacuna on the C people. --- tests/numerics/type_vector_test.h | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/tests/numerics/type_vector_test.h b/tests/numerics/type_vector_test.h index 9ea217a3f52..478e697bdbc 100644 --- a/tests/numerics/type_vector_test.h +++ b/tests/numerics/type_vector_test.h @@ -337,8 +337,10 @@ class TypeVectorTestBase : public CppUnit::TestCase { #endif CPPUNIT_ASSERT(isfinite(avector)); + // numeric_limits isn't specialized for complex, but float + // infinity promotes seamlessly to complex etc. avector(0) = - std::numeric_limits::infinity(); + std::numeric_limits::infinity(); CPPUNIT_ASSERT(!isfinite(avector)); } } From 57999ca41a8df280fd9c53472407c59407cb4985 Mon Sep 17 00:00:00 2001 From: Roy Stogner Date: Tue, 11 Aug 2026 21:55:32 -0500 Subject: [PATCH 7/9] Fix bug in ComplexVectorValueTest logging --- tests/numerics/vector_value_test.C | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/numerics/vector_value_test.C b/tests/numerics/vector_value_test.C index 3dde99f512a..ed0118102b5 100644 --- a/tests/numerics/vector_value_test.C +++ b/tests/numerics/vector_value_test.C @@ -49,7 +49,7 @@ public: this->libmesh_suite_name = "ComplexVectorValueTest"; } - LIBMESH_CPPUNIT_TEST_SUITE( NumberVectorValueTest ); + CPPUNIT_TEST_SUITE( ComplexVectorValueTest ); VECTORVALUETEST From 7c22989f2fd68e983ff5041835eb9bcef857c778 Mon Sep 17 00:00:00 2001 From: Roy Stogner Date: Wed, 12 Aug 2026 13:05:40 -0500 Subject: [PATCH 8/9] LIBMESH_ASSERT_COMPLEX_EQUAL Fixing our failure to test complex instantiations with Number!=complex uncovered the fact that we didn't actually have test macros that work for Number!=complex --- tests/libmesh_cppunit.h | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/tests/libmesh_cppunit.h b/tests/libmesh_cppunit.h index 9c1e34af8a4..80d312a619a 100644 --- a/tests/libmesh_cppunit.h +++ b/tests/libmesh_cppunit.h @@ -19,12 +19,15 @@ CPPUNIT_ASSERT_DOUBLES_EQUAL(expected,actual,tolerance) #endif -#ifdef LIBMESH_USE_COMPLEX_NUMBERS -# define LIBMESH_ASSERT_NUMBERS_EQUAL(expected,actual,tolerance) \ +# define LIBMESH_ASSERT_COMPLEX_EQUAL(expected,actual,tolerance) \ do { \ LIBMESH_ASSERT_FP_EQUAL(libMesh::libmesh_real(expected),libMesh::libmesh_real(actual),tolerance); \ LIBMESH_ASSERT_FP_EQUAL(libMesh::libmesh_imag(expected),libMesh::libmesh_imag(actual),tolerance); \ } while (0) + +#ifdef LIBMESH_USE_COMPLEX_NUMBERS +# define LIBMESH_ASSERT_NUMBERS_EQUAL(expected,actual,tolerance) \ + LIBMESH_ASSERT_COMPLEX_EQUAL(expected,actual,tolerance) #else # define LIBMESH_ASSERT_NUMBERS_EQUAL(expected,actual,tolerance) \ LIBMESH_ASSERT_FP_EQUAL(expected,actual,tolerance) From 694f880b5f6a8608d4069c4d0ae56ec918d93473 Mon Sep 17 00:00:00 2001 From: Roy Stogner Date: Wed, 12 Aug 2026 13:07:34 -0500 Subject: [PATCH 9/9] Use LIBMESH_ASSERT_COMPLEX_EQUAL as appropriate I should probably be less lazy and make something that templates on type and avoids all the "yup, 0i=0i" testing in the real-valued cases, but these are very fast tests anyway. --- tests/numerics/type_vector_test.h | 110 +++++++++++++++--------------- 1 file changed, 55 insertions(+), 55 deletions(-) diff --git a/tests/numerics/type_vector_test.h b/tests/numerics/type_vector_test.h index 478e697bdbc..b6279b3ecf2 100644 --- a/tests/numerics/type_vector_test.h +++ b/tests/numerics/type_vector_test.h @@ -180,12 +180,12 @@ class TypeVectorTestBase : public CppUnit::TestCase { DerivedClass avector = 0; for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , avector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , avector(i) , TOLERANCE*TOLERANCE ); DerivedClass bvector = 2.0; - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , bvector(0) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , bvector(0) , TOLERANCE*TOLERANCE ); for (int i = 1; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , bvector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , bvector(i) , TOLERANCE*TOLERANCE ); } void testScalarMult() @@ -194,9 +194,9 @@ class TypeVectorTestBase : public CppUnit::TestCase { for (int i = 0; i != LIBMESH_DIM; ++i) { - LIBMESH_ASSERT_NUMBERS_EQUAL( 5.0 , ((*m_1_1_1)*5.0)(i) , TOLERANCE*TOLERANCE ); - LIBMESH_ASSERT_NUMBERS_EQUAL( 5.0 , outer_product(*m_1_1_1,5.0)(i) , TOLERANCE*TOLERANCE ); - LIBMESH_ASSERT_NUMBERS_EQUAL( 5.0 , outer_product(5.0,*m_1_1_1)(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 5.0 , ((*m_1_1_1)*5.0)(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 5.0 , outer_product(*m_1_1_1,5.0)(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 5.0 , outer_product(5.0,*m_1_1_1)(i) , TOLERANCE*TOLERANCE ); } } @@ -205,7 +205,7 @@ class TypeVectorTestBase : public CppUnit::TestCase { LOG_UNIT_TEST; for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 1/Real(5) , ((*m_1_1_1)/5.0)(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 1/Real(5) , ((*m_1_1_1)/5.0)(i) , TOLERANCE*TOLERANCE ); } void testScalarMultAssign() @@ -216,7 +216,7 @@ class TypeVectorTestBase : public CppUnit::TestCase { avector*=5.0; for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 5.0 , avector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 5.0 , avector(i) , TOLERANCE*TOLERANCE ); } void testScalarDivAssign() @@ -227,18 +227,18 @@ class TypeVectorTestBase : public CppUnit::TestCase { avector/=5.0; for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 1/Real(5) , avector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 1/Real(5) , avector(i) , TOLERANCE*TOLERANCE ); } void testVectorAdd() { LOG_UNIT_TEST; - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , ((*m_1_1_1)+(*m_n1_1_n1))(0) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , ((*m_1_1_1)+(*m_n1_1_n1))(0) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 1) - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , ((*m_1_1_1)+(*m_n1_1_n1))(1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , ((*m_1_1_1)+(*m_n1_1_n1))(1) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 2) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , ((*m_1_1_1)+(*m_n1_1_n1))(2) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , ((*m_1_1_1)+(*m_n1_1_n1))(2) , TOLERANCE*TOLERANCE ); } void testVectorAddScaled() @@ -249,18 +249,18 @@ class TypeVectorTestBase : public CppUnit::TestCase { avector.add_scaled((*m_1_1_1),0.5); for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 1.5 , avector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 1.5 , avector(i) , TOLERANCE*TOLERANCE ); } void testVectorSub() { LOG_UNIT_TEST; - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , ((*m_1_1_1)-(*m_n1_1_n1))(0) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , ((*m_1_1_1)-(*m_n1_1_n1))(0) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 1) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , ((*m_1_1_1)-(*m_n1_1_n1))(1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , ((*m_1_1_1)-(*m_n1_1_n1))(1) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 2) - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , ((*m_1_1_1)-(*m_n1_1_n1))(2) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , ((*m_1_1_1)-(*m_n1_1_n1))(2) , TOLERANCE*TOLERANCE ); } void testVectorMult() @@ -268,9 +268,9 @@ class TypeVectorTestBase : public CppUnit::TestCase { LOG_UNIT_TEST; if (LIBMESH_DIM == 2) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , (*m_1_1_1)*(*m_n1_1_n1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , (*m_1_1_1)*(*m_n1_1_n1) , TOLERANCE*TOLERANCE ); else - LIBMESH_ASSERT_NUMBERS_EQUAL( -1.0 , (*m_1_1_1)*(*m_n1_1_n1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( -1.0 , (*m_1_1_1)*(*m_n1_1_n1) , TOLERANCE*TOLERANCE ); } void testVectorAddAssign() @@ -281,7 +281,7 @@ class TypeVectorTestBase : public CppUnit::TestCase { avector+=(*m_1_1_1); for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , avector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , avector(i) , TOLERANCE*TOLERANCE ); } void testVectorSubAssign() @@ -291,11 +291,11 @@ class TypeVectorTestBase : public CppUnit::TestCase { DerivedClass avector {*basem_1_1_1}; avector-=(*m_n1_1_n1); - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , avector(0) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , avector(0) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 1) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , avector(1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , avector(1) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 2) - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , avector(2) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , avector(2) , TOLERANCE*TOLERANCE ); } void testValueBase() @@ -377,7 +377,7 @@ class TypeVectorTestBase : public CppUnit::TestCase { LOG_UNIT_TEST; const DerivedClass xvec(1.3); - LIBMESH_ASSERT_NUMBERS_EQUAL( solid_angle(xvec, xvec, xvec), 0, TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( solid_angle(xvec, xvec, xvec), 0, TOLERANCE*TOLERANCE ); #if LIBMESH_DIM > 1 // This is ambiguous with --enable-complex builds? We really need @@ -392,18 +392,18 @@ class TypeVectorTestBase : public CppUnit::TestCase { xypdiag(0) = 0.8; xypdiag(1) = -.8; // Yeah, nothing subtends a non-zero solid angle in 2D either. - LIBMESH_ASSERT_NUMBERS_EQUAL( solid_angle(xvec, xvec, yvec), 0, TOLERANCE*TOLERANCE ); - LIBMESH_ASSERT_NUMBERS_EQUAL( solid_angle(xvec, yvec, xydiag), 0, TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( solid_angle(xvec, xvec, yvec), 0, TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( solid_angle(xvec, yvec, xydiag), 0, TOLERANCE*TOLERANCE ); #endif #if LIBMESH_DIM > 2 const DerivedClass zvec(0.,0.,1.1), xzdiag(0.,.7,.7); - LIBMESH_ASSERT_NUMBERS_EQUAL(solid_angle(xydiag, yvec, zvec), + LIBMESH_ASSERT_COMPLEX_EQUAL(solid_angle(xydiag, yvec, zvec), libMesh::pi/4, TOLERANCE*TOLERANCE); - LIBMESH_ASSERT_NUMBERS_EQUAL(solid_angle(xvec, yvec, xzdiag), + LIBMESH_ASSERT_COMPLEX_EQUAL(solid_angle(xvec, yvec, xzdiag), libMesh::pi/4, TOLERANCE*TOLERANCE); - LIBMESH_ASSERT_NUMBERS_EQUAL(solid_angle(xypdiag, xydiag, zvec), + LIBMESH_ASSERT_COMPLEX_EQUAL(solid_angle(xypdiag, xydiag, zvec), libMesh::pi/2, TOLERANCE*TOLERANCE); // Icosahedron coordinates are a nice analytic test @@ -411,7 +411,7 @@ class TypeVectorTestBase : public CppUnit::TestCase { const DerivedClass icosa1(1., golden, 0.), icosa2(-1., golden, 0.), icosa3(0., 1., golden); - LIBMESH_ASSERT_NUMBERS_EQUAL(solid_angle(icosa1, icosa2, icosa3), + LIBMESH_ASSERT_COMPLEX_EQUAL(solid_angle(icosa1, icosa2, icosa3), libMesh::pi/5, TOLERANCE*TOLERANCE); #endif } @@ -424,20 +424,20 @@ class TypeVectorTestBase : public CppUnit::TestCase { DerivedClass twoone(2); twoone(1) = 1; auto cc1 = circumcenter(origin, e_x, twoone); - LIBMESH_ASSERT_NUMBERS_EQUAL(cc1(0), 0.5, TOLERANCE*TOLERANCE); - LIBMESH_ASSERT_NUMBERS_EQUAL(cc1(1), 1.5, TOLERANCE*TOLERANCE); + LIBMESH_ASSERT_COMPLEX_EQUAL(cc1(0), 0.5, TOLERANCE*TOLERANCE); + LIBMESH_ASSERT_COMPLEX_EQUAL(cc1(1), 1.5, TOLERANCE*TOLERANCE); #if LIBMESH_DIM > 2 - LIBMESH_ASSERT_NUMBERS_EQUAL(cc1(2), 0, TOLERANCE*TOLERANCE); + LIBMESH_ASSERT_COMPLEX_EQUAL(cc1(2), 0, TOLERANCE*TOLERANCE); const DerivedClass twozeroone(2,0,1); auto cc2 = circumcenter(origin, e_x, twozeroone); - LIBMESH_ASSERT_NUMBERS_EQUAL(cc2(0), 0.5, TOLERANCE*TOLERANCE); - LIBMESH_ASSERT_NUMBERS_EQUAL(cc2(1), 0, TOLERANCE*TOLERANCE); - LIBMESH_ASSERT_NUMBERS_EQUAL(cc2(2), 1.5, TOLERANCE*TOLERANCE); + LIBMESH_ASSERT_COMPLEX_EQUAL(cc2(0), 0.5, TOLERANCE*TOLERANCE); + LIBMESH_ASSERT_COMPLEX_EQUAL(cc2(1), 0, TOLERANCE*TOLERANCE); + LIBMESH_ASSERT_COMPLEX_EQUAL(cc2(2), 1.5, TOLERANCE*TOLERANCE); auto cc3 = circumcenter(e_x, twoone, twozeroone); - LIBMESH_ASSERT_NUMBERS_EQUAL(cc3(0), Real(5)/3, TOLERANCE*TOLERANCE); - LIBMESH_ASSERT_NUMBERS_EQUAL(cc3(1), Real(1)/3, TOLERANCE*TOLERANCE); - LIBMESH_ASSERT_NUMBERS_EQUAL(cc3(2), Real(1)/3, TOLERANCE*TOLERANCE); + LIBMESH_ASSERT_COMPLEX_EQUAL(cc3(0), Real(5)/3, TOLERANCE*TOLERANCE); + LIBMESH_ASSERT_COMPLEX_EQUAL(cc3(1), Real(1)/3, TOLERANCE*TOLERANCE); + LIBMESH_ASSERT_COMPLEX_EQUAL(cc3(2), Real(1)/3, TOLERANCE*TOLERANCE); #endif } @@ -488,7 +488,7 @@ class TypeVectorTestBase : public CppUnit::TestCase { LOG_UNIT_TEST; for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 5.0 , ((*basem_1_1_1)*5.0)(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 5.0 , ((*basem_1_1_1)*5.0)(i) , TOLERANCE*TOLERANCE ); } void testScalarDivBase() @@ -496,7 +496,7 @@ class TypeVectorTestBase : public CppUnit::TestCase { LOG_UNIT_TEST; for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 1/Real(5) , ((*basem_1_1_1)/5.0)(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 1/Real(5) , ((*basem_1_1_1)/5.0)(i) , TOLERANCE*TOLERANCE ); } void testScalarMultAssignBase() @@ -507,7 +507,7 @@ class TypeVectorTestBase : public CppUnit::TestCase { avector*=5.0; for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 5.0 , avector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 5.0 , avector(i) , TOLERANCE*TOLERANCE ); } void testScalarDivAssignBase() @@ -518,18 +518,18 @@ class TypeVectorTestBase : public CppUnit::TestCase { avector/=5.0; for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 1/Real(5) , avector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 1/Real(5) , avector(i) , TOLERANCE*TOLERANCE ); } void testVectorAddBase() { LOG_UNIT_TEST; - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , ((*basem_1_1_1)+(*basem_n1_1_n1))(0) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , ((*basem_1_1_1)+(*basem_n1_1_n1))(0) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 1) - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , ((*basem_1_1_1)+(*basem_n1_1_n1))(1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , ((*basem_1_1_1)+(*basem_n1_1_n1))(1) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 2) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , ((*basem_1_1_1)+(*basem_n1_1_n1))(2) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , ((*basem_1_1_1)+(*basem_n1_1_n1))(2) , TOLERANCE*TOLERANCE ); } void testVectorAddScaledBase() @@ -540,18 +540,18 @@ class TypeVectorTestBase : public CppUnit::TestCase { avector.add_scaled((*basem_1_1_1),0.5); for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 1.5 , avector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 1.5 , avector(i) , TOLERANCE*TOLERANCE ); } void testVectorSubBase() { LOG_UNIT_TEST; - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , ((*basem_1_1_1)-(*basem_n1_1_n1))(0) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , ((*basem_1_1_1)-(*basem_n1_1_n1))(0) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 1) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , ((*basem_1_1_1)-(*basem_n1_1_n1))(1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , ((*basem_1_1_1)-(*basem_n1_1_n1))(1) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 2) - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , ((*basem_1_1_1)-(*basem_n1_1_n1))(2) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , ((*basem_1_1_1)-(*basem_n1_1_n1))(2) , TOLERANCE*TOLERANCE ); } void testVectorMultBase() @@ -559,9 +559,9 @@ class TypeVectorTestBase : public CppUnit::TestCase { LOG_UNIT_TEST; if (LIBMESH_DIM == 2) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , (*basem_1_1_1)*(*basem_n1_1_n1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , (*basem_1_1_1)*(*basem_n1_1_n1) , TOLERANCE*TOLERANCE ); else - LIBMESH_ASSERT_NUMBERS_EQUAL( -1.0 , (*basem_1_1_1)*(*basem_n1_1_n1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( -1.0 , (*basem_1_1_1)*(*basem_n1_1_n1) , TOLERANCE*TOLERANCE ); } void testVectorAddAssignBase() @@ -572,7 +572,7 @@ class TypeVectorTestBase : public CppUnit::TestCase { avector+=(*basem_1_1_1); for (int i = 0; i != LIBMESH_DIM; ++i) - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , avector(i) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , avector(i) , TOLERANCE*TOLERANCE ); } void testVectorSubAssignBase() @@ -582,11 +582,11 @@ class TypeVectorTestBase : public CppUnit::TestCase { TypeVector avector(*m_1_1_1); avector-=(*basem_n1_1_n1); - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , avector(0) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , avector(0) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 1) - LIBMESH_ASSERT_NUMBERS_EQUAL( 0.0 , avector(1) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 0.0 , avector(1) , TOLERANCE*TOLERANCE ); if (LIBMESH_DIM > 2) - LIBMESH_ASSERT_NUMBERS_EQUAL( 2.0 , avector(2) , TOLERANCE*TOLERANCE ); + LIBMESH_ASSERT_COMPLEX_EQUAL( 2.0 , avector(2) , TOLERANCE*TOLERANCE ); } void testReplaceAlgebraicType()