Cantilever plate¶

Open In Colab Binder

A plate of dimensions 100mm x 50mm x 0.5mm clamped at one end. It is subjected to concentrated forces of 5N at each corner at the tip.

A reference solution is computed with ABAQUS for the following shell elements:

Type u_max S_max
S3 30.36 423.6
STRI3 30.64 461.7
STRI65 30.90 481.4
In [1]:
import torch

from torchfem import Shell
from torchfem.materials import IsotropicElasticityPlaneStress
from torchfem.mesh import rect_tri

torch.set_default_dtype(torch.float64)

Model setup¶

In [2]:
# Material parameters
material = IsotropicElasticityPlaneStress(E=200000.0, nu=0.3)

# Define nodes and element of the plate
nodes, elements = rect_tri(10, 5, Lx=100.0, Ly=50.0, variant="zigzag")
nodes = torch.hstack([nodes, torch.zeros((nodes.size(0), 1))])

# Define Shell model
cantilever = Shell(nodes, elements, material, thickness=0.5)

# Define masks for boundaries
left = nodes[:, 0] < 0.01
right = nodes[:, 0] > 99.9
front = nodes[:, 1] < 0.01
back = nodes[:, 1] > 49.99
corner_1 = right & front
corner_2 = right & back

# Apply load boundary conditions
cantilever.forces[corner_1, 2] = -5
cantilever.forces[corner_2, 2] = -5

# Apply displacement boundary conditions
cantilever.constraints[left] = True

cantilever.plot(bcs=True)

Solve and extract maxima¶

In [3]:
# Solve
u, f, σ, _, _ = cantilever.solve(aggregate_integration_points=False)

# Extract stresses along thickness integration points
σ_top = σ[0]
σ_middle = σ[1]
σ_bottom = σ[2]

# Compute v. Mises stress at top (z=0.25, integration point 0)
mises_top = torch.sqrt(
    σ_top[:, 0, 0] ** 2
    - σ_top[:, 0, 0] * σ_top[:, 1, 1]
    + σ_top[:, 1, 1] ** 2
    + 3 * σ_top[:, 0, 1] ** 2
)
print(f"u_max: {torch.min(u[:, 0:3]):.2f}, S_max: {torch.max(mises_top):.1f}.")
u_max: -30.72, S_max: 430.6.

Visualize displacement¶

In [4]:
# Displacement on deformed configuration
cantilever.plot(u[:, 0:3], node_property={"u": torch.linalg.norm(u[:, 0:3], dim=1)})

Visualize v. Mises stress¶

In [5]:
# Mises stress on deformed configuration
cantilever.plot(
    u[:, 0:3],
    element_property={"v. Mises stress (top)": mises_top},
    cmap="inferno",
    thickness=True,
)