Skip to content

DrSimulate/gallery

Folders and files

NameName
Last commit message
Last commit date

Latest commit

 

History

59 Commits
 
 
 
 
 
 
 
 
 
 
 
 

Repository files navigation

⚠️ ⚠️ ⚠️ Please allow your browser some time to load the animated GIFs. ⚠️ ⚠️ ⚠️

Banner

Gallery: How to visualize computational mechanics

Welcome to the gallery of visualizations in computational mechanics. These visualizations were created as part of the development of a series of educational video resources on computational mechanics on YouTube (Ref. 2). Below, we highlight selected example animations along with their corresponding source code, generated using the Python libraries Matplotlib and Manim, as well as the 3D animation software Blender. More detailed explanations can be found in the accompanying publication (Ref. 1).

Matplotlib

A straightforward way to create animations is by using the animation functionality of the Python library Matplotlib. Install Matplotlib and import it in Python.

import matplotlib.pyplot as plt
import matplotlib.animation as animation

Example: Stress Tensor Components - Infinitesimal Sphere

Description: We visualize the individual components of the Cauchy stress tensor in a fixed Cartesian basis by illustrating the outward unit normal vectors (in red) and the traction vectors (in white) acting on a spherical infinitesimal volume element. The subroutine animate(frame) is called once per frame and generates a plot of an infinitesimal volume element viewed from a different angle. The function animation.FuncAnimation() combines the static plots from each frame into a continuous animation.

Video link: https://youtu.be/NtTVEzZS3Bg

# !!! full code in matplotlib_animate.py !!!
# sphere
u, v = np.mgrid[0:2 * np.pi:30j, 0:np.pi:20j]
y = radius * np.sin(u) * np.sin(v)
x = radius * np.cos(u) * np.sin(v)
z = radius * np.cos(v)
# loop over stress components
for i in range(6):
    traction = np.matmul(sigma,normal_position_nomalized.T).T
    # animate
    fig = plt.figure(figsize=(_figsize,_figsize),dpi=_dpi)
    ax = plt.axes(projection='3d')
    def animate(frame):
        ax.cla()
        plot_coo(ax,shift=np.zeros(3))
        ax.plot_surface(x,y,z,alpha=.2,color=COLOR0)
        plot_vector_field(ax,normal_position,normal,color=COLOR1)
        plot_vector_field(ax,normal_position,factor*traction,color=COLOR0)
        ax.view_init(elev=20, azim=-60 + frame/2)
        return
    ani = animation.FuncAnimation(fig,animate,frames=FACTOR_FRAMES*FRAMES,interval=1/0.03)
    ani.save(_path + '.mov',codec="png",dpi=_dpi,bitrate=-1,savefig_kwargs={"transparent": True, "facecolor": "none"})

Stress Tensor Components

Example: Stress Tensor Components - Infinitesimal Cube

Description: Similar to the previous example, we visualize the individual components of the Cauchy stress tensor in a fixed Cartesian basis by illustrating the traction vectors acting on a cubic infinitesimal volume element.

Video link: https://youtu.be/NtTVEzZS3Bg

