Skip to content

Shell

Shell

Bases: ShellGeometry, Mechanics

Flat-facet shell model for thin-walled structures.

Triangles follow Krysl, quadrilaterals the MITC4 shear interpolation of Dvorkin and Bathe on the mean plane of their four nodes, so warp is neglected. Each node carries six degrees of freedom (three translations, three rotations). The section is either a homogeneous plane-stress material with a thickness or a layered Laminate.

The drilling rotation is penalized towards the rotation of the membrane field rather than transformed out as in section 3 of Krysl, so drill_penalty remains a tuning parameter on a coarse doubly-curved mesh, and a folded or branched shell, whose nodes carry no unique normal, is not supported.

solve reports stress in the local material frame of each element, and leaves the gradient at the identity, as the formulation never forms one. Its default aggregation averages the through-thickness stations, where a bending stress cancels: pass aggregate_integration_points=False to read one.

Attributes:

  • nodes –

    Nodal coordinates with shape [n_nod, 3].

  • elements –

    Triangle or quadrilateral connectivity with shape [n_elem, 3] or [n_elem, 4].

  • material (Material | None) –

    Vectorized plane-stress material (None for laminate shells).

  • section (Laminate | None) –

    Laminate section (None for homogeneous shells).

  • thickness –

    Element thicknesses with shape [n_elem].

  • offset –

    Reference surface position per element with shape [n_elem], as a fraction of thickness from the mid-plane.

  • orientation –

    Per-element material reference direction with shape [n_elem, 3].

  • forces (Tensor) –

    Applied nodal forces and moments with shape [n_nod, 6].

  • displacements (Tensor) –

    Prescribed nodal displacements and rotations with shape [n_nod, 6].

  • constraints (Tensor) –

    Boolean mask of constrained DOFs with shape [n_nod, 6].

__init__(nodes, elements, material, thickness=1.0, offset=0.0, transverse_nu=0.5, transverse_kappa=5.0 / 6.0, transverse_G=None, drill_penalty=0.001, n_simpson=3, orientation=None)

Initialize the shell FEM problem.

Parameters:

  • nodes (Tensor) –

    Nodal coordinates with shape [n_nod, 3].

  • elements (Tensor) –

    Triangle or quadrilateral connectivity with shape [n_elem, 3] or [n_elem, 4].

  • material (Material | Laminate) –

    Either a single plane-stress Material (homogeneous shell) or a Laminate describing a layered stacking sequence. When a Laminate is passed, thickness and n_simpson are taken from the laminate and the corresponding arguments here are ignored.

  • thickness (Tensor | float, default: 1.0 ) –

    Shell thickness. A float is expanded to all elements, a tensor assigns one thickness per element.

  • offset (Tensor | float, default: 0.0 ) –

    Reference surface position within the section, as a fraction of thickness from the mid-plane along the element normal. +0.5 puts it on the top face.

  • transverse_nu (float, default: 0.5 ) –

    Poisson's ratio used for the shear relaxation of a homogeneous triangle. A quadrilateral needs none.

  • transverse_kappa (float, default: 5.0 / 6.0 ) –

    Shear correction factor, ⅚ for a homogeneous section.

  • transverse_G (list[float] | list[Tensor] | None, default: None ) –

    Pair [G_xz, G_yz] of effective transverse shear moduli, integrated over the thickness. Taken from the material or the laminate when omitted.

  • drill_penalty (float, default: 0.001 ) –

    Stiffness tying the drilling degree of freedom to the in-plane rotation of the membrane field, as a fraction of the shear stiffness of the section.

  • n_simpson (int, default: 3 ) –

    Number of Simpson integration points through the thickness. Must be an odd integer.

  • orientation (Tensor | None, default: None ) –

    Global reference direction from which material/ply angles are measured. It is projected onto each element's surface to define the element's local material 0°-axis. Accepts a single (3,) vector (shared by all elements) or a per-element (n_elem, 3) tensor. Defaults to the global x-direction.

solve(increments=None, max_iter=10, rtol=1e-08, atol=1e-06, stol=1e-10, cutback_factor=0.5, growth_factor=1.1, max_cutbacks=10, verbose=False, method=None, preconditioner=None, device=None, return_intermediate=False, aggregate_integration_points=True, alpha=0.0, differentiable_parameters=None)

Solve the quasi-static finite-element problem by load increments.

