[SOLVED] CLOAD card of .inp file - manually include forces in C3D10 element nodes

About the development of the FEM module/workbench.

Moderator: bernd

EkaitzEsteban
Posts: 100
Joined: Wed Sep 12, 2018 1:31 pm

[SOLVED] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby EkaitzEsteban » Thu Oct 11, 2018 8:52 am

Hello,

I have a question regarding *CLOAD card. I need to manually include forces in C3D10 element nodes correctly.

In writer.py script, one can find the following code, which is the responsible for writing the forces of *CLOAD card.

Code: Select all

  def write_constraints_force(self, f):
        # check shape type of reference shape and get node loads
        self.get_constraints_force_nodeloads()
        # write node loads to file
        f.write('\n***********************************************************\n')
        f.write('** Node loads Constraints\n')
        f.write('** written by {} function\n'.format(sys._getframe().f_code.co_name))
        f.write('*CLOAD\n')
        for femobj in self.force_objects:  # femobj --> dict, FreeCAD document object is femobj['Object']
            f.write('** ' + femobj['Object'].Label + '\n')
            direction_vec = femobj['Object'].DirectionVector
            for ref_shape in femobj['NodeLoadTable']:
                f.write('** ' + ref_shape[0] + '\n')
                for n in sorted(ref_shape[1]):
                    node_load = ref_shape[1][n]
                    if (direction_vec.x != 0.0):
                        v1 = "{:.13E}".format(direction_vec.x * node_load) # Help! Need to find node_load values
                        f.write(str(n) + ',1,' + v1 + '\n')
                    if (direction_vec.y != 0.0):
                        v2 = "{:.13E}".format(direction_vec.y * node_load) # Help! Need to find node_load values
                        f.write(str(n) + ',2,' + v2 + '\n')
                    if (direction_vec.z != 0.0):
                        v3 = "{:.13E}".format(direction_vec.z * node_load) # Help! Need to find node_load values
                        f.write(str(n) + ',3,' + v3 + '\n')
                f.write('\n')
            f.write('\n')
I understand this code. The result is basically a sorted list of 3 components: node number, direction (1=x,y=2,z=3), force magnitude in the node
e.g.

Code: Select all

22,3,1000.000000000000 # the node 22 has a force of 1000 [N] in z direction.
23,3,500.000000000000 # the node 23 has a force of 500 [N] in z direction. 
...


In order to include forces correctly in C3D10 element nodes, I need to define correctly the aforementioned 3 components.

1.- The definition of a manually selected nodes is not a handicap, it is discussed here: https://forum.freecadweb.org/viewtopic. ... 10#p258488
2.- The definition of the direction is not a handicap. x=1,y=2,z=3
3.- My problem is here in the third component: the magnitude of the force for each node.


The last part has a product of two variables: direction_vec.x * node_load direction_vec is not a handicap. The definition of Node_load is my real problem

In order to include forces correctly in C3D10 element nodes, I need to find in FreeCAD folder the python or C++ scripts where the following four functions are defined:

1.- self.get_constraints_force_nodeloads()
2.- self.force_objects
3.- femobj['NodeLoadTable']
4.- node_load = ref_shape[1][n]

I guess that all of them should be in the same script... but I dont know where are defined.

I arrived until here: def FemInputWriter.FemInputWriter.get_constraints_force_nodeloads(self ) but I dont find the python or C++ script.

best regards,

Ekaitz.

More info regarding C3D10 element here: [url]http://web.mit.edu/calculix_v2.7/Calcul ... 3.html[url]
Last edited by EkaitzEsteban on Thu Oct 11, 2018 11:11 am, edited 1 time in total.
EkaitzEsteban
Posts: 100
Joined: Wed Sep 12, 2018 1:31 pm

Re: [Help] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby EkaitzEsteban » Thu Oct 11, 2018 9:23 am

Hello again.

Ok I arrive until here: In writerbase.py, path *\FreeCAD 0.17\Mod\Fem\femsolver\writerbase.py

I find the function which defines the loads for each node:

