JOURNAL OF THEORETICAL

AND APPLIED MECHANICS

42, 3, pp. 539-563, Warsaw 2004

HYBRID, FINITE ELEMENT-ARTIFICIAL NEURAL

NETWORK MODEL FOR COMPOSITE MATERIALS

Marek Lefik

Chair of Geotechnical Engineering and Engineering Structures, Technical University of Łódź

e-mail: emlefik@p.lodz.pl

An application of Artificial Neural Networks for a definition of the ef-
fective constitutive law for a composite is described in the paper. First,
a classical homogenisation procedure is directly interpreted with a use
of this numerical tool. Next, a self-learning Finite Element code (FE
with ANN inside) is used in the case when the effective constitutive
law is deduced from a numerical experiment (substituting here a purely
phenomenological approach). The new contribution to the classical self-
learningprocedure consistsof its adaptationtoacaseof anon-monotonic
loading (non one-to-one load-deformation curve).This new ability of the
method isprincipallydue to the incremental formof the constitutive equ-
ation and the respective scheme of the neural network structure. Also
an organisation of a constitutive data-base containing learning patterns
is suitably modified. It is shown by examples that the training process
is very quick. The error of this method is smaller, comparing to other
schemes of data acquisition.

Keywords:ArtificialNeuralNetworks, compositematerials, self-learning
FEmethod

1. Introduction

The paper is concerned with the presentation of some applications of Ar-
tificial Neural Networks (ANN) in the analysis of composites. Not only a for-
mulation of the constitutive law for a composite but also its representation
inside a finite element code is described in the paper. Themethods presented
here are intended only for composites with a very finemicrostructure.
In the introductory part we explain why the representation of the consti-

tutive law by ANN is useful, and we suggest what kind of problems can be
reasonably treated with it.



540 M.Lefik

The main problem in the analysis of composites is the estimation of the
effective mechanical stiffness, effective thermal conductivity, effective thermal
expansion, effective electrical resistance and so on. Effective – means here –
significant when the structure is considered as a fictious, homogeneous body,
disregarding details of its internal composition. These properties characterise
global (or macro) behaviour, i.e. the behaviour observed at the macrosca-
le. The perturbation of fields of variables of the problem due to the micro-
heterogeneity is insignificant from this perspective. It is obvious that the glo-
bal behaviour is determined at themicrolevel, by geometry of components and
their individual physical properties. Despite this, we want to carry out the es-
timation of global effects without pushing a finite element mesh down, into
the geometry of the microstructure, while the global boundary value problem
under global loading conditions is to be solved. Because of the scale separation
between structural levels a mesh fine enough for the microlevel would result
in a huge number of elements at the macro-structural level, thus, with a huge
number of unknowns in the model. This can be avoided using theory of ho-
mogenisation. The theory of homogenisation (or rather a variety of theories)
can be considered as a fundamental tool of a reasonable analysis of materials
with a microstructure.

Alternatively, a set of experimental observations of the behaviour of a
macro-sample of the heterogeneous material must be performed to define the
effective properties of the composite.Recently thisphenomenological approach
hasbeen frequently substitutedwith a so-called ”numerical experiment” (Ban-
deira et al., 2001; Wriggers and Lohnert, 2002). In the paper we use mainly
such substitution as a source of constitutive data.

We note first that in the frame of the classical homogenisation theory some
regularity of themicrogeometry of the composition is assumed. It is necessary
to simplify the description of a heterogeneous body. There are two principal
idealisations of real structures. Asymptotic homogenisation of periodic media
is applicable to the analysis of structures obtained by translation of an ele-
mentary cell (called the cell of periodicity) while in a self-consistent approach
one assumes that the considered composite is made by multiple, irregular re-
petitions of a scaled pattern. The third possibility is a stochastic description,
which is out of the scope of this paper.

In the case of periodic composites themathematical formalism allows us to
build a qualitative theorywithout any a priori restriction concerning the form
of the resulting effective constitutive law (see, for example Sanchez-Palencia,
1980). However, in the practice of physically non-linear composites the natu-
re of non-linearity must be taken into account at an early stage of the de-



Hybrid, finite element-artificial neural network model... 541

velopment. For example, exponential non-linearity of components requires a
theoretical formulation quite different than that for an elastic-plastic physical
framework (Suquet, 1985). Because of this, many individual, particular the-
ories have been built for particular classes of materials sometimes combined
with particular geometries of the cell of periodicity.

Using asymptotic homogenisation, a solution to the problem of periodic
media typically splits into two boundary value problems. The first one is
formulated at the microlevel, for the repetitive microstructural domain. Its
solution determines effective material characteristics and allows for the con-
struction of stress distribution at the level of microstructure in the phase of
”unsmearing”. The second boundary value problem (the global one) is stated
for the homogeneous, fictious body. The constitutive description of this ma-
terial is obtained at the first stage. The principal difficulties result from the
fact that in the non-linear material or geometrical frame the superposition of
micro- and macro-behaviours is not allowed. Moreover, the coupling between
these two levels is explained by the fact that the solution to themicroBVpro-
blemdependson themicro-geometrywhichchanges during the global solution,
which, in turn, depends on the results of the micro BV problem through the
values of constitutive parameters in the global constitutive law. In the case of
a complicated physical and geometrical composition the mentioned problems
are rather difficult to overcome. This situation is even more difficult from the
point of view of a designer who needs amethod that works well formany geo-
metric andphysical arrangements suitable for various schemes of the structure
that should be considered in the design process. This analysis is simpler if can
be done without revision of the mathematical formulation at each stage of
design and at each level of the structure.