# !!! full code in matplotlib_animate.py !!!
# cube
r = np.linspace(-length/2, length/2, 3)
x, y = np.meshgrid(r, r)
z = np.zeros_like(x)
# loop over stress components
for i in range(6):
    # animate
    fig = plt.figure(figsize=(_figsize,_figsize),dpi=_dpi)
    ax = plt.axes(projection='3d')
    def animate(frame):
        ax.cla()
        plot_coo(ax,shift=np.zeros(3))
        # plot cube faces
        ax.plot_surface(x,y,-length/2.0*np.ones_like(x),alpha=.2,color=COLOR0)
        ax.plot_surface(x,y,length/2.0*np.ones_like(x),alpha=.2,color=COLOR0)
        ax.plot_surface(x,-length/2.0*np.ones_like(x),y,alpha=.2,color=COLOR0)
        ax.plot_surface(x,length/2.0*np.ones_like(x),y,alpha=.2,color=COLOR0)
        ax.plot_surface(-length/2.0*np.ones_like(x),x,y,alpha=.2,color=COLOR0)
        ax.plot_surface(length/2.0*np.ones_like(x),x,y,alpha=.2,color=COLOR0)
        if COMPONENT == 1:
            sigma[0,0] = sigma_given
            traction = np.matmul(sigma,normals_1_plus.T).T
            plot_vector_field(ax,nodes_1_minus,-traction,color=COLOR0)
            plot_vector_field(ax,nodes_1_plus,traction,color=COLOR0)
        if COMPONENT == 2:
            sigma[1,1] = sigma_given
            traction = np.matmul(sigma,normals_2_plus.T).T
            plot_vector_field(ax,nodes_2_minus,-traction,color=COLOR0)
            plot_vector_field(ax,nodes_2_plus,traction,color=COLOR0)
        if COMPONENT == 3:
            sigma[2,2] = sigma_given
            traction = np.matmul(sigma,normals_3_plus.T).T
            plot_vector_field(ax,nodes_3_minus,-traction,color=COLOR0)
            plot_vector_field(ax,nodes_3_plus,traction,color=COLOR0)
        if COMPONENT == 4:
            sigma[0,1] = sigma_given
            sigma[1,0] = sigma_given
            traction1 = np.matmul(sigma,normals_1_plus.T).T
            traction2 = np.matmul(sigma,normals_2_plus.T).T
            plot_vector_field(ax,nodes_1_minus,-traction1,color=COLOR0)
            plot_vector_field(ax,nodes_1_plus,traction1,color=COLOR0)
            plot_vector_field(ax,nodes_2_minus,-traction2,color=COLOR0)
            plot_vector_field(ax,nodes_2_plus,traction2,color=COLOR0)
        if COMPONENT == 5:
            sigma[0,2] = sigma_given
            sigma[2,0] = sigma_given
            traction1 = np.matmul(sigma,normals_1_plus.T).T
            traction3 = np.matmul(sigma,normals_3_plus.T).T
            plot_vector_field(ax,nodes_1_minus,-traction1,color=COLOR0)
            plot_vector_field(ax,nodes_1_plus,traction1,color=COLOR0)
            plot_vector_field(ax,nodes_3_minus,-traction3,color=COLOR0)
            plot_vector_field(ax,nodes_3_plus,traction3,color=COLOR0)
        if COMPONENT == 6:
            sigma[1,2] = sigma_given
            sigma[2,1] = sigma_given
            traction2 = np.matmul(sigma,normals_2_plus.T).T
            traction3 = np.matmul(sigma,normals_3_plus.T).T
            plot_vector_field(ax,nodes_2_minus,-traction2,color=COLOR0)
            plot_vector_field(ax,nodes_2_plus,traction2,color=COLOR0)
            plot_vector_field(ax,nodes_3_minus,-traction3,color=COLOR0)
            plot_vector_field(ax,nodes_3_plus,traction3,color=COLOR0)
    ani = animation.FuncAnimation(fig,animate,frames=FACTOR_FRAMES*FRAMES,interval=1/0.03)
    ani.save(_path + '.mov',codec="png",dpi=_dpi,bitrate=-1,savefig_kwargs={"transparent": True, "facecolor": "none"})
    

Stress Tensor Components

Manim

In the following, we highlight example animations implemented in the Manim Community Edition. To run these examples, install ManimCE by following the installation instructions and import it in Python.

from manim import *

Example: Kinematics 2D

Description: We visualize large deformation kinematics in two dimensions. We define a triangle in the reference configuration. This triangle is then deformed according to the deformation mapping phi(X,t) which depends on time through the ValueTracker.

Video link: https://youtu.be/dT30kLnQNUE

# !!! full code in manim.py !!!
class Kinematics2D(MovingCameraScene):
    def construct(self):
        time = ValueTracker(0.00) # set a scalar parameter that evolves with time
        # reference configuration
        reference = Triangle().stretch(1.3, dim=1).move_to(2*UP + 1.25*RIGHT)
        # deformed configuration
        def phi(X,t):
            x0 = (1 - 0.25*t)*X[0] + t*(X[1]-2)**2 + 2.5*t
            x1 = (1 + 0.25*t)*X[1] + 0.25*t
            x2 = X[2]
            return (x0,x1,x2)
        current = always_redraw(
            lambda: reference.copy().apply_function(lambda X: phi(X,time.get_value()))
        )
        # animate
        self.wait(.25)
        self.play(time.animate.set_value(1.00), rate_func=linear)
        self.wait(.5)
        self.play(time.animate.set_value(0.00), rate_func=linear)
        self.wait(.25)

