DOI:10.14311/APP.2021.30.0047
Acta Polytechnica CTU Proceedings 30:47–52, 2021 © Czech Technical University in Prague, 2021

available online at http://ojs.cvut.cz/ojs/index.php/app

LOCALIZATION ANALYSIS OF DAMAGE FOR
ONE-DIMENSIONAL PERIDYNAMIC MODEL

Karel Mikeša, ∗, Milan Jiráseka, Jan Zemana, Ondřej Rokošb,
Ron H. J. Peerlingsb

a Czech Technical University in Prague, Faculty of Civil Engineering, Department of Mechanics, Thákurova 7,
166 29 Prague, Czech Republic

b Eindhoven University of Technology, Department of Mechanical Engineering, PO Box 513, 5600 MB
Eindhoven, The Netherlands

∗ corresponding author: Karel.Mikes@cvut.cz

Abstract. Peridynamics is a recently developed extension of continuum mechanics, which replaces
the traditional concept of stress by force interactions between material points at a finite distance. The
peridynamic continuum is thus intrinsically nonlocal. In this contribution, a bond-based peridynamic
model with elastic-brittle interactions is considered and the critical strain is defined for each bond
as a function of its length. Various forms of length functions are employed to achieve a variety of
macroscopic responses. A detailed study of three different localization mechanisms is performed for
a one-dimensional periodic unit cell. Furthermore, a convergence study of the adopted finite element
discretization of the peridynamic model is provided and an effective event-driven numerical algorithm
is described.

Keywords: Damage, localization, nonlocal continuum, peridynamics.

1. Bond-based peridynamics

For natural as well as man-made materials, the re-
sulting macroscopic response is dictated by the fine-
scale structure of the material. Individual building
units (such as atoms on the nanoscale, inclusions and
matrix on the microscale, or bricks on the mesoscale)
have a certain characteristic size and their interaction
takes place at finite distance. This can be taken into
account by particle/lattice models, or by various en-
riched continuum theories. In this work, we consider
the bond-based version of the peridynamic model and
study the effect of the specific choice of the critical
strain function on the evolution and localization of
the macroscopic damage.

Peridynamics was introduced in 2000 by
S. A. Silling as a reformulation of the classical
theory of elasticity [1]. In peridynamics, the partial
differential equations of the classical continuum
theory are replaced by integral equations. Since
differentiation of the displacement field is not
needed, peridynamics is especially convenient for
capturing discontinuities such as cracks without the
complications of mathematical singularities.

In a peridynamic model, the material is represented
by infinitely many infinitely small particles that can
interact at a finite distance. Usually, interaction be-
tween two points is considered only if their distance
does not exceed the so-called horizon h, see Figure 1.

The acceleration ü of a particle located at a spatial
position x in the reference configuration at time t can

h

Figure 1. Two-dimensional representation of a dis-
cretized peridynamic model: regular grid of particles
(black dots) and one selected particle (hollow cir-
cle) with interactions (grey lines) within horizon h
(dashed line).

be found from the peridynamic equation of motion

ρü(x, t) =
!

Hx
f (u(x′, t) − u(x, t), x′ − x)dVx′

+ b(x, t)
(1)

ρü(x, t) =
!

Hx
f (u(x′, t)−u(x, t), x′−x)dVx′ + b(x, t)

where Hx is the set of points that interact with the
particle at x, u is the displacement vector field, b

47

http://dx.doi.org/10.14311/APP.2021.30.0047
http://ojs.cvut.cz/ojs/index.php/app


K. Mikeš, M. Jirásek, J. Zeman et al. Acta Polytechnica CTU Proceedings

⇠

h�h

s0(⇠)

⇠

h�h

s0(⇠)

⇠

h�h

s0(⇠)

Figure 2. Various functions describing the depen-
dence of failure strain on the initial bond length:
constant (left), decreasing (middle), and increasing
(right).

is the prescribed body force density field, ρ is the
mass density in the reference configuration, and f is
the so-called pairwise force function. The value of
f is the force vector (per unit volume squared) that
a particle located at the position x′ exerts on the
particle located at x. The pairwise force function
can be written in the form

f (η, ξ) =
ξ + η
|ξ + η|

f (|ξ + η|, ξ) (2)

where

ξ = x′ − x (3)

