In this section, a progressive damage model for composites (PDM) is presented (Curiel Sosa et al., 2008a; Curiel Sosa, 2008b), the results of this model are then compared with the outcome of some stress failure criteria. The formulation concerning the implementation of the damage model is also presented below. In PDM, it is assumed that different damage modes develop simultaneously on the failed composite structure and that they can interact directly influencing each other. The approach presented may be framed within the continuum damage mechanics field. The damage variables represent the state of damage at any stage of deterioration. They are implicitly defined for what they may be considered internal state variables that can provide a quantifiable magnitude of the degradation of the composite material. In PDM, damage variables are referred to damage modes via superposition. Thus, a damage mode is represented for one or more than one damage variables with distinct weights as explained below. The adopted damage modes, generically denoted as y, are modelled by means of a linear combination of growth functions Фт and damage directors vY (Curiel Sosa et al., 2008a). In summary, PDM admits the modelling of different damage modes that can be integrated directly in the description of damage through Equation (6)which is explained on Section (4).
Thus, the damage state may be defined by a series of internal variables w^j, filling the diagonal damage tensor D (see Equation (1)) that represents the state of damage in the composite.
& = D • & (1)
where &T = [a11, a22, v33, a12, (T23, a31 ] is an array formed by the stress components and
& is the socalled effective stress array (Chaboche, 1981). The definition of tensors and properties is conducted in a local system of reference for the lamina. Variables and parameters numeric subscripts refer to the local lamina system of reference. Thus, axis 1 is pointed
Finite Element Analysis of Progressive
Degradation versus Failure Stress Criteria on Composite Damage Mechanics in the longitudinal direction to the fibres whereas the other two axes, i. e. 2 and 3, are in perpendicular direction to fibres.
The damage tensor is built as a diagonal tensor and contains the damage internal variables Wkj, see Equation (2). These are responsible for the degradation of the stiffness components.






The effective stresses a are assumed to fulfil the strain equivalence principle (Lemaitre and Chaboche, 1990) resulting, eventually, in Equation (3).
a = Co • є
Inverting the damage tensor D and substituting Equation (1) into Equation (3), renders the stressstrain constitutive law, Equation (4).
a = D—1 • C0 • є = С(ш) • є (4)
where Co is the stiffness matrix. It should be noticed that the introduction of the degradation internal variables шj yields a nonsymmetric tensor C(d) (see Equation (5)). The matrices A and B, defined in a local system of reference, are introduced in order to read C in a more compact manner.
(1 — Cc? ii) (1 —1/23^32 ) (ЛМи)(уі2+У32Уіз) (}^1і)(уіЗ+Щ2^2з) ЕцЕззА Ед Езз А Ед ЕдА –
(1 0)22)ІУі2+v32^13) (1 —CU22)(1—Уїз^зі) (1 ~(л)д) (у23+v21 УЪ)
ЕцЕззД % ЕцЕззД % ЕцЕггД ,
(1аізз)(Уіз+Уі2У2з) (1ІУзз)(У23+1ДіТДз) (1С^зз) (1У12У21) Е11Е22Д Е11Е22Д Е11Е22Д
where A isa determinant that depends upon the elastic properties,
(1 – u12u21 – 1*231*32 – v31 v13 – 2v21v32v13)
E11E22 E33
(1 – ш12)С12 0 0
0 (1 — Ш23) G23 0
0 0 (1 — ш31)^31
Note that the ‘damaged’ stiffness tensor C Є R6x6 is built as follows,






where O Є R3x3 is a matrix filled with zeros.
2. Definition of damage internal variables
The time variation of damage internal variables is defined as a linear combination of the Ф7 growth functions and the damage directors (Equation (6)).
nmodes
da = £ Ф7у7 (6)
7=1
434 
Advances in 
Composite Materials – Ecodesign and Analysis 
In the above, у denotes a mode of damage and nmodes denotes the total number of failure modes. The growth functions for each damage mode у are computed through Equation (7).
Фу =< V£gY, Є >+ (7)
where V£ is the strain gradient ^ and є the strain rate. < • >+ denotes the nonnegative inner product accounting for the trespassing on the damage surface. The subscript + indicates that the inner product vanishes for negative values. This ensures that there is no growth of damage if the damage surface is not reached. If the strain increment vector is pointing to the interior of the surface (for a generic damage mode у) there is no progression of that particular damage mode. So a simple way to effectively computing this is to perform the nonnegative scalar product as represented in equation (7). In Equation (8), gy are the evolving damage surfaces in the strain space.
gy = tT ■ Gy ■ є – cy (8)
where cy is an empirical parameter defining the damage surface. The variations of these surfaces on the strain space result in Equation (9). It should be noticed that cy are not needed in the numerical scheme as Equation (9) is the one necessary for the computational procedure. In this manner, the number of experimental data, which are difficult, or even impossible with the current techniques to obtain, are sensibly reduced.
V£gy = tT ■ (Gyt + Gy) (9)
After some algebra, Gy secondorder tensors are derived from Equation (3) and from the equivalence of the quadratic forms in stress and strain spaces given by Equation (10) (Curiel Sosa et al., 2008a).
aT ■ Fy ■ a = tT ■ Gy ■ є (10)
Ffl are secondorder tensors are derived from damage surfaces defined on the stress space. The modelling of the unitary damage directors vy is based upon the stiffness components that are degraded when a particular mode of damage occurs. For instance, fibre rupture v(l) affects to the stiffness degradation in (11), (12) and (31) directions,
v(1) = [A<;) 0 0 л12) 0 л311) ]T
The weights Ay may be estimated from experimental observations in a qualitative manner for the corresponding damage mode. This technique is still being researched to provide a more straightforward computational strategy that allows to update vy at every time step of the numerical procedure. Techniques such as Inverse Modelling or Optimization are also possible for a more efficient modelling of vy. However, at present, no attempt of using these techniques is being made.
3. PDM algorithm
The PDM algorithm has been implemented into an inhouse FEM. Additionally, it was coded within AbaqusTM as a vumat subroutine. It is adapted for the majority of the commercial software packages based in the explicit FEM. The computation of stresses, performed by numerical integration, includes the constitutive law expressed by the model described (see the algorithm below). A loop over damage modes is performed for the computation of stress
Finite Element Analysis of Progressive Degradation versus Failure Stress Criteria on Composite Damage Mechanics 
435 
at each quadrature point. This is gathered in step (I) below. The computational algorithm is briefly outlined as follows,
I. Loop over damage modes, for 7 = 1 to 7 = nmodes do:
i. Compute the ‘damaged’ stiffness tensor: C(o).
ii. Generate FY. Note that the 7 superscript may be treated as a third index which may provide clarity to the code.
iii. Calculate: GY = CT ■ FY ■ C. Each 7 gives place to a distinct damage surface in the strain space gY, i. e. one for every damage mode (see Equation (8)). Calculation of gY is not required as GY is the only entity needed for the following steps.
iv. Strain gradient of damage in strain space: Vs gY = єт ■ (GYT + GY).
v. Growth of damage 7: Ф7 =< V£ gY, є > +
vi. Directional damage vector: vY
II. For current gauss point, compute the damage internal variables array as a linear
combination of damage directors and damage mode growth: Co = Ynm0des Ф7vY.