#include #include #include #include #include #include #include "serif/backend/mfem/serialize/configuration.hpp" #include "serif/backend/mfem/serialize/options.hpp" #include "serif/backend/mfem/discretization/element_coordinates.hpp" #include "serif/tests/solver/stellar_fixture.hpp" #include "serif/tests/nonconforming.hpp" namespace { namespace fixture = serif::tests::solver; namespace archive = serif::utils::serialize; namespace backend = serif::backend::mfem::serialize; using Discretization = decltype(fixture::make_discretization(std::declval>())); /** @brief Exercise simplex orientation independently of the moving spherical mesh generator. */ struct FixedResources { ::mfem::Mesh mesh; ::mfem::H1_FECollection scalar_collection{fixture::Density::FunctionalSpace::order, 3}; ::mfem::H1_FECollection displacement_collection{fixture::Displacement::FunctionalSpace::order, 3}; ::mfem::H1_FECollection potential_collection{fixture::Potential::FunctionalSpace::order, 3}; ::mfem::FiniteElementSpace density_space{&mesh, &scalar_collection}; ::mfem::FiniteElementSpace displacement_space{&mesh, &displacement_collection, 3, ::mfem::Ordering::byVDIM}; ::mfem::FiniteElementSpace potential_space{&mesh, &potential_collection}; ::mfem::GridFunction compactification{&density_space}; serif::discretization::backend::mfem::TermBinding density{density_space}; serif::discretization::backend::mfem::TermBinding displacement{displacement_space}; serif::discretization::backend::mfem::TermBinding potential{potential_space}; serif::mapping::DomainMapper mapper{{}, std::make_unique(serif::utils::user_options::KelvinCompactificationOptions{})}; /** @brief Finalize a supported serial mesh before constructing any term binding. */ explicit FixedResources(::mfem::Element::Type geometry) : mesh(make_mesh(geometry)) { compactification = 0.0; } /** @brief Retain nontrivial support and high-order mesh coordinates in the small fixed-domain case. */ static ::mfem::Mesh make_mesh(::mfem::Element::Type geometry) { auto mesh = ::mfem::Mesh::MakeCartesian3D(2, 1, 1, geometry); mesh.GetElement(mesh.GetNE() - 1)->SetAttribute(3); mesh.SetCurvature(2); return mesh; } }; /** @brief Ensure restoration cannot accidentally borrow the application's original resources. */ archive::Node saved_star(backend::SaveOptions options = {}) { const auto resources = fixture::make_resources(1.25, 0, 2); const auto discretization = fixture::make_discretization(resources, 2); return backend::save_discretization(discretization, options); } /** @brief Keep numerical reference choices together with the model and its native resources. */ auto make_configuration() { const auto resources = fixture::make_resources(1.25, 0, 2); const serif::solver::normalization::NormalizationOptions normalization{ .scheme = serif::solver::normalization::LaneEmdenApproximation{ .n = 1.0, .mass = serif::dimensions::MassValue{1.0}, .central_density = serif::dimensions::DensityValue{1.0}} }; return serif::model::make_discretization(serif::model::DiscretizationOptions{ .model = fixture::make_model(), .finite_elements = fixture::make_discretization(resources, 2), .normalization = normalization }); } /** @brief Compare native element coefficients after restoring their exact basis and orientation. */ template void compare_field(const Original& original, const Restored& restored) { const auto& before = original.template binding(); const auto& after = restored.template binding(); REQUIRE(before.size() == after.size()); ::mfem::Vector coefficients(before.size()); for (int dof = 0; dof < coefficients.Size(); ++dof) { coefficients[dof] = 0.125 * static_cast(dof - 7); } coefficients[0] = -0.0; const auto before_elements = before.validated_elements(); const auto after_elements = after.validated_elements(); for (int element = 0; element < before.mesh().GetNE(); ++element) { ::mfem::Vector first, second; before_elements.element_values(element, coefficients, first); after_elements.element_values(element, coefficients, second); REQUIRE(first.Size() == second.Size()); for (int local = 0; local < first.Size(); ++local) { CHECK(std::bit_cast(first[local]) == std::bit_cast(second[local])); } } } } TEST_CASE("Fixed MFEM discretizations retain simplex and tensor element orientation", "[unit][serialize][mfem_archive]") { namespace dependencies = serif::utils::dependencies; for (const ::mfem::Element::Type geometry : {::mfem::Element::HEXAHEDRON, ::mfem::Element::TETRAHEDRON}) { const auto resources = std::make_shared(geometry); const auto bindings = dependencies::make_dependency_context(dependencies::borrow(resources->density), dependencies::borrow(resources->potential), dependencies::borrow(resources->displacement)); const auto original = serif::model::backend::mfem::make_stellar_discretization(resources, bindings, resources->compactification, resources->mapper); const archive::Node record = backend::save_discretization(original); auto restored = backend::restore_discretization(record); INFO((restored ? "restored" : restored.error().path + ": " + restored.error().message)); REQUIRE(restored); CHECK(backend::save_discretization(*restored) == record); compare_field(original, *restored); compare_field(original, *restored); compare_field(original, *restored); } } TEST_CASE("MFEM archives own exact meshes bases support and oriented coefficient interpretations", "[unit][serialize][mfem_archive]") { const auto resources = fixture::make_resources(1.25, 0, 2); const auto original = fixture::make_discretization(resources, 2); const archive::Node record = backend::save_discretization(original); auto restored = backend::restore_discretization(record); INFO((restored ? "restored" : restored.error().path + ": " + restored.error().message)); REQUIRE(restored); CHECK(backend::save_discretization(*restored) == record); CHECK(&restored->template binding().mesh() != resources->mesh.mesh.get()); CHECK(restored->deformation.reference.scalar_space.FEColl() == restored->deformation.reference.displacement_space.FEColl()); CHECK(restored->deformation.reference.displacement_space.GetOrdering() == ::mfem::Ordering::byNODES); CHECK(restored->template binding().space().GetOrdering() == ::mfem::Ordering::byVDIM); CHECK(restored->deformation.model == original.deformation.model); compare_field(original, *restored); compare_field(original, *restored); compare_field(original, *restored); compare_field(original, *restored); compare_field(original, *restored); // Nothing from the writer's resource owner survives this helper call. const archive::Node independent_record = saved_star(); const auto independent = backend::restore_discretization(independent_record); REQUIRE(independent); CHECK(backend::save_discretization(*independent) == independent_record); } TEST_CASE("MFEM archives reject incompatible numbering and resource associations", "[unit][serialize][mfem_archive]") { archive::Node record = saved_star(); auto& resources = record.children.at("resources"); SECTION("Mesh contents disagree with their content identity") { auto& native = resources.children.at("meshes").children.at("0").children.at("native"); std::get(native.value).append("\n"); } SECTION("A matching fingerprint does not make malformed native mesh text valid") { auto& mesh = resources.children.at("meshes").children.at("0"); const std::string malformed = "This is not an MFEM mesh.\n"; mesh.children.at("native") = archive::encode(malformed); mesh.children.at("identity") = archive::encode(archive::fingerprint(malformed)); const auto restored = backend::restore_discretization(record); REQUIRE_FALSE(restored); CHECK(restored.error().message.starts_with("Cannot reconstruct saved MFEM resources:")); return; } SECTION("A collection's basis changes without changing its nominal polynomial order") { const auto space = archive::decode(record.children.at("bindings").children.at("0").children.at("space")); REQUIRE(space); const auto collection_index = archive::decode( resources.children.at("spaces").children.at(std::to_string(*space)).children.at("collection")); REQUIRE(collection_index); auto& collection = resources.children.at("collections").children.at(std::to_string(*collection_index)); const auto family = archive::decode(collection.children.at("family")); REQUIRE(family); REQUIRE(*family == "H1"); collection.children.at("closed_basis") = archive::encode(static_cast(::mfem::BasisType::ClosedUniform)); } SECTION("An element's signed local numbering changes") { auto& numbering = resources.children.at("spaces").children.at("0").children.at("element_numbering"); std::get>(numbering.value)[1] = -100; } SECTION("A term is connected to a missing space") { record.children.at("bindings").children.at("0").children.at("space") = archive::encode(std::uint64_t{900}); } SECTION("A term's reduced map is changed") { auto& reduction = record.children.at("bindings").children.at("0").children.at("coordinates").children.at("full_to_reduced"); std::get>(reduction.value)[0] = 999; } SECTION("A deformation companion is not a displacement space") { record.children.at("reference").children.at("displacement_space") = record.children.at("reference").children.at("scalar_space"); } SECTION("Native MFEM version differs") { resources.children.at("mfem_version") = archive::encode(std::string{"different"}); } CHECK_FALSE(backend::restore_discretization(record)); } TEST_CASE("An omitted MFEM mesh must be resolved by its exact saved identity", "[unit][serialize][mfem_archive]") { const auto resources = fixture::make_resources(1.25, 0, 2); const auto original = fixture::make_discretization(resources, 2); const archive::Node record = backend::save_discretization(original, {.embed_mesh = false}); CHECK_FALSE(backend::restore_discretization(record)); std::map> meshes; for (const ::mfem::Mesh* mesh : {resources->mesh.mesh.get(), resources->mesh.reference_mesh.get()}) { const std::string text = backend::mesh_text(*mesh); meshes.emplace(archive::fingerprint(text), text); } auto restored = backend::restore_discretization(record, [&](std::string_view identity) -> archive::Result { const auto found = meshes.find(identity); if (found == meshes.end()) { return std::unexpected(archive::Error{"mesh", "Unknown mesh identity."}); } return found->second; }); REQUIRE(restored); CHECK(backend::save_discretization(*restored, {.embed_mesh = false}) == record); CHECK_FALSE(backend::restore_discretization(record, [](std::string_view) -> archive::Result { return std::string{"wrong mesh"}; })); } TEST_CASE("MFEM assembly options retain quadrature controls and reject unknown contraction modes", "[unit][serialize][mfem_archive]") { namespace quadrature = serif::discretization::quadrature; using Options = serif::operators::backend::mfem::VolumeOptions; quadrature::RuleSet rules; rules.set_role_control(quadrature::QuadratureRole::discretization, {.fixed_order = 9, .boost = 2}); const Options options{.quadrature_policy = quadrature::Policy{rules}, .contraction_mode = serif::operators::backend::mfem::ContractionMode::dense}; archive::Node record = archive::encode(options); auto restored = archive::decode(record); REQUIRE(restored); CHECK(archive::encode(*restored) == record); record.children.at("contraction") = archive::encode(500); CHECK_FALSE(archive::decode(record)); } TEST_CASE("Configured MFEM models restore their full numerical prescription with independent ownership", "[unit][serialize][mfem_archive]") { using Configuration = decltype(make_configuration()); archive::Node record; std::weak_ptr original_owner; { const auto original = make_configuration(); original_owner = original.finite_elements.owner; record = archive::encode(original); } REQUIRE(original_owner.expired()); REQUIRE(record.type == "model.configured_discretization"); const auto restored = backend::restore_configuration(record); INFO((restored ? "restored" : restored.error().path + ": " + restored.error().message)); REQUIRE(restored); const archive::Node saved_again = archive::encode(*restored); CHECK(saved_again.children.at("finite_elements") == record.children.at("finite_elements")); CHECK(saved_again.children.at("normalization") == record.children.at("normalization")); CHECK(saved_again.children.at("preconditioner") == record.children.at("preconditioner")); CHECK(saved_again.children.at("model").children.at("resolved_equations") == record.children.at("model").children.at("resolved_equations")); SECTION("Unknown configuration versions cannot be restored") { record.version = 2; } SECTION("Missing numerical choices cannot silently become defaults") { record.children.erase("normalization"); } SECTION("Unsupported preconditioners cannot silently become no preconditioning") { record.children.at("preconditioner") = archive::Node{.type = "unsupported", .value = std::string{"Custom preconditioner has no serializer."}}; } CHECK_FALSE(backend::restore_configuration(record)); } TEST_CASE("MFEM archives refuse native parsing when the application selects abort on error", "[unit][serialize][mfem_archive]") { const archive::Node record = saved_star(); struct RestoreErrorAction { ::mfem::ErrorAction previous{::mfem::get_error_action()}; /** @brief Restore the application's prior native error policy even if an assertion stops the test. */ ~RestoreErrorAction() { ::mfem::set_error_action(previous); } } restore_error_action; ::mfem::set_error_action(::mfem::MFEM_ERROR_ABORT); const auto restored = backend::restore_discretization(record); REQUIRE_FALSE(restored); CHECK(restored.error().message.find("requires MFEM_ERROR_THROW") != std::string::npos); CHECK(::mfem::get_error_action() == ::mfem::MFEM_ERROR_ABORT); } TEST_CASE("Nonconforming STROID archives retain independent coordinates after original owners expire", "[unit][serialize][mfem_archive][nonconforming]") { archive::Node record; std::weak_ptr original_owner; const auto restored = [&] { const auto resources = std::make_shared(serif::tests::nonconforming::mesh(2, 0, 2)); original_owner = resources; const auto original = fixture::make_discretization(resources, 2); REQUIRE(resources->potential_space.GetTrueVSize() < resources->potential_space.GetVSize()); record = backend::save_discretization(original); auto result = backend::restore_discretization(record); if (!result) { INFO(result.error().path << ": " << result.error().message); } REQUIRE(result); compare_field(original, *result); compare_field(original, *result); compare_field(original, *result); CHECK(&result->template binding().mesh() != resources->mesh.mesh.get()); return std::move(*result); }(); CHECK(original_owner.expired()); const auto& potential = restored.template binding(); REQUIRE(potential.space().GetTrueVSize() < potential.space().GetVSize()); CHECK(potential.space().GetConformingProlongation() != nullptr); CHECK(restored.deformation.reference.logical_mesh.Nonconforming()); CHECK(backend::save_discretization(restored) == record); mfem::Vector coefficients(potential.size()); for (int index = 0; index < coefficients.Size(); ++index) { coefficients[index] = 0.125 * (index - 7); } const mfem::Vector local = potential.scatter(coefficients); const mfem::Vector independent = potential.gather(local); REQUIRE(independent.Size() == coefficients.Size()); for (int index = 0; index < coefficients.Size(); ++index) { CHECK(independent[index] == coefficients[index]); } }