diff --git a/docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md b/docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md index 174544b..ca81a5f 100644 --- a/docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md +++ b/docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md @@ -809,3 +809,54 @@ | REVIEW1-targeted | `ctest --test-dir .harness/build -C Debug -R SparseAssembly --output-on-failure` | 0 | Exact SparseAssembly suite passed 5/5 | | REVIEW1-full | `ctest --test-dir .harness/build -C Debug --output-on-failure` | 0 | Full accumulated suite passed 47/47 | | REVIEW1-reference | `git diff --exit-code 59da6c6 -- reference/`; `git status --short -- reference/` | 0 | Approved reference artifacts remain unchanged and reference status is empty | + +## Step 19 — essential-constraints + +- task_id: `TASK-19` +- status: `completed` +- changed_files: `include/fesa/constraints/essential_constraints.hpp`, + `src/fesa/constraints/essential_constraints.cpp`, + `tests/unit/constraints/essential_constraints_test.cpp`, + `src/fesa/CMakeLists.txt`, `tests/CMakeLists.txt`, + `docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md`, + `phases/linear-static-3d-euler-beam/index.json` +- requirement_ids: `FESA-REQ-LS3DEB-007`, `FESA-REQ-LS3DEB-022`, + `FESA-REQ-LS3DEB-027`, `FESA-REQ-LS3DEB-034`, + `FESA-REQ-LS3DEB-035` +- test_ids: `T19-CONSTRAINT-001`, `T19-CONSTRAINT-002`, + `T19-CONSTRAINT-003`, `T19-CONSTRAINT-004` + +| stage | exact command | exit_code | expected_or_observed_result | evidence_tail | +| --- | --- | ---: | --- | --- | +| RED-build | `cmake --build .harness/build --config Debug --target fesa_tests` | 1 | Exactly four planned tests were registered before production and the public constraint API was absent | MSVC C1083 reported missing `fesa/constraints/essential_constraints.hpp` from `essential_constraints_test.cpp` | +| GREEN-build | `cmake --build .harness/build --config Debug --target fesa_tests` | 0 | Minimum partition/mapping implementation and all four tests compiled and linked | `essential_constraints.cpp`, its test, `fesa_solver.lib`, and `fesa_unit_tests.exe` built without a FESA warning under `/W4 /WX` | +| GREEN-test | `ctest --test-dir .harness/build -C Debug -R EssentialConstraints --output-on-failure` | 0 | Hand-computed blocks, no/all/mixed constraints, nonzero reconstruction, and invalid dimensions passed | 4/4 exact `EssentialConstraints` tests passed | +| VERIFY-configure | `cmake -S . -B .harness/build -A x64 -DFESA_GTEST_SOURCE_DIR=C:/git/googletest "-DMKL_DIR=C:/Program Files (x86)/Intel/oneAPI/mkl/2026.1/lib/cmake/mkl" "-DTBB_DIR=C:/Program Files (x86)/Intel/oneAPI/tbb/2023.1/lib/cmake/tbb" "-DHDF5_DIR=C:/Program Files/HDF_Group/HDF5/2.1.1/cmake"` | 0 | Approved explicit-dependency MSVC x64 build tree generates | Windows SDK, oneMKL 2026.1, oneTBB, and HDF5 resolved; configure/generate completed | +| VERIFY-build | `cmake --build .harness/build --config Debug` | 0 | Full Debug build passes without a new FESA warning | `fesa_solver.lib` and `fesa_unit_tests.exe` built under `/W4 /WX` | +| VERIFY-targeted | `ctest --test-dir .harness/build -C Debug -R EssentialConstraints --output-on-failure` | 0 | Focused Step 19 suite remains green | 4/4 exact `EssentialConstraints` tests passed | +| VERIFY-discovery | `ctest --test-dir .harness/build -C Debug --show-only=json-v1` | 0 | CTest discovers the accumulated suite and all four exact names | 51 tests discovered, including 4 `EssentialConstraints` tests | +| VERIFY-full | `ctest --test-dir .harness/build -C Debug --output-on-failure` | 0 | Full accumulated C++ suite has zero failures | 51/51 tests passed | +| VERIFY-contract-scans | Public backend and forbidden load/RHS/penalty/MPC scans, exact-test-count scan, `git diff --check`, and reference status/diff checks | 0 | Essential partition remains isolated behind backend-free public types | backend leaks 0; forbidden scope 0; exact tests 4; whitespace clean; reference unchanged | + +- contract_checks: `partition` validates the input CSR, exact square/full-DOF + dimensions, stable increasing free/constrained order, disjoint complete + coverage, free-equation correspondence, and prescribed-vector size before + extracting any block. +- contract_checks: each full CSR slot maps at most once into + `Kff/Kfc/Kcf/Kcc`. Filtering the stable full rows and sorted columns retains + numeric values, local row/column order, and exact structural zeros; every + returned block revalidates through the canonical 0-based CSR constructor. +- contract_checks: no/all/mixed constraint cases include valid zero-sized + rectangular blocks. `gatherFree`, `gatherConstrained`, and + `reconstructFull` reject vector-size mismatch and preserve caller-supplied + constrained values exactly, including nonzero prescribed displacement. +- contract_checks: the module does not assemble loads, form an effective RHS, + implement penalty/MPC behavior, or expose MKL/oneTBB/HDF5/PARDISO types. +- generated_evidence: `.harness/build/src/fesa/Debug/fesa_solver.lib`, + `.harness/build/tests/Debug/fesa_unit_tests.exe` +- reference_diff: unchanged; `git diff --exit-code -- reference/` exit 0 +- handoff: Step 20 may factorize the validated `kff`; Step 21 may use `kfc` + with `DofManager::prescribedValues()` to form the nonzero prescribed + effective RHS without changing this module's mapping responsibility. +- concerns: none; no critical implementation or upstream contract conflict + was found. diff --git a/include/fesa/constraints/essential_constraints.hpp b/include/fesa/constraints/essential_constraints.hpp new file mode 100644 index 0000000..5a5d4c2 --- /dev/null +++ b/include/fesa/constraints/essential_constraints.hpp @@ -0,0 +1,35 @@ +#pragma once + +#include "fesa/core/status.hpp" +#include "fesa/math/sparse_matrix.hpp" +#include "fesa/math/vector.hpp" + +namespace fesa { + +class DofManager; + +struct PartitionedStiffness { + SparseMatrix kff; + SparseMatrix kfc; + SparseMatrix kcf; + SparseMatrix kcc; +}; + +// Applies the DofManager's stable elimination order without owning equation +// numbering, load assembly, or a solver policy. +class EssentialConstraints { +public: + static Result partition( + const SparseMatrix& full, + const DofManager& dofs); + static Vector gatherFree(const Vector& full, const DofManager& dofs); + static Vector gatherConstrained( + const Vector& full, + const DofManager& dofs); + static Vector reconstructFull( + const Vector& freeValues, + const Vector& constrainedValues, + const DofManager& dofs); +}; + +} // namespace fesa diff --git a/src/fesa/CMakeLists.txt b/src/fesa/CMakeLists.txt index 81b6c90..29a765c 100644 --- a/src/fesa/CMakeLists.txt +++ b/src/fesa/CMakeLists.txt @@ -6,6 +6,7 @@ add_library( assembly/parallel_for.cpp assembly/sparse_assembler.cpp build_info.cpp + constraints/essential_constraints.cpp core/diagnostic.cpp core/status.cpp elements/euler_beam_3d.cpp diff --git a/src/fesa/constraints/essential_constraints.cpp b/src/fesa/constraints/essential_constraints.cpp new file mode 100644 index 0000000..79bc934 --- /dev/null +++ b/src/fesa/constraints/essential_constraints.cpp @@ -0,0 +1,261 @@ +#include "fesa/constraints/essential_constraints.hpp" + +#include "fesa/fem/dof_manager.hpp" + +#include +#include +#include +#include +#include +#include + +namespace fesa { +namespace { + +Status constraintFailure( + const std::string& code, + const std::string& identity, + const std::string& message) { + return Status::failure( + FailureCategory::model, + {{Severity::error, + code, + {{}, 0U}, + "ESSENTIAL_CONSTRAINTS", + identity, + message}}); +} + +bool isStrictlyIncreasing(const std::vector& values) { + return std::adjacent_find( + values.begin(), + values.end(), + [](const std::size_t left, const std::size_t right) { + return left >= right; + }) == values.end(); +} + +Status validateDofOrder(const DofManager& dofs) { + const std::size_t fullCount = dofs.fullDofCount(); + const auto& freeDofs = dofs.freeDofs(); + const auto& constrainedDofs = dofs.constrainedDofs(); + if (freeDofs.size() != dofs.freeDofCount() || + constrainedDofs.size() != dofs.constrainedDofCount() || + dofs.prescribedValues().size() != constrainedDofs.size() || + constrainedDofs.size() > fullCount || + freeDofs.size() != fullCount - constrainedDofs.size()) { + return constraintFailure( + "invalid-constraint-dimensions", + std::to_string(fullCount), + "DofManager full, free, constrained, and prescribed dimensions must agree."); + } + if (!isStrictlyIncreasing(freeDofs) || + !isStrictlyIncreasing(constrainedDofs)) { + return constraintFailure( + "invalid-constraint-order", + std::to_string(fullCount), + "Free and constrained DOFs must use stable increasing full-DOF order."); + } + + std::vector ownership(fullCount, 0U); + try { + for (std::size_t equation = 0U; + equation < freeDofs.size(); + ++equation) { + const std::size_t fullDof = freeDofs[equation]; + if (fullDof >= fullCount || ownership[fullDof] != 0U || + dofs.freeEquation(fullDof) != equation) { + return constraintFailure( + "invalid-constraint-order", + std::to_string(fullDof), + "Free equation numbering must match the stable free-DOF order."); + } + ownership[fullDof] = 1U; + } + for (const std::size_t fullDof : constrainedDofs) { + if (fullDof >= fullCount || ownership[fullDof] != 0U || + dofs.freeEquation(fullDof).has_value()) { + return constraintFailure( + "invalid-constraint-order", + std::to_string(fullDof), + "Constrained DOFs must be unique and absent from free equations."); + } + ownership[fullDof] = 2U; + } + } catch (const std::out_of_range&) { + return constraintFailure( + "invalid-constraint-dimensions", + std::to_string(fullCount), + "DofManager equation storage does not cover every full DOF."); + } + if (std::find(ownership.begin(), ownership.end(), 0U) != ownership.end()) { + return constraintFailure( + "invalid-constraint-order", + std::to_string(fullCount), + "Free and constrained DOFs must partition the complete full-DOF range."); + } + return Status::ok(); +} + +Result extractBlock( + const SparseMatrix& full, + const std::vector& rowDofs, + const std::vector& columnDofs) { + const std::size_t absent = (std::numeric_limits::max)(); + std::vector localColumn(full.columns(), absent); + for (std::size_t column = 0U; column < columnDofs.size(); ++column) { + localColumn[columnDofs[column]] = column; + } + + SparsePattern pattern; + pattern.rowOffsets.reserve(rowDofs.size() + 1U); + pattern.rowOffsets.push_back(0U); + std::vector contributions; + contributions.reserve(full.values().size()); + for (std::size_t localRow = 0U; + localRow < rowDofs.size(); + ++localRow) { + const std::size_t fullRow = rowDofs[localRow]; + for (std::size_t position = full.rowOffsets()[fullRow]; + position < full.rowOffsets()[fullRow + 1U]; + ++position) { + const std::size_t column = + localColumn[full.columnIndices()[position]]; + if (column == absent) { + continue; + } + pattern.columnIndices.push_back(column); + // One source CSR entry maps to one block slot, so exact numeric + // values and structural zeros survive without a new reduction. + contributions.push_back({ + localRow, + column, + full.values()[position], + localRow, + position}); + } + pattern.rowOffsets.push_back(pattern.columnIndices.size()); + } + return SparseMatrix::fromCoo( + rowDofs.size(), + columnDofs.size(), + std::move(contributions), + pattern); +} + +void requireDofOrder(const DofManager& dofs) { + if (!validateDofOrder(dofs).isOk()) { + throw std::invalid_argument{ + "DofManager constraint dimensions or order are invalid."}; + } +} + +} // namespace + +Result EssentialConstraints::partition( + const SparseMatrix& full, + const DofManager& dofs) { + const Status matrixStatus = full.validate(); + if (!matrixStatus.isOk()) { + return Result::failure(matrixStatus); + } + if (full.rows() != full.columns() || + full.rows() != dofs.fullDofCount()) { + return Result::failure(constraintFailure( + "invalid-constraint-dimensions", + std::to_string(full.rows()) + "x" + + std::to_string(full.columns()), + "Full stiffness must be square and match the DofManager full dimension.")); + } + const Status dofStatus = validateDofOrder(dofs); + if (!dofStatus.isOk()) { + return Result::failure(dofStatus); + } + + auto kff = extractBlock(full, dofs.freeDofs(), dofs.freeDofs()); + if (!kff.hasValue()) { + return Result::failure(kff.status()); + } + auto kfc = extractBlock(full, dofs.freeDofs(), dofs.constrainedDofs()); + if (!kfc.hasValue()) { + return Result::failure(kfc.status()); + } + auto kcf = extractBlock(full, dofs.constrainedDofs(), dofs.freeDofs()); + if (!kcf.hasValue()) { + return Result::failure(kcf.status()); + } + auto kcc = extractBlock( + full, dofs.constrainedDofs(), dofs.constrainedDofs()); + if (!kcc.hasValue()) { + return Result::failure(kcc.status()); + } + + return Result::success({ + std::move(kff.value()), + std::move(kfc.value()), + std::move(kcf.value()), + std::move(kcc.value())}); +} + +Vector EssentialConstraints::gatherFree( + const Vector& full, + const DofManager& dofs) { + requireDofOrder(dofs); + if (full.size() != dofs.fullDofCount()) { + throw std::invalid_argument{ + "Full vector size must match the DofManager full dimension."}; + } + Vector reduced{dofs.freeDofCount()}; + for (std::size_t equation = 0U; + equation < dofs.freeDofs().size(); + ++equation) { + reduced[equation] = full[dofs.freeDofs()[equation]]; + } + return reduced; +} + +Vector EssentialConstraints::gatherConstrained( + const Vector& full, + const DofManager& dofs) { + requireDofOrder(dofs); + if (full.size() != dofs.fullDofCount()) { + throw std::invalid_argument{ + "Full vector size must match the DofManager full dimension."}; + } + Vector reduced{dofs.constrainedDofCount()}; + for (std::size_t index = 0U; + index < dofs.constrainedDofs().size(); + ++index) { + reduced[index] = full[dofs.constrainedDofs()[index]]; + } + return reduced; +} + +Vector EssentialConstraints::reconstructFull( + const Vector& freeValues, + const Vector& constrainedValues, + const DofManager& dofs) { + requireDofOrder(dofs); + if (freeValues.size() != dofs.freeDofCount() || + constrainedValues.size() != dofs.constrainedDofCount()) { + throw std::invalid_argument{ + "Reduced vector sizes must match the DofManager order."}; + } + + Vector full{dofs.fullDofCount()}; + for (std::size_t equation = 0U; + equation < dofs.freeDofs().size(); + ++equation) { + full[dofs.freeDofs()[equation]] = freeValues[equation]; + } + // Preserve caller-supplied dc exactly; nonzero prescribed displacement is + // never replaced with an implicit homogeneous constraint. + for (std::size_t index = 0U; + index < dofs.constrainedDofs().size(); + ++index) { + full[dofs.constrainedDofs()[index]] = constrainedValues[index]; + } + return full; +} + +} // namespace fesa diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 1e6659b..256ab6b 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -7,6 +7,7 @@ add_executable( unit/analysis/analysis_state_test.cpp unit/assembly/parallel_for_test.cpp unit/assembly/sparse_assembler_test.cpp + unit/constraints/essential_constraints_test.cpp unit/core/diagnostic_test.cpp unit/core/source_identity_test.cpp unit/core/status_test.cpp diff --git a/tests/unit/constraints/essential_constraints_test.cpp b/tests/unit/constraints/essential_constraints_test.cpp new file mode 100644 index 0000000..71b7438 --- /dev/null +++ b/tests/unit/constraints/essential_constraints_test.cpp @@ -0,0 +1,243 @@ +#include "fesa/constraints/essential_constraints.hpp" + +#include "fesa/analysis/analysis_model.hpp" +#include "fesa/fem/dof_manager.hpp" +#include "fesa/model/domain.hpp" + +#include + +#include +#include +#include +#include +#include +#include + +namespace { + +fesa::DofManager makeDofs( + std::vector boundaries) { + const std::filesystem::path source{"models/essential-constraints.inp"}; + fesa::ModelDefinition definition{}; + definition.sourcePath = source; + definition.sourceContentIdentity = "fnv1a64:1234567890abcdef"; + definition.nodes = { + {{"Beam-1", 1, "1"}, {0.0, 0.0, 0.0}, {source, 2U}}}; + definition.steps = {{ + "Step-1", + std::move(boundaries), + {}, + 0.1, + 1.0, + 0.01, + 1.0, + {source, 10U}}}; + + auto domain = fesa::Domain::create(std::move(definition)); + EXPECT_TRUE(domain.hasValue()); + auto model = fesa::AnalysisModel::create(domain.value()); + EXPECT_TRUE(model.hasValue()); + auto dofs = fesa::DofManager::create(model.value()); + EXPECT_TRUE(dofs.hasValue()); + return std::move(dofs.value()); +} + +fesa::SparseMatrix makeMatrix( + const std::size_t rows, + const std::size_t columns, + const std::vector& denseValues) { + EXPECT_EQ(denseValues.size(), rows * columns); + fesa::SparsePattern pattern; + std::vector contributions; + pattern.rowOffsets.reserve(rows + 1U); + pattern.rowOffsets.push_back(0U); + for (std::size_t row = 0U; row < rows; ++row) { + for (std::size_t column = 0U; column < columns; ++column) { + pattern.columnIndices.push_back(column); + contributions.push_back({ + row, + column, + denseValues[row * columns + column], + row, + column}); + } + pattern.rowOffsets.push_back(pattern.columnIndices.size()); + } + + auto matrix = fesa::SparseMatrix::fromCoo( + rows, columns, std::move(contributions), pattern); + EXPECT_TRUE(matrix.hasValue()); + return std::move(matrix.value()); +} + +std::vector sequentialDense(const std::size_t size) { + std::vector values(size * size); + for (std::size_t row = 0U; row < size; ++row) { + for (std::size_t column = 0U; column < size; ++column) { + values[row * size + column] = + static_cast(row * 10U + column + 1U); + } + } + return values; +} + +void expectShape( + const fesa::SparseMatrix& matrix, + const std::size_t rows, + const std::size_t columns) { + EXPECT_EQ(matrix.rows(), rows); + EXPECT_EQ(matrix.columns(), columns); + EXPECT_TRUE(matrix.validate().isOk()); +} + +} // namespace + +TEST(EssentialConstraints, ExtractsHandComputedBlocksInStableOrder) { + const auto dofs = makeDofs({ + {"1", 2, 2, 2.5, {{}, 12U}}, + {"1", 5, 5, -3.25, {{}, 13U}}}); + auto fullValues = sequentialDense(6U); + fullValues[2U * 6U + 4U] = 0.0; + const auto full = makeMatrix(6U, 6U, fullValues); + + auto result = fesa::EssentialConstraints::partition(full, dofs); + ASSERT_TRUE(result.hasValue()); + const auto& blocks = result.value(); + + EXPECT_EQ(blocks.kff.rowOffsets(), (std::vector{0U, 4U, 8U, 12U, 16U})); + EXPECT_EQ( + blocks.kff.columnIndices(), + (std::vector{ + 0U, 1U, 2U, 3U, 0U, 1U, 2U, 3U, + 0U, 1U, 2U, 3U, 0U, 1U, 2U, 3U})); + EXPECT_EQ( + blocks.kff.values(), + (std::vector{ + 1.0, 3.0, 4.0, 6.0, + 21.0, 23.0, 24.0, 26.0, + 31.0, 33.0, 34.0, 36.0, + 51.0, 53.0, 54.0, 56.0})); + EXPECT_EQ(blocks.kfc.rowOffsets(), (std::vector{0U, 2U, 4U, 6U, 8U})); + EXPECT_EQ( + blocks.kfc.columnIndices(), + (std::vector{0U, 1U, 0U, 1U, 0U, 1U, 0U, 1U})); + EXPECT_EQ( + blocks.kfc.values(), + (std::vector{2.0, 5.0, 22.0, 0.0, 32.0, 35.0, 52.0, 55.0})); + EXPECT_EQ(blocks.kcf.rowOffsets(), (std::vector{0U, 4U, 8U})); + EXPECT_EQ( + blocks.kcf.columnIndices(), + (std::vector{0U, 1U, 2U, 3U, 0U, 1U, 2U, 3U})); + EXPECT_EQ( + blocks.kcf.values(), + (std::vector{11.0, 13.0, 14.0, 16.0, 41.0, 43.0, 44.0, 46.0})); + EXPECT_EQ(blocks.kcc.rowOffsets(), (std::vector{0U, 2U, 4U})); + EXPECT_EQ(blocks.kcc.columnIndices(), (std::vector{0U, 1U, 0U, 1U})); + EXPECT_EQ(blocks.kcc.values(), (std::vector{12.0, 15.0, 42.0, 45.0})); + + EXPECT_EQ(blocks.kfc.values()[3U], 0.0); + EXPECT_TRUE(blocks.kff.validate().isOk()); + EXPECT_TRUE(blocks.kfc.validate().isOk()); + EXPECT_TRUE(blocks.kcf.validate().isOk()); + EXPECT_TRUE(blocks.kcc.validate().isOk()); +} + +TEST(EssentialConstraints, HandlesNoAllAndMixedConstraints) { + const auto full = makeMatrix(6U, 6U, sequentialDense(6U)); + + const auto noConstraints = makeDofs({}); + auto none = fesa::EssentialConstraints::partition(full, noConstraints); + ASSERT_TRUE(none.hasValue()); + expectShape(none.value().kff, 6U, 6U); + expectShape(none.value().kfc, 6U, 0U); + expectShape(none.value().kcf, 0U, 6U); + expectShape(none.value().kcc, 0U, 0U); + EXPECT_EQ(none.value().kff.values(), full.values()); + + const auto allConstraints = makeDofs({{"1", 1, 6, 1.0, {{}, 12U}}}); + auto all = fesa::EssentialConstraints::partition(full, allConstraints); + ASSERT_TRUE(all.hasValue()); + expectShape(all.value().kff, 0U, 0U); + expectShape(all.value().kfc, 0U, 6U); + expectShape(all.value().kcf, 6U, 0U); + expectShape(all.value().kcc, 6U, 6U); + EXPECT_EQ(all.value().kcc.values(), full.values()); + + const auto mixedConstraints = makeDofs({{"1", 3, 4, 0.0, {{}, 12U}}}); + auto mixed = fesa::EssentialConstraints::partition(full, mixedConstraints); + ASSERT_TRUE(mixed.hasValue()); + expectShape(mixed.value().kff, 4U, 4U); + expectShape(mixed.value().kfc, 4U, 2U); + expectShape(mixed.value().kcf, 2U, 4U); + expectShape(mixed.value().kcc, 2U, 2U); +} + +TEST(EssentialConstraints, ReconstructsNonzeroPrescribedValues) { + const auto dofs = makeDofs({ + {"1", 2, 2, 2.5, {{}, 12U}}, + {"1", 5, 5, -3.25, {{}, 13U}}}); + fesa::Vector full{6U}; + full[0U] = 10.0; + full[1U] = 2.5; + full[2U] = 20.0; + full[3U] = 30.0; + full[4U] = -3.25; + full[5U] = 40.0; + + const auto free = fesa::EssentialConstraints::gatherFree(full, dofs); + const auto constrained = + fesa::EssentialConstraints::gatherConstrained(full, dofs); + EXPECT_EQ(free.size(), 4U); + EXPECT_DOUBLE_EQ(free[0U], 10.0); + EXPECT_DOUBLE_EQ(free[1U], 20.0); + EXPECT_DOUBLE_EQ(free[2U], 30.0); + EXPECT_DOUBLE_EQ(free[3U], 40.0); + EXPECT_EQ(constrained.size(), 2U); + EXPECT_DOUBLE_EQ(constrained[0U], 2.5); + EXPECT_DOUBLE_EQ(constrained[1U], -3.25); + EXPECT_EQ(constrained[0U], dofs.prescribedValues()[0U]); + EXPECT_EQ(constrained[1U], dofs.prescribedValues()[1U]); + + const auto reconstructed = fesa::EssentialConstraints::reconstructFull( + free, dofs.prescribedValues(), dofs); + ASSERT_EQ(reconstructed.size(), full.size()); + for (std::size_t index = 0U; index < full.size(); ++index) { + EXPECT_DOUBLE_EQ(reconstructed[index], full[index]); + } +} + +TEST(EssentialConstraints, RejectsDimensionOrOrderMismatch) { + const auto dofs = makeDofs({{"1", 2, 2, 1.0, {{}, 12U}}}); + const auto wrongSquare = makeMatrix(5U, 5U, sequentialDense(5U)); + auto wrongDimension = + fesa::EssentialConstraints::partition(wrongSquare, dofs); + ASSERT_FALSE(wrongDimension.hasValue()); + EXPECT_EQ( + wrongDimension.status().failureCategory(), + fesa::FailureCategory::model); + ASSERT_EQ(wrongDimension.status().diagnostics().size(), 1U); + EXPECT_EQ( + wrongDimension.status().diagnostics()[0U].code, + "invalid-constraint-dimensions"); + + const auto rectangular = makeMatrix( + 6U, 5U, std::vector(30U, 0.0)); + auto wrongOrder = fesa::EssentialConstraints::partition(rectangular, dofs); + ASSERT_FALSE(wrongOrder.hasValue()); + EXPECT_EQ( + wrongOrder.status().diagnostics()[0U].code, + "invalid-constraint-dimensions"); + + EXPECT_THROW( + static_cast(fesa::EssentialConstraints::gatherFree( + fesa::Vector{5U}, dofs)), + std::invalid_argument); + EXPECT_THROW( + static_cast(fesa::EssentialConstraints::gatherConstrained( + fesa::Vector{7U}, dofs)), + std::invalid_argument); + EXPECT_THROW( + static_cast(fesa::EssentialConstraints::reconstructFull( + fesa::Vector{4U}, fesa::Vector{2U}, dofs)), + std::invalid_argument); +}