diff --git a/src/fesa/analysis/eigen_reduction.cpp b/src/fesa/analysis/eigen_reduction.cpp index b7e98fc..35b2ccb 100644 --- a/src/fesa/analysis/eigen_reduction.cpp +++ b/src/fesa/analysis/eigen_reduction.cpp @@ -118,6 +118,22 @@ Result MakeLift(const Matrix& mass, if (positive == 0U) return Result::Success(std::move(result)); const Matrix scaled = ScaledMass(mass, length); + // Raw frame tangents can enter this internal boundary independently of + // Domain validation. Check the actual operands before any dense dispatch. + for (std::size_t row = 0; row < 6U; ++row) { + for (std::size_t j = 0; j < basis.Columns(); ++j) { + if (!std::isfinite(basis(row, j))) + return Result::Failure( + Failure("invalid-eigen-mass-space", node, + "Invalid positive mass basis.")); + } + for (std::size_t j = 0; j < 6U; ++j) { + if (!std::isfinite(scaled(row, j))) + return Result::Failure( + Failure("invalid-eigen-mass-space", node, + "Scaled mass is nonfinite.")); + } + } Matrix mp = dense_blas_internal::MultiplyTransposeLeft(basis, scaled) .Multiply(basis); double scale = 0.0; diff --git a/tests/unit/analysis/eigen_reduction_test.cpp b/tests/unit/analysis/eigen_reduction_test.cpp index bc31dd3..a7f9108 100644 --- a/tests/unit/analysis/eigen_reduction_test.cpp +++ b/tests/unit/analysis/eigen_reduction_test.cpp @@ -6,6 +6,7 @@ #include #include #include +#include #include #include @@ -92,6 +93,38 @@ TEST(EigenReduction, RemovesOnlyActualFreeDirectorNull) { .HasValue()); } +TEST(EigenReduction, RejectsNonfiniteUsedBasisAndPreservesPriorCandidate) { + const fesa::ShellNodeInitialFrame frame{0U, {0, 0, 1}, {1, 0, 0}, {0, 1, 0}}; + const auto mass = NodalMass(frame.director); + fesa::Matrix stiffness(6U, 6U); + for (std::size_t i = 0; i < 5U; ++i) stiffness(i, i) = mass(i, i) * 7.0; + stiffness(5U, 5U) = .25; + const auto matrix = Sparse(stiffness); + const std::vector mask(6U, 0U); + const auto prior = fesa::eigen_internal::BuildEigenReduction( + matrix, {mass}, {frame}, mask, {.25}, 2.0); + ASSERT_TRUE(prior.HasValue()); + const auto prior_values = prior.Value().whitened_operator.Values(); + for (const double invalid : {std::numeric_limits::quiet_NaN(), + std::numeric_limits::infinity()}) { + for (const bool first_tangent : {false, true}) { + auto bad_frame = frame; + (first_tangent ? bad_frame.tangent_a : bad_frame.tangent_b)[0] = invalid; + const auto rejected = fesa::eigen_internal::BuildEigenReduction( + matrix, {mass}, {bad_frame}, mask, {.25}, 2.0); + ASSERT_FALSE(rejected.HasValue()); + EXPECT_EQ(rejected.GetStatus().Category(), fesa::FailureCategory::kSolver); + EXPECT_EQ(rejected.GetStatus().Diagnostics()[0].code, + "invalid-eigen-mass-space"); + EXPECT_EQ(prior.Value().whitened_operator.Values(), prior_values); + } + } + const auto rerun = fesa::eigen_internal::BuildEigenReduction( + matrix, {mass}, {frame}, mask, {.25}, 2.0); + ASSERT_TRUE(rerun.HasValue()); + EXPECT_EQ(rerun.Value().whitened_operator.Values(), prior_values); +} + TEST(EigenReduction, RejectsPositiveButUnresolvableMassWithoutDeletingDofs) { const double p = 1e-8; const std::array d{2 * p / (1 + p * p), 0,