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 ====================
Isotropic:
E = 1.50e+10, v = 0.25
planeStress = False
thickness = 3.60e+02
solver:scipy
============= Boundary Conditions =============
Unspecified.
=================== Results ===================
W def = 620258634.76
Svm max = 2788891.49
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: 10.746 ms
Boundary Conditions: 21.696 µs
Matrix: 4.604 ms
Solver: 2.056 ms
PostProcessing: 489.235 µs
13 import matplotlib.pyplot as plt
14 import numpy as np
15
16 from EasyFEA import Display, Models, ElemType, Simulations
17 from EasyFEA.Geoms import Points
18
19 if __name__ == "__main__":
20 Display.Clear()
21
22 # ----------------------------------------------
23 # Configuration
24 # ----------------------------------------------
25
26 # geom
27 dim = 2
28 h = 180 # m (thickness)
29 thickness = 2 * h
30
31 # model
32 coef = 1e6
33 E = 15000 * coef # Pa (Young's modulus)
34 v = 0.25 # Poisson's ratio
35
36 # load
37 g = 9.81 # m/s^2 (acceleration due to gravity)
38 ro = 2400 # kg/m^3 (density)
39 w = 1000 # kg/m^3 (density)
40
41 # ----------------------------------------------
42 # Mesh
43 # ----------------------------------------------
44
45 contour = Points([(0, 0), (h, 0), (0, h)], h / 10)
46
47 if dim == 2:
48 mesh = contour.Mesh_2D([], ElemType.TRI6)
49 print(f"err area = {np.abs(mesh.area - h**2 / 2) / mesh.area:.3e}")
50 elif dim == 3:
51 mesh = contour.Mesh_Extrude([], [0, 0, -thickness], [3], ElemType.PRISM15)
52 print(
53 f"error volume = {np.abs(mesh.volume - h**2 / 2 * thickness) / mesh.volume:.3e}"
54 )
55
56 nodes_x0 = mesh.Nodes_Conditions(lambda x, y, z: x == 0)
57 nodes_y0 = mesh.Nodes_Conditions(lambda x, y, z: y == 0)
58
59 # ----------------------------------------------
60 # Simulation
61 # ----------------------------------------------
62
63 material = Models.Elastic.Isotropic(
64 dim, E, v, planeStress=False, thickness=thickness
65 )
66 simu = Simulations.Elastic(mesh, material)
67
68 simu.add_dirichlet(nodes_y0, [0] * dim, simu.Get_unknowns())
69 simu.add_surfLoad(
70 nodes_x0, [lambda x, y, z: w * g * (h - y)], ["x"], description="[w*g*(h-y)]"
71 )
72 simu.add_volumeLoad(mesh.nodes, [-ro * g], ["y"], description="[-ro*g]")
73
74 sol = simu.Solve()
75 simu.Save_Iter()
76
77 # ----------------------------------------------
78 # Results
79 # ----------------------------------------------
80 print(simu)
81
82 Display.Plot_Mesh(simu, h / 10 / np.abs(sol.max()))
83 Display.Plot_BoundaryConditions(simu)
84 Display.Plot_Result(simu, "Svm", nodeValues=True, coef=1 / coef, ncolors=20)
85
86 plt.show()
Total running time of the script: (0 minutes 0.251 seconds)


