2. Composite material modeling fundamentals

Before any progressive failure analysis can be trusted, the constitutive behavior of a single lamina must be defined correctly. A unidirectional composite ply is not isotropic: it is stiff and strong along the fiber direction and far more compliant and weaker transverse to it. This directional dependence is the foundation of everything that follows — ply stresses, failure indices, delamination onset, and final structural collapse all inherit the assumptions made here. This section establishes the orthotropic elasticity framework, the engineering constants that populate it, and the coordinate transformations that let a single lamina definition be reused at any fiber angle.

From Isotropy to Orthotropy

An isotropic material needs only two independent elastic constants, typically Young's modulus E and Poisson's ratio ν. A unidirectional lamina requires four independent in-plane constants because it has two perpendicular planes of material symmetry: one normal to the fiber axis and one normal to the transverse in-plane axis. The standard engineering constants are:

• E₁ — longitudinal modulus, measured along the fiber direction (typically 120–180 GPa for carbon/epoxy).
• E₂ — transverse modulus, measured perpendicular to the fibers in the ply plane (typically 8–12 GPa).
• ν₁₂ — major Poisson's ratio, the transverse contraction under longitudinal extension (typically 0.25–0.35).
• G₁₂ — in-plane shear modulus (typically 4–6 GPa).

The ratio E₁/E₂ is often 15–20, which is precisely why ignoring orthotropy produces meaningless results. For a plane-stress lamina (thin ply, no through-thickness stress), the reduced stiffness matrix Q in the material (1–2) coordinate system is:

Q11 = E1 / (1 - ν12·ν21)
Q22 = E2 / (1 - ν12·ν21)
Q12 = ν12·E2 / (1 - ν12·ν21)
Q66 = G12
where ν21 = ν12 · E2 / E1   (reciprocity: ν12·E2 = ν21·E1)

The stress–strain relation is then {σ} = [Q]{ε} with σ = (σ₁, σ₂, τ₁₂) and ε = (ε₁, ε₂, γ₁₂). Note that Q₁₂ couples the two normal strains through Poisson's effect, while Q₆₆ is fully decoupled from the normal terms — a direct consequence of orthotropy.

Why Coordinate Transformation Is Unavoidable

Real laminates stack plies at angles such as 0°, ±45°, and 90°. The material axes (1, 2) of each ply therefore do not align with the structural axes (x, y) used to assemble the global stiffness. A transformation is required to rotate [Q] from material coordinates into structural coordinates, producing the transformed reduced stiffness matrix [Q̄].

Let θ be the angle from the structural x-axis to the fiber 1-axis, measured counterclockwise positive. The transformation uses m = cos θ and n = sin θ:

Q̄11 = Q11·m⁴ + 2(Q12 + 2Q66)·m²n² + Q22·n⁴
Q̄22 = Q11·n⁴ + 2(Q12 + 2Q66)·m²n² + Q22·m⁴
Q̄12 = (Q11 + Q22 - 4Q66)·m²n² + Q12·(m⁴ + n⁴)
Q̄16 = (Q11 - Q12 - 2Q66)·m³n - (Q22 - Q12 - 2Q66)·mn³
Q̄26 = (Q11 - Q12 - 2Q66)·mn³ - (Q22 - Q12 - 2Q66)·m³n
Q̄66 = (Q11 + Q22 - 2Q12 - 2Q66)·m²n² + Q66·(m⁴ + n⁴)

The critical observation is that Q̄16 and Q̄26 are nonzero for any angle other than 0° or 90°. These terms represent extension–shear coupling: pulling an off-axis ply along x also generates shear strain, and shearing it generates normal strain. This coupling is the physical origin of bend–twist coupling in unsymmetric laminates and of the complex stress states that drive matrix cracking in angle plies. A simulation that omits Q̄16 and Q̄26 will underpredict shear stresses and misplace failure initiation.

Worked Example: Transforming a Carbon/Epoxy Ply

Consider a standard carbon/epoxy lamina with E₁ = 140 GPa, E₂ = 10 GPa, ν₁₂ = 0.30, G₁₂ = 5.0 GPa. First compute ν₂₁ = 0.30 × 10/140 = 0.02143, then the denominator 1 − ν₁₂ν₂₁ = 1 − 0.006429 = 0.99357.

Q11 = 140 / 0.99357 = 140.9 GPa
Q22 = 10  / 0.99357 = 10.06 GPa
Q12 = 0.30 × 10 / 0.99357 = 3.019 GPa
Q66 = 5.0 GPa

Now rotate to θ = 45° (m = n = 0.7071, m² = n² = 0.5, m⁴ = n⁴ = 0.25, m²n² = 0.25):

Q̄11 = 140.9(0.25) + 2(3.019 + 10.0)(0.25) + 10.06(0.25)
     = 35.23 + 6.51 + 2.515 = 44.25 GPa
Q̄12 = (140.9 + 10.06 - 20.0)(0.25) + 3.019(0.5) = 32.74 + 1.510 = 34.25 GPa
Q̄16 = (140.9 - 3.019 - 10.0)(0.3536) - (10.06 - 3.019 - 10.0)(0.3536)
     = 45.22 - (-1.046) = 46.27 GPa
Q̄66 = (140.9 + 10.06 - 6.038 - 10.0)(0.25) + 5.0(0.5) = 33.73 + 2.50 = 36.23 GPa

The large Q̄16 = 46.27 GPa confirms strong extension–shear coupling at 45°. If a ±45° laminate is loaded in tension, each ply experiences significant shear stress even though the applied load is purely axial — a result that isotropic intuition would miss entirely.

From Ply to Laminate: Classical Lamination Theory

Once each ply has been transformed into structural coordinates, the laminate stiffness is assembled through Classical Lamination Theory (CLT). The extensional stiffness matrix A, coupling matrix B, and bending stiffness matrix D are obtained by integrating through the thickness:

A_ij = Σ Q̄_ij^(k) · (z_k − z_{k−1})
B_ij = ½ Σ Q̄_ij^(k) · (z_k² − z_{k−1}²)
D_ij = ⅓ Σ Q̄_ij^(k) · (z_k³ − z_{k−1}³)

Here z_k are the through-thickness coordinates of each ply interface, measured from the laminate mid-plane. For a symmetric laminate, B = 0 and there is no extension–bending coupling. For a balanced laminate (equal +θ and −θ plies), the A₁₆ and A₂₆ terms cancel, eliminating in-plane extension–shear coupling at the laminate level even though each individual ply retains nonzero Q̄16 and Q̄26. This distinction between ply-level and laminate-level coupling is essential: a simulation may show no global shear coupling while individual plies still carry significant shear stress.

Practical Modeling Checklist

When setting up a composite model in FEA or a CLT script, verify the following before running any failure analysis:

The orthotropic elasticity framework, the four independent lamina constants, and the rotation to structural coordinates form the complete input needed for CLT and for any progressive failure model. Every subsequent chapter — stress recovery, failure criteria, degradation rules — operates on the [Q̄] matrices and A/B/D stiffnesses defined here. Getting these right is not optional; it is the difference between a simulation that predicts failure and one that merely produces numbers.