is the relative position of the interacting particles in
the reference configuration,

η = u(x′, t) − u(x, t) (4)

is their relative displacement, and f is the scalar bond
force governed by a constitutive relation.

In this work, bonds are considered as perfectly
elastic-brittle, in line with the concept of microelastic
brittle material, adopted from [2]. The corresponding
bond force is defined as

f (|ξ + η|, ξ) = c s(|ξ + η|, |ξ|) µ(t, ξ) (5)

where c is a spring constant called the bond elastic
micromodulus,

s(|ξ + η|, |ξ|) =
|ξ + η| − |ξ|

|ξ|
(6)

is the bond stretch and

µ(t, ξ) =
"

1 if s(t′, ξ) < s0 for all 0 ≤ t′ ≤ t,
0 otherwise

(7)

is a history-dependent damage function that jumps
from 1 to 0 when the bond breaks, which is as-
sumed to happen when its stretch attains the critical
level, s0.

In this work, a fixed value of the elastic micromod-
ulus c is considered for all bonds, whereas the critical
stretch for bond failure, s0, is evaluated for each link
as a function of its length. Three different forms of
length functions, namely constant, linear decreasing,
and linear increasing, are employed, as illustrated in
Figure 2.

l

ux

h = 0.3lx 0.1l

⇠ = 0.1l

⇠ = 0.2l

⇠ = 0.3l

Figure 3. Geometry of an illustrative example of
the one-dimensional finite element model of a peri-
odic unit cell of length l with horizon h = 0.3l: pri-
mary nodes (filled circles), slave nodes (hollow cir-
cles), virtual node with prescribed displacement (gray
circle), and an expanded illustration of the considered
interactions with different lengths ξ (gray solid lines).
Dashed arrows denote connections of slave nodes with
primary nodes.

2. Numerical implementation
For numerical implementation, we use an analogy
with the finite element analysis, adopted from [3],
to represent a bond-based peridynamic system by a
finite element truss model. Thus, individual bonds
on the microscopic level are modelled by truss ele-
ments, so-called links, and the material parameters
of the bonds, c and s0, are related to the macroscopic
Young modulus E and critical strain ε0, respectively.

2.1. One-dimensional model
Damage localization and evolution of failure mecha-
nisms in a one-dimensional periodic lattice submitted
to positive strain is investigated for different forms of
the critical strain length function. The geometry of
an illustrative example of the finite element model
of a one-dimensional periodic unit cell is depicted in
Figure 3. All nodes in the periodic unit cell (except
the last one) are modeled as primary nodes that carry
independent degrees of freedom (horizontal displace-
ments). To enforce periodicity, additional slave nodes
are constructed as periodic images of primary nodes
on one side of the periodic cell. The number of pe-
riodic slave nodes is given by the size of the hori-
zon. Loading in the form of prescribed displacement
is applied to the cell through one virtual node which
is linked to all slave nodes simultaneously. The dis-
placement of each slave node is evaluated as the sum
of the unknown displacement of its master (one of
the primary nodes) and the prescribed displacement
of the virtual node.

Bonds, in the form of one-dimensional truss ele-
ments, are placed between all pairs of primary nodes
as well as between primary nodes and slave nodes
within the given horizon h; see the solid grey lines in
Figure 3. The cross-sectional area of a truss element
is the same for all links and it is scaled to provide a
unit sum of overlapping link areas representing the
macroscopic cross-sectional area.

48



vol. 30/2021 One-dimensional peridynamic model damage analysis

During the loading process, the evolution of dam-
age in a macroscopic cell is modelled as a sequence
of failures (breakages) of individual links followed by
stress redistribution. The macroscopic damage ω can
be evaluated for each point x of the unit cell as

ω(x) =
Nf (x)
N (x)

(8)

where N is the number of all links passing through
the given point x and Nf is the number of broken
links (i.e., those with µ = 0) passing through x.

Theoretically, it is possible that multiple links can
reach the critical strain simultaneously. Then the
subsequent solution is not unique and additional
choices must be made to distinguish one link or a
group of links preferred for failure in the given step.
For an ideal lattice geometry without imperfections,
such a situation occurs quite often, not only at the
beginning of the loading, but also almost every time
when broken links form a symmetric pattern. To
avoid this phenomenon and to more closely mimic the
behavior of real materials, the initial position of the
primary nodes is randomized by a small perturbation
with a uniform distribution and a maximum value
l × 10−6. The macroscopic response of the random-
ized model is evaluated statistically based on multiple
random realizations.

