How can a certified admissible constitutive operator model phase-field fracture in composite shells?

A replaceable operator replaces one fixed stored energy

Conventional phase-field models of fracture in fibre-reinforced shells are written against a single closed-form stored energy. The stress, the consistent tangent and the plane-stress condensation are then derived and verified for that expression alone, so exchanging the energy means re-deriving and re-verifying each of them. The paper treats that obstacle as structural rather than computational and moves the constitutive law into a replaceable operator.

The operator sits at stage 4 of an eight-stage formulation. Stages 1-3 carry the geometrically exact Reissner-Mindlin shell kinematics to the Green-Lagrange strain, to its Cholesky pull-back, and to the invariant basis. Stages 5-8 handle the crack driving force, channel degradation, assembly of the regularised functional and the solve. Those downstream stages are independent of the representation seated in stage 4, so a closed form and a learned representation pass the same admission gate.

An exact closure identity makes the channels a change of basis

The strain is normalised by the Cholesky factor of the reference metric and expressed on four mutually orthogonal invariants whose closure identity is exact. Because the identity is exact, separating the stored energy into a fibre channel, an inter-fibre channel and a fibre-transverse interaction channel is a change of basis rather than a modelling assumption. The three channel potentials become the only material input, and everything downstream is written in those potentials and their derivatives.

On a curved midsurface the metric, its Cholesky factor and the fibre direction vary through the thickness through the shifter I - zeta b, and the invariants inherit that dependence. The basis is therefore defined on curved as well as flat midsurfaces, and its chart invariance is certified in the numerical studies. The fibre-transverse interaction energy is stored once and degraded by the product of the fields whose mechanisms it couples, so the coupling appears in the residual as well as in the Hessian.

Six conditions define which operators are admissible

Admissibility is stated as six conditions in Definition A.1. Within them, the two linear invariants are proved convex in the deformation gradient directly. The corresponding statement for the two quadratic invariants is left open and discharged instead by a pointwise certificate, so convexity of the full operator is not claimed analytically.

Two further properties are proved to survive an exchange of representation. The geometric part of the tangent is nonzero and independent of the material tangent, so K = K_mat + K_geo still requires both terms even when the operator is learned. And the plane-stress condensation stays regular wherever the residual stiffness is positive. Because plane-stress enforcement is a constitutive operation performed at every quadrature point, its regularity is a property of the energy rather than of the element.

Three tiers separate algebra, solver state and training data

Tier A certifies the algebra. It checks the residual and tangent against first and second variations of the functional across all blocks, using a control that omits the geometric tangent. Finite-difference probes are flagged rather than scored as failures whenever the perturbation crosses one of the two switching surfaces introduced by the Macaulay split.

Tier B certifies the solver state. The accepted-equilibrium residual is measured on the complementarity object minimised by the active-set method rather than on plain stationarity, the active set must be unchanged after a classify-reset-solve cycle, and definiteness of the assembled tangent is decided at the accepted state by a shifted-inertia count. Because that count is a local, state-dependent signal that may legitimately fail near a limit state, the integer count is reported rather than the verdict alone.

Tier C certifies learning. A constitutive state enters a training set only if it passes both Tier A and Tier B at that state. This is a checkable admissibility criterion for data, orthogonal to whether the underlying law is analytic or learned.

Notched IM7/8552 shells verify the discretisation before fracture is interpreted

The numerical studies are ordered so that each relies only on what precedes it. The single-edge-notched tension specimen is 36 x 72 mm with a shell thickness of 1 mm and a notch a0 = 9 mm at mid-height; the bottom row of nodes fixes all three displacement components and leaves the rotations free, and the top edge carries the prescribed displacement. The fibre direction is swept over theta = 0, 30, 60 and 90 degrees measured from the horizontal, so theta = 90 degrees places fibres parallel to the loading direction. A clamped bending specimen shares the plate dimensions, carries its notch 6 mm from the clamped edge, and is loaded out of plane at the tip.

The material is IM7/8552, and the property table separates the measured ply constants from the invariant coefficients derived from them, since the latter are what the operator consumes. The relations that produce them are printed with the table so the mapping from a datasheet to the invariant basis can be reproduced. The fibre toughness exceeds the inter-fibre toughness by a factor of 1630, and this ratio governs the orientation dependence reported in the studies.

Sections 7.3 and 7.4 verify the discretisation itself, in the thin regime and on curved geometry, before any quantity computed with it is interpreted. The thin-limit study examines Pi/t^3 and its spread under refinement, with rotations in the per-cell basis as the control. The curved-geometry chain is checked against a flat limit and chart invariance, with a spoofed geometry and a planted metric as controls. Only after those checks do the later studies use the certified model to examine how fibre orientation, the interface field and curvature control fracture.

