European Conference on Spacecraft Structures, Materials & Environmental
Testing, 1-4 April 2014, Braunschweig, Germany
[Year]
DAMAGE MODELING OF COMPOSITES AND RELIABILITY ANALYSIS
Michael Bruyneel (1,2), Samih Zein (1)
(1)
SAMTECH (A Siemens Industry Software Company), Liège Science Park, Belgium, Michael.bruyneel@lmsintl.com
(2)
Aerospace & Mechanical Engineering Department, University of Liège, Belgium:
ABSTRACT
In this paper, the intra-laminar damage of laminated
composites is modelled and analysed with the SAMCEF
finite element code. A continuum damage mechanics
model is used, and its parameters are identified based on
physical tests conducted at the coupon level.
Considering the variability from the tests and the
sensitivity of the computed mechanical response
depending on the values of the model parameters, a
reliability analysis is considered to evaluate the most
sensitive parameters in the model.
1. INTRODUCTION
In order to propose predictive simulation tools, it is
important to use material models able to represent the
different modes of degradation of the plies forming the
laminated composite structure. Sometimes, damage at
the interface between the plies, that is delamination,
must also be taken into account in the model. In this
paper, laminated composites made of unidirectional
plies are considered. Damage inside the plies alone, that
is intra-laminar damage, is considered in this paper.
The intra-laminar damage model used here is based on
the continuum damage mechanics and was initially
developed by the Ladevèze’s team in Cachan [1].
Damage variables impacting the stiffness of the ply are
associated to the different failure modes, representing
the fibre breaking, the matrix cracking and the decohesion between fibres and matrix. Plasticity in the
matrix is also taken into account. Such an advanced
damage model includes lots of parameters, which must
be identified based on test results at the coupon level.
Test results present some variability, even when they
are conducted on coupons coming from the same plate,
because of the presence of small defects arising from
the manufacturing process. It is therefore important to
take into account such dispersions of the material
properties and the influence they may have on the
mechanical response of the composite.
In this paper, problems at the coupon level are
addressed for laminates made up of unidirectional plies.
The intra-laminar damage model available in the
SAMCEF finite element code is first presented. The set
of parameters is provided and the basics of the
parameter identification procedure of such material
models are briefly explained. The damage model is then
used in a reliability analysis based on the polynomial
chaos expansion approach [2,3], in order to determine
the most sensitive parameters as well as the probability
of failure depending on the parameters variability. The
polynomial chaos expansion (PCE) is an efficient
numerical method for performing a reliability analysis.
It relates the output of a nonlinear system with the
uncertainty in its input parameters using a
multidimensional polynomial approximation (the socalled PCE). Numerically, such an approximation can
be obtained by using a regression method with a
suitable design of experiments. The cost of this
approximation depends on the size of the design of
experiments. If the design of experiments is large and
the system is modelled with a computationally
expensive FEA model, the PCE approximation becomes
infeasible. In papers [2,3], an algorithm is proposed to
generate efficiently a design of experiments of a size
defined by the user, in order to make the PCE
approximation feasible in computational time. It is an
optimization algorithm which seeks to find the best
design of experiments in the D-optimal sense for the
PCE. This algorithm is a coupling between genetic
algorithms and the Fedorov exchange algorithm.
2. DAMAGE MODEL OF UD PLIES FOR
LAMINATES
The ply damage model relies on the development
proposed in Ladeveze and Le Dantec [1]. For intralaminar damage, the following potential (1) with
damage, named ed, is used, where d11, d22 and d12 are the
damages related to the fibres, the transverse and the
shear directions, respectively.
ed =
+
2
σ 11
2(1 − d11 ) E10
σ 22
2
+
2(1 − d 22 ) E20
+
−
0
ν 12
E10
σ 22
2 E20
σ 11σ 22
2
−
+
2
σ 12
(1)
0
2(1 − d12 )G12
These damage variables allow considering damage
associated to the fibre direction, cracks in the transverse
direction and de-cohesion between fibres and matrix.
The thermodynamic forces represent the effect of the
loading in the corresponding mode. These
thermodynamic forces are derived from the potential
and manage the evolution of the damages via relations
of the form d11 = g11 (Y11), d22 = g22 (Y12,Y22) and d12 =
g12 (Y12,Y22). For instance, the thermodynamic force
associated to shear is given in (2). Finally the model can
be coupled to plasticity with isotropic hardening. The
yield criterion is defined as a function of the effective
European Conference on Spacecraft Structures, Materials & Environmental
Testing, 1-4 April 2014, Braunschweig, Germany
stresses impacting the matrix behaviour (3).
Y12 =
2
σ12
0
2(1 − d12 )2 G12
2
2
f (σ~, p ) = σ~12
+ a 2σ~22
− R0 − R ( p ) ≤ 0
[Year]
3. DAMAGE MODEL PARAMETERS
IDENTIFICATION
(2)
(3)
σ~ij = σ ij /(1 − d ij )
R ( p ) = Kpγ
In (3), R0 is the initial yield stress. The non-linear
behaviours taken into account in this model are
illustrated in Figure 1: non linearity in the fibre
direction, non-linearity including plasticity in the
matrix.
The value of some parameters must be provided in order
to feed the damage model. These parameters are the
elastic properties (E01, E02, G012, ν12), the coefficients of
the plastic law (R0, K, γ, a), as well as the parameters
associated to the damage law (e.g. b, Y012, YC12, YS12;
other parameters are also used for the behaviour along
the fibre direction). When only traction is considered
(and not compression), the parameter identification
procedure is based on physical tests conducted on
coupons made up of 3 different stacking sequences. One
of these sequences is [±45]ns. Quasi-static cyclic loading
is conducted (Figure 1): the initial stiffness is
determined, as well as the evolution of the damage
variable linked to the decrease in the stiffness during the
cyclic loading. The procedure is explained in [1,4].
From the test results, it is clear that some variability
exists in the mechanical response of the laminates. This
dispersion can be reproduced by varying the value of
the material parameters in the finite element model.
4. INFLUENCE OF THE PARAMETERS ON
THE DAMAGE MODEL RESPONSE
The parameters of the damage model are determined at
the coupon level (Figure 2), based on test results. It is
clear that some variability appears in the tests results,
and that most of the time some “mean” values are used
to feed the model.
Figure 1. Non-linear behaviour of the damage model
In this paper, the evolution of the shear damage variable
d12 as a function of the thermodynamic force is given by
the following equation:
d12 =
0
Y − Y12
C
Y12
−
0
Y12
if
S
Y ≤ Y12
d12 = 1 otherwise
The equivalent thermodynamic force Y makes a link
between the transverse and shear effects, via:
Y = Y12 + bY22
where b is a coupling coefficient.
Figure 2. Model of the coupon used for the parameters
identification
In Figure 3, the mechanical response of a coupon made
up of a [±45]2s laminate is illustrated. It is clear that the
global mechanical behavior is non-linear, that
permanent deformation is observed after an unloading,
together with damage associated to the decrease in the
material stiffness. The reference values (from one test)
are superposed to the results of the simulation, for the
identified parameters of the damage model. A very good
agreement is observed. It is clear that a variability exists
at the physical test level, even when the tested coupons
are machined from the same composite plate. This is
also concluded when the values of the identified
parameters are modified, as illustrated in Figure 4. In
that case, the simulation and test results are no more in
good agreement.
European Conference on Spacecraft Structures, Materials & Environmental
Testing, 1-4 April 2014, Braunschweig, Germany
Figure 3. Mechanical response of the [±45]2s
[Year]
Figure 6. Influence of the damage evolution law on the
mechanical response
Figure 7. Influence of the initial shear stiffness on the
mechanical response
Figure 4. Mechanical response of the [±45] 2s for
modified values of the identified model parameters
In Figure 5, the (reference) values of the parameters are
given by: G012 = 4500 MPa, K = 400 MPa, YC12 = 6 MPa
and θ = 45° in the [±θ]2s laminate. In Figures 6 to 9, the
values of these parameters are changed a little bit, and
their influence on the mechanical response of the
laminate can be estimated. The evolutions of d12 and of
the hardening law are given in Figures 10 and 11,
respectively.
Figure 5. Influence of the parameters variation on the
non-linear mechanical response of the coupon
Figure 8. Influence of the hardening law on the
mechanical response
Figure 9. Influence of the fibre orientation on the
mechanical response
European Conference on Spacecraft Structures, Materials & Environmental
Testing, 1-4 April 2014, Braunschweig, Germany
[Year]
approximated using the Polynomial Chaos Expansion as
follows (PCE) up to degree p:
~
S (ζ ) =
n
ai Ψi (ζ )
i =0
Figure 10. Damage law for different values of the YC12
parameter
~
where S is the approximation of S about its nominal
value, Ψi are the multivariate Hermite polynomial,
ai are the PCE coefficients and ζ is the vector of normal
Gaussian random vairables of . The number n of terms
in the PCE is: n=(N+p)!/N!p! where N is the number of
uncertain parameters. The PCE coefficients are
computed using the proposed method in [2]. A Doptimal design of experiments is generated and the
regression method is used to compute the PCE
coefficients. According to the orthogonal property of the
PCE, the mean value and variance of S are
approximated by
n
ai2
a0 and
i =1
Figure 11. Hardening law for different values of the K
parameter
The following questions then arise: 1/ how can we
quantify the validity of the identified model parameters
knowing that some variability exists at the tests level?;
2/what is the most sensitive parameter of the model,
which should be identified in a very accurate way? 3/
how can we be sure that the predicted strength is in an
acceptable range?
5. RELIABILITY ANALYSIS USING THE
POLYNOMIAL CHAOS EXPANSION
Some variability exists in the mechanical properties of
composites. The design paramerter xi is defined as an
independent gaussian random variable defined as
follows:
xi = xi + σ iξi ,
where xi is its mean value, σi is its standard deviation
and ξi a normal Gaussian variable. Note that xi can be
any of the mechanical parameters cited above. This
uncertainty on the parameters is propagated through the
mechanical system and its response, for example the
stress, also becomes a random variable. The reliability
analysis consists in computing the mean value and the
variance of the response of the mechanical system. It
also consists in computing the probability that the
response does not exceed a certain threshold. This
probability is referred to as probability of failure. Let S
be the response of the mechanical system. S is
respectively. One can also easily deduce the Sobol’s
indices from these coefficients (see [3]). The variance
of S is the sum of variances due to the perturbation of
each parameter and the joint perturbations of
parameters. The Sobol’s index is the part the
perturbation of a parameter or a set of parameters with
respect to the total variance. The probability of failure
defined by Pr[S>q], the probability that S exceed a
threshold q is computed by using the PCE
approximation of S and the Monte Carlo Method. Mt
samples of ζ are drawn randomly and the corresponding
~
values S (ζ ) are computed. If Mq samples correspond
~
to S (ζ ) > q ,
Pr[ S > q ] ≈
Mq
Mt
.
6. APPLICATIONS
6.1 [±45]2s laminate
The problem is run with the following values of the
parameters:
Table 1. The parameters and their variability
Parameter
G0
K
θ
YC
Mean value
5500 MPa
400 MPa
45°
60 MPa
Standard deviation
150 MPa
15 MPa
1.5°
0.33 MPa
The ratios between the mean value and the standard
European Conference on Spacecraft Structures, Materials & Environmental
Testing, 1-4 April 2014, Braunschweig, Germany
deviation are quite similar, so there is no unbalanced
variable uncertainty in the problem, which would bias
the sensitivity analysis using the Sobol’s indices. A
DOE of 30 experiments is conducted, in order to
initialize the population. Two functions are considered
in the problem: the first one concerns the value of the
shear strain in the ply, calculated for an applied force
equal to 2000N, and the second one is related to the
maximum force the composite can sustain. The relative
sensitivities of these functions with respect to the
parameter uncertainties are expressed in terms of
Sobol’s indices. From Table 2, it is concluded that the
angle deviation is the most sensitive parameter when
both functions are considered.
Table 2. Results of the first application
Sensitive
variables
G012 alone
G012 and K
G012 and θ
G012 and Y
K alone
K and θ
K and Y
θ alone
θ and Y
Y
Sum of all the
contributions
Sobol’s index of
ε12
0.14533
7.86508e-06
0.000954293
5.37728e-05
0.060854
0.00139441
3.43177e-05
0.744552
0.00100739
0.0458122
Sobol’s index of
Fmax
0.107495
0.000291474
0.000405614
0.00183905
0.000871712
7.77772e-05
0.000610869
0.421117
0.00250115
0.464791
1
1
6.2 [±67.5]2s laminate
The methodology is applied on a plate made up of a
[±67.5]2s laminate. The same parameters are used, with
the same variability as in Table 1. In Figure 12, the
mechanical responses of the plate are plotted, for the
reference values and for some perturbations. It is clear
that the variability in the fiber orientation is here the
most important parameter. This result is confirmed by
the robust analysis, as reported in Table 3.
Figure 12. Mechanical response (Force/strain) of the
plate
[Year]
Table 3. Results of the second application
Sensitive
variables
G012 alone
G012 and K
G012 and θ
G012 and Y
K alone
K and θ
K and Y
θ alone
θ and Y
Y
Sum of all the
contributions
Sobol’s index of
Fmax
0.0549
5.37e-5
0.004
1.48e-7
0.00139
3.29e-5
8.19e-5
0.923045
0.00081
0.01504
1
7. CONCLUSIONS
In this paper, the intra-laminar damage model for
laminated composites available in SAMCEF was
presented. It is a continuum damage mechanics model,
whose parameters are identified based on physical tests
conducted at the coupon level. Considering the
variability from the tests and the sensitivity of the
computed mechanical response depending on the values
of the model parameters, a reliability analysis was
considered to evaluate the most sensitive parameters in
the model.
8. REFERENCES
1. Ladeveze, P. & Le Dantec, S. (1992). Damage
modelling of the elementary ply for laminated
composites, Composites Science and Technology,
43, 257-267.
2. Zein, S., Colson, B. & Glineur, F. (2013). An
Efficient Sampling Method for Regression Based
Polynomial Chaos Expansion, Communications in
computational physics, 13(4), 1173-1188.
3. B. Sudret (2013). Global sensitivity analysis using
polynomial chaos expansions, Reliab. Eng. Syst.
Safe., 93, pp. 964-979.
4. Bruyneel, M., Delsemme, J.P., Goupil, A.C.,
Jetteur, P., Lequesne, C., Naito, T., Urushiyama, Y.
(2014). Damage modelling of laminated
composites: validation of the intra-laminar law of
SAMCEF at the coupon level for UD plies.
ECCM16, Seville, Spain, 22-26 June, 2014.