2.2. Solution algorithm
To capture the exact sequence of breaking links, it
is important to investigate not only the beginning
of failure but also the full failure process of individ-
ual links and take into account that a new link can
start failing before another (already failing) link be-
comes fully stress-free. For elastic-brittle links, fail-
ure means breaking, i.e., decreasing link stress under
constant link strain. In this case, neglecting the pos-
sibility of partial unloading usually has only a mi-
nor effect on the macroscopic response. However,
for elastic-damaging links with softening (which are
not considered in this contribution), the possibility of
partial unloading is essential.

It is obvious that the size of the loading increment
to reach link failure or unloading is different in each
step. Therefore, conventional solvers for an incre-
mental nonlinear structural analysis are inconvenient
for the present problem because they may miss the
correct failure mechanism if the prescribed increment
is too large. To solve the present problem exactly, an
event-driven algorithm with step size control intro-
duced in [4] has been employed. Since the bond-level
constitutive law is piecewise linear, the step size can
be selected such that the tangential stiffness matrix
does not change during the step. Therefore, no equi-
librium iteration is necessary. In one solution step,
either a positive or a negative force increment is con-
sidered to achieve the onset of failure or complete
failure of one link. Furthermore, an effective method
of inelastic forces can be used to avoid the necessity

of decomposing the tangential stiffness matrix in each
step; for more details see section 4.2 in [4].

3. Results
In this section, results are presented for the afore-
mentioned finite element representation of the peridy-
namic model of a one-dimensional periodic cell with
elastic-brittle interactions subjected to positive (av-
erage) strain. The modelled cell is considered to have
a length of l = 1 m, with the horizon h = 0.1l and
cross-sectional area A = l2. The links have a constant
Young modulus E = 100 N m−2 and critical strain
function with a maximum value ε0 = 1 × 10−2 and
various forms, namely constant, linearly decreasing,
and linearly increasing with link length.

3.1. Convergence study
First of all, a convergence study of the peridynamic
model with decreasing computational grid spacing is
performed. For that purpose, five simulations with
the same macroscopic parameters but with 10, 20,
50, 100, and 200 primary nodes have been performed.
The simplest case of constant critical strain for dif-
ferent link lengths has been considered. The results
are shown in Figure 4 in terms of normalized force-
displacement diagrams, in which the elastic part of
the displacement, F/Kel, is subtracted from the to-
tal macroscopic displacement u and only the inelastic
part is plotted. F is the macroscopic reaction force
and Kel denotes the initial elastic stiffness. For the
purpose of clarity, only one representative random re-
alization is shown for each discretization.

Figure 4. Convergence of normalized force-
displacement diagrams (elastic part of displacements
has been subtracted) for one-dimensional periodic lat-
tice with fixed size l, horizon h = 0.1l, for discretiza-
tions with different numbers of nodes. For each dis-
cretization, only one representative random realiza-
tion is shown.

The maximum force corresponding to the initial
failure of the first link as well as the maximum dis-
placement for the final damage mechanism are pre-
dicted exactly for all refinements. However, for a
small number of nodes, significant oscillation can be

49



K. Mikeš, M. Jirásek, J. Zeman et al. Acta Polytechnica CTU Proceedings

Figure 5. Results for constant critical strain: final damage profile (top left), positions and lengths of broken
links (bottom left), and normalized force-displacement diagram where the elastic part of the displacement has been
subtracted (right).

Figure 6. Results for critical strain linearly decreasing with link length: final damage profile (top left), positions
and lengths of broken links (bottom left), and normalized force-displacement diagram where the elastic part of the
displacement has been subtracted (right).

observed, while with an increasing number of nodes,
the macroscopic response tends to a smooth curve.
Discretization with 100 nodes turned out to be suffi-
cient to capture the important characteristics of the
macroscopic response.

Note that for 10 primary nodes, there is only one
link between each pair of neighbouring nodes and
therefore damage localizes in one step by failure of the
first link and the corresponding response is a straight
line, see the gray line in Figure 4.