Code: Select all

 def get_constraints_force_nodeloads(self):
        # check shape type of reference shape
        for femobj in self.force_objects:  # femobj --> dict, FreeCAD document object is femobj['Object']
            FreeCAD.Console.PrintMessage("Constraint force: " + femobj['Object'].Name + '\n')
            frc_obj = femobj['Object']
            if femobj['RefShapeType'] == 'Vertex':
                # print("load on vertices --> we do not need the femelement_table and femnodes_mesh for node load calculation")
                pass
            elif femobj['RefShapeType'] == 'Face' and FemMeshTools.is_solid_femmesh(self.femmesh) and not FemMeshTools.has_no_face_data(self.femmesh):
                # print("solid_mesh with face data --> we do not need the femelement_table but we need the femnodes_mesh for node load calculation")
                if not self.femnodes_mesh:
                    self.femnodes_mesh = self.femmesh.Nodes
            else:
                # print("mesh without needed data --> we need the femelement_table and femnodes_mesh for node load calculation")
                if not self.femnodes_mesh:
                    self.femnodes_mesh = self.femmesh.Nodes
                if not self.femelement_table:
                    self.femelement_table = FemMeshTools.get_femelement_table(self.femmesh)
        # get node loads
        for femobj in self.force_objects:  # femobj --> dict, FreeCAD document object is femobj['Object']
            frc_obj = femobj['Object']
            if frc_obj.Force == 0:
                FreeCAD.Console.PrintMessage('  Warning --> Force = 0\n')
            if femobj['RefShapeType'] == 'Vertex':  # point load on vertieces
                femobj['NodeLoadTable'] = FemMeshTools.get_force_obj_vertex_nodeload_table(self.femmesh, frc_obj)
            elif femobj['RefShapeType'] == 'Edge':  # line load on edges
                femobj['NodeLoadTable'] = FemMeshTools.get_force_obj_edge_nodeload_table(self.femmesh, self.femelement_table, self.femnodes_mesh, frc_obj)
            elif femobj['RefShapeType'] == 'Face':  # area load on faces
                femobj['NodeLoadTable'] = FemMeshTools.get_force_obj_face_nodeload_table(self.femmesh, self.femelement_table, self.femnodes_mesh, frc_obj) 


Someone can help me with the calculus of the load magnitude?

best regards,
Ekaitz.
User avatar
bernd
Posts: 8329
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: [Help] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby bernd » Thu Oct 11, 2018 10:38 am

I'm on the way. I have not read all your posts. For face loads these two links might help you ...

https://github.com/FreeCAD/FreeCAD/blob ... se.py#L185

https://github.com/FreeCAD/FreeCAD/blob ... 1016-L1150 line loads are in the same module

We do calculate the area and multiplicate by the factors from the element definitions according the book mentioned in source code

hope that helps bernd
EkaitzEsteban
Posts: 100
Joined: Wed Sep 12, 2018 1:31 pm

Re: [Help] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby EkaitzEsteban » Thu Oct 11, 2018 11:03 am

bernd wrote:
Thu Oct 11, 2018 10:38 am
I'm on the way. I have not read all your posts. For face loads these two links might help you ...

https://github.com/FreeCAD/FreeCAD/blob ... se.py#L185

https://github.com/FreeCAD/FreeCAD/blob ... 1016-L1150 line loads are in the same module

We do calculate the area and multiplicate by the factors from the element definitions according the book mentioned in source code

hope that helps bernd
Hi bernd,

You are so fast! ATM I was just looking Meshtools.py

Line 715

Code: Select all

def get_force_obj_face_nodeload_table(femmesh, femelement_table, femnodes_mesh, frc_obj):
in order to understand the force calculus....

Ok I think your tip in Line 1016 would help me a lot.

Code: Select all

def get_ref_facenodes_areas(femnodes_mesh, face_table):


best regards,
Ekaitz.
EkaitzEsteban
Posts: 100
Joined: Wed Sep 12, 2018 1:31 pm

Re: [SOLVED] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby EkaitzEsteban » Thu Oct 11, 2018 11:21 am

Hello,

Now I understand much better the force application of C3D10 elements.

In the *CLOAD card I check that using a single force in Z direction as shown in the following code:

Code: Select all

***********************************************************
** Node loads Constraints
** written by write_constraints_force function
*CLOAD
** FemConstraintForce
** node loads on shape: Cone:Face3
2,3,-0.0000000000000E+00 # Corner node
11,3,-0.0000000000000E+00 # Corner node
12,3,-0.0000000000000E+00 # Corner node
13,3,-0.0000000000000E+00 # Corner node
14,3,-7.5026359679759E-01
15,3,-7.5026359679759E-01
16,3,-7.5026359679759E-01
17,3,-7.5026359679759E-01
58,3,-0.0000000000000E+00 # Corner node
59,3,-1.5005271935952E+00
60,3,-1.5005271935952E+00
61,3,-1.5005271935952E+00
62,3,-1.5005271935952E+00
All nodes placed at the corner of the element were zero. It seemed strange to me.

The explanation can be found in page 208 from the book G. Lakshmi Narasaiah, Finite Element Analysis:
ExplanationForces.PNG
ExplanationForces.PNG (26.73 KiB) Viewed 244 times
The force magnitude is applied equally in the mid-side nodes.


EDIT: Calculix assumes this hypothesis?
Last edited by EkaitzEsteban on Thu Oct 11, 2018 11:54 am, edited 1 time in total.
EkaitzEsteban
Posts: 100
Joined: Wed Sep 12, 2018 1:31 pm

Re: [SOLVED] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby EkaitzEsteban » Thu Oct 11, 2018 11:42 am

Hello again,

