I think the post below is better for the finish as a program, so please refer to it.
I made a linear elastic two-dimensional FEM stress analysis program by Python.
The element is an isoparametric element with 4 nodes per element, and each element has 4 Gauss integral points.
The loads that can be handled are as follows.
$ \ {\ boldsymbol {\ epsilon_0} \} $ is distortion due to temperature change
$ \ {\ boldsymbol {\ epsilon} \} $ is the distortion of any point in the element
$ \ {\ boldsymbol {v} \} $ is the displacement of any point in the element
Since Python is an interpreter, I was worried whether it would be useful in terms of calculation speed, but the calculation time for the model with the number of nodes: 3417, the number of elements: 3257, and the degree of freedom: 6834 is 7.9 seconds, which is a usable level. It was confirmed that.
The program is long, so I put a link to what I pasted in Gist. In addition, contour drawing only fills the elements with different colors according to the stress level. For my purpose, it is not necessary to insert it separately and draw a contour line.
npoin nele nsec npfix nlod NSTR #Basic amount
t E po alpha gamma gkh gkv #Material property
..... (1~nsec) .....
node1 node2 node3 node4 isec #element-Nodal relationship, material property number
..... (1~nele) .....
x y deltaT #Nodal coordinates, nodal temperature changes
..... (1~npoin) .....
node kox koy rdisx rdisy #Displacement constraint conditions
..... (1~npfix) .....
node fx fy #External force
..... (1~nlod) .....
npoin, nele, nsec td> | Number of nodes, number of elements, number of material characteristics td> tr> |
npfix, nlod td> | Number of constrained nodes, number of loaded nodes td> tr> |
NSTR td> | Stress state (plane stress: 0, plane stress: 1) td> tr> |
t, E, po, alpha td> | Plate thickness, elastic modulus, Poisson's ratio, coefficient of linear expansion td> tr> |
gamma, gkh, gkv td> | Unit volume weight, horizontal and vertical acceleration (ratio of g) td> tr> |
x, y, deltaT td> | Node x coordinate, node y coordinate, node temperature change td> tr> |
node, kox, koy td> | Constrain node number, presence / absence of x and y direction constraints (constraint: 1, freedom: 0) td> tr> |
rdisx, rdisy td> | Displacement in x and y directions (enter 0 even if unconstrained) td> tr> |
node, fx, fy td> | Load node number, x-direction load, y-direction load td> tr> |
npoin nele nsec npfix nlod NSTR
4 1 1 2 2 1
sec t E po alpha gamma gkh gkv
1 1.0000000e+00 1.0000000e+03 0.0000000e+00 1.0000000e-05 2.3000000e+00 0.000 0.000
node x y fx fy deltaT kox koy
1 0.0000000e+00 0.0000000e+00 0.0000000e+00 0.0000000e+00 0.0000000e+00 1 1
2 1.0000000e+00 0.0000000e+00 0.0000000e+00 0.0000000e+00 0.0000000e+00 0 1
3 1.0000000e+00 1.0000000e+00 0.0000000e+00 1.0000000e+01 0.0000000e+00 0 0
4 0.0000000e+00 1.0000000e+00 0.0000000e+00 1.0000000e+01 0.0000000e+00 0 0
node kox koy rdis_x rdis_y
1 1 1 0.0000000e+00 0.0000000e+00
2 0 1 0.0000000e+00 0.0000000e+00
elem i j k l sec
1 1 2 3 4 1
node dis-x dis-y
1 0.0000000e+00 0.0000000e+00
2 8.8817842e-19 0.0000000e+00
3 1.7763568e-18 2.0000000e-02
4 4.1448326e-18 2.0000000e-02
elem sig_x sig_y tau_xy p1 p2 ang
1 -7.4014868e-16 2.0000000e+01 1.1333527e-16 2.0000000e+01 0.0000000e+00 9.0000000e+01
n=8 time=0.044 sec
node, dis-x, dis-y td> | Node number, displacement in x direction, displacement in y direction td> tr> |
elme, sig_x, sig_y, tau_xy td> | Element number, x-direction direct stress, y-direction direct stress, shear stress td> tr> |
p1, p2, ang td> | First principal stress, second principal stress, direction of first principal stress td> tr> |
n, time td> | total degrees of freedom, calculation time td> tr> |
that's all
Recommended Posts