3.2. Damage mechanisms
In this subsection, the influence of the form of the
critical strain function is analyzed. The discretiza-
tion with 100 nodes is used based on the previous
study. The results are shown in Figures 5, 6, and 7

for the constant, decreasing, and increasing functions
of the critical strain, respectively. For each function,
10 independent random realizations have been eval-
uated to analyze the influence of randomness to the
mechanisms of failure. However, for the purpose of
clarity, only the results of one representative random
realization are depicted.

The resulting mechanisms of failure are illustrated
by the final damage profiles (top left in each figure),
evaluated according to Eq. (8), and by the positions
of all broken links (bottom left); here, each blue dot
represent one broken link. The length of the broken
link ξ/l is plotted against the position of its center.
The values are normalized by the cell size l. More-
over, the corresponding inelastic force-displacement
diagram (right) is shown to illustrate the macroscopic

50



vol. 30/2021 One-dimensional peridynamic model damage analysis

Figure 7. Results for critical strain linearly increasing with link length: final damage profile (top left), positions
and lengths of broken links (bottom left), and normalized force-displacement diagram where the elastic part of the
displacement has been subtracted (right).

behavior during damage evolution.

3.2.1. Constant critical strain
In the initial (undamaged) configuration, the strain of
the links is uniform. Consequently, for the first case
of constant critical strain, there is no preferred length
of links for failure in the first step. In our simulations,
one of the shortest links is made weaker to initialize
damage localization. Once the shortest link breaks,
the stress is redistributed, two links in the second
shortest layer spanning the broken link become the
most loaded links and one of them fails. Then, each
subsequent broken link is situated in the next longer
layer spanning all previously broken links until the
first link in the longest layer is reached. Once this se-
quence of propagating broken links from the shortest
to the longest layer is finished, similar sequences fol-
low starting from the second shortest to the longest,
the third shortest to the longest, etc., until the final
triangular pattern visible in Figure 5 (bottom left) is
formed.

Note that if the shortest layer of links is not artifi-
cially favored at the beginning, the initial failure can
occur for a different link length due to the random-
ization. However, in such a case, the closest shortest
links break in the next step and the same failure se-
quences can be observed with only the exception of
the first link.

For constant critical strain, the position of max-
imum damage (ω = 1), i.e., the position of macro-
scopic failure, is dictated by the first broken link. The
final damage mechanism is simple. The resulting fail-
ure pattern and associated failure properties such as
fracture energy can be evaluated from the geometry
of the cell. Also the force-displacement diagram does
not exhibit any significant irregularities.

3.2.2. Decreasing critical strain
In the case of a critical strain decreasing with link
length, the longest layer of links is the weakest
one and the first broken link appears in this layer.
Subsequently, other links in the weakest layer are
breaking one by one until the weakest layer is com-
pletely broken. This process is repeated for each
next longest layer of links until a critical localization
layer is reached. The gradual breaking of individ-
ual layers corresponds to horizontal parts of individ-
ual steps (jumps) in the force-displacement diagram.
Particularly, four steps can be observed at the lev-
els of normalized force about 0.01×10−2, 0.08×10−2,
0.12×10−2, and 0.135×10−2, see Fig. 6 (right). From
the macroscopic point of view, this process can be
understood as a uniform growth of damage.

Once the critical localization layer is reached, bro-
ken links localize in a triangular pattern similar to
the previous case of constant critical strain; see Fig. 6
(bottom left). In the force-displacement diagram, lo-
calization corresponds to the final unloading part,
where both macroscopic force and macroscopic dis-
placement are decreasing.

The position of maximum damage is not affected
by the first broken link but it is determined during the
localization of the damage in the critical localization
layer.

3.2.3. Increasing critical strain
For the last case of increasing critical strain, the
shortest layer of links is the weakest one and the first
broken link appears here. However, in this case, dam-
age of the first link is followed neither by damage of
the whole weakest layer nor by localization. On the
contrary, links continue breaking in different layers
without a specific pattern. Consequently, the macro-
scopic damage is distributed randomly and no clear

51



K. Mikeš, M. Jirásek, J. Zeman et al. Acta Polytechnica CTU Proceedings

