Note
Go to the end to download the full example code.
Elas3#
Hydraulic dam subjected to water pressure and its own weight.
err area = 0.000e+00
==================== Mesh ====================
Element type: TRI6
Ne = 133, Nn = 302
==================== Model ====================
ElasIsot:
E = 1.50e+10, v = 0.25
planeStress = False
thickness = 3.60e+02
solver : pypardiso
============= Boundary Conditions =============
Unspecified.
=================== Results ===================
W def = 620258634.76
Svm max = 2788622.22
Evm max = 0.02 %
Ux max = 2.06e-02
Ux min = 0.00e+00
Uy max = 4.64e-04
Uy min = -9.77e-03
=================== TicTac ===================
Mesh : 334.186 ms
Boundary Conditions : 809.669 µs
Matrix : 2.215 s
Solver : 5.392 s
Display : 471.464 ms
PostProcessing : 245.131 ms
Resolution hyperelastic : 6.573 s
PyVista_Interface : 17.072 s
12 from EasyFEA import Display, Models, plt, np, ElemType, Simulations
13 from EasyFEA.Geoms import Points
14
15 if __name__ == "__main__":
16 Display.Clear()
17
18 # ----------------------------------------------
19 # Configuration
20 # ----------------------------------------------
21
22 # geom
23 dim = 2
24 h = 180 # m (thickness)
25 thickness = 2 * h
26
27 # model
28 coef = 1e6
29 E = 15000 * coef # Pa (Young's modulus)
30 v = 0.25 # Poisson's ratio
31
32 # load
33 g = 9.81 # m/s^2 (acceleration due to gravity)
34 ro = 2400 # kg/m^3 (density)
35 w = 1000 # kg/m^3 (density)
36
37 # ----------------------------------------------
38 # Mesh
39 # ----------------------------------------------
40
41 contour = Points([(0, 0), (h, 0), (0, h)], h / 10)
42
43 if dim == 2:
44 mesh = contour.Mesh_2D([], ElemType.TRI6)
45 print(f"err area = {np.abs(mesh.area - h**2 / 2) / mesh.area:.3e}")
46 elif dim == 3:
47 mesh = contour.Mesh_Extrude([], [0, 0, -thickness], [3], ElemType.PRISM15)
48 print(
49 f"error volume = {np.abs(mesh.volume - h**2 / 2 * thickness) / mesh.volume:.3e}"
50 )
51
52 nodes_x0 = mesh.Nodes_Conditions(lambda x, y, z: x == 0)
53 nodes_y0 = mesh.Nodes_Conditions(lambda x, y, z: y == 0)
54
55 # ----------------------------------------------
56 # Simulation
57 # ----------------------------------------------
58
59 material = Models.ElasIsot(dim, E, v, planeStress=False, thickness=thickness)
60 simu = Simulations.ElasticSimu(mesh, material)
61
62 simu.add_dirichlet(nodes_y0, [0] * dim, simu.Get_unknowns())
63 simu.add_surfLoad(
64 nodes_x0, [lambda x, y, z: w * g * (h - y)], ["x"], description="[w*g*(h-y)]"
65 )
66 simu.add_volumeLoad(mesh.nodes, [-ro * g], ["y"], description="[-ro*g]")
67
68 sol = simu.Solve()
69 simu.Save_Iter()
70
71 # ----------------------------------------------
72 # Results
73 # ----------------------------------------------
74 print(simu)
75
76 Display.Plot_Mesh(simu, h / 10 / np.abs(sol.max()))
77 Display.Plot_BoundaryConditions(simu)
78 Display.Plot_Result(simu, "Svm", nodeValues=True, coef=1 / coef, ncolors=20)
79
80 plt.show()
Total running time of the script: (0 minutes 0.221 seconds)