Wriggers and others (Wriggers and Lohnert, 2002; Bandeira et al., 2001)
have made an important step toward such a solution. The proposed method
consists of some numerical experiments using FEmethod at the meso- or mi-
crolevel. The proposed approach is consistentwith fundamental findings of the
homogenisation theory: as in the theory of homogenisation the global behavio-
ur is deduced froma small sample of a compositematerial but now themethod
is fully numerical. It assumes the application of pre-existing numerical proce-
dures in the modelling of contact, friction and other effects at the microlevel
(inter-particularmechanics is considered for granular assemblies for example).
Another important direction in the numerical homogenisation is represented
by Kouznetzova, Brekelmans and co-workers (see Kouznetzova et al., 2001).
Also in this approach an expected form of the resulting, effective constitutive
law should be imposed a priori.



542 M.Lefik

A self-consistent approach can be considered as themost suitable theoreti-
cal frame resulting with a closed-form expression giving the effective material
characteristics as a function of the volume fraction of the forcing material di-
spersed in thematrix. This approach, (see Zaoui, 1997; Kroner, 1978; Boutin,
2000; Hashin, 1968) is less sensitive on the morphology of the sample of the
composite. Boutin, Zaoui andmany other authors have developed expressions
for different effective properties in the case of three separate levels of porosity:
micro, meso andmacro (Boutin et al., 1998; Zaoui, 1997). Unfortunately, also
thismethod assumes some qualitative a priori decisions concerning the global
constitutive law. Finally, it seems to be very difficult to take into account all
possible non-linear, unilateral microfeatures defining ”generic heterogeneity”
for this method. Also the geometry of the scalable cell should be very simple
in this approach (concentric ellipsoids).

Our method of solving the problem of the effective constitutive law for a
materialwithamicrostructure is entirely numerical. It is basedon the assump-
tions that the Cauchy stress (or its increment) in the fictitious, homogenised
body is an unknown, non-linear function of mean strains, strain increment,
physical parameters defining themicrophysics of composition, and the geome-
try of the representative volume. The history of loading-unloading can affect
all this data. We assume a form of this function representing it by an Artifi-
cial Neural Network with hidden layers. At the input of this network we have
all independent variables, mentioned in the above sentence. At the output
we obtain a state of stress in the point of the body at the macrolevel. Since
the components of the strain tensor are present among the input variables,
the Artificial Neural Network acts as a numerical representation of the effec-
tive constitutive law. It is well known that the ANN can be considered as a
universal approximation of a function, functional or operator (a formal proof
of this statement can be found in Chen and Chen (1995)). Because of this,
this representation does not constrain the generality of the global, effective
constitutive law we are looking for. To define the ANN we must determine
all weights of the links between neurones belonging to various layers of the
network. The coefficients in this very rich and flexible representation are defi-
ned during an iterative process of ”training”. The source of knowledge for the
learning can be a real or numerical experiment. If a numerical experiment is
considered, an auxiliary FE modelling process is similar to that proposed by
Wriggers and Lohnert (2002) with the exception that in our case the form of
the constitutive equation is defined by the ANN thus, it does not restrict the
possible effective behaviour. In this case, like in the periodic homogenisation,
we discover rather than impose the global constitutive law, but now – also



Hybrid, finite element-artificial neural network model... 543

for a possibly non-linear effective behaviour. The FE analysis of a sample of
the composite is a source of ”examples” for the training process. These exam-
ples form a knowledge base containing the constitutive data. The considered
sample of the composite that defines the microstructure consists of only one
repetitive cell or a representative volume. In examples presented in this pa-
per the representative volume will be replaced, for simplicity, by an array of
cells. The equivalence between the composite structure and its homogeneous
counterpart is defined in the paper in an original manner, consistent with the
proposed numerical procedure (different from that, described in the frame of
the theory of homogenisation).

The use of the ANN in the role of the constitutive operator is not new.
Classical application of Artificial Neural Networks for constitutive modelling
of concrete was originally proposed by Ghaboussi et al. (1991). An improved
technique of ANN approximation for a problem of mechanical behaviour of
drained and undrained sand is presented in Ghaboussi and Sidarta (1998). In
the latest surveys Waszczyszyn (1997, 1998), Yagawa and Okuda (1996) the
role of neural computing in the constitutive modelling is clearly pointed out.
A similar approach is employed in Garcia et al. (2000), Gawin et al. (2001),
Kortesis and Panagiotopoulos (1993), Mucha (1997), Penumadu and Zhao
(1999) andmany other papers.

The essence of the ANN technique is to construct an application that
attributes a given set of output vectors to a given set of input vectors. When
applied to the constitutive description, the physical nature of these input-
output data is clearly determined by measured quantities: strains-stresses or
displacements-forces. The neural ”black box” operator, replacing the existing
symbolic description, cannot be directly used in the symbolic development of
the formulae and is useless for a closed-form solution. It can be very efficient
as an element of a FE code, as it is shown in Gawin et al. (2001), Mucha
(1997). A hybrid FE-ANN code is described also in Shin and Pande (2000,
2001). The authors show that the insertion of the constitutive law presented
in the formof a neural operator leads to some qualitative improvements in the
application ofFE in engineeringpractice.Namely, theANNrepresentation can
bemodified to reduce the error of a FE numerical experiment with respect to
the real experiment. Our representation of the constitutive lawwith the ANN
is slightly different. It is incremental while in Shin andPande (2000, 2001) the
ε−σ functions are directly approximated.

Obviously, there is no need to use the presented method in the case of a
simple micro-geometry of the composite and when the behaviour of compo-
nents is linear. In this case the classical description of the effective behaviour



544 M.Lefik

