eu DESIGN OF ONE-DIMENSIONAL STRUCTURES WITH PERIODIC INELASTIC MATERIALS

This work shows the design of one-dimensional periodic structures composed by inelastic materials. These structures are composed by multilayer unit cells with different materials (linear and non-linear). In this work, we design three types of unit cells. Newmark scheme for integration and finite element method for space discretization are used in the algorithms. The unit cell design proposed in this work can be extended to multidimensional problems. AMS Subject Classification: 74B05, 74C05, 74J05, 74S05, 74H15


Introduction
A composite materials can be modeled with a regular repetition of a cell in space.A periodic structures are composed for this class of material.A detailed study of materials with computational algorithms, constitutive equations was developed for several authors.Simo and Hughes [1] describe inelastic materials; Hughes [2] analyses elastic materials using finite element method; Cook [3] present elastic and plastic cases for dynamic problems; Ju [4] developed a damage theory coupling for elastoplastic materials; Wood [5] shows integration schemes for transient treatment; Kushwaha et al. [6] detail periodic structures composed by composite materials formed by unit cells; Lee and Yang [7] study waves propagation in structures with composite periodic materials; Hussein et.al. [8] study wave propagation in elastic periodic materials through numerical one-dimensional simulation using finite element method; Hussein et al. [9] proposed a design methodology for the characterization of periodic material using binary numbers; Fabro et.al. [10] study structural wave propagation in waveguides with slowly randomly varying material and geometrical properties along the axis of propagation; among others.Navarro and Souza-Braun [11] analyzed wave propagation in one-dimensional structure of periodic inelastic composites.Based on investigations by Hussein et al. ([8], [9]) and our own previous works ( [11] and [12]), we applied a methodology for design of periodic structure and we applied it for three types of unit cells.

