hermit.fenics.ops

The FEniCSx custom operations. Called only from hermit._solve / hermit.outputs / hermit.transfer.

ShellSolveOp – implicit; assemble R/J, linear (or Newton) solve, adjoint. ShellScalarFormsOp (scalar output forms) – explicit, dense Jacobians. ShellFieldFormsOp (field outputs – one UFL expression per FieldSpec into its target space, by project / interpolate / average / midpoint) – explicit, reverse-mode only (one cotangent solve or matvec per field, never a dense Jacobian). project / average / midpoint keep mesh-coordinate derivatives; interpolate does not.

Port of StateOperation / OutputOperation / OutputFieldOperation (femo dev_coupling), self-contained (no FEA). The shell tangent is symmetric, so one MUMPS factorization of dR/dw serves both the forward solve and the adjoint.

Classes

`FieldSpec`	One field output: a UFL expression + how to represent it as FE coefficients.
`ShellFieldFormsOp`	Field outputs -- several UFL expressions -> FE coefficients, in one op.
`ShellScalarFormsOp`	Assemble several scalar UFL output forms + their partials in one op.
`ShellSolveOp`	`orientation`, given, is `(name, space)` with `name` in

Functions

`orientation_cos_sin_spec`(pde, orientation_name, ...[, ...])	`FieldSpec` for `(cos theta, sin theta)` of the fibre orientation angle.
`strain_fields`(pde[, space, method, frame])	`name -> FieldSpec` for the membrane strain / bending curvature / transverse

Module Contents

class hermit.fenics.ops.FieldSpec

One field output: a UFL expression + how to represent it as FE coefficients.

space is (family, degree) (blocked by n_components when the op builds the FunctionSpace). method:

"project" (default) – L2 projection M c = b. Any space; one mass factorization; full derivatives including mesh coordinates.
"interpolate" – collocate expr at the space’s nodal points. Any Lagrange space; cheaper (no solve); exact where expr already lies in the space. Has no mesh-coordinate (shape) derivative – use it only when the geometry is fixed.
"average" / "midpoint" – (1/vol_K) * integral_K(expr), i.e. L2 projection onto ("DG", 0) (diagonal mass -> a divide, no solve); "average" uses full quadrature, "midpoint" a single centroid point. DG0 only; full derivatives including mesh coordinates.

args: tuple = ('disp_solid',)

expr: object

frame: str = 'local'

kind: str = 'scalar'

method: str = 'project'

n_components: int

space: tuple = ('DG', 2)

class hermit.fenics.ops.ShellFieldFormsOp(pde, fields: dict, differentiable_geometry: bool = False, coefficients=None, quadrature_degree=4)

Bases: csdl_alpha.CustomExplicitOperation

Field outputs – several UFL expressions -> FE coefficients, in one op.

fields maps output name -> FieldSpec. Each output is (n_scalar_dofs, n_components) in that field’s target space, FE dof order (wrap in hermit._field.Field for user ordering, point evaluation and frame changes).

Reverse-mode only (compute_jacvec_product), like femo’s OutputFieldOperation: one cotangent solve / matvec per field, never a dense Jacobian. method="interpolate" has no mesh_nodes derivative – the op raises if differentiable_geometry is set with an interpolated field.

compute(input_vals, output_vals)

compute_jacvec_product(inputs, outputs, d_inputs, d_outputs, mode)

evaluate(fe)

arg_names

fields

gdim

nel

pde

tdim

class hermit.fenics.ops.ShellScalarFormsOp(pde, forms: dict, differentiable_geometry: bool = False, orientation=None, coefficients=None)

Bases: csdl_alpha.CustomExplicitOperation

Assemble several scalar UFL output forms + their partials in one op.

Port of femo’s OutputOperation (multi-form). forms maps output name -> (ufl_form, tuple of arg names). Function args are written onto the persistent ShellPDE functions; mesh_nodes writes the mesh geometry (scratch buffer, restored after each call) and, when differentiable, contributes a d(form)/d(SpatialCoordinate) partial scattered back to user node ordering.

