diff --git a/include/numsim-materials/materials/isotropic_damage.h b/include/numsim-materials/materials/isotropic_damage.h index 889da36..d061a93 100644 --- a/include/numsim-materials/materials/isotropic_damage.h +++ b/include/numsim-materials/materials/isotropic_damage.h @@ -24,6 +24,8 @@ namespace numsim::materials { /// damage_source::damage, damage_source::d_damage — from propagation law /// state_source::d_equivalent_strain — from state function /// yield_source::is_yielding — from yield function +/// The tangent may come from a different material than the stress: set the +/// optional "tangent_source". Absent, both come from "elastic_source". template class isotropic_damage final : public material_base, Traits> { @@ -50,8 +52,13 @@ class isotropic_damage final // inputs m_stress(base::template add_input( m_elastic_source, "stress", EdgeKind::Global)), + // Absent: same material as the stress. contains() rather than an + // empty-string test, so "" is an error and not a fallback. m_tangent(base::template add_input( - m_elastic_source, "tangent", EdgeKind::Global)), + base::m_parameter_handler.contains("tangent_source") + ? base::template get_parameter("tangent_source") + : m_elastic_source, + "tangent", EdgeKind::Global)), m_damage(base::template add_input( m_damage_source, "damage", EdgeKind::Global)), m_d_damage(base::template add_input( @@ -65,6 +72,8 @@ class isotropic_damage final static input_parameter_controller parameters() { input_parameter_controller para{base::parameters()}; para.template insert("elastic_source").template add(); + // No check == optional; declared so the JSON schema still knows the key. + para.template insert("tangent_source"); para.template insert("damage_source").template add(); para.template insert("state_source").template add(); para.template insert("yield_source").template add(); diff --git a/include/numsim-materials/materials/weighted_sum.h b/include/numsim-materials/materials/weighted_sum.h index 3c19b29..210eafd 100644 --- a/include/numsim-materials/materials/weighted_sum.h +++ b/include/numsim-materials/materials/weighted_sum.h @@ -1,6 +1,7 @@ #ifndef WEIGHTED_SUM_H #define WEIGHTED_SUM_H +#include #include #include "numsim-materials/core/material_base.h" @@ -18,6 +19,13 @@ namespace numsim::materials { /// Parameters: /// "name": material name /// "terms": vector> — (weight_name, constituent_name) pairs +/// "tangent_sources": optional vector — one per term, naming a +/// different material for that term's tangent; "" keeps the term's own. +/// +/// One entry per term buys ALIGNMENT, not correctness: an entry naming the +/// wrong material wires and runs, and shows only as a wrong summed tangent and +/// slow Newton convergence. Uncheckable, since decoupling stress from tangent +/// is the point. template class weighted_sum final : public material_base, Traits> { @@ -38,20 +46,43 @@ class weighted_sum final m_terms_param(base::template get_parameter("terms")), m_weight_property(base::template get_parameter("weight_property")), m_stress_property(base::template get_parameter("stress_property")), - m_tangent_property(base::template get_parameter("tangent_property")) + m_tangent_property(base::template get_parameter("tangent_property")), + m_tangent_sources(base::m_parameter_handler.contains("tangent_sources") + ? base::template get_parameter>("tangent_sources") + : std::vector{}) { + // Positional: a shorter list shifts every override onto the wrong + // constituent, and both names still resolve. One entry per term, with "" + // meaning "keeps its own" so a non-leading term can be overridden alone. + if (!m_tangent_sources.empty() && + m_tangent_sources.size() != m_terms_param.size()) + throw std::runtime_error( + "weighted_sum '" + base::name() + "': tangent_sources has " + + std::to_string(m_tangent_sources.size()) + " entries but there are " + + std::to_string(m_terms_param.size()) + + " terms — supply one per term (\"\" keeps a term's own tangent) or " + "omit it entirely"); + // Dynamically create inputs for each term + std::size_t i = 0; for (const auto& [weight_name, mat_name] : m_terms_param) { + const bool overridden = + i < m_tangent_sources.size() && !m_tangent_sources[i].empty(); + const std::string& tangent_owner = + overridden ? m_tangent_sources[i] : mat_name; auto& w = base::template add_input(weight_name, m_weight_property, EdgeKind::Global); auto& s = base::template add_input(mat_name, m_stress_property, EdgeKind::Global); - auto& c = base::template add_input(mat_name, m_tangent_property, EdgeKind::Global); + auto& c = base::template add_input(tangent_owner, m_tangent_property, EdgeKind::Global); m_terms.push_back({&w, &s, &c}); + ++i; } - } static input_parameter_controller parameters() { input_parameter_controller para{base::parameters()}; + // Optional (no check): absent means every term uses its stress material. + // If given, one entry per term; "" keeps that term's own tangent. + para.template insert>("tangent_sources"); para.template insert("terms").template add(); para.template insert("weight_property") .template add(std::string{"value"}); @@ -62,10 +93,9 @@ class weighted_sum final return para; } - void update() override { - update_stress(); - update_tangent(); - } + // No update() override: the engine drives each output through its add_output + // callback, so one would never run and would strand any output wired only + // into it. (backward_euler/vector_newton route update() themselves.) void update_stress() noexcept { m_stress = tensor2{}; @@ -94,6 +124,7 @@ class weighted_sum final const std::string& m_weight_property; const std::string& m_stress_property; const std::string& m_tangent_property; + const std::vector m_tangent_sources; std::vector m_terms; }; diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index d03b294..1f94db2 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -24,6 +24,7 @@ add_numsim_test(test_material_point_evaluator test_material_point_evaluator.cpp) add_numsim_test(test_plane_stress_evaluator test_plane_stress_evaluator.cpp) add_numsim_test(test_umat_interface test_umat_interface.cpp) add_numsim_test(test_tangent_generator test_tangent_generator.cpp) +add_numsim_test(test_weighted_sum test_weighted_sum.cpp) target_link_libraries(test_umat_interface PRIVATE Threads::Threads) # Data dumper for plotting (not a test — standalone executable) diff --git a/tests/test_weighted_sum.cpp b/tests/test_weighted_sum.cpp new file mode 100644 index 0000000..fbc194d --- /dev/null +++ b/tests/test_weighted_sum.cpp @@ -0,0 +1,213 @@ +#include +#include +#include +#include +#include +#include "numsim-materials/core/material_context.h" +#include "numsim-materials/materials/constant_scalar.h" +#include "numsim-materials/materials/isotropic_tangent.h" +#include "numsim-materials/materials/linear_elasticity.h" +#include "numsim-materials/materials/linear_stress.h" +#include "numsim-materials/materials/scalar_identity_weight.h" +#include "numsim-materials/materials/scalar_stepper.h" +#include "numsim-materials/materials/weighted_sum.h" +#include "numsim-materials/umat/external_state_source.h" + +namespace { + +namespace nm = numsim::materials; + +using policy = nm::material_policy_default; +using T = policy::value_type; +using ctx_type = nm::material_context; +using param_type = policy::ParameterHandler; +using tensor2 = tmech::tensor; +using tensor4 = tmech::tensor; +using terms_type = std::vector>; + +constexpr T KA = 100.0, GA = 40.0; +constexpr T KB = 300.0, GB = 140.0; + +/// A weight in [0,1] published as "value", via scalar_identity_weight over a +/// scalar_stepper's history. +void add_weight(ctx_type& ctx, const std::string& name, T increment) { + param_type p; + p.insert("name", name + "_drv"); + p.insert("increment", increment); + ctx.create>(p); + p.clear(); + p.insert("name", name); + p.insert("source", name + "_drv::state"); + ctx.create>(p); +} + +void add_monolithic(ctx_type& ctx, const std::string& name, T k, T g) { + param_type p; + p.insert("name", name); + p.insert("strain_producer_name", "strain_in"); + p.insert("K", k); + p.insert("G", g); + ctx.create>(p); +} + +/// A constituent whose stiffness lives in its OWN material, so its stress and +/// tangent come from two different names. +void add_decomposed(ctx_type& ctx, const std::string& name, T k, T g) { + param_type p; + p.insert("name", name + "_K"); + p.insert("value", k); + ctx.create>(p); + p.clear(); + p.insert("name", name + "_G"); + p.insert("value", g); + ctx.create>(p); + p.clear(); + p.insert("name", name + "_stiff"); + p.insert("K_source", name + "_K"); + p.insert("G_source", name + "_G"); + ctx.create>(p); + p.clear(); + p.insert("name", name); + p.insert("tangent_source", name + "_stiff"); + p.insert("strain_source", "strain_in"); + ctx.create>(p); +} + +nm::external_strain_source& add_strain(ctx_type& ctx) { + param_type p; + p.insert("name", "strain_in"); + return ctx.create>(p); +} + +tensor2 uniaxial(T v) { + tensor2 e; + e.fill(0.0); + e(0, 0) = v; + return e; +} + +// --------------------------------------------------------------------------- +// Baseline: the mixture rule itself +// --------------------------------------------------------------------------- + +/// Two monolithic constituents at weights 0.25 and 0.5. Both stress and tangent +/// must be the weighted sum. +TEST(WeightedSum, SumsStressAndTangentByWeight) { + ctx_type ctx; + auto& src = add_strain(ctx); + add_weight(ctx, "wA", 0.25); + add_weight(ctx, "wB", 0.5); + add_monolithic(ctx, "matA", KA, GA); + add_monolithic(ctx, "matB", KB, GB); + + param_type p; + p.insert("name", "mix"); + p.insert("terms", {{"wA", "matA"}, {"wB", "matB"}}); + ctx.create>(p); + ctx.finalize(); + + src.bind(uniaxial(0.001), uniaxial(0.001)); + ctx.update(); + + const T cA = KA + 4.0 * GA / 3.0; + const T cB = KB + 4.0 * GB / 3.0; + EXPECT_NEAR(ctx.get("mix", "stress")(0, 0), + (0.25 * cA + 0.5 * cB) * 0.001, 1e-10); + EXPECT_NEAR(ctx.get("mix", "tangent")(0, 0, 0, 0), + 0.25 * cA + 0.5 * cB, 1e-9); +} + +// --------------------------------------------------------------------------- +// tangent_sources +// --------------------------------------------------------------------------- + +/// Absent: every term takes its tangent from the material producing its stress. +TEST(WeightedSum, AbsentTangentSourcesUsesEachTermsOwnMaterial) { + ctx_type ctx; + auto& src = add_strain(ctx); + add_weight(ctx, "wA", 1.0); + add_monolithic(ctx, "matA", KA, GA); + + param_type p; + p.insert("name", "mix"); + p.insert("terms", {{"wA", "matA"}}); + ctx.create>(p); + ctx.finalize(); + + src.bind(uniaxial(0.001), uniaxial(0.001)); + ctx.update(); + EXPECT_NEAR(ctx.get("mix", "tangent")(0, 0, 0, 0), + KA + 4.0 * GA / 3.0, 1e-9); +} + +/// What tangent_sources exists for: a constituent whose stress and tangent have +/// different owners. +TEST(WeightedSum, OverridesOneTermsTangentOwner) { + ctx_type ctx; + auto& src = add_strain(ctx); + add_weight(ctx, "wA", 1.0); + add_decomposed(ctx, "matA", KA, GA); + + param_type p; + p.insert("name", "mix"); + p.insert("terms", {{"wA", "matA"}}); + p.insert>("tangent_sources", {"matA_stiff"}); + ctx.create>(p); + ctx.finalize(); + + src.bind(uniaxial(0.001), uniaxial(0.001)); + ctx.update(); + EXPECT_NEAR(ctx.get("mix", "tangent")(0, 0, 0, 0), + KA + 4.0 * GA / 3.0, 1e-9); +} + +/// An empty entry keeps that term's own tangent, so a NON-LEADING term can be +/// overridden alone; without it a positional list could only override a prefix. +TEST(WeightedSum, EmptyEntryKeepsATermsOwnTangent) { + ctx_type ctx; + auto& src = add_strain(ctx); + add_weight(ctx, "wA", 0.5); + add_weight(ctx, "wB", 0.5); + add_monolithic(ctx, "matA", KA, GA); // keeps its own tangent + add_decomposed(ctx, "matB", KB, GB); // tangent lives elsewhere + + param_type p; + p.insert("name", "mix"); + p.insert("terms", {{"wA", "matA"}, {"wB", "matB"}}); + p.insert>("tangent_sources", {"", "matB_stiff"}); + ctx.create>(p); + ctx.finalize(); + + src.bind(uniaxial(0.001), uniaxial(0.001)); + ctx.update(); + + const T cA = KA + 4.0 * GA / 3.0; + const T cB = KB + 4.0 * GB / 3.0; + EXPECT_NEAR(ctx.get("mix", "tangent")(0, 0, 0, 0), + 0.5 * cA + 0.5 * cB, 1e-9); +} + +/// A shorter list is rejected: positional matching would apply the override to +/// the WRONG constituent, and both names still resolve. The only symptom would +/// be a wrong summed tangent — slow Newton, correct stresses. +TEST(WeightedSum, RejectsATangentSourcesListThatDoesNotMatchTheTermCount) { + auto build = [](std::vector sources) { + ctx_type ctx; + add_strain(ctx); + add_weight(ctx, "wA", 0.5); + add_weight(ctx, "wB", 0.5); + add_monolithic(ctx, "matA", KA, GA); + add_decomposed(ctx, "matB", KB, GB); + param_type p; + p.insert("name", "mix"); + p.insert("terms", {{"wA", "matA"}, {"wB", "matB"}}); + p.insert>("tangent_sources", std::move(sources)); + ctx.create>(p); + }; + + EXPECT_THROW(build({"matB_stiff"}), std::runtime_error); // too short + EXPECT_THROW(build({"", "matB_stiff", "extra"}), std::runtime_error); // long + EXPECT_NO_THROW(build({"", "matB_stiff"})); // exact +} + +} // namespace