From ebdb2519dc86eb1e279cfbddebaab0734634cd7a Mon Sep 17 00:00:00 2001 From: NINI Date: Mon, 4 May 2026 23:05:01 +0900 Subject: [PATCH] test: add mitc4 patch benchmark coverage --- PLAN.md | 8 +- PROGRESS.md | 28 ++- docs/VERIFICATION_PLAN.md | 13 ++ .../index.json | 2 +- tests/test_main.cpp | 207 ++++++++++++++++++ 5 files changed, 252 insertions(+), 6 deletions(-) diff --git a/PLAN.md b/PLAN.md index 9003ecf..edef15c 100644 --- a/PLAN.md +++ b/PLAN.md @@ -13,7 +13,7 @@ Every new agent session must read this file together with `PROGRESS.md` before p - If an item becomes obsolete, move it to `PROGRESS.md` with a short reason instead of silently deleting it. ## Current Objective -Continue the new Phase 1 rebaseline plan in `phases/1-linear-static-mitc4-rebaseline`, starting with P1R-11 MITC4 patch, locking-sensitivity, and benchmark tests. The old `phases/1-linear-static-mitc4` path is historical and superseded by the paper-based MITC4 formulation reset. +Continue the new Phase 1 rebaseline plan in `phases/1-linear-static-mitc4-rebaseline`, starting with P1R-12 assembly, solver-adapter boundary, constrained solve, and full-vector RF recovery revalidation. The old `phases/1-linear-static-mitc4` path is historical and superseded by the paper-based MITC4 formulation reset. ## Required Reading For New Agents 1. `AGENTS.md` @@ -36,7 +36,7 @@ Continue the new Phase 1 rebaseline plan in `phases/1-linear-static-mitc4-rebase ## Active Phase Files - Active phase directory: `phases/1-linear-static-mitc4-rebaseline` - Execute with: `python scripts/execute.py 1-linear-static-mitc4-rebaseline` -- Step numbering is zero-based. `step0.md` is complete and recorded in `phases/1-linear-static-mitc4-rebaseline/step0-audit.md`; `step1.md` is complete and created the `quad_02_phase1.inp` normalized reference path; `step2.md` is complete and revalidated core harness guardrails; `step3.md` is complete and revalidated the Phase 1 parser/domain subset; `step4.md` is complete and strengthened validation/singular diagnostics; `step5.md` is complete and revalidated the DofManager/reaction foundation; `step6.md` is complete and revalidated the minimum result model plus displacement CSV comparator; `step7.md` is complete and revalidated MITC4 natural coordinates, tying points, center directors, and integration bases; `step8.md` is complete and revalidated degenerated-continuum displacement, direct covariant strain rows, and MITC shear tying rows; `step9.md` is complete and revalidated plane-stress material, convected-to-local transform, and `2 x 2 x 2` material integration scaffolding; `step10.md` is complete and revalidated MITC4 stiffness, internal force, six-DOF transform, and drilling stabilization; `step15.md` is the independent evaluator closeout. +- Step numbering is zero-based. `step0.md` is complete and recorded in `phases/1-linear-static-mitc4-rebaseline/step0-audit.md`; `step1.md` is complete and created the `quad_02_phase1.inp` normalized reference path; `step2.md` is complete and revalidated core harness guardrails; `step3.md` is complete and revalidated the Phase 1 parser/domain subset; `step4.md` is complete and strengthened validation/singular diagnostics; `step5.md` is complete and revalidated the DofManager/reaction foundation; `step6.md` is complete and revalidated the minimum result model plus displacement CSV comparator; `step7.md` is complete and revalidated MITC4 natural coordinates, tying points, center directors, and integration bases; `step8.md` is complete and revalidated degenerated-continuum displacement, direct covariant strain rows, and MITC shear tying rows; `step9.md` is complete and revalidated plane-stress material, convected-to-local transform, and `2 x 2 x 2` material integration scaffolding; `step10.md` is complete and revalidated MITC4 stiffness, internal force, six-DOF transform, and drilling stabilization; `step11.md` is complete and added MITC4 membrane, bending, shear, twist, drilling-sensitivity, and thin-cantilever locking-sensitivity tests; `step15.md` is the independent evaluator closeout. - Every step file contains a sprint contract with objective, required reading, scope, allowed files, explicit non-goals, tests to write first, reference artifacts, acceptance command, evaluator checklist, and handoff requirements. - Historical phase directory: `phases/1-linear-static-mitc4` - Historical phase status: blocked/superseded. Do not resume the old P1-15/P1-16 path unless the user explicitly requests recovery of that exact phase. @@ -76,7 +76,7 @@ Each gate should be satisfied before moving to the next implementation band unle | G1 - Build and validation | satisfied | Build system, test framework, and `scripts/validate_workspace.py` run real checks. | Validation command output | | G2 - Parser and domain | satisfied | Parser subset revalidated in step 3; validation and singular diagnostics revalidated in step 4. | Parser acceptance/rejection tests, validation negative tests, and validation output | | G3 - DOF/math/results infrastructure | partial | Core aliases, DOF mapping, validation harness, model diagnostic context, DofManager, sparse-connectivity inputs, full-vector reaction formula, result model metadata, and displacement CSV comparator were revalidated in steps 2, 5, and 6; assembly remains for step 12. | P1R-02, P1R-05, and P1R-06 validation output | -| G4 - MITC4 element readiness | partial | MITC4 formulation was rewritten from local papers; Steps 7 through 10 rebuilt geometry/director/local-basis scaffolding, displacement interpolation, direct covariant strain rows, MITC shear tying rows, plane-stress material, convected-to-local transform, `2 x 2 x 2` material integration scaffolding, stiffness/internal force, six-DOF transform, and drilling stabilization, while patch/benchmark coverage remains for step 11. | P1R-07 through P1R-10 validation output; future element tests | +| G4 - MITC4 element readiness | satisfied | MITC4 formulation was rewritten from local papers; Steps 7 through 11 rebuilt geometry/director/local-basis scaffolding, displacement interpolation, direct covariant strain rows, MITC shear tying rows, plane-stress material, convected-to-local transform, `2 x 2 x 2` material integration scaffolding, stiffness/internal force, six-DOF transform, drilling stabilization, and patch/locking-sensitivity tests. | P1R-07 through P1R-11 validation output | | G5 - End-to-end solver | reopened | Linear static path must be revalidated through steps 13 and 14 after the MITC4 rebuild and `quad_02` compatibility path. | Future integration/reference regression output | ## Phase 1 Implementation Milestones @@ -92,7 +92,7 @@ All milestones are intended to become one or more self-contained sprint contract | P1R-08 | completed | MITC4 generator | Implement degenerated-continuum displacement, covariant strain rows, and MITC shear tying. | P1R-07 | Finite-difference and tying interpolation tests | | P1R-09 | completed | MITC4 generator | Implement material matrix, transform, and `2 x 2 x 2` integration scaffolding. | P1R-08 | Material/integration tests | | P1R-10 | completed | MITC4 generator | Assemble MITC4 stiffness/internal force with six-DOF transform and drilling stabilization. | P1R-09, P1R-05 | Symmetry, rigid body, drilling sensitivity tests | -| P1R-11 | pending | verification generator | Add MITC4 patch, locking-sensitivity, and benchmark tests. | P1R-10 | Membrane/bending/shear/twist/locking tests | +| P1R-11 | completed | verification generator | Add MITC4 patch, locking-sensitivity, and benchmark tests. | P1R-10 | Membrane/bending/shear/twist/locking tests | | P1R-12 | pending | assembly generator | Rebuild assembly, solver adapter boundary, constrained solve, and full-vector RF recovery. | P1R-05, P1R-10 | Assembly and full-vector reaction tests | | P1R-13 | pending | analysis generator | Rebuild linear static workflow from input to U/RF result fields. | P1R-03, P1R-04, P1R-06, P1R-12 | End-to-end linear static tests | | P1R-14 | pending | reference generator | Run stored reference displacement regression using accepted Phase 1-compatible S4 cases. | P1R-13 | At least one automated CSV displacement regression | diff --git a/PROGRESS.md b/PROGRESS.md index 0c30401..ddbac24 100644 --- a/PROGRESS.md +++ b/PROGRESS.md @@ -13,10 +13,36 @@ Every new agent session must read this file together with `PLAN.md` before plann - Do not remove history unless the user explicitly asks for archival cleanup. ## Current Status -Phase 1 has a new rebaseline phase definition in `phases/1-linear-static-mitc4-rebaseline`. Steps 0 through 10 are complete. `quad_02_phase1.inp` is now the normalized Phase 1-compatible input path for the stored `quad_02` S4 reference pair, while the original `quad_02.inp` remains preserved unsupported provenance. Core numeric aliases, DOF mapping, validation harness, model diagnostic context, the Phase 1 parser/domain subset, validation/singular diagnostics, DofManager/reaction foundation, minimum result model metadata, displacement CSV comparator foundation, MITC4 geometry/director scaffolding, MITC4 displacement/strain/tying row scaffolding, MITC4 material/transform/integration scaffolding, and MITC4 stiffness/drilling/internal-force scaffolding have been revalidated. The old `phases/1-linear-static-mitc4` path is historical and superseded after the MITC4 formulation reset. +Phase 1 has a new rebaseline phase definition in `phases/1-linear-static-mitc4-rebaseline`. Steps 0 through 11 are complete. `quad_02_phase1.inp` is now the normalized Phase 1-compatible input path for the stored `quad_02` S4 reference pair, while the original `quad_02.inp` remains preserved unsupported provenance. Core numeric aliases, DOF mapping, validation harness, model diagnostic context, the Phase 1 parser/domain subset, validation/singular diagnostics, DofManager/reaction foundation, minimum result model metadata, displacement CSV comparator foundation, MITC4 geometry/director scaffolding, MITC4 displacement/strain/tying row scaffolding, MITC4 material/transform/integration scaffolding, MITC4 stiffness/drilling/internal-force scaffolding, and MITC4 patch/locking-sensitivity tests have been revalidated. The old `phases/1-linear-static-mitc4` path is historical and superseded after the MITC4 formulation reset. ## Completed Work +### 2026-05-04 - P1R-11 MITC4 patch and benchmark tests completed +Author: Codex + +Changed files: +- `tests/test_main.cpp` +- `docs/VERIFICATION_PLAN.md` +- `phases/1-linear-static-mitc4-rebaseline/index.json` +- `PLAN.md` +- `PROGRESS.md` + +Summary: +- Added reusable MITC4 test helpers for analytic displacement/rotation fields, local strain sampling through the tied strain-row path, and a one-element cantilever smoke model. +- Added constant membrane patch tests that verify uniform local `[eps11, eps22, gamma12]` at all `2 x 2 x 2` samples and positive physical energy. +- Added pure bending, pure transverse shear, and pure twist patch tests to exercise through-thickness bending behavior, MITC transverse shear tying, and bilinear twist/shear interpolation. +- Added a drilling-sensitivity negative check proving membrane patch energy is not hidden by the artificial drilling scale when no drilling DOF participates. +- Added a thin one-element cantilever strip test that verifies tip displacement grows materially as thickness decreases, providing early shear-locking sensitivity coverage before stored Abaqus reference regression. +- Documented the Step 11 analytic patch layer and kept Scordelis-Lo as a deferred benchmark scaffold until pressure-load support and a stored curved-shell reference artifact are defined. + +Verification: +- `python scripts/validate_workspace.py` configured CMake, built `fesa_core` and `fesa_tests`, and ran CTest successfully. +- CTest result: 1 test executable passed. + +Follow-up: +- Continue with P1R-12 assembly, solver-adapter boundary, constrained solve, and full-vector RF recovery revalidation. +- Step 11 intentionally did not execute Abaqus, add pressure loads, add stored CSV reference regression, or tune benchmark tolerances against Abaqus artifacts. + ### 2026-05-04 - P1R-10 MITC4 stiffness and drilling completed Author: Codex diff --git a/docs/VERIFICATION_PLAN.md b/docs/VERIFICATION_PLAN.md index 5cb1805..a478a8f 100644 --- a/docs/VERIFICATION_PLAN.md +++ b/docs/VERIFICATION_PLAN.md @@ -205,6 +205,19 @@ Final benchmark tolerances must be stored with each reference case. | Pinched cylinder | Curved shell bending/membrane behavior | loaded-point displacement | | Twisted beam | Warped shell and drilling sensitivity | tip displacement/rotation | +### Step 11 Analytic Patch And Benchmark Scaffold +The first MITC4 patch-test layer is kept local and analytic so it does not depend on Abaqus execution or stored CSV artifacts. + +Current in-repository coverage: +- Constant membrane patch: verifies uniform local `[eps11, eps22, gamma12]` strain at all `2 x 2 x 2` samples and positive physical energy. +- Pure bending patch: verifies membrane-free midsurface bending behavior, through-thickness antisymmetry, and positive physical bending energy. +- Pure transverse shear patch: verifies constant `gamma13` reproduction through the MITC tying path. +- Pure twist patch: verifies bilinear transverse-shear reproduction for `w = twist * x * y`, covering the locking-sensitive twist/shear interpolation path before curved benchmarks. +- Drilling sensitivity negative check: verifies membrane patch energy is unchanged by drilling scale changes when no drilling DOF participates. +- Thin one-element cantilever strip: verifies tip displacement increases materially as thickness decreases, providing an early shear-locking sensitivity smoke test rather than a final convergence benchmark. + +Scordelis-Lo roof remains a named benchmark scaffold for Phase 1 planning, but it is not executable in the current Step 11 test layer because Phase 1 does not yet support pressure loads and does not yet have an accepted curved-shell Abaqus reference artifact. The first executable Scordelis-Lo sprint must define load representation, normalized input compatibility, stored displacement reference data, and scale-aware displacement tolerances before it is used as an automated acceptance case. + ## Minimum Phase 1 Acceptance Before MITC4 Phase 1 is considered credible: - All parser smoke tests pass. diff --git a/phases/1-linear-static-mitc4-rebaseline/index.json b/phases/1-linear-static-mitc4-rebaseline/index.json index 594ff69..396fc85 100644 --- a/phases/1-linear-static-mitc4-rebaseline/index.json +++ b/phases/1-linear-static-mitc4-rebaseline/index.json @@ -13,7 +13,7 @@ { "step": 8, "name": "mitc4-covariant-strain-tying", "status": "completed" }, { "step": 9, "name": "mitc4-material-integration", "status": "completed" }, { "step": 10, "name": "mitc4-stiffness-drilling", "status": "completed" }, - { "step": 11, "name": "mitc4-patch-benchmark-tests", "status": "pending" }, + { "step": 11, "name": "mitc4-patch-benchmark-tests", "status": "completed" }, { "step": 12, "name": "assembly-sparse-solver-path", "status": "pending" }, { "step": 13, "name": "linear-static-workflow", "status": "pending" }, { "step": 14, "name": "stored-reference-regression", "status": "pending" }, diff --git a/tests/test_main.cpp b/tests/test_main.cpp index 09c992a..047fef2 100644 --- a/tests/test_main.cpp +++ b/tests/test_main.cpp @@ -160,6 +160,76 @@ fesa::Real quadraticEnergy(const fesa::DenseMatrix& stiffness, const std::vector return energy; } +std::array flatUnitSquareCoordinates() { + return {fesa::Vec3{0.0, 0.0, 0.0}, + fesa::Vec3{1.0, 0.0, 0.0}, + fesa::Vec3{1.0, 1.0, 0.0}, + fesa::Vec3{0.0, 1.0, 0.0}}; +} + +using MITC4DofField = std::function; + +fesa::MITC4ElementDofVector elementDofsFromFields(const std::array& coords, + MITC4DofField translation, + MITC4DofField rotation) { + fesa::MITC4ElementDofVector values{}; + for (std::size_t node = 0; node < coords.size(); ++node) { + const auto displacement = translation(coords[node]); + const auto rotational_dof = rotation(coords[node]); + values[node * 6 + 0] = displacement.x; + values[node * 6 + 1] = displacement.y; + values[node * 6 + 2] = displacement.z; + values[node * 6 + 3] = rotational_dof.x; + values[node * 6 + 4] = rotational_dof.y; + values[node * 6 + 5] = rotational_dof.z; + } + return values; +} + +std::vector toVector(const fesa::MITC4ElementDofVector& values) { + return {values.begin(), values.end()}; +} + +fesa::MITC4StrainVector localStrainAt(const fesa::MITC4Geometry& geometry, + const fesa::MITC4ElementDofVector& values, + fesa::Real xi, + fesa::Real eta, + fesa::Real zeta) { + const auto rows = fesa::mitc4TiedCovariantStrainRows(geometry, xi, eta, zeta); + FESA_CHECK(rows.ok()); + const auto convected_strain = fesa::evaluateMITC4StrainRows(rows, values); + FESA_CHECK(convected_strain.ok()); + const auto basis = fesa::computeMITC4IntegrationBasis(geometry, xi, eta, zeta); + FESA_CHECK(basis.ok()); + const auto transform = fesa::mitc4CovariantToLocalStrainTransform(basis); + FESA_CHECK(transform.ok()); + return fesa::multiplyMITC4MaterialMatrix(transform.matrix, convected_strain.values); +} + +fesa::Real singleElementCantileverTipUz(fesa::Real thickness) { + fesa::Domain domain; + domain.nodes[1] = {1, {0.0, 0.0, 0.0}}; + domain.nodes[2] = {2, {4.0, 0.0, 0.0}}; + domain.nodes[3] = {3, {4.0, 1.0, 0.0}}; + domain.nodes[4] = {4, {0.0, 1.0, 0.0}}; + domain.elements[1] = {1, fesa::ElementType::MITC4, {1, 2, 3, 4}, "SHELLS"}; + domain.node_sets["clamped"] = {"CLAMPED", {1, 4}}; + domain.element_sets["shells"] = {"SHELLS", {1}}; + domain.materials["shell_steel"] = {"SHELL_STEEL", 210000.0, 0.30}; + domain.shell_sections.push_back({"SHELLS", "SHELL_STEEL", thickness}); + domain.boundary_conditions.push_back({"CLAMPED", 1, 6, 0.0}); + domain.loads.push_back({"2", 3, -0.5}); + domain.loads.push_back({"3", 3, -0.5}); + + fesa::LinearStaticAnalysis analysis; + const auto result = analysis.run(domain); + FESA_CHECK(result.ok()); + fesa::DofManager dofs(domain); + const auto node2_uz = result.state.u_full[static_cast(dofs.fullIndex(2, fesa::Dof::UZ))]; + const auto node3_uz = result.state.u_full[static_cast(dofs.fullIndex(3, fesa::Dof::UZ))]; + return 0.5 * (node2_uz + node3_uz); +} + fesa::Domain singleElementValidationDomain() { fesa::Domain domain; domain.nodes[1] = {1, {0, 0, 0}}; @@ -1316,6 +1386,143 @@ FESA_TEST(mitc4_rigid_translation_energy_is_zero_and_drilling_energy_is_document FESA_CHECK_NEAR(quadraticEnergy(result.global, drilling), 4.0 * result.drilling_stiffness, 1.0e-8); } +FESA_TEST(mitc4_constant_membrane_patch_strain_is_uniform) { + const auto coords = flatUnitSquareCoordinates(); + const auto geometry = fesa::buildMITC4Geometry(coords, 0.2); + FESA_CHECK(geometry.ok()); + const fesa::Real eps_x = 0.010; + const fesa::Real eps_y = -0.004; + const fesa::Real gamma_xy = 0.006; + const auto dofs = elementDofsFromFields( + coords, + [&](const fesa::Vec3& p) { + return fesa::Vec3{eps_x * p.x + 0.5 * gamma_xy * p.y, + eps_y * p.y + 0.5 * gamma_xy * p.x, + 0.0}; + }, + [](const fesa::Vec3&) { return fesa::Vec3{}; }); + + for (const auto& sample : fesa::mitc4GaussQuadrature2x2x2()) { + const auto strain = localStrainAt(geometry, dofs, sample.xi, sample.eta, sample.zeta); + FESA_CHECK_NEAR(strain[0], eps_x, 1.0e-12); + FESA_CHECK_NEAR(strain[1], eps_y, 1.0e-12); + FESA_CHECK_NEAR(strain[2], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[3], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[4], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[5], gamma_xy, 1.0e-12); + } + + const auto stiffness = fesa::mitc4ElementStiffness(coords, 210000.0, 0.30, 0.2); + FESA_CHECK(stiffness.ok()); + const auto energy = quadraticEnergy(stiffness.local_without_drilling, toVector(dofs)); + FESA_CHECK(energy > 0.0); +} + +FESA_TEST(mitc4_pure_bending_patch_is_membrane_free_and_thickness_antisymmetric) { + const auto coords = flatUnitSquareCoordinates(); + const fesa::Real thickness = 0.2; + const auto geometry = fesa::buildMITC4Geometry(coords, thickness); + FESA_CHECK(geometry.ok()); + const fesa::Real curvature_x = 0.025; + const auto dofs = elementDofsFromFields( + coords, + [](const fesa::Vec3&) { return fesa::Vec3{}; }, + [&](const fesa::Vec3& p) { return fesa::Vec3{0.0, curvature_x * p.x, 0.0}; }); + + const fesa::Real xi = -0.30; + const fesa::Real eta = 0.45; + const fesa::Real zeta = 1.0 / std::sqrt(3.0); + const auto top = localStrainAt(geometry, dofs, xi, eta, zeta); + const auto mid = localStrainAt(geometry, dofs, xi, eta, 0.0); + const auto bottom = localStrainAt(geometry, dofs, xi, eta, -zeta); + + FESA_CHECK(std::fabs(top[0]) > 1.0e-6); + FESA_CHECK_NEAR(top[0] + bottom[0], 0.0, 1.0e-12); + FESA_CHECK_NEAR(mid[0], 0.0, 1.0e-12); + FESA_CHECK_NEAR(top[1] + bottom[1], 0.0, 1.0e-12); + FESA_CHECK_NEAR(mid[1], 0.0, 1.0e-12); + + const auto stiffness = fesa::mitc4ElementStiffness(coords, 210000.0, 0.30, thickness); + FESA_CHECK(stiffness.ok()); + const auto energy = quadraticEnergy(stiffness.local_without_drilling, toVector(dofs)); + FESA_CHECK(energy > 0.0); +} + +FESA_TEST(mitc4_pure_transverse_shear_patch_is_constant) { + const auto coords = flatUnitSquareCoordinates(); + const auto geometry = fesa::buildMITC4Geometry(coords, 0.2); + FESA_CHECK(geometry.ok()); + const fesa::Real gamma_xz = 0.012; + const auto dofs = elementDofsFromFields( + coords, + [&](const fesa::Vec3& p) { return fesa::Vec3{0.0, 0.0, gamma_xz * p.x}; }, + [](const fesa::Vec3&) { return fesa::Vec3{}; }); + + for (const auto& sample : fesa::mitc4GaussQuadrature2x2x2()) { + const auto strain = localStrainAt(geometry, dofs, sample.xi, sample.eta, sample.zeta); + FESA_CHECK_NEAR(strain[0], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[1], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[2], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[3], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[4], gamma_xz, 1.0e-12); + FESA_CHECK_NEAR(strain[5], 0.0, 1.0e-12); + } +} + +FESA_TEST(mitc4_pure_twist_patch_reproduces_bilinear_transverse_shear) { + const auto coords = flatUnitSquareCoordinates(); + const auto geometry = fesa::buildMITC4Geometry(coords, 0.2); + FESA_CHECK(geometry.ok()); + const fesa::Real twist = 0.018; + const auto dofs = elementDofsFromFields( + coords, + [&](const fesa::Vec3& p) { return fesa::Vec3{0.0, 0.0, twist * p.x * p.y}; }, + [](const fesa::Vec3&) { return fesa::Vec3{}; }); + + for (const auto& sample : fesa::mitc4GaussQuadrature2x2x2()) { + const auto x = 0.5 * (sample.xi + 1.0); + const auto y = 0.5 * (sample.eta + 1.0); + const auto strain = localStrainAt(geometry, dofs, sample.xi, sample.eta, sample.zeta); + FESA_CHECK_NEAR(strain[0], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[1], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[2], 0.0, 1.0e-12); + FESA_CHECK_NEAR(strain[3], twist * x, 1.0e-12); + FESA_CHECK_NEAR(strain[4], twist * y, 1.0e-12); + FESA_CHECK_NEAR(strain[5], 0.0, 1.0e-12); + } +} + +FESA_TEST(mitc4_membrane_patch_energy_is_not_hidden_by_drilling_scale) { + const auto coords = flatUnitSquareCoordinates(); + const auto dofs = elementDofsFromFields( + coords, + [](const fesa::Vec3& p) { return fesa::Vec3{0.003 * p.x, -0.002 * p.y, 0.0}; }, + [](const fesa::Vec3&) { return fesa::Vec3{}; }); + + fesa::ElementStiffnessOptions no_drilling; + no_drilling.drilling_stiffness_scale = 0.0; + fesa::ElementStiffnessOptions strong_drilling; + strong_drilling.drilling_stiffness_scale = 1.0e-1; + const auto physical = fesa::mitc4ElementStiffness(coords, 210000.0, 0.30, 0.2, no_drilling); + const auto stabilized = fesa::mitc4ElementStiffness(coords, 210000.0, 0.30, 0.2, strong_drilling); + FESA_CHECK(physical.ok()); + FESA_CHECK(stabilized.ok()); + + const auto physical_energy = quadraticEnergy(physical.local_with_drilling, toVector(dofs)); + const auto stabilized_energy = quadraticEnergy(stabilized.local_with_drilling, toVector(dofs)); + FESA_CHECK(physical_energy > 0.0); + FESA_CHECK_NEAR(stabilized_energy, physical_energy, physical_energy * 1.0e-12); +} + +FESA_TEST(mitc4_single_element_cantilever_is_thickness_sensitive) { + const auto thick_tip = singleElementCantileverTipUz(0.20); + const auto thin_tip = singleElementCantileverTipUz(0.10); + FESA_CHECK(thick_tip < 0.0); + FESA_CHECK(thin_tip < thick_tip); + const auto response_ratio = std::fabs(thin_tip / thick_tip); + FESA_CHECK(response_ratio > 2.0); +} + FESA_TEST(mitc4_internal_force_is_stiffness_times_displacement_for_linear_phase1) { const std::array coords = {{{0, 0, 0}, {1, 0, 0}, {1, 1, 0}, {0, 1, 0}}}; const auto result = fesa::mitc4ElementStiffness(coords, 1000.0, 0.25, 0.2);