Fem: Load heat flux from CalculiX results - fixes #12117

src/Mod/Fem/App/FemVTKTools.cpp

@@ -659,6 +659,7 @@ std::map<std::string, std::string> _getFreeCADMechResultVectorProperties()
     resFCVecProp["PS1Vector"] = "Major Principal Stress Vector";
     resFCVecProp["PS2Vector"] = "Intermediate Principal Stress Vector";
     resFCVecProp["PS3Vector"] = "Minor Principal Stress Vector";
+    resFCVecProp["HeatFlux"] = "Heat Flux";
 
     return resFCVecProp;
 }

src/Mod/Fem/feminout/importCcxFrdResults.py

@@ -311,6 +311,7 @@ def read_frd_result(
     mode_strain = {}
     mode_peeq = {}
     mode_temp = {}
+    mode_heatflux = {}
     mode_massflow = {}
     mode_networkpressure = {}
 
@@ -322,6 +323,7 @@ def read_frd_result(
     mode_strain_found = False
     mode_peeq_found = False
     mode_temp_found = False
+    mode_heatflux_found = False
     mode_massflow_found = False
     mode_networkpressure_found = False
     end_of_section_found = False
@@ -641,6 +643,19 @@ def read_frd_result(
             temperature = float(line[13:25])
             mode_temp[elem] = (temperature)
 
+        # Check if we found heat flux section
+        if line[5:9] == "FLUX":
+            mode_heatflux_found = True
+        if mode_heatflux_found and (line[1:3] == "-1"):
+            # we found a heat_flux line
+            elem = int(line[4:13])
+            mode_heatflux_x = float(line[13:25])
+            mode_heatflux_y = float(line[25:37])
+            mode_heatflux_z = float(line[37:49])
+            mode_heatflux[elem] = FreeCAD.Vector(mode_heatflux_x, mode_heatflux_y, mode_heatflux_z)
+
+
+
         # Check if we found a mass flow section
         if line[5:11] == "MAFLOW":
             mode_massflow_found = True
@@ -713,6 +728,12 @@ def read_frd_result(
                 mode_temp_found = False
                 node_element_section = False
 
+            if mode_heatflux_found:
+                mode_results["heatflux"] = mode_heatflux
+                mode_heatflux = {}
+                mode_heatflux_found = False
+                node_element_section = False
+
             if mode_massflow_found:
                 mode_results["mflow"] = mode_massflow
                 mode_massflow = {}

src/Mod/Fem/feminout/importToolsFem.py

@@ -455,6 +455,11 @@ def fill_femresult_mechanical(
                 res_obj.Temperature = list(map((lambda x: x), Temperature.values()))
             res_obj.Time = step_time
 
+    if "heatflux" in result_set:
+        HeatFlux = result_set["heatflux"]
+        if HeatFlux:
+            res_obj.HeatFlux = list(map((lambda x: x), HeatFlux.values()))
+
     # fill res_obj.MassFlow
     if "mflow" in result_set:
         MassFlow = result_set["mflow"]

src/Mod/Fem/femobjects/result_mechanical.py

@@ -230,6 +230,16 @@ class ResultMechanical(base_fempythonobject.BaseFemPythonObject):
             "Temperature field",
             True
         )
+        obj.addProperty(
+            "App::PropertyVectorList",
+            "HeatFlux",
+            "NodeData",
+            "List of heat flux vectors",
+            True
+        )
+        obj.setPropertyStatus("HeatFlux", "LockDynamic")
+
+
         obj.setPropertyStatus("Temperature", "LockDynamic")
         obj.addProperty(
             "App::PropertyFloatList",

src/Mod/Fem/femsolver/calculix/write_step_output.py

@@ -51,10 +51,13 @@ def write_step_output(f, ccxwriter):
         f.write("U\n")
     if not ccxwriter.member.geos_fluidsection:
         f.write("*EL FILE\n")
+        variables = "S, E"
+        if ccxwriter.analysis_type == "thermomech":
+            variables += ", HFL"
         if ccxwriter.solver_obj.MaterialNonlinearity == "nonlinear":
-            f.write("S, E, PEEQ\n")
-        else:
-            f.write("S, E\n")
+            variables += ", PEEQ"
+
+        f.write(variables + "\n")
 
         # dat file
         # reaction forces: freecad.org/tracker/view.php?id=2934

src/Mod/Fem/femtest/data/calculix/thermomech_bimetall.inp

@@ -7111,7 +7111,7 @@ ConstraintTemperatureNormal,11,11,273.0
 *NODE FILE
 U, NT
 *EL FILE
-S, E
+S, E, HFL
 ** outputs --> dat file
 ** reaction forces for Constraint fixed
 *NODE PRINT, NSET=ConstraintFixed, TOTALS=ONLY