Parameters:

  • increments (Tensor | None, default: None ) –

    Load scale factors, typically [0, 1]. They may rise and fall, so a load cycle is expressed as a sequence like [0, 1, 0]. Results are always returned at exactly these values. If a Newton solve does not converge, the increment is subdivided internally and retried, and the substep is grown again after each success.

  • max_iter (int, default: 10 ) –

    Maximum Newton iterations before an increment is cut back.

  • rtol (float, default: 1e-08 ) –

    Relative residual tolerance for Newton convergence.

  • atol (float, default: 1e-06 ) –

    Absolute residual tolerance for Newton convergence.

  • stol (float, default: 1e-10 ) –

    Tolerance used by iterative linear solvers.

  • cutback_factor (float, default: 0.5 ) –

    Factor applied to the substep size after a Newton solve failed to converge.

  • growth_factor (float, default: 1.1 ) –

    Factor applied to the substep size after a Newton solve converged, capped at the requested increment.

  • max_cutbacks (int, default: 10 ) –

    Number of successive cutbacks accepted within an increment before the solve is given up.

  • verbose (bool, default: False ) –

    If True, reports the solver configuration and a table of per-increment progress, updated in place inside notebooks.

  • method (Literal['direct', 'cg', 'bicgstab'] | None, default: None ) –

    Linear solver method, chosen by size and tangent symmetry when omitted.

  • preconditioner (Literal['amg', 'jacobi', 'none'] | None, default: None ) –

    Preconditioner for an iterative method, chosen by device and available backends when omitted.

  • device (str | None, default: None ) –

    Optional device hint for the linear solver backend.

  • return_intermediate (bool, default: False ) –

    If True, returns values for all increments.

  • aggregate_integration_points (bool, default: True ) –

    If True, averages flux, gradient, and state over integration points.

  • alpha (float, default: 0.0 ) –

    Damping factor for viscous stabilization. Dissipated energy is accumulated in self.stabilization_energy.

  • differentiable_parameters (Tensor | Iterable[Tensor] | None, default: None ) –

    Explicit parameter(s) to differentiate through implicit Newton/sparse solves. Accepts either a single tensor or an iterable of tensors.

Returns:

  • tuple[Tensor, Tensor, Tensor, Tensor, Tensor] –

    Tuple of displacement, internal force, flux, gradient, and material state. If return_intermediate is True, each tensor includes an increment dimension as the leading axis.

integrate_surface_load(mask, load)

Consistent nodal loads from a load per unit area, e.g. a pressure.

A shell element is its own surface, so the loaded surface is made up of the elements whose nodes all lie in mask.

Parameters:

  • mask (Tensor) –

    Boolean nodal mask with shape [n_nod] selecting the surface.

  • load (float | Tensor) –

    Load per unit area. A float is a pressure acting along the element normal, while shape [3] or [n_elem, 3] is a traction in global coordinates.

Returns:

  • Tensor –

    Nodal loads with shape [n_nod, k], to be added to forces[:, 0:3]

  • Tensor –

    or heat_flux.

plot(u=0.0, node_property=None, element_property=None, orientations=None, thickness=True, mirror=(False, False, False), show_undeformed=False, axes=False, bcs=False, plotter=None, camera=None, **kwargs)

Plot the shell mesh with PyVista, optionally with results.

Parameters:

  • u (float | Tensor, default: 0.0 ) –

    Nodal displacements added to the positions, e.g. to plot the deformed configuration. Defaults to 0.0 (undeformed).

  • node_property (Tensor | dict[str, Tensor] | None, default: None ) –

    Nodal field, optionally keyed by its color bar title, e.g. {"u": u[:, :3]}.

  • element_property (Tensor | dict[str, Tensor] | None, default: None ) –

    Element field, keyed like node_property.

  • orientations (Tensor | None, default: None ) –

    Per-element direction vectors with shape [n_elem, k, 3] with k <= 3, e.g. the local frames self.t, drawn on the unmirrored mesh as red, green, and blue arrows of the mean element size.

  • thickness (bool, default: True ) –

    If True, extrudes elements by their thickness.

  • mirror (tuple[bool, bool, bool], default: (False, False, False) ) –

    Mirrors the mesh about the (x, y, z) planes, e.g. to visualize symmetric halves. Warns if the nodes on a mirrored plane are not constrained to enforce that symmetry.

  • show_undeformed (bool, default: False ) –

    If True, draws the undeformed mesh as a grey wireframe.

  • axes (bool, default: False ) –

    If True, shows labeled coordinate axes around the mesh.

  • bcs (bool, default: False ) –

    If True, renders boundary conditions on the unmirrored mesh: arrows for forces and prescribed displacements, spheres at displacement tips, and a cone per constrained DOF. Rotational DOFs use doubled heads, the usual convention for moments. Constraints enforcing the symmetry of a mirrored plane are skipped, since the mirrored copy shows that symmetry already.

  • plotter (Plotter | None, default: None ) –

    PyVista plotter. Defaults to None.

  • camera (CameraPositionOptions | None, default: None ) –

    Camera position, either a plane ("xy", "xz", "yz"), "iso", or an explicit position, focal point and view up. Defaults to None.

  • **kwargs –

    Forwarded to pyvista.Plotter.add_mesh.