So, If I want to simulate an uniform pressure normal to all surfaces of a single tetrahedron (element C3D10), I need to:
1.- Find mid-side nodes of the element.
2.- Calculate the area of each surface of a tetrathedon.
3.- Calculate the direction of the vector for each force (normal to all surfaces).
4.- Apply the formula of the book to the mid-side nodes.

Finally, write everything in the .inp file. I am right?

Is it necessary to write the values of the nodes placed at the corner in the .inp file?

best regards,
EKaitz.
fandaL
Posts: 338
Joined: Thu Jul 24, 2014 8:29 am

Re: [SOLVED] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby fandaL » Thu Oct 11, 2018 1:43 pm

EkaitzEsteban wrote:
Thu Oct 11, 2018 11:21 am
The explanation can be found in page 208 from the book G. Lakshmi Narasaiah, Finite Element Analysis:

ExplanationForces.PNG

The force magnitude is applied equally in the mid-side nodes.

EDIT: Calculix assumes this hypothesis?
yes

EkaitzEsteban wrote:
Thu Oct 11, 2018 11:42 am
Is it necessary to write the values of the nodes placed at the corner in the .inp file?
If there is 0 load, I don’t see a reason to write it.

EkaitzEsteban wrote:
Thu Oct 11, 2018 11:42 am
So, If I want to simulate an uniform pressure normal to all surfaces of a single tetrahedron (element C3D10), I need to:
1.- Find mid-side nodes of the element.
2.- Calculate the area of each surface of a tetrathedon.
3.- Calculate the direction of the vector for each force (normal to all surfaces).
4.- Apply the formula of the book to the mid-side nodes.

Finally, write everything in the .inp file. I am right?
See also *DLOAD card which is used to apply pressure on the element face (shell edge), or gravity load or centrifugal load according to the element density.
http://www.feacluster.com/CalculiX/ccx_ ... de185.html
http://www.feacluster.com/CalculiX/ccx_ ... de187.html
http://www.feacluster.com/CalculiX/ccx_ ... de241.html
User avatar
bernd
Posts: 8329
Joined: Sun Sep 08, 2013 8:07 pm
Location: Zürich, Switzerland

Re: [SOLVED] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby bernd » Thu Oct 11, 2018 1:51 pm

EkaitzEsteban wrote:
Thu Oct 11, 2018 11:21 am
EDIT: Calculix assumes this hypothesis?
To be honest I'm not that deep at FEM in theory. I used the factors from the bood for a 2nd order element, as it is what CalculiX has. I developed the factors and made a test If you use the Calculix cantilever example from StartWB and load it with a tension force the front face keeps in direction of the beam. If you use different factors the face drifts to one side ... There are some topics on the forum in this regard which show some pics in this regard.

Would be great if someone with good FEM background knowledge could proof if we gone do the right thing.

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

Re: [SOLVED] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby bernd » Thu Oct 11, 2018 1:53 pm

EkaitzEsteban wrote:
Thu Oct 11, 2018 11:42 am
So, If I want to simulate an uniform pressure normal to all surfaces of a single tetrahedron (element C3D10), I need to:
use DLOAD as fandal has written before.
EkaitzEsteban
Posts: 100
Joined: Wed Sep 12, 2018 1:31 pm

Re: [SOLVED] CLOAD card of .inp file - manually include forces in C3D10 element nodes

Postby EkaitzEsteban » Thu Oct 11, 2018 2:14 pm

fandaL wrote:
Thu Oct 11, 2018 1:43 pm
EkaitzEsteban wrote:
Thu Oct 11, 2018 11:21 am
The explanation can be found in page 208 from the book G. Lakshmi Narasaiah, Finite Element Analysis:

ExplanationForces.PNG

The force magnitude is applied equally in the mid-side nodes.

EDIT: Calculix assumes this hypothesis?
yes

EkaitzEsteban wrote:
Thu Oct 11, 2018 11:42 am
Is it necessary to write the values of the nodes placed at the corner in the .inp file?
If there is 0 load, I don’t see a reason to write it.

EkaitzEsteban wrote:
Thu Oct 11, 2018 11:42 am
So, If I want to simulate an uniform pressure normal to all surfaces of a single tetrahedron (element C3D10), I need to:
1.- Find mid-side nodes of the element.
2.- Calculate the area of each surface of a tetrathedon.
3.- Calculate the direction of the vector for each force (normal to all surfaces).
4.- Apply the formula of the book to the mid-side nodes.

Finally, write everything in the .inp file. I am right?
See also *DLOAD card which is used to apply pressure on the element face (shell edge), or gravity load or centrifugal load according to the element density.
http://www.feacluster.com/CalculiX/ccx_ ... de185.html
http://www.feacluster.com/CalculiX/ccx_ ... de187.html
http://www.feacluster.com/CalculiX/ccx_ ... de241.html
Thank you FandaL, awesome!