Note

Go to the end to download the full example code.

Elas9#

Wave propagation.

Elas9
0.000  %
5.000  %
10.000  %
15.000  %
20.000  %
25.000  %
30.000  %
35.000  %
40.000  %
45.000  %
50.000  %
55.000  %
60.000  %
65.000  %
70.000  %
75.000  %
80.000  %
85.000  %
90.000  %
95.000  %
100.000  %
Generate movie 01/21 (4.76 %) 2.40 s
Generate movie 02/21 (9.52 %) 2.18 s
Generate movie 03/21 (14.29 %) 2.19 s
Generate movie 04/21 (19.05 %) 1.95 s
Generate movie 05/21 (23.81 %) 1.86 s
Generate movie 06/21 (28.57 %) 1.78 s
Generate movie 07/21 (33.33 %) 1.60 s
Generate movie 08/21 (38.10 %) 1.46 s
Generate movie 09/21 (42.86 %) 1.39 s
Generate movie 10/21 (47.62 %) 1.25 s
Generate movie 11/21 (52.38 %) 1.16 s
Generate movie 12/21 (57.14 %) 1.06 s
Generate movie 13/21 (61.90 %) 916.33 ms
Generate movie 14/21 (66.67 %) 800.67 ms
Generate movie 15/21 (71.43 %) 700.98 ms
Generate movie 16/21 (76.19 %) 572.55 ms
Generate movie 17/21 (80.95 %) 457.47 ms
Generate movie 18/21 (85.71 %) 350.06 ms
Generate movie 19/21 (90.48 %) 229.15 ms
Generate movie 20/21 (95.24 %) 115.17 ms
Generate movie 21/21 (100.00 %) 0.00 µs


 11 # TODO: Compare results with analytical values.
 12
 13 from EasyFEA import Folder, Display, Models, Tic, plt, ElemType, Simulations, PyVista
 14 from EasyFEA.Geoms import Domain, Circle, Line
 15
 16 if __name__ == "__main__":
 17     Display.Clear()
 18
 19     # ----------------------------------------------
 20     # Configuration
 21     # ----------------------------------------------
 22
 23     # outputs
 24     folder = Folder.Join(Folder.RESULTS_DIR, "LinearizedElasticity", "Dynamic2")
 25     plotModel = False
 26     plotIter = False
 27     makeMovie = True
 28     result = "speed_norm"
 29
 30     # Define geometric parameters
 31     a = 1
 32     meshSize = a / 50
 33     diam = a / 10
 34     r = diam / 2
 35
 36     # Time parameters
 37     tMax = 1e-6
 38     Nt = 20
 39     dt = tMax / Nt
 40
 41     # Load parameters
 42     load = 1e-3
 43     f0 = 2
 44     a0 = 1
 45     t0 = dt * 4
 46
 47     # ----------------------------------------------
 48     # Mesh
 49     # ----------------------------------------------
 50
 51     # Define the domain and create the mesh
 52     domain = Domain((-a / 2, -a / 2), (a / 2, a / 2), meshSize)
 53     circle = Circle((0, 0), diam, meshSize, isHollow=False)
 54     line = Line((0, 0), (diam / 4, 0))
 55     mesh = domain.Mesh_2D([circle], ElemType.TRI6, cracks=[line])
 56
 57     # Plot the model if specified
 58     if plotModel:
 59         Display.Plot_Tags(mesh)
 60         plt.show()
 61
 62     # Get nodes for boundary conditions and loading
 63     nodesBorders = mesh.Nodes_Tags(["L0", "L1", "L2", "L3"])
 64     nodesLoad = mesh.Nodes_Point((0, 0))
 65
 66     # ----------------------------------------------
 67     # Simulation
 68     # ----------------------------------------------
 69
 70     # Define material properties
 71     material = Models.ElasIsot(2, E=210000e6, v=0.3, planeStress=False, thickness=1)
 72     lmbda = material.get_lambda()
 73     mu = material.get_mu()
 74
 75     # Create the simulation object
 76     simu = Simulations.ElasticSimu(mesh, material)
 77     simu.Set_Rayleigh_Damping_Coefs(1e-10, 1e-10)
 78     simu.Solver_Set_Hyperbolic_Algorithm(beta=1 / 4, gamma=1 / 2, dt=dt)
 79
 80     # Plot the result at the initial iteration if specified
 81     if plotIter:
 82         ax = Display.Plot_Result(simu, result, nodeValues=True, title=result)
 83
 84     # Create a timer object
 85     tic = Tic()
 86
 87     # Time loop
 88     t = 0
 89     while t <= tMax:
 90
 91         simu.Bc_Init()
 92
 93         # Set displacement boundary conditions
 94         simu.add_dirichlet(nodesBorders, [0, 0], ["x", "y"], description="[0,0]")
 95
 96         # Set Neumann boundary conditions (loading) at t = t0
 97         if t == t0:
 98             simu.add_neumann(nodesLoad, [load], ["x"])
 99
100         # Solve the simulation
101         simu.Solve()
102
103         # Save the iteration results
104         simu.Save_Iter()
105
106         # Print the progress
107         print(f"{t / tMax * 100:.3f}  %", end="\r")
108
109         # Update the plot at each iteration if specified
110         if plotIter:
111             ax = Display.Plot_Result(simu, result, nodeValues=True, ax=ax, title=result)
112             plt.pause(1e-12)
113
114         t += dt
115
116     # ----------------------------------------------
117     # Results
118     # ----------------------------------------------
119
120     if makeMovie:
121         PyVista.Movie_simu(
122             simu,
123             f"{result}",
124             folder,
125             f"{result}.gif",
126         )
127
128     plt.show()

Total running time of the script: (0 minutes 18.002 seconds)

Gallery generated by Sphinx-Gallery