Failure and scaling properties of a softening interface connected to an elastic block

The damage growth in a softening interface connected to an elastic block is analysed. The elastic block, assumed to be infinite, is modelled as a two-dimensional continuum and the interface is one-dimensional with a constitutive response which follows a scalar damage model. The solution technique is based on the equilibrium of the interfacial forces resulting from the deformation of the elastic block and from the interface constitutive response. The interface failure process is compared to that of a hierarchical model which was obtained analytically (Delaplace et al., 1998). The two are found to be similar, without an internal length scaling the distribution of damage at the inception of macro-cracking. Finally, scale effects on the occurrence of bifurcation and instability are considered. It is shown that bifurcation may occur prior to or after the limit point under displacement control, depending on the elastic block height or stiffness.


Introduction
One of the important characteristics of quasi-brittle heterogeneous materials in tension is the existence of a strain softening regime. This peculiarity of the stress-strain response induces strain localisation which can be understood as a sudden increase of the strains in a very narrow region. In tension, the inception of strain localisation is often seen as a bifurcation from an homogeneous deformation regime to a non-homogeneous, continuous or discontinuous one, occurring in a band of material oriented perpendicularly to the applied tensile load. If the band is small with respect to the specimen size, it can be considered as a softening interface which is a precursor of the macro-crack observed at complete failure.
Among the popular techniques for describing this type of material and structural response is the approach where the material response in the softening band, collapsed onto a line, is expressed as a function of the stress versus the relative displacement across the band (the cohesive crack model). This function is, of course central in the constitutive model and there have been many proposals, assuming that the stress versus relative displacement relation is linear, bi-linear, or exponential (see for a recent review (Bažant and Planas, 1998)). Once a crack is initiated, that is when the material stress has been exhausted in the course of loading in a region of the band, another issue is the fracture process zone which occurs ahead of the crack due to progressive failure. The size of this process zone, here the length of the interface ahead of the crack where softening occurs, induces size effect which is typical of quasi-brittle structures. Hence it is critical to understand how the fracture process zone forms in the course of loading.
In the present contribution, we consider the response of a strain softening solid in tension from the viewpoint of an interface problem. Hence the width of the softening zone is assumed to be small compared to the specimen length. We will use a one-dimensional continuous damage model for describing the interface behaviour. Our intent is to understand how failure is initiated, what is the process of damage localisation in the interface, and what scaling parameters can be expected as outcomes of the redistribution of forces in the interface and of the interaction between the interface and the rest of the specimen.
In order to better understand how damage develops in the interface, a first study was made, based on a hierarchical model which allows to extract the interface behaviour from the global system response (Delaplace et al., 1998). The interface was composed of a set of elastic brittle fibres, and heterogeneity was included through a random variable, the threshold where the fibres break. With this model, we have shown that at the beginning of the loading, damage develops homogeneously all over the interface. Then damage becomes localised in the interface, first at a macro scale, that is the half system size, and further in a narrower and narrower region. In the limit of an infinite number of fibres corresponding to a continuous description of the interface, this process was understood as a cascade of bifurcations. At the initiation of the macro-crack, the damage profile in the interface was computed analytically, and damage was shown to spread all over the interface. In other words, the fracture process zone resulting from the progressive failure process considered did not exhibit a scale other than the entire length of the interface.
Although the hierarchical description is equivalent to the classical elastic one and provides a convenient setting for deriving analytical results and is more efficient for computations, it has two drawbacks: (i) because the interaction between two points on the interface occurs through hierarchical levels, correlation lengths are not straightforward to interpret physically; (ii) during the failure process, the discrete nature of the hierarchical description becomes more and more prominent and can result into unrealistic responses such as strength recoveries. In order to confirm the above results and to arrive at an easier interpretation, we will use here a more classical representation where two-dimensional elasticity is used for the description of the elastic block instead of a hierarchical model. The tensile specimen will be modelled by an elastic block connected to an interface with a continuous damage law for its behaviour, without any imperfection at the interface. In order to detect possible bifurcation points along the specimen response, we will implement classical criteria (see for instance (Pijaudier-Cabot and Huerta 1991;De Borst, 1987) in order to treat loss of uniqueness or loss of stability.
In the first part of the paper, we present the model geometry and the discretisation technique which yields to the basic equilibrium equations, paying specific attention to the boundary conditions. The bases of the stability and uniqueness criteria are recalled in section three. Then, the interface failure process is compared to the hierarchical model. The two are found to be similar. Finally scale effects on the occurrence of bifurcation and instability are considered.

Discretisation of the problem
We consider in this contribution the problem of an elastic block connected to a rigid substrate through an interface. The elastic block is a two-dimensional body, infinite in direction x of height h (Figure 1). At y = h, the boundary of the elastic block is constrained to follow a fixed displacement u c . The material properties of the block are its Poisson ratio ν, and its Young modulus E. Periodic loading conditions are used on the boundary y = 0 in the x direction. Hence, the periodic cell is an elastic body connected to an interface at x = 0. L and h are the width and the height of the cell respectively. The interface is one-dimensional which means that the shear stiffness and tension-shear coupling terms are neglected. The interface is a row of n identical uniaxial elements of width l = L/n (its behaviour will be detailed in Section 2.4). As we will see later L/l will be the discretisation parameter in the solution technique. Note that because the interface is one-dimensional, the rigid substrate on which the interface is connected can be seen as an axis of symmetry. Therefore, the problem considered is also equivalent to that of an infinite elastic band of width 2h, with a softening interface in the middle, which is subjected to fixed displacements on its upper and lower boundaries.
In order to solve the equilibrium equations in the block and at the interface, one can use the finite element method. Another technique, which will be implemented here uses, to some extent, the boundary element method for computing the displacement and stresses in the elastic block and to write the equilibrium of each uniaxial element of the interface. In fact, one needs only to know the relationship between the boundary tractions and the boundary displacements in the upper (y = h) and lower (y = 0) contours of the periodic cell. In the following, we are going to derive from known Green functions these relationships, starting from a single elastic block and then implementing the periodicity conditions. Then, the equilibrium at the interface level will be introduced along with the nonlinear constitutive response of the interface.

THE ELASTIC BLOCK
We first consider the elastic block without periodic conditions. The vertical displacement of each point of the block is called u(x, y). x varies between −∞ and +∞, as y varies between 0 and h. Note that for y = h, we have (1) Because we are just focusing on the displacement along the interface, we will call in the following u(x) the displacement u(x, 0) in order to simplify the notation. For a perfectly elastic block, a point load P applied at the lower surface (y = 0) at the position x 0 induces a displacement on this boundary which is (see for instance (Timoshenko, 1948)) where the symbol |x| refers to the absolute value of x. In our case, we will deal with distributed loads which are due to the interface response. Because of the discretisation of the interface, these distributed forces will be piecewise constant over the length l corresponding to the discretisation parameter. Hence, we consider now a constant pressure distribution q over a distance l, which starts at position x 0 . At points (with the coordinate y = 0) located outside the pressure distribution (Figure 2a), and for At points (with the coordinate y = 0) located outside the pressure distribution, and for x > For the points located inside the pressure distribution, the displacement becomes ( Figure 2b) Note that the two expressions in (3, 5) are obtained under the condition that the displacement on the upper boundary of the block at y = h satisfies the prescribed fixed displacement u c .

PERIODICITY OF THE PRESSURE DISTRIBUTION
Let us come back to our problem with periodic loads. The periodicity condition for the displacements is u(x) = u(x + L). The periodicity of the pressure distribution entering in (3, 5) needs also to be taken into account. This can be performed with a superposition technique where the influences of pressure distributions spaced out by L are added (see appendix for the complete derivation). An other way to find this result is to consider a Fourier transform version of the Green function. Then, the mode k = 0 can be treated in real space as the mere uniaxial compression, like that introduced in the appendix. From (32, 36), we observe that for a large number m of left and right elastic blocks (cells) in the superposition scheme, the displacement behaves as  where u ∞ is the exact value of the displacement, and A a constant. Note that this result could be probably demonstrated analytically by considering the limit of the series defined in the appendix. A numerical verification of (6) has been performed by plotting the displacement versus the inverse of number m of blocks ( Figure 3). The intersection with the vertical axis provides the exact displacement for an infinite periodic block resulting from the superposition of an infinite number of pressure distributions.

DISCRETISATION OF THE ELASTIC BLOCK RESPONSE
The last step necessary to compute the displacement is to write the previous (32, 36) in a discrete form. We choose to compute the lower surface displacement in the middle of each part of width l. For a distributed pressure q j on the element j and a measured displacement u i on element i, the continuous variable x (the distance between the applied load and the measured displacement) is changed into (32) x = l/2 inside the distributed load (36).
According to the equations derived in the appendix, the relation between u i and q j becomes We introduce now the distributed load vector q = {q 1 , q 2 , . . . , q n } and the corresponding displacement vector u = {u 1 , u 2 , . . . , u n }. The discretised equations can be cast in the form where G is the Green matrix of the problem. Note that G is symmetric and is independent of u. The displacement of the block surface at the level y = 0 of the interface can be computed for any pressure distribution and with the fixed displacement u=u c at y = h. Figure 4 shows an example of such a computation where a unit distributed load (tension) is locally applied between x = 7 and x = 8 and the block length is L/ l = 32.

THE INTERFACE BEHAVIOUR
The constitutive relation of the interface is assumed to derive from a damage model. Because the interface is one-dimensional, damage is a scalar variable, called D, that evolves between 0 (undamaged material) and 1 (broken material). For any interface element, the force-displacement relation is where C is the stiffness of the undamaged material. We introduce the simplest form for the damage assuming it is simply proportional to the displacement This relation is in fact, the asymptotic behaviour of a fibre bundle model with an infinite number of fibres of unit stiffness, and with a uniform distribution between 0 and 1 for the thresholds where the fibres break irreversibly (Daniels 1945;Hemmer and Hansen, 1992). For simplicity, we introduce a dimensionless displacementũ = u/u 0 . With these notations, we have simply D =ũ.
In the following, geometrical distances will be made dimensionless by scaling them by the size of one element l. Forces will be scaled by Cu 0 / l 2 . Thus we introduce the reduced variableẼ = El 2 /Cu 0 , andf = f l/(Cu 0 ). Note thatẼ is the only the pertinent parameter in our problem. Equation (10) becomes For an element i, the equivalent distributed load q i over the element width l is Finally, if the element is undergoing loading, the relation betweenq i andũ i is In the situation of unloading, we havẽ where the secant modulus K sc i = 1 − max(ũ i ) over the loading history. The corresponding rate relations are From now on, we will omit the ∼ superscript for reduced variables.

EQUILIBRIUM OF THE INTERFACE
The set of (9) has been derived from the equilibrium of the elastic block. We have now to write the equilibrium of the interface, that is to write that the distributed forces applied on the lower surface of the block result from the interface constitutive relations. Because the interface has a nonlinear constitutive response, we write this condition in a rate form. Equation (9) whereu = {u 1 ,u 2 , . . . ,u n } the vector of nodal velocities andq = {q 1 ,q 2 , . . . ,q n } the vector of loading rate are introduced. In the same way, the interface response is written from equation (16) where K(u) is just a diagonal matrix. In these two last equations, the vectors of the rate of distributed forces should be equal in order to ensure equilibrium of the interface. We will define now Hence, the global rate equation that governs the system is Or where M = K + G is a symmetric matrix. Note that this equation has been derived assuming that the displacement at the upper surface of the elastic block is controlled with a rateu c . Load control conditions would provide a slightly different system of equations.

LOSS OF UNIQUENESS
The system has an unique solution unless M is singular, i.e.
Because M is real and symmetric, therefore it can be diagonalised, this condition corresponds to the existence of at least one vanishing eigenvalue.

LOSS OF STABILITY
For our problem, the condition that guarantees the stability of the system is (Hill, 1959) u t .M.u > 0 ∀u (23) The condition of critical equilibrium (u t .M.u = 0) is det M = 0 oru is orthogonal to Mu.

POST-BIFURCATION RESPONSE
When the condition det M = 0 is reached, i.e. when at least one eigenvalue, called here λ 1 , vanishes, the solution of the rate equilibrium equations is no longer unique and a classical perturbation technique can be implemented in order to search for another solution in addition to the fundamental oneu (De Borst, 1987;1988). We call v 1 the eigenvector associated to λ 1 . One can easily show that the following displacement rate vectoru * also satisfies the equilibrium equation for any arbitrary scalar β We have also to determine if the singularity of M corresponds to a limit point (under displacement controlled conditions) or to a bifurcation point. This can be performed following the steps in (De Borst, 1987;1988). The solution verifies u could be written in the M-eigenvectors v i basė The equilibrium equation (21) becomes All the eigenvectors are independent, and then for all i we have For i = 1, it leads to v t i .u g = 0. Then, a limit point corresponds tou g = 0 (loss of stability). A bifurcation point corresponds to v t i ⊥u g (loss of uniqueness).
For our system, the first eigenvalue that vanishes is unique for a limit point, and is at least double for a bifurcation point because of the periodic boundary conditions. In these equations, however, the eigenvalues depend on M, and then on K. But from (16), K is not a single valued matrix. It depends on the strain rate field in the interface, i.e. if damage grows or not in each element. It follows that a complete bifurcation and stability analysis would require the investigation of all the possible loading/unloading combinations in the rate constitutive relations of the interface. This can be hardly performed, except for a very coarse discretisation of the interface. Computations for coarse discretisations show that it may be sufficient to consider the tangent stiffness matrix computed under the assumption that all the interface elements which were loading (growth of damage) during the previous displacement step are still undergoing loading.

BIFURCATION
The post bifurcation analysis is applied to the system with a number L of elements in the interface that varies between 4 and 32. The first singularity encountered is as just mentioned either a limit point (i.e. the maximum global displacement) or a bifurcation point. At the bifurcation point and according to the considerations in the previous section, we look at the case where the vanishing eigenvalue is at least double. Then, two independent eigenvectors are found. They have a L-periodic 'sine' form, and are shifted forward with a phase of L/4. Note that because of the periodic boundary conditions, all eigenvectors of this form could be associated to the vanishing eigenvalues. One of these eigenvectors is then added to the perturbationu in accord with (25), and the loading could be continued according to the procedure defined in (De Borst, 1987;Pijaudier-Cabot and Huerta, 1991). In the all computations only one bifurcation point was found during the course of loading. After bifurcation occurred, the equilibrium equations had a single solution in the post-bifurcation regime.
The response of a 32-element system is shown in Figure 5. It is similar to the hierarchical model response (Delaplace et al., 1998). This result is also obtained if we introduce artificially an initial defect (initial nonzero damage for instance) in the interface. In this case, the damage profile in the interface is determined by the defect location and bifurcation does not occur, same as in the hierarchical model with initial disorder in the interface. In the presence of an initial defect, the global behaviour is only slightly affected.

DAMAGE EVOLUTION
One of the main motivations of this work is to study the evolution of damage in the interface. This evolution is plotted in Figure 6. From the very beginning of the loading to the first bifurcation point, damage D grows homogeneously in the interface. It condenses progressively into a narrower and narrower region down to just one element, which breaks (i.e. D = 1 in this element). Then, a macro-crack is initiated in this element and propagates in the interface. This last stage is not analysed in this paper. This damage evolution is very close to what was encountered with the hierarchical model. An interesting analytical result was obtained with the hierarchical model about the damage profile at the initiation of the first macro-crack. This unusual profile was where x is the distance, in the hierarchical sense, between the macro-crack location and the location at which damage is computed. This analytical result could not be found with the present approach, but, as shown in Figure 6, the damage profile at the initiation of the first crack is well fitted with such an inverse law with the distance. Note that the damage spreads all over the interface, and that there is no characteristic length which scales the damage distribution at the inception of macro-cracking. This is contrary to what is expected in the Hillerborg approach (Hillerborg et al., 1976), at least in the regime of inception of cracking. It is important to note the present conclusion does not pertain to the macro-crack propagation regime which is not analysed here.

SCALE EFFECT
Two scaling parameters are relevant in this problem. The first one is the block height, called h, and the second is the ratio between the block stiffness E and the interface stiffness C.
If we change these two parameters independently, we observe the same effect: the initial slope in the force displacement response changes (Figure 7). It can be easily understood by considering an equivalent system: a spring of stiffness E and of length h which is coupled in series with one damageable element influences the initial slope of the response, same as in the one-dimensional analysis devised in (Bažant and Cedolin, 1991). A more relevant variable is the so called variable g 0 (Section 2.5). This variable depends of course on E and h, but its influence is directly related to the global behaviour of the system through (21). Hence, a constant g 0 gives the same initial stiffness and the same characteristic if the interface displacement is homogeneous. Under this condition, the influence of h and E can be conveniently analysed: The main effect of h and E is to shift the position of the first  bifurcation point. More specifically, bifurcation can be encountered before or after the limit point. For a fixed g 0 , if h and E increase, the first bifurcation is delayed during the loading. On the lower left zoom of the Figure 8, the first post-bifurcation branch (dashed-bold) is in the stable regime. The solid becomes unstable after the maximum displacement is reached (fine-dashed). Upon an increase of h and E (lower right zoom) the response of the system becomes unstable before bifurcation occurs (fine continuous curve).

Conclusion
We have analysed the growth of damage in a system that is composed of a softening interface coupled in series to an elastic block. This study follows a first one where a hierarchical model was used in order to obtain main results analytically. We have used here a more classical representation (two-dimensional continuum), instead of a hierarchical decomposition. We have focused attention on the post bifurcation response, where damage localises progressively. The following conclusions can be drawn -Damage grows first homogeneously in the interface. At the first bifurcation point, damage localises continuously from a macroscopic scale to a narrower and narrower region. Then, just one element fails and a macro-crack is initiated.
-The damage profile is computed at the initiation of the macro-crack. For the chosen interface behaviour, an inverse power law with the distance was obtained for the hierarchical model and the numerical result is very close the analytical one.
-The scaling parameters in the problem are the block height and the relative stiffnesses of the elastic block and of the interface. These two parameters influence the relative location of the bifurcation point and of the limit point on the response of the solid. There is no length scale in the damage profile at the inception of macro-cracking.
Finally, it should be noted that obtaining experimentally the response of such a softening interface should be quite difficult, due to the loss of stability and most of all due to the bifurcation of equilibrium states. The same problem is encountered when the constitutive response of a material is sought from a tensile test.

Appendix: Displacement due to periodic pressure distribution
We derive here the displacement corresponding to a periodic distribution of pressure using a superposition. Let us recall the elastic block properties h is the block height; -L is the width of one cell; -ν is the Poisson ratio; -E is the Young modulus; -u c is the displacement of the upper surface (constant all along the block).
The lower surface displacement verifies the relation: u(x) = u(x + L). For one cell, the pressure is spread over a distance l, and starts at the position x 0 . In the following equations, x is assumed to vary between 0 and L, as x 0 varies between 0 and L − l. From (3), the displacement is A m is the left load contribution, B refers to the pressure distribution applied on the block where the displacement is measured, and C m is the right load contribution. In the same way, if the measured displacement is inside the pressure distribution, the displacement is with Figure 9. Superposition of the two problems.
D m is the left pressure contribution, E refers to the pressure distribution applied on the block where the displacement is measured, and F m is the right pressure contribution. This superposition scheme leads to an infinite displacement, due to the constant expression (1 − ν)/(πE)ql. Of course, this unrealistic physical feature is due to the periodic conditions. In order to avoid this divergence, one solution is to superpose two problems (Figure 9) -the first problem is the infinite elastic block, loaded by a pressure distribution q over a length l, that starts at position x 0 , and an uniform pressure distribution −ql/L all over the surface. Hence, the load resultant along the surface is zero, and the divergent term vanishes in (8). -The second problem is the same elastic block, loaded by just an uniform pressure distribution ql/L. In this case, the induced surface displacement is finite and well known for all x Finally, the addition of these two displacements provides a finite displacement as expected on physical grounds.