Kinematics 2D

Example: Kinematics 2D - Deformation Gradient Tensor Field

Description: We visualize large deformation kinematics in two dimensions. To visualize the deformation gradient tensor field, we animate the deformation of infinitesimal elements at selected points.

Video link: https://youtu.be/1GwAEnegaRs

# !!! full code in manim.py !!!
class Kinematics2D(MovingCameraScene):
    def construct(self):
        time = ValueTracker(0.00) # set a scalar parameter that evolves with time
        reference = get_reference()
        current = always_redraw(
            lambda: reference.copy().apply_function(lambda X: my_phi(X,time.get_value()))
        )
        # animate
        self.wait(.25)
        self.play(time.animate.set_value(1.00), rate_func=linear)
        self.wait(.5)
        self.play(time.animate.set_value(0.00), rate_func=linear)
        self.wait(.25)

Kinematics 2D - Deformation Gradient Tensor Field

Example: Different Stress Tensors

Description: Different stress tensors are visualized by animating undeformed and deformed area elements and tractions acting on them for a given state of deformation and stress at a material point. The Cauchy stress gives the current traction if multiplied by the normal of the deformed area element. Multiplying this traction by the area of the deformed element yields the current force. If multiplied by the normal of the undeformed area element, the first Piola-Kirchhoff stress gives an auxiliary current traction, which after multiplication by the area of the undeformed element yields the current force. Finally, the second Piola-Kirchhoff stress gives an auxiliary reference traction if multiplied by the normal of the undeformed area element. Multiplying this reference traction by the area of the undeformed element yields the reference force which is related to the current force by the inverse of the deformation gradient. Animating a rotation of the area elements, while updating the tractions accordingly, illustrates that stress tensors encode the tractions acting on arbitrarily oriented surfaces.

Video link: https://youtu.be/5G0rzB1X70k

# !!! full code in manim.py !!!
class AreaElementRotationCauchy(ThreeDScene):
    def construct(self):
        N0 = normalize(np.array([0,0,1]))
        t_tracker = ValueTracker(0)
        def N_path(t):
            angle = t  # one rotation if t ∈ [0, 2π]
            R_y = np.array([
                [ np.cos(angle), 0.0, np.sin(angle)],
                [ 0.0,           1.0, 0.0          ],
                [-np.sin(angle), 0.0, np.cos(angle)]
            ])
            N = R_y @ N0
            return normalize(N)
        area_element = always_redraw(
            lambda: get_area_element_normal(
                normal=get_N()
            ).copy().apply_matrix(MY_F)
        )
        arrow_t_i = VGroup()
        arrow_t_i.add(always_redraw(lambda: get_arrow_t_i(0)))
        arrow_n_i = VGroup()
        arrow_n_i.add(always_redraw(lambda: get_arrow_n_i(0)))
        # animate
        self.wait(0.25)
        self.play(t_tracker.animate.set_value(2 * PI),run_time=3)
        self.wait(0.25)

Stress Tensors

Example: Tensor Components

Description: We visualize how the deformation gradient tensor transforms a set of unit vectors by illustrating their deformed images. In addition, we demonstrate how the tensor components change under a rotation of the basis. This example highlights that the physical transformation described by the tensor remains unchanged, even though its component representation depends on the chosen basis.

Video link: https://youtu.be/1GwAEnegaRs

# !!! full code in manim.py !!!
class TensorComponents(MovingCameraScene):
    def construct(self):
        F_np = np.array([[2.0, 1.0],[0.5, 1.5]])
        F_3d = np.array([[2.0, 1.0, 0.0],[0.5, 1.5, 0.0],[0.0, 0.0, 1.0]])
        alpha = ValueTracker(0.00) # set a scalar parameter that varies with the angle of rotation
        # rotating grid
        grid_rot = always_redraw(lambda:
            GRID.copy().apply_matrix(Q(alpha.get_value()))
        )
        # tensor components
        F_rot_text = always_redraw(lambda:
            MathTex(r"\boldsymbol{F} = ", matrix2text(rotate_matrix(F_np,-alpha.get_value()))
                            ).move_to(LEFT * 11)
            )
        # tensor visualization
        scale = 2
        Line_coords = scale * get_Line_coords()
        lines = get_lines(ORIGIN,Line_coords,F_3d)
        Volume = Circle(color=WHITE,stroke_width=6).scale(scale)
        volume = Volume.copy().apply_function(lambda X: F_3d @ X)
        self.add(grid_rot,F_rot_text,lines,volume)
        # animate
        self.wait(0.125)
        self.play(alpha.animate.set_value(30 * DEGREES), run_time=3)
        self.wait(0.25)
        self.play(alpha.animate.set_value(0 * DEGREES), run_time=4)
        self.wait(0.125)