solve_modes(n_modes)

Compute the natural frequencies and mode shapes.

Solves the generalized eigenvalue problem

\[\mathbf{K}\boldsymbol{\phi} = \omega^2 \mathbf{M}\boldsymbol{\phi}\]

Parameters:

  • n_modes (int) –

    Number of eigenpairs to compute.

Returns:

  • tuple[Tensor, Tensor] –

    Tuple (omega_sq, modes) where omega_sq has shape [n_modes] (squared angular frequencies, differentiable) and modes has shape [n_modes, n_nod, n_dof_per_node] (detached).

ShellHeat

Bases: ShellGeometry, Heat

Heat conduction model for thin-walled structures.

Uses the same flat facets, local frames and plotting as Shell, with one temperature per node, so the section conducts in-plane only.

solve reports flux and gradient in the local material frame, as Shell reports its stress.

Attributes:

  • nodes –

    Nodal coordinates with shape [n_nod, 3].

  • elements –

    Triangle or quadrilateral connectivity with shape [n_elem, 3] or [n_elem, 4].

  • material (Material | None) –

    Vectorized plane thermal material.

  • thickness –

    Element thicknesses with shape [n_elem].

  • orientation –

    Per-element material reference direction with shape [n_elem, 3].

  • heat_flux (Tensor) –

    Applied nodal heat sources with shape [n_nod, 1].

  • temperatures (Tensor) –

    Prescribed nodal temperatures with shape [n_nod, 1].

  • constraints (Tensor) –

    Boolean mask of constrained DOFs with shape [n_nod, 1].

__init__(nodes, elements, material, thickness=1.0, offset=0.0, orientation=None)

Initialize the shell FEM problem.

Parameters:

  • nodes (Tensor) –

    Nodal coordinates with shape [n_nod, 3].

  • elements (Tensor) –

    Triangle or quadrilateral connectivity with shape [n_elem, 3] or [n_elem, 4].

  • material (Material | None) –

    Pointwise material model, or None for a shell whose section carries the materials.

  • thickness (Tensor | float, default: 1.0 ) –

    Shell thickness. A float is expanded to all elements, a tensor assigns one thickness per element.

  • offset (Tensor | float, default: 0.0 ) –

    Reference surface position within the section, as a fraction of thickness from the mid-plane along the element normal. +0.5 puts it on the top face.

  • orientation (Tensor | None, default: None ) –

    Global reference direction from which material/ply angles are measured. It is projected onto each element's surface to define the element's local material 0°-axis. Accepts a single (3,) vector (shared by all elements) or a per-element (n_elem, 3) tensor. Defaults to the global x-direction.

solve(increments=None, max_iter=10, rtol=1e-08, atol=1e-06, stol=1e-10, cutback_factor=0.5, growth_factor=1.1, max_cutbacks=10, verbose=False, method=None, preconditioner=None, device=None, return_intermediate=False, aggregate_integration_points=True, alpha=0.0, differentiable_parameters=None)

Solve the quasi-static finite-element problem by load increments.