localization is observed. In the end, the majority of
broken links are situated in several weakest layers.
However, these weak layers (except the last one) are
not completely broken. On the other hand, in several
regions, damage significantly propagates also into lay-
ers of stronger links, see Fig. 7 (bottom left).

In the central region, the broken links form an
unusual mechanism of failure where no macroscopic
cross-section exhibits complete damage (ω = 1), even
when the macroscopic reaction force completely van-
ishes. Such a mechanism is formed by overlapping
sliding parts. As a result of this, the final macroscopic
strain is localized not only into one but into three seg-
ments between neighbouring nodes where two of them
exhibit positive strain increments and the remaining
one a negative strain increment of equal magnitude.

Unlike the previous cases of constant and decreas-
ing forms of critical strain function, in this case, the
order of broken links is unstable and the position and
shape of the final damage pattern as well as the evo-
lution of the damage before localization are affected
by the randomization of the initial lattice geometry.
As a result, significant oscillations can be observed in
the force-displacement diagram, see Fig. 7 (right).

4. Conclusions
In this contribution, a finite element truss model rep-
resenting a peridynamic system with elastic-brittle
interactions has been introduced, together with an ac-
curate and efficient event-driven numerical algorithm
which is able to capture exactly not only the failure
but also the full unloading of individual links. The
convergence study of the presented model has shown
that upon refinement the system tends to a smooth
macroscopic response.

The influence of various forms of the critical strain
function on the type of macroscopic response has been
investigated for a one-dimensional periodic unit cell.
Constant, decreasing, and increasing critical strain
functions have been investigated and the resulting
responses with different types of damage evolution
have been analyzed, namely direct damage localiza-
tion, uniform damage evolution preceding localiza-
tion, and random damage distribution, respectively.

The presented results demonstrate that the con-
cept of critical strain values depending on link length
can be used to reach various types of macroscopic
responses even for such a simple geometry (1D) and
one of the simplest nonlinear material models (elastic-
brittle). The significant advantage of this concept is
that the presented elastic-brittle model can be fit-
ted to represent various macroscopic behaviors, but
at the same time the micro-level constitutive law re-
mains piecewise linear and therefore the exact and
efficient event-driven algorithm can be used.

In addition, the spectrum of different response
types can be further extended by considering vari-
ous micromodulus functions or by adopting elastic-
damaging links with softening.

List of symbols
h Material horizon [m]
Hx Neighborhood of a particle
Vx Volume of a particle neighborhood
x Particle coordinates [m]
u Displacement vector [m]
ü Acceleration vector [m s−2]
ρ Mass density [kg m−3]
b Body force density [N m−3]
t Time [s]
ξ Relative position of two particles [m]
η Relative displacement of two particles [m]
f pairwise force function [N m−6]
f scalar bond force [N m−6]
s Bond stretch [–]
x Horizontal coordinate [m]
l Cell size [m]
A Cross-sectional area of links [m2]
ux Prescribed horizontal displacement [m]
u Total lattice displacement [m]
F Reaction force [N]
Kel Initial elastic stiffness [N m−1]
c Spring constant [N m−2]
s0 Critical stretch for bond failure [–]
E Young modulus [N m−2]
ε0 Critical strain [–]
µ Damage function [–]
ω Damage [–]
N Number of links [–]
Nf Number of broken links [–]

Acknowledgements
This research has been performed in the Center of Ad-
vanced Applied Sciences (CAAS), financially supported
by the European Regional Development Fund (project
No. CZ.02.1.01/0.0/0.0/16_019/0000778). Financial
support received from the Czech Technical University in
Prague under project SGS20/038/OHK1/1T/11 is grate-
fully acknowledged.

References
[1] S. A. Silling. Reformulation of elasticity theory for

discontinuities and long-range forces. Journal of the
Mechanics and Physics of Solids 48(1):175–209, 2000.

[2] S. A. Silling, E. Askari. A meshfree method based on
the peridynamic model of solid mechanics. Computers
& structures 83(17-18):1526–1535, 2005.

[3] R. W. Macek, S. A. Silling. Peridynamics via finite
element analysis. Finite Elements in Analysis and
Design 43(15):1169–1178, 2007.

[4] M. Jirásek, Z. P. Bažant. Macroscopic fracture
characteristics of random particle systems.
International Journal of Fracture 69(3):201–228, 1994.

52