orientation, given, is (name, space) (name in ("fiber_angle", "fiber_direction")) – needed whenever one of forms’ arg-name tuples includes it (elastic_energy depends on theta for an oriented composite; a strain-field form does not). Unlike the fixed-space _ARG_TO_PDE args, an orientation coefficient’s space is caller-chosen, so it is resolved through pde.coefficient(name, space) (see that method) rather than a fixed pde attribute – the same seam ShellSolveOp uses, so both ops reference the exact same persistent Function the form the caller built was differentiated against.

coefficients, given, is a dict {name: space} generalizing that same seam to any other arg name a caller-built form needs on a caller-chosen space (Phase F’s per-load-term / arbitrary-constitutive_space forms – e.g. hm.solve’s ShellPDE.compliance_form(loads=...), whose args are per-term names like "traction_0" with no fixed pde attribute at all). None (default, every pre-Phase-F call site) leaves every non-orientation arg on the fixed _ARG_TO_PDE route, unchanged.

compute(input_vals, output_vals)

compute_derivatives(input_vals, output_vals, derivatives)

evaluate(fe)

arg_names

forms

gdim

pde

class hermit.fenics.ops.ShellSolveOp(pde, bc, arg_names, options, form_cache=None, orientation=None, material=None, loads=None)

Bases: csdl_alpha.experimental.CustomImplicitOperation

orientation, given, is (name, space) with name in ("fiber_angle", "fiber_direction") – the in-form Teps(theta) rotation (see hermit.fenics.elastic_model.ElasticModel). name must also appear in arg_names for it to be a declared (differentiable) CSDL input.

material, given, is a dict {"A": space, "B": space, "D": space, "As": space} (all four): each coefficient is routed through pde.coefficient(name, space) instead of the PDE’s fixed VABD/VAs space, so an arbitrary caller-chosen constitutive_space reaches the solve. None (default) keeps the fixed-space Functions.

loads, given, is a list of (argname, kind, space) triples: kind in "traction" / "moment" / "pressure", each routed through pde.coefficient(argname, space) and contributing its own residual/compliance term (see ElasticModel.weak_residual / ShellPDE.compliance_form’s load_terms), rather than being reduced onto the fixed f/m VF coefficients. Every argname must also appear in arg_names. None (default) keeps the fixed-VF f/m mechanism; an empty list selects the per-term path with zero distributed terms (e.g. a point-load-only Loads).

The residual is not structurally fixed once orientation, material spaces, load composition, and multi-region penalty terms can all vary – form_cache is keyed by a structural signature of all four, not by a single fixed string, so two structurally different compositions on the same ShellDomain never share a (wrong) compiled residual.

apply_inverse_jacobian(inputs, outputs, d_outputs, d_residuals, mode)

compute_jacvec_product(inputs, outputs, d_inputs, d_outputs, d_residuals, mode)

evaluate(fe)

solve_residual_equations(inputs, outputs)

Vc

arg_names

bc

dR_dw

gdim

ndof

options

orientation = None

pde

residual

state_name = 'disp_solid'

hermit.fenics.ops.orientation_cos_sin_spec(pde, orientation_name, coeff_space, target_space=('DG', 0), method='average')

FieldSpec for (cos theta, sin theta) of the fibre orientation angle.

The angle between a global fiber_direction and each element’s in-plane frame is a function of the mesh: moving nodes rotates the frame. Building the cos/sin in UFL – the same expression the shell form itself uses – keeps that dependence in the differentiated path, so the geometry derivative reaches the Tsai-Wu outputs. Evaluating the angle in numpy instead freezes it at the reference geometry and silently drops a term from the shape gradient.

fiber_angle is element-relative and carries no geometry dependence, but it is accepted here too so both orientation kinds take one path.

hermit.fenics.ops.strain_fields(pde, space=None, method=None, frame=None)

name -> FieldSpec for the membrane strain / bending curvature / transverse shear (Voigt / engineering convention). Natural space ("DG", 2), natural method "project"; override with space / method.

frame selects the components:

"local" – the element in-plane frame: strain / curvature as a length-3 in-plane Voigt vector, shear as a length-2 vector. Only well defined on a discontinuous space (adjacent elements have different in-plane frames, so a continuous space would average frame-inconsistent components at shared nodes).
"global" – global Cartesian: strain / curvature as a length-6 symmetric tensor [xx, yy, zz, 2yz, 2xz, 2xy] (engineering shear), shear as a length-3 vector. Frame-consistent -> works on CG spaces.

frame=None picks "local" for a DG space and "global" for a CG space.