is much better. The ANN representation can be useful in a case, when the
physical or geometrical non-linearity dominates themicro-description and the
FE analysis at the level of the representative volume (or a direct experiment)
is the only possible tool for its effective analysis.

The paper contains academic examples illustrating the assumed strategy
of the modelling.We refer the interested reader to our previous papers (Lefik
and Schrefler, 2002a,b) in which some aspects of practical applications of the
presentedmethodto theanalysis ofmechanical behaviourofa superconducting
cable is shown.

2. The effective constitutive law for acomposite and its

representation by an Artificial Neural Network with hidden

layers

2.1. Effective properties of the composite

For usual analysis of a two scale composite we define two sets of co-
ordinates: local y related to the single, repetitive cell of periodicity Y and
the global x for the entire body Ω (composed of a finite number of the cells,
Fig.1).

Fig. 1. The cell of periodicity and two systems of coordinates



Hybrid, finite element-artificial neural network model... 545

Thedimensionless characteristic length ε of the cell Y is treated as a small
parameter tending to zero (the ”length” of the cell Y related to the ”length”
of Ω). Let the displacement field uε(x,y) be an ε-depending solution to the
problem

find uεi ∈V and σ
ε
ij ∈L

2 such that:
(2.1)

∀vi(x)∈V

∫

Ω

σεkl(u)vk,l dΩ=

∫

Sf

Fivi dSf

where the stress tensor σ is computed via (2.2)1 and remains inside the ad-
missible set P (2.2)2 defined by any usual yield criterion (an associated flow
rule completes the constitutive description)

σεij(u
ε)(x,y)= aijkl(y)ekl(u

ε(x,y)) σεij ∈P(y) (2.2)

V is the usual space of kinematically admissible displacements, e(v) is the
linearized strain tensor computed at v, F denotes a vector of external load.
Components of the fourth order elasticity tensor a are piece-wise constant
functions of ywith discontinuities along regular surfaces, and satisfy the clas-
sical conditions of symmetry, ellipticity and positivity.

Let u0, σ0 be the solution to the ”homogenised” problem, i.e. problem
(2.1) in which the variable material coefficients a(y) are replaced with some
unknown constant values ah. We suppose that the periodicity of material
characteristics imposes an analogous periodical perturbation on the quantities
describing themechanical behaviour of the body.Hence, for the displacements
we have

uε(x)≡u0(x)+εu1(x,y) (2.3)

u1 Y -periodic.

The same can be obtained for strains and the stress tensor via a simple
total differential formula with respect to x

eε(x)= e0(x,y)+εe1(x,y)+ ...+0(ε2)
(2.4)

σε(x)=σ0(x,y)+εσ1(x,y)+ ...+0(ε2)

Thanks to (2.4), problem (2.1) splits into a sequence of problems of the order
ε0, ε1 and ε2.
The problemof the order ε0 leads to equation (2.5) inwhich the perturba-

tion of the displacement is computed according to (2.6) (this last expression



546 M.Lefik

is not an assumption, it follows from the detailed analysis of the problem, see
for example Sanchez-Palencia (1980) and references given there)

find χ
pq
i ∈VY such that: (2.5)

∀vi ∈VY

∫

Y