Tensor Components

Example: Finite Element Method

Description: An unknown function has infinitely many degrees of freedom prior to discretization. Finite element discretization reduces this to a finite number, which can be visualized by adjusting the associated parameters.

Video link: https://youtu.be/xZpESocdvn4

# !!! full code in manim.py !!!
class FiniteElements(MovingCameraScene):
    def construct(self):
        variable1 = Variable(0, MathTex(r"u_0"), num_decimal_places=2).scale(0.5)
        variable2 = Variable(0, MathTex(r"u_1"), num_decimal_places=2).scale(0.5)
        variable3 = Variable(0, MathTex(r"u_2"), num_decimal_places=2).scale(0.5)
        variable4 = Variable(0, MathTex(r"u_3"), num_decimal_places=2).scale(0.5)
        variable5 = Variable(0, MathTex(r"u_4"), num_decimal_places=2).scale(0.5)
        variable6 = Variable(0, MathTex(r"u_5"), num_decimal_places=2).scale(0.5)
        ax_u = Axes(
            x_range=[0, 1.05, 0.2],
            y_range=[-0.5, 0.5, 1],
            x_length=5.5,
            y_length=1.75,
            axis_config={"include_tip": False},
        )
        labels_u = VGroup(
            MathTex(r"x").next_to(ax_u.coords_to_point(1.05, 0),RIGHT),
            MathTex(r"u(x)").next_to(ax_u.coords_to_point(0, 0.5),UP,buff=0).shift(UP * 0.1),
            )
        graph_uN = ax_u.plot(lambda x : func_uN(x,np.array([variable1.tracker.get_value(),
                                                            variable2.tracker.get_value(),
                                                            variable3.tracker.get_value(),
                                                            variable4.tracker.get_value(),
                                                            variable5.tracker.get_value(),
                                                            variable6.tracker.get_value()])), x_range=[0,1], use_smoothing=False, color=YELLOW)
        
        graph_uN.add_updater(lambda m : m.become(
            ax_u.plot(lambda x : func_uN(x,np.array([variable1.tracker.get_value(),
                                                     variable2.tracker.get_value(),
                                                     variable3.tracker.get_value(),
                                                     variable4.tracker.get_value(),
                                                     variable5.tracker.get_value(),
                                                     variable6.tracker.get_value()])), x_range=[0,1], use_smoothing=False, color=YELLOW)
            ))
        # animate
        self.wait(0.125)
        x_val = np.linspace(0,1,6)
        u_val = func_u_sin(x_val)
        self.wait()
        self.play(variable1.tracker.animate.set_value(u_val[0]),
                  variable2.tracker.animate.set_value(u_val[1]),
                  variable3.tracker.animate.set_value(u_val[2]),
                  variable4.tracker.animate.set_value(u_val[3]),
                  variable5.tracker.animate.set_value(u_val[4]),
                  variable6.tracker.animate.set_value(u_val[5]))

Finite Elements

Example: Finite Element Method - Mesh Refinement

Description: Mesh refinement is visualized by sucessively increasing the number of finite elements.

Video link: https://youtu.be/1wSE6iQiScg

