Note

Go to the end to download the full example code.

Dyna2#

Wave propagation.

Dyna2
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  %
Movie_func 0/20
Movie_func 1/20 (5.00 %) 2.95 s
Movie_func 2/20 (10.00 %) 2.82 s
Movie_func 3/20 (15.00 %) 2.61 s
Movie_func 4/20 (20.00 %) 2.49 s
Movie_func 5/20 (25.00 %) 2.39 s
Movie_func 6/20 (30.00 %) 2.17 s
Movie_func 7/20 (35.00 %) 1.99 s
Movie_func 8/20 (40.00 %) 1.92 s
Movie_func 9/20 (45.00 %) 1.75 s
Movie_func 10/20 (50.00 %) 1.54 s
Movie_func 11/20 (55.00 %) 1.39 s
Movie_func 12/20 (60.00 %) 1.23 s
Movie_func 13/20 (65.00 %) 1.09 s
Movie_func 14/20 (70.00 %) 924.24 ms
Movie_func 15/20 (75.00 %) 771.54 ms
Movie_func 16/20 (80.00 %) 618.13 ms
Movie_func 17/20 (85.00 %) 467.40 ms
Movie_func 18/20 (90.00 %) 307.31 ms
Movie_func 19/20 (95.00 %) 154.68 ms
Movie_func 20/20 (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, "Dynamics", "Dyna2")
 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 5.260 seconds)

Gallery generated by Sphinx-Gallery