Parameters:

  • increments (Tensor | None, default: None ) –

    Load scale factors, typically [0, 1]. They may rise and fall, so a load cycle is expressed as a sequence like [0, 1, 0]. Results are always returned at exactly these values. If a Newton solve does not converge, the increment is subdivided internally and retried, and the substep is grown again after each success.

  • max_iter (int, default: 10 ) –

    Maximum Newton iterations before an increment is cut back.

  • rtol (float, default: 1e-08 ) –

    Relative residual tolerance for Newton convergence.

  • atol (float, default: 1e-06 ) –

    Absolute residual tolerance for Newton convergence.

  • stol (float, default: 1e-10 ) –

    Tolerance used by iterative linear solvers.

  • cutback_factor (float, default: 0.5 ) –

    Factor applied to the substep size after a Newton solve failed to converge.

  • growth_factor (float, default: 1.1 ) –

    Factor applied to the substep size after a Newton solve converged, capped at the requested increment.

  • max_cutbacks (int, default: 10 ) –

    Number of successive cutbacks accepted within an increment before the solve is given up.

  • verbose (bool, default: False ) –

    If True, reports the solver configuration and a table of per-increment progress, updated in place inside notebooks.

  • method (Literal['direct', 'cg', 'bicgstab'] | None, default: None ) –

    Linear solver method, chosen by size and tangent symmetry when omitted.

  • preconditioner (Literal['amg', 'jacobi', 'none'] | None, default: None ) –

    Preconditioner for an iterative method, chosen by device and available backends when omitted.

  • device (str | None, default: None ) –

    Optional device hint for the linear solver backend.

  • return_intermediate (bool, default: False ) –

    If True, returns values for all increments.

  • aggregate_integration_points (bool, default: True ) –

    If True, averages flux, gradient, and state over integration points.

  • alpha (float, default: 0.0 ) –

    Damping factor for viscous stabilization. Dissipated energy is accumulated in self.stabilization_energy.

  • differentiable_parameters (Tensor | Iterable[Tensor] | None, default: None ) –

    Explicit parameter(s) to differentiate through implicit Newton/sparse solves. Accepts either a single tensor or an iterable of tensors.

Returns:

  • tuple[Tensor, Tensor, Tensor, Tensor, Tensor] –

    Tuple of displacement, internal force, flux, gradient, and material state. If return_intermediate is True, each tensor includes an increment dimension as the leading axis.

integrate_surface_load(mask, load)

Consistent nodal loads from a load per unit area, e.g. a pressure.

A shell element is its own surface, so the loaded surface is made up of the elements whose nodes all lie in mask.

Parameters:

  • mask (Tensor) –

    Boolean nodal mask with shape [n_nod] selecting the surface.

  • load (float | Tensor) –

    Load per unit area. A float is a pressure acting along the element normal, while shape [3] or [n_elem, 3] is a traction in global coordinates.

Returns:

  • Tensor –

    Nodal loads with shape [n_nod, k], to be added to forces[:, 0:3]

  • Tensor –

    or heat_flux.

plot(u=0.0, node_property=None, element_property=None, orientations=None, thickness=True, mirror=(False, False, False), show_undeformed=False, axes=False, bcs=False, plotter=None, camera=None, **kwargs)

Plot the shell mesh with PyVista, optionally with results.

Parameters:

  • u (float | Tensor, default: 0.0 ) –

    Nodal displacements added to the positions, e.g. to plot the deformed configuration. Defaults to 0.0 (undeformed).

  • node_property (Tensor | dict[str, Tensor] | None, default: None ) –

    Nodal field, optionally keyed by its color bar title, e.g. {"u": u[:, :3]}.

  • element_property (Tensor | dict[str, Tensor] | None, default: None ) –

    Element field, keyed like node_property.

  • orientations (Tensor | None, default: None ) –

    Per-element direction vectors with shape [n_elem, k, 3] with k <= 3, e.g. the local frames self.t, drawn on the unmirrored mesh as red, green, and blue arrows of the mean element size.

  • thickness (bool, default: True ) –

    If True, extrudes elements by their thickness.

  • mirror (tuple[bool, bool, bool], default: (False, False, False) ) –

    Mirrors the mesh about the (x, y, z) planes, e.g. to visualize symmetric halves. Warns if the nodes on a mirrored plane are not constrained to enforce that symmetry.

  • show_undeformed (bool, default: False ) –

    If True, draws the undeformed mesh as a grey wireframe.

  • axes (bool, default: False ) –

    If True, shows labeled coordinate axes around the mesh.

  • bcs (bool, default: False ) –

    If True, renders boundary conditions on the unmirrored mesh: arrows for forces and prescribed displacements, spheres at displacement tips, and a cone per constrained DOF. Rotational DOFs use doubled heads, the usual convention for moments. Constraints enforcing the symmetry of a mirrored plane are skipped, since the mirrored copy shows that symmetry already.

  • plotter (Plotter | None, default: None ) –

    PyVista plotter. Defaults to None.

  • camera (CameraPositionOptions | None, default: None ) –

    Camera position, either a plane ("xy", "xz", "yz"), "iso", or an explicit position, focal point and view up. Defaults to None.

  • **kwargs –

    Forwarded to pyvista.Plotter.add_mesh.