FEM CalculiX 2D mechanical analyses (plane stress, plane strain and a…

…xisymmetric) (FreeCAD#12562)
FEA-eng · Mar 2, 2024 · 7eae061 · 7eae061
1 parent 6c8ee59
commit 7eae061
Show file tree

Hide file tree

Showing 4 changed files with 92 additions and 3 deletions.
diff --git a/src/Mod/Fem/femmesh/meshtools.py b/src/Mod/Fem/femmesh/meshtools.py
@@ -2519,6 +2519,61 @@ def beam_reduced_integration(
     f.truncate()
     f.close()
 
+def plane_stress(
+    fileName
+):
+    # replace shell elements with plane stress elements
+    f = open(fileName, "r+")
+    lines = f.readlines()
+    f.seek(0)
+    for line in lines:
+        if line.find("S3") != -1:
+            line = line.replace("S3", "CPS3")
+        if line.find("S6") != -1:
+            line = line.replace("S6", "CPS6")
+        if line.find("S4") != -1:
+            line = line.replace("S4", "CPS4")     
+        f.write(line)
+
+    f.truncate()
+    f.close()
+def plane_strain(
+    fileName
+):
+    # replace shell elements with plane strain elements
+    f = open(fileName, "r+")
+    lines = f.readlines()
+    f.seek(0)
+    for line in lines:
+        if line.find("S3") != -1:
+            line = line.replace("S3", "CPE3")
+        if line.find("S6") != -1:
+            line = line.replace("S6", "CPE6")
+        if line.find("S4") != -1:
+            line = line.replace("S4", "CPE4")      
+        f.write(line)
+
+    f.truncate()
+    f.close()
+def axisymmetric(
+    fileName
+):
+    # replace shell elements with axisymmetric elements
+    f = open(fileName, "r+")
+    lines = f.readlines()
+    f.seek(0)
+    for line in lines:
+        if line.find("S3") != -1:
+            line = line.replace("S3", "CAX3")
+        if line.find("S6") != -1:
+            line = line.replace("S6", "CAX6")
+        if line.find("S4") != -1:
+            line = line.replace("S4", "CAX4")    
+        f.write(line)
+
+    f.truncate()
+    f.close()
+
 # ************************************************************************************************
 def sub_shape_at_global_placement(obj, sub_name):
     sub_sh = obj.getSubObject(sub_name)

diff --git a/src/Mod/Fem/femsolver/calculix/solver.py b/src/Mod/Fem/femsolver/calculix/solver.py
@@ -381,6 +381,21 @@ def add_attributes(obj, ccx_prefs):
         )
         obj.BeamReducedIntegration = True
 
+    if not hasattr(obj, "ModelSpace"):
+        model_space_types = [
+            "3D",
+            "plane stress",
+            "plane strain",
+            "axisymmetric"
+        ]
+        obj.addProperty(
+            "App::PropertyEnumeration",
+            "ModelSpace",
+            "Fem",
+            "Type of model space"
+        )
+        obj.ModelSpace = model_space_types
+
 """
 Should there be some equation object for Calculix too?
 

diff --git a/src/Mod/Fem/femsolver/calculix/write_femelement_geometry.py b/src/Mod/Fem/femsolver/calculix/write_femelement_geometry.py
@@ -109,11 +109,14 @@ def write_femelement_geometry(f, ccxwriter):
                 f.write(section_geo)
             elif "shellthickness_obj" in matgeoset:  # shell mesh
                 shellth_obj = matgeoset["shellthickness_obj"]
-                section_def = "*SHELL SECTION, {}{}\n".format(elsetdef, material)
+                if ccxwriter.solver_obj.ModelSpace == "3D":
+                    section_def = "*SHELL SECTION, {}{}\n".format(elsetdef, material)
+                else:
+                    section_def = "*SOLID SECTION, {}{}\n".format(elsetdef, material)
                 thickness = shellth_obj.Thickness.getValueAs("mm").Value
                 section_geo = "{:.13G}\n".format(thickness)
                 f.write(section_def)
-                f.write(section_geo)
+                f.write(section_geo)    
             else:  # solid mesh
                 section_def = "*SOLID SECTION, {}{}\n".format(elsetdef, material)
                 f.write(section_def)

diff --git a/src/Mod/Fem/femsolver/calculix/write_mesh.py b/src/Mod/Fem/femsolver/calculix/write_mesh.py
@@ -54,7 +54,15 @@ def write_mesh(ccxwriter):
         # Use reduced integration beam elements if this option is enabled in ccx solver settings
         if ccxwriter.solver_obj.BeamReducedIntegration:
             meshtools.beam_reduced_integration(ccxwriter.femmesh_file)
-
+
+        # Use 2D elements if model space is not set to 3D
+        if ccxwriter.solver_obj.ModelSpace == "plane stress":
+            meshtools.plane_stress(ccxwriter.femmesh_file)
+        if ccxwriter.solver_obj.ModelSpace == "plane strain":
+            meshtools.plane_strain(ccxwriter.femmesh_file)
+        if ccxwriter.solver_obj.ModelSpace == "axisymmetric":
+            meshtools.axisymmetric(ccxwriter.femmesh_file)
+
         inpfile = codecs.open(ccxwriter.file_name, "w", encoding="utf-8")
         inpfile.write("{}\n".format(59 * "*"))
         inpfile.write("** {}\n".format(write_name))
@@ -77,6 +85,14 @@ def write_mesh(ccxwriter):
         if ccxwriter.solver_obj.BeamReducedIntegration:
             meshtools.beam_reduced_integration(ccxwriter.femmesh_file)
 
+        # Use 2D elements if model space is not set to 3D
+        if ccxwriter.solver_obj.ModelSpace == "plane stress":
+            meshtools.plane_stress(ccxwriter.femmesh_file)
+        if ccxwriter.solver_obj.ModelSpace == "plane strain":
+            meshtools.plane_strain(ccxwriter.femmesh_file)
+        if ccxwriter.solver_obj.ModelSpace == "axisymmetric":
+            meshtools.axisymmetric(ccxwriter.femmesh_file)
+
         # reopen file with "append" to add all the rest
         inpfile = codecs.open(ccxwriter.femmesh_file, "a", encoding="utf-8")
         inpfile.write("\n\n")