Fibre orientation sets the damage envelope and load capacity

The notched-shell studies show that fibre orientation governs the damage envelope and the load capacity. The orientation sweep is paired with a matched peak-load state as the control, so the comparison is made at a common load level rather than at each solve's own end point. The mechanism-competition study examines load capacity and driving shares against a seeded notch baseline.

This matters because a model carrying a single scalar damage variable can reproduce a load-displacement curve but cannot report which mechanism consumed the energy. Fibre rupture, inter-fibre cracking of the matrix and interfacial separation each carry their own strength, their own toughness and their own dependence on the angle between fibre and load. The channel decomposition is what lets the formulation report which mechanism is loaded rather than only how much damage accumulated.

Curvature drives a through-thickness asymmetry one field cannot hold

On the curved strip, surface-resolved damage is plotted on the upper and lower surfaces at R = 25 mm under circumferential loading at the last converged state, on one shared scale fixed to 0-1, with the seeded notch common to both surfaces. The difference between the two surface fields is plotted on a diverging scale with the zero contour drawn. The same difference for the flat pair is identically zero at every step, so the field measures curvature alone and no residual asymmetry of numerical origin can be mistaken for the effect.

The separation between the two loading configurations is read at four radii at one common damage level rather than at each solve's own end point, since the legs terminate at different states. The flat value is -1.9 x 10^-14, so the curved values stand about thirteen orders above it. The separation is monotone in 1/R and follows a sub-linear law with exponent 0.375; the best line through the origin departs from the data by 51 per cent. An O(h/R) argument applied to the strain asymmetry alone would predict a linear law, so damage does not inherit the scaling of the field that drives it, which is consistent with a threshold process in which the driving force must exceed a critical value before any of it is converted into crack surface. The exponent is reported at one mesh and one reading level.

Two limitations attach to this result. The thin limit is not shown, because pairing curvature with thickness would require a surface-resolved thickness ladder, whereas the thickness legs available here are mechanism-resolved with damage uniform through the thickness and no route for the observable to vanish in h. And the study demonstrates a formulation capability rather than a calibrated structural prediction, since the interfacial toughness is uncalibrated and the geometry is a benchmark.

The learned operator is certified but not deployed

A corpus of 65,000 constitutive states was collected from eleven solves, each state admitted only if it passed Tier A and Tier B at the point where it was recorded, which is the Tier C criterion applied literally. In the inter-fibre channel plane, none of the 6000 held-out states falls outside the training box, so the held-out error measures interpolation rather than extrapolation.

The distilled network reaches 0.013 per cent relative error on the training set and 0.165 per cent on the held-out solve for the fibre channel, and 0.337 per cent and 0.605 per cent for the inter-fibre channel. The inter-fibre channel is the harder of the two, as expected, since it depends on two invariants rather than one and its state space is strongly non-uniform, with most samples concentrated near one invariant value and a thin tail reaching the corners of the support.

Two statements must be kept apart. The first concerns capacity: an admissible container built from convex activations on the quadratic slots and unconstrained weights on the linear slots reproduces the closed form to better than one per cent at interpolation scale while remaining admissible, since none of the six conditions is relaxed to obtain it. The second concerns deployment: the operator that produced every other result in this work runs with the network identically zero, so its error against the closed form is exactly zero by construction. Predictive quality on states outside the corpus is neither established nor claimed, and predictive learned constitutive modelling on this socket is left to a companion study.

Watch the channel decomposition, not the load curve alone

For readers comparing composite-shell fracture models, the decision-relevant feature of this paper is not the load-displacement response, which a single scalar damage variable can already reproduce. It is whether the formulation can report which mechanism consumed the energy and whether that reporting survives a change of constitutive representation. The exact closure identity is what makes the fibre, inter-fibre and interaction channels a change of basis rather than a modelling assumption, and the six admissibility conditions plus the three-tier protocol are what let a closed form and a learned operator pass the same gate.

The clearest scenario split in the evidence is between flat and curved geometry. On the flat pair the upper-lower surface damage difference is identically zero at every step, while on the curved strip it stands about thirteen orders above the flat value and grows sub-linearly in 1/R with exponent 0.375. A single midsurface damage field cannot represent that asymmetry, so curvature is the condition under which the two-surface treatment earns its cost. The fibre-orientation sweep, in turn, is where the 1630 toughness ratio between fibre and inter-fibre channels becomes visible in the damage envelope and load capacity.

The priority for anyone reusing this work is therefore to check the channel decomposition and the admissibility gate before trusting any downstream fracture quantity, and to treat the curvature asymmetry as a formulation capability demonstrated on a benchmark rather than a calibrated structural prediction, since the interfacial toughness is uncalibrated and the thin limit is not shown.

Sources used in this answer