Cantilever plate

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	PEEQ_max
STRI65	46.33	304.9	0.4846%

import torch

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

torch.set_default_dtype(torch.float64)

Material model

# Material properties
E = 200000.0
nu = 0.3
sigma_y = 250.0
k = 10000.0


# Hardening function
def sigma_f(q):
    return sigma_y + k * q


# Derivative of the hardening function
def sigma_f_prime(q):
    return k


material = IsotropicPlasticityPlaneStress(E, nu, sigma_f, sigma_f_prime)

Model setup

# 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, n_simpson=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

# Solve
increments = torch.linspace(0, 1.0, 21)
u, f, σ, F, ɑ = cantilever.solve(
    increments=increments,
    return_intermediate=True,
    aggregate_integration_points=False,
)

# Extract final plastic strains along thickness integration points
ɑ_top = ɑ[-1, 0, :, 0]
ɑ_middle = ɑ[-1, 2, :, 0]
ɑ_bottom = ɑ[-1, 4, :, 0]

# Extract final stresses along thickness integration points
σ_top = σ[-1, 0, :, :]
σ_middle = σ[-1, 2, :, :]
σ_bottom = σ[-1, 4, :, :]

mises_top = torch.sqrt(
    σ_top[:, 0, 0] ** 2
    - σ_top[:, 0, 0] * σ_top[:, 1, 1]
    + σ_top[:, 1, 1] ** 2
    + 3 * σ_top[:, 0, 1] ** 2
)

mises_middle = torch.sqrt(
    σ_middle[:, 0, 0] ** 2
    - σ_middle[:, 0, 0] * σ_middle[:, 1, 1]
    + σ_middle[:, 1, 1] ** 2
    + 3 * σ_middle[:, 0, 1] ** 2
)

mises_bottom = torch.sqrt(
    σ_bottom[:, 0, 0] ** 2
    - σ_bottom[:, 0, 0] * σ_bottom[:, 1, 1]
    + σ_bottom[:, 1, 1] ** 2
    + 3 * σ_bottom[:, 0, 1] ** 2
)

print(f"Maximum displacement: {torch.linalg.norm(u[-1, :, 0:3], dim=1).max():.4f}")
print(f"Maximum plastic strain at top: {ɑ_top.max():.4f}")
print(f"Maximum plastic strain at middle: {ɑ_middle.max():.4f}")
print(f"Maximum plastic strain at bottom: {ɑ_bottom.max():.4f}")
print(f"Maximum von Mises stress at top: {mises_top.max():.4f}")
print(f"Maximum von Mises stress at middle: {mises_middle.max():.4f}")
print(f"Maximum von Mises stress at bottom: {mises_bottom.max():.4f}")

Maximum displacement: 45.6149
Maximum plastic strain at top: 0.0038
Maximum plastic strain at middle: 0.0000
Maximum plastic strain at bottom: 0.0038
Maximum von Mises stress at top: 288.1621
Maximum von Mises stress at middle: 0.0000
Maximum von Mises stress at bottom: 288.1621

Visualize results on deformed geometry

# Displacement on deformed configuration
cantilever.plot(
    u[-1, :, 0:3],
    element_property={"von Mises Stress (top)": mises_top},
    cmap="inferno",
)

# Displacement on deformed configuration
cantilever.plot(
    u[-1, :, 0:3], element_property={"PEEQ (top) in % ": ɑ_top * 100}, cmap="Blues"
)