fix: enforce thermal boundary dof endpoints

include/fesa/analysis/analysis_model.h

@@ -20,8 +20,10 @@ class AnalysisModel {
   /// @brief Returns the non-owned Domain backing this view.
   const Domain& GetDomain() const noexcept;
 
-  /// @brief Returns the sole active static step.
+  /// @brief Returns the sole active structural step.
-  const StepDefinition& Step() const noexcept;
+  /// @throws std::logic_error when this is a thermal view.
+  /// @note Thermal callers must use ThermalStep().
+  const StepDefinition& Step() const;
 
   /// @brief Reports whether this view selects the uncoupled thermal step.
   bool IsThermal() const noexcept;
@@ -56,6 +58,12 @@ class AnalysisModel {
   /// @brief Returns connected thermal node indices in stable Domain order.
   const std::vector<EntityIndex>& ActiveThermalNodes() const noexcept;
 
+  /// @brief Returns thermal-material collection positions in source order.
+  const std::vector<EntityIndex>& ActiveThermalMaterials() const noexcept;
+
+  /// @brief Returns thermal-section collection positions in source order.
+  const std::vector<EntityIndex>& ActiveThermalSections() const noexcept;
+
   /// @brief Returns active surface pressure records without copying Domain data.
   const std::vector<SurfacePressureLoadDefinition>& ActiveSurfacePressureLoads()
       const noexcept;
@@ -74,6 +82,8 @@ class AnalysisModel {
   std::vector<EntityIndex> active_loads_;
   std::vector<EntityIndex> active_thermal_elements_;
   std::vector<EntityIndex> active_thermal_nodes_;
+  std::vector<EntityIndex> active_thermal_materials_;
+  std::vector<EntityIndex> active_thermal_sections_;
   bool is_thermal_{false};
 };
 

include/fesa/elements/ds4_thermal.h

@@ -2,6 +2,7 @@
 #define FESA_ELEMENTS_DS4_THERMAL_H_
 
 #include <array>
+#include <limits>
 #include <utility>
 #include <vector>
 
@@ -14,6 +15,9 @@ namespace fesa {
 /// state.
 class Ds4ThermalDefinition final : public ElementDefinition {
  public:
+  /// @brief Marks that DS4 has no entry in Domain::Properties().
+  static constexpr EntityIndex kNoPropertyIndex =
+      (std::numeric_limits<EntityIndex>::max)();
   Ds4ThermalDefinition(SourceEntityId source_id,
                        std::array<EntityIndex, 4> node_indices,
                        EntityIndex section_index, SourceLocation location);
@@ -74,7 +78,7 @@ inline EntityIndex Ds4ThermalDefinition::MaterialIndex() const noexcept {
 }
 
 inline EntityIndex Ds4ThermalDefinition::PropertyIndex() const noexcept {
-  return section_index_;
+  return kNoPropertyIndex;
 }
 
 inline EntityIndex Ds4ThermalDefinition::ThermalSectionIndex() const noexcept {

include/fesa/model/model_types.h

@@ -239,6 +239,8 @@ struct ThermalMaterialDefinition {
 };
 
 /// @brief Stores one centered homogeneous thermal shell assignment.
+/// @note material_index addresses ModelDefinition::thermal_materials and,
+/// after Domain construction, Domain::ThermalMaterials(), never Properties().
 struct ThermalShellSectionDefinition {
   SourceEntityId source_id;
   std::string name;
@@ -268,6 +270,16 @@ struct ThermalInitialConditionDefinition {
   SourceLocation location;
 };
 
+/// @brief Identifies the fixed bottom-to-top thermal station source DOFs.
+/// @note This is separate from the six-component structural DofComponent.
+enum class ThermalDofComponent : int {
+  kBottom = 11,
+  kLower = 12,
+  kMiddle = 13,
+  kUpper = 14,
+  kTop = 15
+};
+
 /// @brief Stores one prescribed thermal source DOF value.
 struct ThermalBoundaryConditionDefinition {
   EntityIndex node_index;

src/fesa/analysis/analysis_model.cpp

@@ -1,6 +1,8 @@
 #include "fesa/analysis/analysis_model.h"
 
+#include <algorithm>
 #include <cstddef>
+#include <stdexcept>
 #include <string>
 #include <vector>
 
@@ -40,7 +42,10 @@ Result<AnalysisModel> AnalysisModel::Create(const Domain& domain) {
 
 const Domain& AnalysisModel::GetDomain() const noexcept { return *domain_; }
 
-const StepDefinition& AnalysisModel::Step() const noexcept {
+const StepDefinition& AnalysisModel::Step() const {
+  if (is_thermal_) {
+    throw std::logic_error{"Thermal AnalysisModel requires ThermalStep()."};
+  }
   return domain_->Steps()[0U];
 }
 
@@ -92,6 +97,16 @@ const std::vector<EntityIndex>& AnalysisModel::ActiveThermalNodes()
   return active_thermal_nodes_;
 }
 
+const std::vector<EntityIndex>& AnalysisModel::ActiveThermalMaterials()
+    const noexcept {
+  return active_thermal_materials_;
+}
+
+const std::vector<EntityIndex>& AnalysisModel::ActiveThermalSections()
+    const noexcept {
+  return active_thermal_sections_;
+}
+
 const std::vector<SurfacePressureLoadDefinition>&
 AnalysisModel::ActiveSurfacePressureLoads() const noexcept {
   static const std::vector<SurfacePressureLoadDefinition> kNoSurfacePressures;
@@ -109,6 +124,12 @@ AnalysisModel::AnalysisModel(const Domain& domain) : domain_{&domain} {
          ++index) {
       active_thermal_elements_.push_back(
           static_cast<EntityIndex>(first_element + index));
+      active_thermal_sections_.push_back(
+          domain.ThermalElements()[index].ThermalSectionIndex());
+      active_thermal_materials_.push_back(
+          domain.ThermalSections()[domain.ThermalElements()[index]
+                                       .ThermalSectionIndex()]
+              .material_index);
       for (const EntityIndex node :
            domain.ThermalElements()[index].NodeIndices()) {
         active_nodes[node] = true;
@@ -119,6 +140,17 @@ AnalysisModel::AnalysisModel(const Domain& domain) : domain_{&domain} {
         active_thermal_nodes_.push_back(static_cast<EntityIndex>(index));
       }
     }
+    std::sort(active_thermal_sections_.begin(), active_thermal_sections_.end());
+    active_thermal_sections_.erase(
+        std::unique(active_thermal_sections_.begin(),
+                    active_thermal_sections_.end()),
+        active_thermal_sections_.end());
+    std::sort(active_thermal_materials_.begin(),
+              active_thermal_materials_.end());
+    active_thermal_materials_.erase(
+        std::unique(active_thermal_materials_.begin(),
+                    active_thermal_materials_.end()),
+        active_thermal_materials_.end());
     return;
   }
   std::vector<bool> reachable_materials(domain.Materials().Size(), false);

src/fesa/model/domain.cpp

@@ -20,7 +20,7 @@ Status InvalidThermalDefinition(const SourceLocation& location,
         "THERMAL", identity, message}});
 }
 
-bool IsThermalDof(const int dof) { return dof >= 11 && dof <= 15; }
+bool IsThermalDof(const int dof) { return dof == 11 || dof == 15; }
 
 Status ValidateThermalTargets(
     const std::vector<ThermalBoundaryConditionDefinition>& boundaries,

tests/unit/model/thermal_semantics_test.cpp

@@ -4,6 +4,7 @@
 
 #include <array>
 #include <limits>
+#include <stdexcept>
 #include <type_traits>
 #include <utility>
 
@@ -41,6 +42,14 @@ fesa::ModelDefinition ValidDefinition() {
 
 }  // namespace
 
+TEST(ThermalDomain, StructuralDefinitionHasNoThermalStep) {
+  fesa::ModelDefinition definition{};
+  definition.steps = {{"Static", {}, {}, 0.1, 1.0, 0.01, 1.0, Location()}};
+  const auto result = fesa::Domain::Create(std::move(definition));
+  ASSERT_TRUE(result.HasValue());
+  EXPECT_EQ(result.Value().ThermalStep(), nullptr);
+}
+
 TEST(IsotropicThermalMaterial, RequiresFinitePositivePropertiesAndIdentity) {
   const auto valid = fesa::IsotropicThermalMaterial::Create(
       {"", 7, "Thermal-Steel"}, "Thermal-Steel", 45.0, 7800.0, 500.0,
@@ -56,6 +65,15 @@ TEST(IsotropicThermalMaterial, RequiresFinitePositivePropertiesAndIdentity) {
       std::numeric_limits<double>::infinity(), 7800.0, 500.0, Location());
   ASSERT_FALSE(invalid.HasValue());
   EXPECT_EQ(invalid.GetStatus().Category(), fesa::FailureCategory::kModel);
+
+  EXPECT_FALSE(fesa::IsotropicThermalMaterial::Create(
+                   {"", 7, "Thermal-Steel"}, "Thermal-Steel", 45.0, 0.0,
+                   500.0, Location())
+                   .HasValue());
+  EXPECT_FALSE(fesa::IsotropicThermalMaterial::Create(
+                   {"", 7, "Thermal-Steel"}, "Thermal-Steel", 45.0, 7800.0,
+                   std::numeric_limits<double>::quiet_NaN(), Location())
+                   .HasValue());
 }
 
 TEST(ThermalDomain, OwnsFiveStationDefinitionsWithoutEquationState) {
@@ -70,12 +88,16 @@ TEST(ThermalDomain, OwnsFiveStationDefinitionsWithoutEquationState) {
   EXPECT_EQ(domain.ThermalElements().Size(), 1U);
   EXPECT_EQ(domain.ThermalElements()[0U].Kind(),
             fesa::ElementDefinitionKind::kDs4Thermal);
+  EXPECT_EQ(domain.ThermalElements()[0U].PropertyIndex(),
+            fesa::Ds4ThermalDefinition::kNoPropertyIndex);
   EXPECT_EQ(domain.ThermalElements()[0U].NodeIndices(),
             (std::vector<fesa::EntityIndex>{0U, 1U, 2U, 3U}));
   EXPECT_EQ(domain.ThermalInitialConditions()[0U].station_values,
             (std::array<double, 5>{12.5, 12.5, 12.5, 12.5, 12.5}));
   EXPECT_EQ(domain.ThermalStep()->Boundaries()[0U].value, 20.0);
   EXPECT_EQ(domain.ThermalStep()->Loads()[0U].source_order, 0U);
+  EXPECT_EQ(static_cast<int>(fesa::ThermalDofComponent::kBottom), 11);
+  EXPECT_EQ(static_cast<int>(fesa::ThermalDofComponent::kTop), 15);
 }
 
 TEST(ThermalDomain, RejectsInactiveTargetsConflictsAndMixedElements) {
@@ -97,6 +119,13 @@ TEST(ThermalDomain, RejectsInactiveTargetsConflictsAndMixedElements) {
   const auto mixed = fesa::Domain::Create(std::move(invalid));
   ASSERT_FALSE(mixed.HasValue());
   EXPECT_EQ(mixed.GetStatus().Category(), fesa::FailureCategory::kModel);
+
+  invalid = ValidDefinition();
+  invalid.thermal_sections[0U].thickness = 0.0;
+  const auto invalid_thickness = fesa::Domain::Create(std::move(invalid));
+  ASSERT_FALSE(invalid_thickness.HasValue());
+  EXPECT_EQ(invalid_thickness.GetStatus().Category(),
+            fesa::FailureCategory::kModel);
 }
 
 TEST(ThermalDomain, AllowsEqualInitialTargetsAndStepBoundaryOverride) {
@@ -118,6 +147,87 @@ TEST(ThermalDomain, AllowsEqualInitialTargetsAndStepBoundaryOverride) {
   EXPECT_TRUE(overriding_step.HasValue());
 }
 
+TEST(ThermalDomain, RejectsInteriorThermalDofsAndInactiveFluxTargets) {
+  for (const int dof : {12, 13, 14}) {
+    auto definition = ValidDefinition();
+    definition.thermal_step->boundaries = {{0U, dof, 1.0, Location()}};
+    const auto boundary = fesa::Domain::Create(std::move(definition));
+    EXPECT_FALSE(boundary.HasValue());
+
+    definition = ValidDefinition();
+    definition.thermal_step->loads = {{0U, dof, 1.0, Location()}};
+    const auto flux = fesa::Domain::Create(std::move(definition));
+    EXPECT_FALSE(flux.HasValue());
+  }
+
+  auto definition = ValidDefinition();
+  definition.thermal_step->loads = {{4U, 11, 1.0, Location()}};
+  const auto inactive_flux = fesa::Domain::Create(std::move(definition));
+  EXPECT_FALSE(inactive_flux.HasValue());
+
+  definition = ValidDefinition();
+  definition.thermal_model_boundaries = {{4U, 11, 1.0, Location()}};
+  const auto inactive_boundary = fesa::Domain::Create(std::move(definition));
+  EXPECT_FALSE(inactive_boundary.HasValue());
+}
+
+TEST(ThermalDomain, PreservesAdditiveFluxSourceOrder) {
+  auto definition = ValidDefinition();
+  definition.thermal_step->loads = {{0U, 11, -5.0, Location()},
+                                    {1U, 11, 3.0, Location()}};
+  const auto result = fesa::Domain::Create(std::move(definition));
+  ASSERT_TRUE(result.HasValue());
+  const auto& loads = result.Value().ThermalStep()->Loads();
+  ASSERT_EQ(loads.size(), 2U);
+  EXPECT_EQ(loads[0U].source_order, 0U);
+  EXPECT_EQ(loads[1U].source_order, 1U);
+  EXPECT_DOUBLE_EQ(loads[0U].power, -5.0);
+  EXPECT_DOUBLE_EQ(loads[1U].power, 3.0);
+}
+
+TEST(ThermalAnalysisModel,
+     KeepsThermalSectionMaterialIndicesSeparateFromDomainMaterials) {
+  auto definition = ValidDefinition();
+  definition.materials = {{"Unused-Structural", 210.0e9, 0.3, Location()}};
+  definition.thermal_materials.push_back(
+      {{"", 2, "Thermal-Aluminum"}, "Thermal-Aluminum", 120.0, 2700.0,
+       900.0, Location()});
+  definition.thermal_sections = {
+      {{"", 1, "Section-A"}, "Section-A", 0.01, 1U, Location()},
+      {{"", 2, "Section-B"}, "Section-B", 0.02, 0U, Location()}};
+  definition.nodes.push_back(
+      {{"Plate", 6, "6"}, {3.0, 0.0, 0.0}, Location()});
+  definition.nodes.push_back(
+      {{"Plate", 7, "7"}, {4.0, 0.0, 0.0}, Location()});
+  definition.nodes.push_back(
+      {{"Plate", 8, "8"}, {4.0, 1.0, 0.0}, Location()});
+  definition.nodes.push_back(
+      {{"Plate", 9, "9"}, {3.0, 1.0, 0.0}, Location()});
+  definition.thermal_elements = {
+      {{"Plate", 8, "8"}, {0U, 1U, 2U, 3U}, 0U, Location()},
+      {{"Plate", 9, "9"}, {5U, 6U, 7U, 8U}, 1U, Location()}};
+
+  const auto domain_result = fesa::Domain::Create(std::move(definition));
+  ASSERT_TRUE(domain_result.HasValue());
+  const fesa::Domain& domain = domain_result.Value();
+  ASSERT_EQ(domain.Materials().Size(), 3U);
+  ASSERT_EQ(domain.ThermalMaterials().Size(), 2U);
+  EXPECT_EQ(domain.ThermalSections()[0U].material_index, 1U);
+  EXPECT_EQ(domain.ThermalSections()[1U].material_index, 0U);
+  EXPECT_EQ(domain.ThermalElements()[0U].MaterialIndex(), 2U);
+  EXPECT_EQ(domain.ThermalElements()[1U].MaterialIndex(), 1U);
+  EXPECT_EQ(domain.ThermalElements()[0U].PropertyIndex(),
+            fesa::Ds4ThermalDefinition::kNoPropertyIndex);
+
+  const auto model_result = fesa::AnalysisModel::Create(domain);
+  ASSERT_TRUE(model_result.HasValue());
+  EXPECT_EQ(model_result.Value().ActiveThermalMaterials(),
+            (std::vector<fesa::EntityIndex>{0U, 1U}));
+  EXPECT_EQ(model_result.Value().ActiveThermalSections(),
+            (std::vector<fesa::EntityIndex>{0U, 1U}));
+  EXPECT_TRUE(model_result.Value().ActiveProperties().empty());
+}
+
 TEST(ThermalAnalysisModel, ViewsConnectedDs4EntitiesWithoutOwningDomain) {
   auto domain_result = fesa::Domain::Create(ValidDefinition());
   ASSERT_TRUE(domain_result.HasValue());
@@ -128,9 +238,24 @@ TEST(ThermalAnalysisModel, ViewsConnectedDs4EntitiesWithoutOwningDomain) {
   const fesa::AnalysisModel& model = model_result.Value();
   EXPECT_TRUE(model.IsThermal());
   EXPECT_EQ(model.ThermalStep(), domain.ThermalStep());
+  EXPECT_THROW(model.Step(), std::logic_error);
   EXPECT_EQ(model.ActiveThermalElements(),
             (std::vector<fesa::EntityIndex>{0U}));
   EXPECT_EQ(model.ActiveThermalNodes(),
             (std::vector<fesa::EntityIndex>{0U, 1U, 2U, 3U}));
+  EXPECT_EQ(model.ActiveThermalMaterials(),
+            (std::vector<fesa::EntityIndex>{0U}));
+  EXPECT_EQ(model.ActiveThermalSections(),
+            (std::vector<fesa::EntityIndex>{0U}));
+  EXPECT_TRUE(model.ActiveProperties().empty());
+  EXPECT_EQ(domain.ThermalSections()[0U].material_index, 0U);
+  EXPECT_EQ(domain.ThermalElements()[0U].MaterialIndex(), 0U);
   EXPECT_EQ(&model.GetDomain(), &domain);
+
+  fesa::Domain moved_domain = std::move(domain_result.Value());
+  const auto moved_model = fesa::AnalysisModel::Create(moved_domain);
+  ASSERT_TRUE(moved_model.HasValue());
+  EXPECT_EQ(&moved_model.Value().GetDomain(), &moved_domain);
+  EXPECT_DOUBLE_EQ(moved_domain.ThermalMaterials()[0U].Conductivity(), 45.0);
+  EXPECT_EQ(moved_domain.ThermalElements()[0U].SourceId().source_label, 8);
 }