# !!! full code in manim.py !!!
class FiniteElementsMeshRefinement(MovingCameraScene):
    def construct(self):
        ax_u = Axes(
            x_range=[0, 1.05, 0.1],
            y_range=[-0.5, 0.5, 1],
            x_length=5.5,
            y_length=1.75,
            axis_config={"include_tip": False},
        )
        labels_u = VGroup(
            MathTex(r"x").next_to(ax_u.coords_to_point(1.05, 0),RIGHT),
            MathTex(r"u(x)").next_to(ax_u.coords_to_point(0, 0.5),UP,buff=0).shift(UP * 0.1),
            )
        graph_u_analytical = ax_u.plot(lambda x : 0.5*f*x**2 - 0.5*f*L*x, x_range=[0,1], use_smoothing=False, color=WHITE)
        _, u2 = fem_1d_poisson(2, L=L, f=f)
        graph_uN2 = ax_u.plot(lambda x : func_uN2(x,u2), x_range=[0,1], use_smoothing=False, color=YELLOW)
        _, u3 = fem_1d_poisson(3, L=L, f=f)
        graph_uN3 = ax_u.plot(lambda x : func_uN3(x,u3), x_range=[0,1], use_smoothing=False, color=YELLOW)
        _, u4 = fem_1d_poisson(4, L=L, f=f)
        graph_uN4 = ax_u.plot(lambda x : func_uN4(x,u4), x_range=[0,1], use_smoothing=False, color=YELLOW)
        _, u5 = fem_1d_poisson(5, L=L, f=f)
        graph_uN5 = ax_u.plot(lambda x : func_uN5(x,u5), x_range=[0,1], use_smoothing=False, color=YELLOW)
        # animate
        self.wait(0.25)
        self.play(
            FadeIn(v_lines_2, graph_uN2),
            )
        self.wait(0.5)
        self.play(
            FadeOut(v_lines_2, graph_uN2),
            FadeIn(v_lines_3, graph_uN3),
            )
        self.wait(0.5)
        self.play(
            FadeOut(v_lines_3, graph_uN3),
            FadeIn(v_lines_4, graph_uN4),
            )
        self.wait(0.5)
        self.play(
            FadeOut(v_lines_4, graph_uN4),
            FadeIn(v_lines_5, graph_uN5),
            )
        self.wait(0.5)
        self.play(
            FadeOut(v_lines_5, graph_uN5),
            )
        self.wait(0.25)

Finite Elements

Blender

Objects animated using Matplotlib or Manim generally do not support realistic surfaces and textures. Such features can be created using Blender. Deformations of the object can be imported from simulation results or hard-coded directly into the mesh.

Example: Finite Element Simulation Results in Blender

Description: The results of finite element simulations can be visualized in Blender. The example below shows a rubber specimen whose deformation was computed using the finite element method and subsequently exported from ParaView to Blender with the Stop-motion-OBJ add-on.

Video link: https://youtu.be/svIs3-0t2LY

Dogbone Simulation Rubber

Example: Imposed Object Deformation in Blender

Description: Blender’s Python scripting interface allows direct modification and deformation of meshes. For example, a prescribed deformation mapping phi(X,t) can be applied to a predefined object, such as the water drop shown below.

Video link: https://youtu.be/1GwAEnegaRs

# !!! full code in blender.py !!!
# material
mat = bpy.data.materials.new(name="Water")
# define deformation
def phi(XYZ,t=0):
    X = XYZ[0]
    Y = XYZ[1]
    Z = XYZ[2]
    x = (1 + 0.75*t)*X
    y = Y
    z = (1 - 0.25*t)*Z + 0.35*t*X**2 - 1.5*t
    return np.array([x, y, z])
# animate
for frame in range(1,int(TOTAL_FRAMES)+1):
    t = frame / TOTAL_FRAMES
    bpy.ops.mesh.primitive_uv_sphere_add(radius=4)
    obj = bpy.context.active_object
    mesh = obj.data
    reference = np.array([v.co[:] for v in mesh.vertices])
    for i, v in enumerate(mesh.vertices):
        v.co = Vector(phi(reference[i],t))
    for face in mesh.polygons:
        face.use_smooth = True
    visible_on_frame(obj,frame)

Water Drop

Software

All presented animations were created using Matplotlib (version 3.10.1), ManimCE (version 0.18.1), and Blender (version 4.5).

References

  1. Flaschel, Moritz. How to visualize computational mechanics: Animating finite elements, continuum mechanics, and tensor calculus. Computer Applications in Engineering Education, 2026. DOI: 10.1002/cae.70189

  2. Flaschel, Moritz. Educational video resources on computational mechanics. YouTube, 2024. URL: https://www.youtube.com/@DrSimulate

About

Gallery of visualizations in computational mechanics.

Resources

Stars

Watchers

Forks

Releases

No releases published

Packages

 
 
 

Contributors

Languages