aijkl(y)(δipδjq+χ
pq
i,j(y)
(y)vk,l(y)(y) dΩ=0

u1i(x,y)= epq(x)(u
0(x))χ

pq
i (y)+Ci(x) (2.6)

We call χpq(y) the homogenisation functions for displacements.
For ah uniquely defined via the classical homogenisation procedure uε

converges weakly to u0 as ε tends to zero. In this case, the tensor ah defines
an effective constitutive relationship between e(u0) and σ0 (a new admissible
set in the stress space Phmakes apartof this effective constitutivedescription)

ahijpq = |Y |
−1
∫

Y

aijkl(y)(δkpδlq+χ
pq
k,l(y)
) dY (2.7)

σ0ij(x)= a
h
ijklekl(u

0) (2.8)

The heterogeneous structure can now be studied as a homogeneous one with
effective material coefficients given by (2.7). The global displacements, strains
and average stresses can be computed at that moment.
In this paper, instead of looking for ah and Ph we are going to use a

numerical representation of the effective constitutive law by a suitably trained
ANN N. The input-output vector of N contains components of the strains
e(u0) and stress tensor σ0 (in the vector-like notation)

σ0 =N@e(u0) σi =
∑

q

w
(2
qi)f(q)(epw

(1)
p(q)
+θ
(1)
(q)
)+θ

(2)
i (2.9)

Symbol@denotes a result of an action of an operator at the specified operand.
In this case it is thevalue of theneural operator at the inputvector.Repetition
of the subscript indicates summation over all its range (unless it is enclosed
by paranthesis). Superscripts denote the number of a layer of neurones.
The parameters w, θ and f in (2.9)2 will be interpreted as weights, biases

and transfer functions of the Artificial Neural Network defined and discussed
in the next section. The finite element code with N included as a material
description subroutine can be used to solve numerically problem (2.10) for u0

find u0i ∈V and σ
0
ij ∈L

2 such that:
(2.10)

∀vi(x)∈V

∫

Ωh

σ0kl(u)vk,l dΩ=

∫

Sf

Fivi dSf σ
0 =Nh@e(u0)



Hybrid, finite element-artificial neural network model... 547

The correct, effective constitutive law in (2.10), representedby theANNcalled
Nh assures the following requirement

∀F ‖uε−u0‖L2(Ω) ¬ τ τ ¬ ε‖u
1‖L2(Y ) as ε→ 0 (2.11)

Estimation (2.11)2 follows from a direct application of Cauchy’s inequality.
Expression (2.11) defines a correspondence between the heterogeneous the ho-
mogeneous (homogenised) fictitious body inamanner suitable for ourpurpose.
In the numerical practice, comparison (2.11) of the effective solution with the
exact one will be checked only in a few strategic points. Condition (2.11)
formulated for any F is verifiable in practice only for a finite number of pre-
scribed loads. Searching for the Nh, we employ an iterative procedure called
the ”training of Artificial Neural Network” proposed in all textbooks devoted
to ANNs, as for example Hertz et al. (1991), Osowski (1996), Tadeusiewicz
(1993) and shortly described in the next section.

3. Artificial Neural Network for constitutive description and a

hybrid ANN-FE code

An inspiration for the application of Artificial Neural Networks (ANNs) to
different branches of engineering sciences comes from the analysis of transmis-
sion and transformation of signals in human or animal neural systems. The
importance of themethod significantly increased during last years. Nowadays,
theANNs are successfully used in the computer-aidedmanagement,modelling
of different physical dynamic processes depending onmany fuzzy variables. In
this application, the ANN will be included into the classical Finite Element
procedure as a constitutive subroutine, as it is described in this section.

3.1. Artificial neural network with hidden layer for incremental repre-

sentation of the constitutive law

TheANN can be considered as a universal, non-linear operator that trans-
forms a set of suitably interpreted discrete values into another set of nume-
rical data. Considering the structure of this operator, the ANN appears as a
collection of some simple processing units (called neurones or nodes) that are
mutually connected by linkswith adjustable strengths (seeFig.1). This logical
system, with nodes organized in layers, transforms the input signal presented
at its ”input” nodes into the output signal produced at the ”output layer”.



548 M.Lefik

The sequence of signals presented at the input and that expected to be obta-
ined at the output layer, are called the input and target patterns, respectively.
The activity of the network consists of transformation of each input pattern
into the output signal, and is defined by expression (3.1). After each cycle
(forward transmission of the input signal and back propagation of the output
error), the weights of connection aremodified to reduce the error between the
network response and the given output pattern

oi = f
(

f
(

w(2)qr f(ipw
(1)
pq )+ b

(1)
r

)

w
(3)
ri + b

(2)
i

)

σi = oi (3.1)

In this representation, f is a given sigmoid function of one variable x, Eq.
(3.2),while w and bareweights andbiasesbeingconstant duringtransmission
of the signal

fi(x)=
pi[1− exp(−λix)]

1+exp(−λix)
(3.2)

All independent variables are in i. It can be proved that this form appro-
ximates all continuous functions of several variables. A graph related to this
formula is commonly known as a neural-like network. In Fig.2 a special form
of such a network, useful for our purpose, is illustrated.

In Fig.2 we show a canonical scheme of the Artificial Neural Network for
incremental approximation of the constitutive law.

Fig. 2. A scheme of the ANN that is used throughout the paper. The light arrows
show spontaneous activity of the net along a given path in the space of
deformations. ρ is an internal parameter that allows for many suitable

interpretations, depending on the field of application



Hybrid, finite element-artificial neural network model... 549

For such a operational scheme there exist some well developed methods
for finding the coefficients: weights and biases.

Onecanobserve that the formof (2.9) for theapproximationof the effective
constitutive law is identical with expression (3.1), suitably truncated.

3.2. An incremental form of the constitutive equation

Inour formerpapers (Lefik, 1997, 2001; LefikandSchrefler, 2002, 2003),we
showed that the incremental formof the constitutive equation ismore suitable
for approximation by anArtificial NeuralNetwork.We shortly repeat here the
arguments that support this idea.

We consider only a special form of the incremental constitutive law, na-
mely (3.3)

σ̇= g(σ,F, Ḟ ,ρ) (3.3)

In this representation F denotes the gradient of deformation, g is a function.
It is well known that g in the representation (3.3) can be chosen to fulfill all
conditions defining the admissible subspace P in the space of stresses. In this
way, the solution σ to differential equation (3.3) can be statically admissible.
This formulation, on one hand involves only unknown function g that can be
well approximated by the ANN, on the other hand does not introduce any
”geometric” object in the stress space like a yield surface or other potential
usually defined in the case of non associated plastic flow.These are fundamen-
tal arguments supporting this representation when anANN approximation is
used. The third important argument is following. The approximation of the
incremental form is simpler since the ANN learns only local and short seg-
ments of the stress-strain graph. In the hidden layers some neurons specialize
in switching betweendifferent ”current states of stress”.Thus, thenetwork can
be small and learns quickly. This is shown by examples by Lefik and Schrefler
(2003). Some graphs taken from that paper are presented in Fig.3. In this
Figure we see an approximation of the experimental graph of the mechanical
histeresis of a superconducting cable, described byNijuhuis et al. (1998). The
graph in this Fig.3a shows an autonomous activity of the trained ANN along
the given ”path”, as it is expressed by equation (3.4), below. The approxi-
mation is realistic and the ANN – very small. Such a result is impossible to
obtain by a non-incremental network. In that case a typical result is shown in
Fig.3b.

The incremental approximation of the constitutive law by the ANN can
be defined as follows.



550 M.Lefik

Fig. 3. Experimental graph (dashed line) taken fromNijuhuis et al. (1998),
approximated by the incremental and non-incremental ANN (solid line).

Autonomous generation of the graph by the networks for the points never shown in
the training, according to Lefik and Schrefler (2002b): (a) incremental approach,

(b) non-incremental one

The trainedANN is an operator that performs the the following operation

(σt,F t,ρt,∆F t)→ (σt+1,F t+1,ρt+1) ∆F t = Ḟ
t
∆t (3.4)

3.3. Insertion of the ANN into a Finite Element code

We show next that representations (2.9), (3.1) work well inside a fini-
te element code. Let us suppose that the basic equation for the standard;
displacement-based finite element method can be written in the form

∫

V

B
>

N : τ dV
0 =

∫

S

N
>

Nt dS+

∫

V

N
>

Nf dV

(3.5)

ε(du)=BNdu

where V is the volume of the body in the reference configuration, BN is the
matrix that operating on the vector of admissible variation of independent



Hybrid, finite element-artificial neural network model... 551

variables to give the strainmeasure ε(du) coupledwith thematerial stress τ .
In what follows in the paper, we will consider small transformations, thus an
infinitesimal strain tensor. The index N means that B is constructed on the
basis of the approximation of u by a set of interpolation functions N(x) on
appropriate finite elements. On the right-hand side of (2.8) t is a stress vector
given on the part S of the boundary, while f are body forces acting on the
elementary volume dV .

Since the considered material behaviour is non-linear, the Newton algori-
thm will be applied to solve the system of equations (2.8). The Jacobian of
the left-hand side of (2.8) can be written as follows

J =

∫

V

[dτ : ε(du)+τ : dε(du)] dV 0 (3.6)

The first term in integral (2.9) can be computed using a usual constitutive
assumption

dτ =D : dε (3.7)

where Dij = ∂τi/∂εj.

We can rewrite the above equations, taking variation with respect to the
independent variables of the problem u and obtaining (by definition) the
stiffness matrix K

KmainMN =

∫

V 0

BM :D :BN dV
0 (3.8)

The second term in integral (3.6) represents the initial stress matrix.

Using the assumed representation of the constitutive law by the ANN,
instead of (3.7) we obtain

dτ =Nd,σ@dε (3.9)

The index d denotes that the network quality is best for some given value of
the increment d, σmeans that the stress increment is computed at the current
value of τ =σ. It is clear thatwemust replace the neural operator in (3.1) by
the matrix D or simply, construct this matrix using the given representation
of the constitutive law. This will be done by trial incrementation of ε. Let us
suppose that both tensors dτ and dε are represented by column vectors

dσ= [dτ1,dτ2,dτ1] = [Nd,σ@dε1,Nd,σ@dε2,Nd,σ@dε1]

[dεt] = [dε1,dε2,dε3] (3.10)

D= [dσ][dεt]−1



552 M.Lefik

Thematrix of trial vectors dεt is always proportional to the strains in the last
equilibratedpoint (σ,ε) during theNewton iteration process (preceding step).
The trial vectors cannot be arbitrary because Nd,σ@dε 6=−Nd,σ@(−dε), and
in fact two different tangent stiffnessmatrices can be defined in any point: one
for loading and the other for unloading. It is supposed then that the loading
(unloading) is continued during current increment in the Newton iterations.
Formulae (3.10)1,2 are used here instead of computing the derivatives of the
neural network with respect to input values.

The stress in the second term in integral (3.6) is computedusing the neural
network in the recall mode for a given constant step dε, until the strain ε
at the trial solution in the current step is reached. The ANN acts here in
the autonomous activity mode as it is defined in Section 3. This process cor-
responds to classical integration of the incremental constitutive equation for
updating σ. It always starts in the last equilibrated point, and the increment
dε is proportional to the one defined for this step (loading or unloading).

4. Classical homogenisation procedure directly interpreted by the

Artificial Neural Network

It can be shown that the vector of homogenisation functions, Eq. (2.5),
appears as a perturbationfield in the numerical experiment inwhich the cell of
periodicity is loadedwith a uniformdisplacement of points on theborder.This
imposes a given average strain state equal to the homogeneous one, computed
from formula (2.9) in which u is to be applied on the border

ui = e
ave
ij xj (4.1)

The corresponding stress can be obtained according to expression (4.2) with
t used for the stress vector on the boundary of the cell

σaveij =

∫

∂Y

xi⊗ tj ds (4.2)

Since Eq. (2.9) is a universal approximator, the introduced assumption does
not restrict the possible effective behaviour. Like in the periodic homogeni-
sation, we discover rather than impose the global constitutive law. All the
constitutive information is in the data set furnished to the network.



Hybrid, finite element-artificial neural network model... 553

4.1. Constitutive data-base for the training

Aconstitutive data-base contains pairs of input vectors andoutput vectors
for the training of the network. The input and output is accorded with the
structure of the ANN (Fig.2) and is always of the form:

{(input vector), (output vector)}

For the case of the incremental law, for the k-th pattern, we have

{[εi,σ(εi),ηi,∆εi],∆σαβi}k (4.3)

Additionally, in the below shown example, we must satisfy the following
condition imposedby isotropyof the approximated constitutive law (@denotes
the action of the neural network operator on the input data)

∀i if N@



















εi
σi
η
∆εi



















= {∆σi} then ∀Θ Θ
>
Θ=1 we have

(4.4)

N@



















Θ>εiΘ

Θ>σiΘ

η

Θ>∆εiΘ



















= {Θ>∆σiΘ}

To satisfy condition (4.4), wemust train thenetworkwith some supplementary
data: the new subset of patterns is of the form

{

[Θ>εiΘ,Θ
>σiΘ,ηi,Θ

>∆εiΘ],Θ
>∆σαβiΘ

}

(4.5)

The subscript k refers to the pattern obtained from the ith experimental
pointby transformationvia the kth rotationmatrixΘ.The total number K of
these additional terms in thematrix of patternsdependson thenumberof trial
rotation needed to train the network up to a satisfactory level of tolerance. An
artificial construction of the experimental data has been necessary to perform
a 2D numerical experiment.
In this case, the source of constitutive data is a numerical experiment de-

fined with expressions (4.1) and (4.2) executed on a single cell of periodicity
of the composite. The results of FE computations for different (assumed and
imposed over the cell) values of a small deformation tensor or tensor of defor-
mation gradient must be completed by the manipulation prescribed by (4.4).



554 M.Lefik

4.2. Example

As an example, we consider a hyperelastic composite, the repetitive cell
of which consists of a neohookean material with a circular hole of the radius
equal to the distance between the neighbouring, circular pore.

In Fig.4 this situation is shown. In Fig.5 some of deformed configurations
are presented. The number of such numerical experimentswas 1260, including
incrementation.The full constitutive data base containing all rotated datawas
much bigger. As far as the increments with the negative sign were concerned
it was about 15000 patterns. The training was organized by epochs that con-
tained about 1250 patterns. After the first decisive reduction of the error of
the training, most of the patterns used as the test data revealed very good
coherence with the training results.

Fig. 4. A single cell of periodicity of the composite

5. The self-learning Finite Element code for deducing effective

constitutive relationships from a numerical experiment

The self-learning finite element procedure is based on the standard FE
code including the Newton-Raphson type iterative method. The neural con-
stitutivemodel, Eq. (3.1), can be applied for the solution to general structural
problems in the samemanner as the conventional one.All details can be found



Hybrid, finite element-artificial neural network model... 555

Fig. 5. Examples of trial deformations of the cell of periodicity – a source of the
constitutive data

in Shin and Pande (2000, 2001). The main advantage of this substitution is
the possibility tomodify the constitutive model simply by retraining the neu-
ral network. This feature results in a self-learning code, which evaluates the
suppliedANNmodel andadjusts itself to a desirable response of the structure.
This strategy permits the constitutive relationship to be built directly from
experimental data or, alternatively, from the exact FE solution. In our case,
it is solution of (2.1) for a structure composed of many repetitive cells of the
given internal structure (may be very complicated). If the considered number
of cells is large enoughwe are able to capture its effective (global) properties.



556 M.Lefik

The procedure starts with calculation of the trial FE solution using any
initial constitutive law (e.g. linear constitutive matrix). Pairs (eI,σI) at each
Gauss point on any considered load level are saved for future training of the
ANN. The displacements at some strategic points of the structure are moni-
tored δI and compared with the known data describing the behaviour of the
compositematerial δD. Those vectors constitute, respectively, the output and
the target of the whole procedure.
The current adjustment ismeasuredbymeans of thenormalized root-mean

square error calculated as follows (instead of (2.6)1)

Err=
1

δmax
D
− δmin
D

√

√

√

√

1

MN

N
∑

n=1

M
∑

m=1

(δm
II,n
)2 <τ (5.1)

In (5.1), N is the number of load increments, M – number of monitored
displacements, the vector δII = δI − δD represents the errors made by the
ANNat themonitored points and δmaxD , δ

min
D are themaximumandminimum

values of |δD|, respectively.
In the following, the displacement differences δII are treated as the load

for the new solution. The obtained strain-stress pairs (eII,σII) are used to
correct the current guess (eI,σI)

etrain = eI σ
train =σI +σII (5.2)

Suchpreparedtrainingdata (5.2) areused for retraining theneuralnetwork
that forms the constitutivemodel. The newweights of theANNare saved and
the next step of the self-learning procedure starts with the corrected ANN
constitutivemodel. If the value of E becomes lower than the tolerance level τ
? or the admissible number of self-learning steps is reached, the process is
finished.According to Shin andPande (2000, 2001), this procedure converges.

5.1. Example

Weconsider the same example as inSection 4.3.This time, however, a sam-
ple containing several cells is considered. The non-deformed cell of periodicity
formsa squarewith ahole of diameter equal to half of the edge. For acquisition
of the constitutive data we perform three different numerical FE experiments
(because of (2.6)). All these experiments are qualitatively different than those
proposedbyShinandPande.We imposeuniformdisplacements of borders of a
portion of the composite. The nature of these kinematic loadings is easily seen
inFig.8: two different tests of tension (free andconstrained tension) anda she-
ar test. The controlled quantities are now values of reactions in the boundary



Hybrid, finite element-artificial neural network model... 557

nodes. A test like this was never proposed by the authors of the self-learning
method since its original application was acquisition of the constitutive data
from real and in situ observations. It is obvious that the point-wise reaction
is not an observable quantity. In the context of the numerical experiment this
strategyworks quitewell, even in the analysis of compositematerials forwhich
the stress jumps in the small scale of the body. In the approximation of the
constitutive law an incremental constitutive description and the corresponding
architecture of the ANN has been used (described by Lefik (2001)).

Fig. 6. Three numerical tests on a sample of the composite used simultaneously for
the preparation of the data-base for ANN training

5.2. Discussion of the results

In Fig.7 two testing finite element computations are shown. These com-
putations have been executed for a loading never used in the training and
with the FE code including the ANN inside. A relatively small ANN with 10
and 7 neurones in hidden layers have been trained with the constitutive data
collected from each first Gauss point of the homogeneous problem (the black,
rectangularmesh in Fig.7). The performance of the trainedANNwas checked
with the loading scheme never used in the training. This scheme is shown in
Fig.7. The coarse, rectangular mesh (dark line) represents the displacements
of the ”effective”, homogeneous body under action of the horizontal shearing



558 M.Lefik

Fig. 7. Two numerical tests on a sample of the composite that have never been used
in the training. The coarsmesh – results obtained for the effective model, the light

deformedmesh – ”exact” solution



Hybrid, finite element-artificial neural network model... 559

stress vector applied to the upper edge of the rectangular sample. This compu-
tationwas carried out byFE-ANN code (the FE code publishedbyBonet and
Wood (1997) wasmodified and used). The triangularmesh (light lines) repre-
sents the ”true” deformed configuration. The coarse mesh is defined for the
homogenised body. The qualitative accordance is good, themaximum error in
displacements does not exceed 5%.

In Fig.8, the difference of the return curve and the one obtained with the
loading is used as the measure of the error. It is known that for hyperela-
stic composites these two curves should coincide. What is shown in Fig.8 is
obtained for the example trained within the direct approach. In Fig.8b the
curves coincide much better that in Fig.8a. The self-learning method of data
acquisition is more efficient. In this case the ANN is simply ”tailored” for the
Finite Element code in which will be included as a constitutive procedure.

Fig. 8. Comparison of loading and unloading stress-strain curves used as themeasure
of the error of the ANN representation of the constitutive law inside the FE code

In the example, only three cases of load are reported. As far as the co-
nvergence to zero is concerned, see Bandeira et al. (2001), a relatively small
trial portion of the composite is sufficient to predict the overall behaviour of
the non-homogeneous material. In our example we consider only 20 cells of
periodicity but even in this case we notice that the perturbation of the mean
displacements is negligible.

6. Conclusions

We conclude that the effective constitutive law for a composite can by
approximated by the Artificial Neural Network with hidden layers. This ap-
proximation does not constrain the form of unknown effective relationships.



560 M.Lefik

The following condition must be fulfilled to assure successful approximation:

• The constitutive law must be formulated in the incremental form

• The ANN must be included as a subroutine in the Finite Element pro-
cedure.

Moreover:

• The self-learning procedure seems to be a very effective tool for the
identification of the global behaviour of composites. Themicrostructure
of thematerial (position of themonitoredpoint on the cell of periodicity)
does not influence the identification abilities of the program.

• The ANN representation of any constitutive law is a flexible tool for
the representation of the effective (global) behaviour of materials with a
complex internal microstructure.

• This representation is very suitable for analysis of composites since it is
”automatic” in the sense that it does not require any a priori choice or
adaptation of the existing constitutive theory for the description of the
observed material behaviour.

• The convergence is surprisingly fast. Four steps are enough to obtain a
qualitatively goodmodel.

• A representation of the effective constitutive law can be very simple (a
network of the architecture 3-6-3 is sufficient in the first example).

• The usefulness of the hybrid FE-ANN code has been confirmed since it
opens up new possibilities in comparison with the standard FE codes
in the sense that the constitutive models can be easily modified and its
simplicity accelerates work of the program.

These conclusions have been supported by examples, not by a theoretical
proof.

References

1. Bandeira A.A., Wriggers P., Pimenta P.M., 2001, Homogenizationme-
thod leading to interface laws of contactmechanics – a finite element approach
for large 3D deformation using ALM, IJMN

2. Bonet J., Wood R.D., 1997, Nonlinear Continuum Mechanics for Finite
Element Analysis, Cambridge University Press

3. Boutin C., 2000, Study of permeability by periodic and self-consistent homo-
genization,Eur. J. Mech., A, 19, 603-632



Hybrid, finite element-artificial neural network model... 561

4. Boutin C., Royer P., Auriault J.L., 1998, Acoustic absorption of porous
surfacing with double porosity, Int. J. Solid Structures, 35, (34-35), 4709-4737

5. Chen T., Chen H., 1995, Universal approximation to non-linear operators
by neural networks with arbitrary activation functions and its application to
dynamical systems, IEEE Trans. on Neural Networks, 6, 4, 911-917

6. Garcia S., Romo M.P., Taboada-Urtuzuastegui V., 2000, Knowledge-
basedmodellingof sandbehaviour,Proceedings of ECCOMAS2000, Barcelona,
11-14

7. GawinD., LefikM., SchreflerB.A., 2001,ANNapproach to sorption hy-
steresiswithin a coupled hygro- thermo-mechanicalFEanalysis, Int. J. Numer.
Meth. Engng., 50, 299-323

8. Ghaboussi J., Garrett J.H., Wu X., 1991, Knowledge-basedmodelling of
material behaviour with neural networks, Journal of Engineering Mechanics,
117, 132-151

9. Ghaboussi J., Sidarta D.E., 1998, New nested adaptive neural networks
(NANN) for constitutive modelling,Computers and Geotechnics, 22, 1, 29-52

10. Hashin Z., 1968,Assessment of self-consistent scheme approximation: conduc-
tivity of particulate composites, J. Comp. Mater., 2, 284-304

11. Hertz J., Krogh A., Palmer G.R., 1991, Introduction to the Theory of
Neural Computation, Lecture Notes Vol. I, Santa Fe Institute Studies in the
Sciences of Complexity, Addison-Wesley,

12. Hu Y.H., Hwang J-N., Eds., 2002,Handbook of Neural Network Signal Pro-
cessing, CRCPRESS

13. Kortesis S., Panagiotopoulos P.D., 1993,Neural networks for computing
in structural analysis: Methods and prospects of applications, Int. J. Num.
Meth. Eng., 36, 2305-2318

14. Kouznetsova V., Brekelmans W., Baaijens F., 2001, An approach to
micro-macromedelling of heterogeneousmaterials,Comp. Mech., 27, 37-48

15. Kroner E., 1978, Self-consistent scheme and graded disorder in poly crystal
elasticity, J. Phys., 8, 22-61

16. Lefik M., 1997,Use of artificial neural network to define a non-linear effective
constitutive law for a composite,Proc. 13th Polish Conf. ComputerMethods in
Mechanics PCCMM’97, Poznań, 725-732

17. Lefik M., 2001,Modified BP artificial neural network as an incremental non-
linear constitutivemodel,Proceedings of European Conference on Computatio-
nal Mechanics, ECCM-2001, on CD

18. Lefik M., Schrefler B.A., 2002a, Artificial neural network for parameter
identifications foran elasto-plasticmodel of super-conductingcableunder cyclic
loading,Computers and Structures, 80, 22, 1699-1713



562 M.Lefik

19. Lefik M., Schrefler B.A., 2002b, One-dimensional model of cable-in-
conduit superconductors under cyclic loading using artificial neural networks,
Fusion Engineering and Design, 60, 2, 105-117

20. LefikM., SchreflerB.A., 2003,Artificial neural network as an incremental
non-linear constitutive model for a finite element code, Computer Methods in
Applied Mechanics and Engineering, Elsevier, Volume/Issue 192/28-30, 3265-
3283

21. Mucha G., Waszczyszyn Z., 1997, Hybrid neural-network/computational
program for bending analysis of elastoplastic beams, Proc. of the XIII Polish
Conf. on Computer Methods in Mechanics, 949-956

22. Nijuhuis A., Noordman N.H.W., Ten Kate H.H.J., 1998, Mechanical
andElectrical testing of an ITERCS1Model Coil Conductor underTransverse
Loading in a Cryogenic Press, PreliminaryReport, University of Twente

23. Osowski S., 1996, Sieci neuronowe w ujęciu algorytmicznym, Wydawnictwo
NaukowoTechniczne,Warszawa

24. Penumadu D., Zhao R., 1999, Triaxial compression behaviour of sand and
gravel using artificial neural networks (ANN),Computers and Geotechnics, 24,
207-230

25. Pichler B.,MangH., 2001,Parameter identification for sophisticatedmate-
rialmodels bymeans of iterativebackanalyseson soft computing,ECCM-2001,
Cracow, Poland

26. Shin H.S., Pande G.N., 2000, On self-learning finite element codes based on
monitored response of structues,Computers and Geotechnics, 27, 161-178

27. Shin H.S., Pande G.N., 2001, ”Intelligent” Finite Elements, in:Computatio-
nalMechanics –NewFrontiers forNewMillenium, S.ValliappanandN.Khalili
(Edit.), Elsevier Science

28. Sanchez-PalenciaE., 1980,Non-HomogeneousMedia andVibrationTheory,
Springer V.

29. Sikora Z., Ossowski R., Ichikawa Y., Tkacz K., 1998, Neural networks
as a tool for constitutive modelling, in: Oka and Yashima, Edit., Localization
and Bifurcation Theory for Soils and Rocks, Balkema, Rotterdam

30. Suquet P.M., 1985, Elements of Homogenisation for Inelastic Solid Mecha-
nics, in: Lecture Notes in Physics, 272, Springer Verlag

31. Tadeusiewicz R., 1993, Sieci Neuronowe, Akademicka OficynaWydawnicza

32. Waszczyszyn Z., 1998, Some new results in applications of backpropagation
neural networks in structural and civil engineering, in:Advances in Engineering
Computational Technology, Civil-Comp Press, Edinburgh, 173-187



Hybrid, finite element-artificial neural network model... 563

33. Waszczyszyn Z., 2000, Neural networks in plasticity: some new results and
prospects of applications, European Congress on Computational Methods in
Applied Sciences and Engineering ECCOMAS 2000, on CD

34. Wriggers P., Lohnert S., 2002, Virtual testing of composite materials,
WCCM V, Vienna, Austria

35. Yagawa G., Okuda H., 1996, Neural networks in computational mechanics,
Archives of Computational Methods in Engineering, 3, 4, 435-512

36. Zaoui A., 1997, Structural morphology and constitutive behaviour of hetero-
geneousmaterials, in:ContinuumMicromechanics, P. Suquet (Edit.), Springer,
291-347

Zastosowanie sztucznych sieci neuronowych w modelowaniu

numerycznym kompozytów przy pomocy metody elementów skończonych

Streszczenie

Wartykuleopisanozastosowaniesztucznych sieci neuronowychdookreśleniaefek-
tywnego związku konstytutywnego dla kompozytów. To narzędzie numeryczne użyte
zostało dwojako: do bezpośredniego zapisu wyników otrzymanychw ramach klasycz-
nejmetody homogenizacji oraz downioskowania owłasnościach efektywnych na pod-
stawie eksperymentu numerycznego (zastępującego eksperyment rzeczywisty) wyko-
nanegonamałej, lecz reprezentatywnej próbce kompozytu.Wtymdrugimprzypadku
zastosowano schemat „samouczącego się” programumetody elementów skończonych,
wktórymzwiązek konstytutywnyopisany jest siecią neuronową.Schemat ten zaadap-
towano tak, żemoże być użytywprzypadku obciążeń niemonotonicznych orazwtedy,
gdy zależność: miara odkształcenia–miara naprężenia nie jest wzajemnie jednoznacz-
na. Te nowe możliwości uzyskane zostały dzięki przedstawieniu związku konstytu-
tywnego w formie przyrostowej oraz opracowania odpowiedniej do tego budowy sieci
neuronowej. Schemat „samouczącego się” programuMES charakteryzuje się tym, że
proces formułowania nieznanego związku konstytutywnego jest szybki, a zgodność
modelu numerycznego z eksperymentemwiększa niż dla innychmetod.

Manuscript received April 13, 2004; accepted for print April 26, 2004