Stress Result in Shells

About the development of the FEM module/workbench.

Moderator: bernd

User avatar
bernd
Posts: 7534
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Stress Result in Shells

Postby bernd » Wed Mar 14, 2018 1:13 pm

some code to test FreeCAD FEM mesh export of all element types to all file types with ... No crash for no element for vtk. I only need to add the missing ones ... You need to have FreeCAD 0.17.13414 or highter

Code: Select all

import Fem
doc = App.ActiveDocument
path = '/home/hugo/Desktop/'
filetype = '.vtk'  # possibilities are: *.unv *.med *.dat *.stl *.inp *.vtk *.vtu


#################################################
# Beams
# 2 node line --> seg2 ##########################
seg2 = Fem.FemMesh()
seg2.addNode( 0,  0, 0, 1)
seg2.addNode(10,  0, 0, 2)
seg2.addEdge(1, 2)
seg2
obj = doc.addObject("Fem::FemMeshObject","seg2")
obj.FemMesh = seg2
obj.Placement.Base = (0,110,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
seg2.write(path + 'seg2' + filetype)


# 3 node line --> seg3 ##########################
seg3 = Fem.FemMesh()
seg3.addNode( 0,  0, 0, 1)
seg3.addNode(10,  0, 0, 2)
seg3.addNode( 5,  0, 0, 3)
seg3.addEdge([1, 2, 3])
seg3
obj = doc.addObject("Fem::FemMeshObject","seg3")
obj.FemMesh = seg3
obj.Placement.Base = (30,110,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
seg3.write(path + 'seg3' + filetype)


#################################################
# Shells
# 3 node triangle --> tria3 #####################
tria3 = Fem.FemMesh()
tria3.addNode( 0,  0, 0, 1)
tria3.addNode( 6, 12, 0, 2)
tria3.addNode(12,  0, 0, 3)
tria3.addFace([1,  2,  3])
tria3
obj = doc.addObject("Fem::FemMeshObject","tria3")
obj.FemMesh = tria3
obj.Placement.Base = (0,80,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
obj.ViewObject.BackfaceCulling = False
tria3.write(path + 'tria3' + filetype)


# 6 node triangle --> tria 6 ####################
tria6 = Fem.FemMesh()
tria6.addNode( 0,  0, 0, 1)
tria6.addNode( 6, 12, 0, 2)
tria6.addNode(12,  0, 0, 3)
tria6.addNode( 3,  6, 0, 4)
tria6.addNode( 9,  6, 0, 5)
tria6.addNode( 6,  0, 0, 6)
tria6.addFace([1,  2, 3,  4, 5, 6])
tria6
obj = doc.addObject("Fem::FemMeshObject","tria6")
obj.FemMesh = tria6
obj.Placement.Base = (30,80,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
obj.ViewObject.BackfaceCulling = False
tria6.write(path + 'tria6' + filetype)


# 4 node quad --> quad4 #########################
quad4 = Fem.FemMesh()
quad4.addNode( 0, 10, 0, 1)
quad4.addNode(10, 10, 0, 2)
quad4.addNode(10,  0, 0, 3)
quad4.addNode( 0,  0, 0, 4)
quad4.addFace([1,  2,  3, 4])
quad4
obj = doc.addObject("Fem::FemMeshObject","quad4")
obj.FemMesh = quad4
obj.Placement.Base = (60,80,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
obj.ViewObject.BackfaceCulling = False
quad4.write(path + 'quad4' + filetype)


# 8 node quad --> quad8 #########################
quad8 = Fem.FemMesh()
quad8.addNode( 0, 10, 0, 1)
quad8.addNode(10, 10, 0, 2)
quad8.addNode(10,  0, 0, 3)
quad8.addNode( 0,  0, 0, 4)
quad8.addNode( 5, 10, 0, 5)
quad8.addNode(10,  5, 0, 6)
quad8.addNode( 5,  0, 0, 7)
quad8.addNode( 0,  5, 0, 8)
quad8.addFace([1,  2,  3, 4, 5, 6, 7, 8])
quad8
obj = doc.addObject("Fem::FemMeshObject","quad8")
obj.FemMesh = quad8
obj.ViewObject.BackfaceCulling = False
quad8.write(path + 'quad8' + filetype)


#################################################
# Volumes
# 4 node tetrahedron --> tetra4 #################
tetra4 = Fem.FemMesh()
tetra4.addNode( 6, 12, 18, 1)
tetra4.addNode( 0,  0, 18, 2)
tetra4.addNode(12,  0, 18, 3)
tetra4.addNode( 6,  6,  0, 4)
tetra4.addVolume([1,  2,  3, 4])
tetra4
obj = doc.addObject("Fem::FemMeshObject","tetra4")
obj.FemMesh = tetra4
obj.Placement.Base = (0,50,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
tetra4.write(path + 'tetra4' + filetype)


# 10 node tetrahedron --> tetra10 ###############
tetra10 = Fem.FemMesh()
tetra10.addNode( 6, 12, 18, 1)
tetra10.addNode( 0,  0, 18, 2)
tetra10.addNode(12,  0, 18, 3)
tetra10.addNode( 6,  6,  0, 4)

tetra10.addNode( 3,  6, 18, 5)
tetra10.addNode( 6,  0, 18, 6)
tetra10.addNode( 9,  6, 18, 7)

tetra10.addNode( 6,  9,  9, 8)
tetra10.addNode( 3,  3,  9, 9)
tetra10.addNode( 9,  3,  9,10)
tetra10.addVolume([1, 2, 3, 4, 5, 6, 7, 8, 9, 10])
tetra10
obj = doc.addObject("Fem::FemMeshObject","tetra10")
obj.FemMesh = tetra10
obj.Placement.Base = (30,50,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
tetra10.write(path + 'tetra10' + filetype)

# 8 node Hexahedron --> hexa8 ###################
hexa8 = Fem.FemMesh()
hexa8.addNode( 0, 10, 10, 1)
hexa8.addNode( 0,  0, 10, 2)
hexa8.addNode(10,  0, 10, 3)
hexa8.addNode(10, 10, 10, 4)
hexa8.addNode( 0, 10,  0, 5)
hexa8.addNode( 0,  0,  0, 6)
hexa8.addNode(10,  0,  0, 7)
hexa8.addNode(10, 10,  0, 8)
hexa8.addVolume([1, 2, 3, 4, 5, 6, 7, 8])
hexa8
obj = doc.addObject("Fem::FemMeshObject","hexa8")
obj.FemMesh = hexa8
obj.Placement.Base = (60,50,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
hexa8.write(path + 'hexa8' + filetype)


# 20 node Hexahedron --> hexa20 #################
hexa20 = Fem.FemMesh()
hexa20.addNode( 0, 10, 10,  1)
hexa20.addNode( 0,  0, 10,  2)
hexa20.addNode(10,  0, 10,  3)
hexa20.addNode(10, 10, 10,  4)
hexa20.addNode( 0, 10,  0,  5)
hexa20.addNode( 0,  0,  0,  6)
hexa20.addNode(10,  0,  0,  7)
hexa20.addNode(10, 10,  0,  8)

hexa20.addNode( 0,  5, 10,  9)
hexa20.addNode( 5,  0, 10, 10)
hexa20.addNode(10,  5, 10, 11)
hexa20.addNode( 5, 10, 10, 12)

hexa20.addNode( 0,  5,  0, 13)
hexa20.addNode( 5,  0,  0, 14)
hexa20.addNode(10,  5,  0, 15)
hexa20.addNode( 5, 10,  0, 16)

hexa20.addNode( 0, 10,  5, 17)
hexa20.addNode( 0,  0,  5, 18)
hexa20.addNode(10,  0,  5, 19)
hexa20.addNode(10, 10,  5, 20)
hexa20.addVolume([1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20])
hexa20
obj = doc.addObject("Fem::FemMeshObject","hexa20")
obj.FemMesh = hexa20
obj.Placement.Base = (90,50,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
hexa20.write(path + 'hexa20' + filetype)


# 6 node pentahedron --> penta6 #################
penta6 = Fem.FemMesh()
penta6.addNode(10,10,10, 1)
penta6.addNode( 0, 0,10, 2)
penta6.addNode(20, 0,10, 3)
penta6.addNode(10,10, 0, 4)
penta6.addNode( 0, 0, 0, 5)
penta6.addNode(20, 0, 0, 6)
penta6.addVolume([1, 2, 3, 4, 5, 6])
penta6
obj = doc.addObject("Fem::FemMeshObject","penta6")
obj.FemMesh = penta6
obj.Placement.Base = (0,0,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
penta6.write(path + 'penta6' + filetype)


# 15 node pentahedron --> penta15 ###############
penta15 = Fem.FemMesh()
penta15.addNode(10,10,10, 1)
penta15.addNode( 0, 0,10, 2)
penta15.addNode(20, 0,10, 3)
penta15.addNode(10,10, 0, 4)
penta15.addNode( 0, 0, 0, 5)
penta15.addNode(20, 0, 0, 6)

penta15.addNode( 5, 5,10, 7)
penta15.addNode(10, 0,10, 8)
penta15.addNode(15, 5,10, 9)

penta15.addNode( 5, 5, 0,10)
penta15.addNode(10, 0, 0,11)
penta15.addNode(15, 5, 0,12)

penta15.addNode(10,10, 5,13)
penta15.addNode( 0, 0, 5,14)
penta15.addNode(20, 0, 5,15)
penta15.addVolume([1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15])
penta15
obj = doc.addObject("Fem::FemMeshObject","penta15")
obj.FemMesh = penta15
obj.Placement.Base = (40,0,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
penta15.write(path + 'penta15' + filetype)


# 5 node pyramid --> pyra5 ######################
pyra5 = Fem.FemMesh()
pyra5.addNode( 0,20, 0, 1)
pyra5.addNode(20,20, 0, 2)
pyra5.addNode(20, 0, 0, 3)
pyra5.addNode( 0, 0, 0, 4)
pyra5.addNode(10,10,10, 5)
pyra5.addVolume([1, 2, 3, 4, 5])
pyra5
obj = doc.addObject("Fem::FemMeshObject","pyra5")
obj.FemMesh = pyra5
obj.Placement.Base = (80,0,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
pyra5.write(path + 'pyra5' + filetype)


# 13 node pyramid --> pyra13 ####################
pyra13 = Fem.FemMesh()
pyra13.addNode( 0,20, 0, 1)
pyra13.addNode(20,20, 0, 2)
pyra13.addNode(20, 0, 0, 3)
pyra13.addNode( 0, 0, 0, 4)
pyra13.addNode(10,10,10, 5)

pyra13.addNode(10,20, 0, 6)
pyra13.addNode(20,10, 0, 7)
pyra13.addNode(10, 0, 0, 8)
pyra13.addNode( 0,10, 0, 9)

pyra13.addNode( 5,15, 5,10)
pyra13.addNode(15,15, 5,11)
pyra13.addNode(15, 5, 5,12)
pyra5 = Fem.FemMesh()
pyra13.addNode( 5, 5, 5,13)
pyra13.addVolume([1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13])
pyra13
obj = doc.addObject("Fem::FemMeshObject","pyra13")
obj.FemMesh = pyra13
obj.Placement.Base = (120,0,0)
obj.ViewObject.DisplayMode = "Faces, Wireframe & Nodes"
pyra13.write(path + 'pyra13' + filetype)
User avatar
bernd
Posts: 7534
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Stress Result in Shells

Postby bernd » Mon Mar 19, 2018 12:17 pm

User avatar
bernd
Posts: 7534
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Stress Result in Shells

Postby bernd » Fri Mar 23, 2018 1:17 pm

@HarryvL
Would you test all kind of shell post processing? The known issues should be fixed.
User avatar
HarryvL
Posts: 958
Joined: Sat Jan 06, 2018 7:38 pm

Re: Stress Result in Shells

Postby HarryvL » Fri Mar 23, 2018 1:52 pm

Hi Bernd, I will !! Nice one to try on the RHS knee joint. Harry

PS: I have been busy to develop some nice applications/examples for the modal analysis option. Hope to report on that soon.
User avatar
bernd
Posts: 7534
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Stress Result in Shells

Postby bernd » Fri Mar 23, 2018 8:01 pm

HarryvL wrote:
Fri Mar 23, 2018 1:52 pm
Hi Bernd, I will !! Nice one to try on the RHS knee joint. Harry

PS: I have been busy to develop some nice applications/examples for the modal analysis option. Hope to report on that soon.
:D I'm curious
User avatar
HarryvL
Posts: 958
Joined: Sat Jan 06, 2018 7:38 pm

Re: Stress Result in Shells

Postby HarryvL » Sat Mar 24, 2018 5:37 am

I updated my FC-daily version through PPA, ran the original test case for a strip under pure bending (Note: I did not re-build the model, is that a problem?) and obtained the following promising result for the maximum principal stress (210MPa tension at the underside of the strip) :

Test56_PS_Max_Pure_Bending_2.png
Test56_PS_Max_Pure_Bending_2.png (309.16 KiB) Viewed 219 times

However, when zooming in, an anomaly is visible. The stresses at the center nodes of the 15-node prism should all be equal (i.e. zero), but I get strange fluctuations. The CGX plot in my original post shows that Calculix does give the correct results, so perhaps something is still not quite right between importing the FRD file and displaying the results for a 3D Shell? Could it perhaps be the node numbering at the center? Printing the stress values at the center nodes (SXX should be zero for all nodes) may give a clue. Unfortunately I don't know how to find them through Python.

Test56_PS_Max_Pure_Bending_Detail.png
Test56_PS_Max_Pure_Bending_Detail.png (339.54 KiB) Viewed 219 times

When I try to create a VTK pipeline from the CalculiX_static_results object FC crashes. Perhaps updating through PPA does not capture the latest updates?

OS: Ubuntu 16.04.4 LTS
Word size of OS: 64-bit
Word size of FreeCAD: 64-bit
Version: 0.17.13452 (Git)
Build type: None
Branch: master
Hash: 34633c144de9133c1f9aeb7da783f369cae9bfaf
Python version: 2.7.12
Qt version: 4.8.7
Coin version: 4.0.0a
OCC version: 7.1.0
Locale: English/UnitedStates (en_US)
User avatar
HarryvL
Posts: 958
Joined: Sat Jan 06, 2018 7:38 pm

Re: Stress Result in Shells

Postby HarryvL » Sat Mar 24, 2018 5:48 am

and here is the input file:
Test56.fcstd
(185.85 KiB) Downloaded 21 times
User avatar
HarryvL
Posts: 958
Joined: Sat Jan 06, 2018 7:38 pm

Re: Stress Result in Shells

Postby HarryvL » Sat Mar 24, 2018 5:54 am

PS: For those who want to reproduce my results to within 5 decimal places ;) please note that I manually set Poisson's ratio to zero in the INP file after writing by FC, but before running Calculix. The reason for this I explained in my original post.
User avatar
bernd
Posts: 7534
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Stress Result in Shells

Postby bernd » Sat Mar 24, 2018 8:30 am

I will have a look.
In the regardof poisons ratio, you just can overwrite material parameter in material task panel. The changed values will even be saved in the material object.
User avatar
bernd
Posts: 7534
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: Stress Result in Shells

Postby bernd » Sat Mar 24, 2018 8:30 am

I will have a look.
In the regardof poisons ratio, you just can overwrite material parameter in material task panel. The changed values will even be saved in the material object.