2x2 square with axial load - Example from SolidsPy

[adtzlr]

This example solves the 2x2 linear-elastic plane-stress example from the docs of SolidsPy in FElupe. Results are identical.

import felupe as fem
import numpy as np

mesh = fem.Rectangle(a=(0, 0), b=(2, 2), n=3)
region = fem.RegionQuad(mesh)
field = fem.FieldContainer([fem.Field(region, dim=2)])

bottom_corners = np.logical_and(mesh.y == 0, np.logical_or(mesh.x == 0, mesh.x == 2))
boundaries = fem.BoundaryDict(
    bottom_corners=fem.Boundary(field[0], mask=bottom_corners, skip=(1, 0)),
    bottom_center=fem.Boundary(field[0], fy=0, fx=1, mode="and"),
)
umat = fem.constitution.LinearElasticPlaneStress(E=1.0, nu=0.3)
solid = fem.SolidBody(umat=umat, field=field)

top_corners = np.logical_and(mesh.y == 2, np.logical_or(mesh.x == 0, mesh.x == 2))
top_center = np.logical_and(mesh.y == 2, mesh.x == 1)
loads = [
    fem.PointLoad(field, points=np.where(top_corners)[0], values=[0.0, 1.0]),
    fem.PointLoad(field, points=np.where(top_center)[0], values=[0.0, 2.0]),
]

loads[1].plot(plotter=loads[0].plot(scale=0.1), scale=0.2).show()

step = fem.Step(items=[solid, *loads], boundaries=boundaries)
job = fem.Job(steps=[step]).evaluate()

field.plot("Displacement", component=0).show()
field.plot("Displacement", component=1).show()