Mathematical Model
The governing equations obtained from the study of dynamic response of an elastoplastic-damage problem are represented by the following mathematical model.
Let Ω a body with boundary ∂Ω in a one-dimensional space ℜ.The space variable is denoted by x ∈ Ω, the time variable by t ∈ I =]0, T [ and u(x, t) ∈ ℜ represent the displacement of a point x in the body.The linearized strain is represented by ǫ = ∂u(x,t) ∂x and the Cauchy stress σ is related to ǫ through a constitutive relation.The body force is written as f : Ω×I → ℜ and considering that it has axial motion, the problem is governed by Eqn.(1), [2] where Ω =]0, L[.The boundary and initial conditions are specified as where σ n is the normal stress on the body.The sets ∂Ω g and ∂Ω h are Dirichlet and Neumann types.We have ∂Ω g ∩ ∂Ω z = ∅ and ∂Ω g ∪ ∂Ω z = ∂Ω.
The normal stress depends on the kind of the materials.The constitutive equations for the materials have a linear dependence for elasticity problems and nonlinear for elastoplasticity problems (see [1]).For the linear range, the elastic constitutive equation relates linearly stress with strain as where the material parameter E is the Young's Modulus which is assumed to be constant over domain Ω .For materials with nonlinear characteristics (plasticity), the stress is a nonlinear function (possibly history-dependent) of ǫ.The total strain ǫ splits into two parts ǫ = ǫ e + ǫ p , where the superscript e referred to the elastic part and p the plastic part.The stress is written as In many metals subjected to cyclic loading, it is experimentally observed that the center of the yield surface experiences a motion in the direction of the plastic flow.A phenomenological model that captures this kinematic effect is constructed by introducing an additional internal variable, denotes by q and called back stress, which defines the location of the center of the yield surface.
The one-dimensional model of rate-independent plasticity with isotropic and kinematic hardening is examined.The hardening is isotropic in the sense that at any state of loading, the admissible relative stresses are constrained to lie in the closed interval [−(σ Y + Kα), σ Y + Kα] ⊂ ℜ.This leads to a yield criterion of the form where σ Y > 0 (flow stress) and K ≥ 0 (plastic modulus) are given constants.The variable α : [0, T ] → ℜ is a nonnegative function of the amount of plastic flow, called an internal hardening variable.If K < 0, one speaks of a strainsoftening response.The hardening is linear in the amount of plastic flow ( α = | ǫp |).The irreversible mechanism that governs the evolution of plastic flow, which is defined by the flow rule, remains unchanged.Then, we assume that ǫp = γ sign(σ − q) where γ ≥ 0 is the rate at which slip takes place and q is the back stress, which defines the location of the center of the yield surface.The evolution of the back stress q is defined by Zieglers rule where H is the kinematic hardening modulus.The irreversible nature of plastic flow is again captured by means of the Karush-Kuhn-Tucker loading/unloading conditions [4], which in the present context read γ ≥ 0, F (σ, q, α) ≤ 0, γF (σ, q, α) = 0, (11) where γ ≥ 0 is determined by the consistency condition γ Ḟ (σ, q, α) = 0, ( if The model is now completed with a damage dissipation mechanism (see [4]).The stress σ is related with an effective stress σ associated to the surface region that effectively resists the load where D is a variable bounded between 0 and 1, where D = 0 means undamaged material and D = 1 a fully broken one.Damage phenomenology has a strain equivalence principle, which states that any constitutive equation for a damaged material may be derived in the same way as for a virgin material, except that the stress is replaced by the effective stress.The total potential energy function of the virgin material can be defined as The stress can be calculated as The state of damage in the material at current time t is characterized by means of a damage criterion where r t is the current damage threshold defined by r t = max{r 0 , ω s }, s ∈ [0, t] with r 0 denoting the initial damage threshold before any loading is applied.The equations of evolution for D t and r t are expressed as where μ ≥ 0 is a damage consistency parameter which defines damage Kuhn-Tucker conditions The Lemaitre damage material model was developed by Lemaitre [13] for an isotropic linear elastic virgin material with stress-strain law as follows where D represents the extent of damage with the damage evolution law where A and B are the characteristic parameters of the material, ǫ D 0 is the initial damage threshold strain and ǭ is the effective strain defined by ǭ = with H is the Heaviside function [13].
The value of μ is determined by the damage consistency condition that reads μ = ωt [4].

Results
This section we analyze a structural problem composed of an elastic periodic material.The problem simulated is a one-dimensional uniform bar with prescribed load (F (0, t). Figure 1(a) shows the analyzed structured.On the left is applied force and on the right the end is fixed boundary condition with zero displacement.There are many possible designs for periodic structures.At this time, we are assuming the structure is divided by unit cells.The design of the unit cells of banded materials following the works of Hussein et al.[ [8], [9]].The unit cell is composed of two layer of fiber ('f ') and one of matrix ('m') materials (see Fig. 1(b)).The fiber is represent by the number '1' and the matrix with '0'.A periodic structure design with a set of N Cell unit cells.The governing equation of motion is expressed by Eqn.(1).The bar area is A = 0.000628m 2 and length is L = 0.5m.It is initially at rest and is modeled by 80 equal-length linear-displacement elements, so that, L e = 0.00625m.
In the case of periodic composite material.we study the problem for different values of N Cell , i.e., N Cell = 1, 2, 3. Note that, d = 10L e , and as we used 80 finite elements, we have a maximum of 8 unit cells, so when we consider for example N Cell = 1 at the left end, the rest part of the structure is homogeneous.Each unit cell has a different design with ten finite-element.In   2(a)) we can see a unit cell with four fiber, two matrix, and four fiber-layer, namely type A (T A).In Fig. 2(b))and (c) we define two different design for periodic material, namely T B and T C, respectively.The homogeneous material is represent by the number '2'.So, a structure with one unit cell of the type A and 80 finite-elements is represent by vector {1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 2, 2, 2, ..., 2, 2, 2}.The computational procedure for integration is given by the central-difference method.

Elastic Constitutive Model
We consider linear elastic material (6), with properties: density ρ = 7, 900 kg/m 3 and Young modulus E = 210.10 9 Pa.We investigate a homogeneous material with a constant Elastic modulus E and a periodic composite material with Elastic modulus E(x).
The dimensions of the unit cell are d f /d = δ, where d is the length of the unit cell and d f is the dimension of the layers of fiber.The material property ratios are ρ f /ρ m = 3 and E f /E m = 12 for density and elastic modulus, respectively.As stated by Hussein et al., the dispersive characteristics of a periodic composite is determined by the ratios of material properties, and not the absolute properties of the individual constituents.Consider the uniform linear-displacement element, with lumped mass, the highest frequency of this element is (ω max ) e = 2c/L e where c = E/ρ is the acoustic wave speed.Using the central-difference method for integration, we have that the scheme is conditionally stable and requires ∆t ≤ 2/(ω max ) e .So, this method require a CFL condition ∆t ≤ L e /c.Namely a stable time-step as ∆t s = L e /c, we have for this problem ∆t s = 1.212.10−6 s.The wave speeds for fiber and matrix are computed by c f = E f /ρ f and c m = E m /ρ m , respectively.A prescribed constant force is applied at the left end, F = F max .The maximum load is F max = 500N .
The equivalent homogeneous structure have same total mass and effective static stiffness computed.So, the homogeneous medium has equivalent averaged elastic modulus, Ē, and the density, ρ, of the heterogeneous unit cell.The values of Ē and ρ are computed using the rule of mixtures, by Eqns.( 22) and (23), respectively.
For unit cell of the type A (T A) we have δ = 0.8 and for types B and C (T B,T C) δ = 0.5.For T A, the speeds are c f = 1.848c and c m = 0.924c.So, we assume the fiber wave speed for the stability limit, and the time-step stable is ∆t f s = 0.5411∆t s .Then, a time-step of ∆t = 0.5∆t s = 6.06.10 −7 s is used.For this case, we have a Courant number C n = ∆t/∆t f s = 0.92 based on fiber layer.The fiber and matrix properties are: ρ f = 9, 115.4 kg/m 3 , ρ m = 3, 038.5 kg/m 3 , E f = 6.72. 10 11 Pa, E m = 5.6.10 10 Pa.
For T B and T B, the speeds are c f = 2.082c and c m = 1.041c.So, we assume the fiber wave speed for the stability limit, and the time-step stable is ∆t f s = 0.480∆t s .Then, a time-step of ∆t = 0.45∆t s = 5.454.10 −7 s is used.For this case, we have a Courant number C n = ∆t/∆t f s = 0.9375 based on fiber layer.The fiber and matrix properties are: ρ f = 11, 850 kg/m 3 , ρ m = 3, 950 kg/m 3 , E f = 1.365.10 12 Pa, E m = 1.1375.10 11Pa.
Figure 3 shows the displacement time history at x = 0 m from time t = 0 s to time t = 2, 000∆t = 1.21.10 −3 s for a homogeneous structure and for periodic structures composed by unit cells of the type T A, T B, and T C .You should be able to verify that the exact solution is a saw-tooth with amplitude 3.75.10−6 m and period approximately 0.4.10 −3 s (the time for a plane wave to make three trips from one end of the bar to the other), due to wave refraction from the built-in end for constant tip loading.The homogeneous case (Fig. 3(a)), we observe a linear variation in the autonomous conservative system.For periodic materials the waves are not linear with dissipative characteristics.
Figure 4 displays the stress time history at x = 0.253125 m (element 41).The stress at this point is zero until sufficient time has elapsed (t ∼= 0.05.10 −3 s) for the stress wave to propagate from the instantaneously loaded tip to the bar midpoint.The increase in mean compressive stress from ∼ −8.10 5 Pa to sim − 17.10 5 Pa at ∼ 0.15.10 −3 s is due to wave refraction from the built-in end.This phenomena increasing and decreasing the stress repeat periodically.Figures 4(a)-(d) present stress time-history at x = 0.253125 m for homogeneous structure, periodic composite structure with 8 cells (type T A, T B, and T C),respectively, from time t = 0 s to time t = 2, 000∆t = 1.21.10 −3 s. Figure 5 show the evolution of stress at x = 0 m for homogeneous (a) and periodic composite structure with 8 cells until t = 2, 000 time-steps, for unit cells of type T A (b), T B (c), and T C (d), respectively.
Figure 6 displays maximum displacement,u[m], versus position x[m] for a homogeneous structure.For each position x the displacement is plotted from zero to the maximum values during the simulation time 1.21.10 −3 s.We can see a linear variation from the loaded to fixed end.
Figure 7 shows maximum displacement,u[m], versus position x[m] for a periodic composite structure with 1 cell, 2 cells, 3 cells, and 8 cells, of the type T A. We observed a switch between fiber and matrix in the periodic structure.The first unit cell effects the behavior of the wave at the initial part of the structure.The second unit cell effects the subsequent region, and so on for other units cells.The homogeneous material is non-dispersive and is therefore not frequency-banded.When the elastic periodic structure is excited by a harmonic function there is a dispersive phenomena limiting the frequency.We observed a switch between fiber and matrix in the periodic structure.The first unit cell effects the behavior of the wave at the initial part of the structure.The second unit cell effects the subsequent region, and so on for other units cells.The homogeneous material is non-dispersive and is therefore not frequency-banded.terial these relations are lines.We can see for homogeneous material one line with inclination represents by E. For composite structure there are two more lines for inclination of the fiber and matrix of the unit cells (cases, (b), (c), and (d), unit cells of the type T A).In case (d), there line for homogeneous material is not present because the whole rod is compose by unit cells.

Elastic-Plastic Constitutive Model
Here, the uniform bar has a instantaneous tip loading (360, 000N ).The material is an elastic-plastic isotropic kinematic hardening model.The flow stress is σ Y = 27E7P a, for plastic modulus is K = E/3 = 7E10P a, and the kinematic hardening modulus is H = 21E9P a. Similarly the case of elastic linear model, the material property ratios for fiber and matrix are , and H f /H m = 12.As before, a time-step of ∆t = 6.06.10 −7 s is used for preserves stability.The number of temporal shows the displacement timehistory, respectively.We can see in Figure 11(a) two separate displacement waves arrive at different times.These correspond to an elastic wave traveling at speed c, which arrives first, followed by a a plastic wave traveling at speed c p = E t /ρ = 0.55c, where E t = E ×(K +H)/(E +(K +H)).Thus, at x = 0m, an elastic wave arrives at t = 0.194.10−3 s (see also Fig. 3) and a plastic wave arrives at t = 0.353.10−3 s (Fig. 11(a)), where at this time the pattern of the response changes by the composition of the elastic and plastic waves.Similarly to Figures 7 -9, Figures 12, 13, and 14 plot maximum displacement versus position during the simulation time.These figures are divided in four parts, each one for a type of periodic composite structure with 1 cell, 2 cells, 3 cells, and 8 cells, respectively.The regions for fiber and matrix in theses cells can be observed.After, these cells a homogeneous part complete the structure.Figures 11,12        Comparing these figures with homogeneous case, we can see that damage took longer to appear.We also observed that stress reaches higher values.

Conclusion
This study involves the analysis of the influence of three periodic material represented by unit cells.It is analyzed the elastic-plastic-damage dynamics of one-dimensional structures.The problem is modeled with a generic non-linear material.The Finite Element Method is applied for space discretization.Timestepping algorithm for non-linear dynamics is also presented.We present simulations for homogeneous and periodic structures.We observed that the periodic materials are more resistant than the homogeneous so take longer to start the process of damage.

Acknowledgment
The second author thanks the following Brazilian agency for financial supports CNPq -"Conselho Nacional de Desenvolvimento Científico e Tecnológico"

Figure 1 :
Figure 1: (a) Structure problem; (b) Isolate unit cell (black and white represent fiber and matrix, respectively).

Figures 8 and 9
Figures 8 and 9 show maximum displacement,u[m], versus position x[m] for unit cells of the type T B and T C, respectively.Each figure is composed by periodic composite structure with 1 cell, 2 cells, 3 cells, and 8 cells, items (a), (b), (c), and (d), respectively.

Figure 10 Figure 4 :
Figure 10  plots relations between stress and strain.For linear elastic ma-

Figure 11 (
b) shows the maximum displacement versus position during simulation time.Figure11(c) presents the phase portrait including the transient period.After this period we can see a periodic movement.Figure11(d)shows the stress versus strain in the structure during the simulation time.
, and 13 present maximum displacement for unit cells of the type T A, T B and T C, respectively.We can see different behaviors