Analysis of Elastic-Plastic Waves in a Thin-Walled Tube By a Novel Lie-Group Differential Algebraic Equations Method

Chein-Shan Liuand Satya N.Atluri

1 Introduction

The combined axial and torsional testing of thin-walled tubes is ideal for the study of constitutive equations of metals;see,for example,Nadai(1950)and Hill(1950).The thin-walled tubular specimen is usually subjected to a combination of axial loadP(t)and torqueT(t).With an appropriate feedback arrangement,the lengthZ(t)and the relative twist angle Θ(t),as well asP(t)andT(t)can serve as control variables.Thus the axial-torsional testing of a thin-walled tube may have(P,T),(Z,Θ),(P,Θ),and(Z,T)as control input pairs.Under the assumption of uniform deformation and stress distribution in the main parallel segment of the thin wall of the specimen,the four control pairs can be correspondingly related to(σxx,σxθ),(εxx,εxθ),(σxx,εxθ),and(εxx,σxθ).They are,respectively,pure stress control,strain control,and mixed controls;see,for example,Klisinski,Mroz and Runesson(1992).In a small deformation theory,the stress and strain paths for the considered tests are such that their rates are of the forms:

where the superimposed dot denotes a differentiation with respect to timet.Notice that˙σ21=˙σ12,˙ε21=˙ε12,and˙ε33=˙ε22.In cylindrical coordinates(x,θ,r),ε11=εxxis the axial strain,ε22= εθθis the hoop(or circumferential)strain,ε33= εrris the radial strain,and ε12= εxθis the shear strain in the thin wall of the tube,whereas σ11= σxxis the axial stress and σ12= σxθis the shear stress in the thin wall of the tube.

The experiments involving combined longitudinal and torsional plastic waves in a thin-walled tube of an aluminum alloy have been reported by Lipkin and Clifton(1970a).Two types of plastic waves were observed,which involved coupled longitudinal and torsional motions.In the course of an impulse loading,the metallic material behaves plastically when stress is over the yield point,and the phenomenon of stress wave propagation in elastic-plastic solids has been studied for a long time[Cristescu(1967);Nowacki(1978)].In an earlier analytical work,Clifton(1966)has considered the combined longitudinal and torsional plastic wave propagation in isotropically-hardening materials,and these two waves are called fast and slow plastic waves depending on whether the plastic wave speed is grater or less than the elastic shear wave speed.The analysis was extended by Lipkin and Clifton(1970b)to the general class of elastic-plastic materials characterized by a smooth loading surface.Goel and Malvern(1970)presented solutions to the same problem for the case of combined kinematic and isotropic hardening materials.Most theoretical results and computational techniques are based on the method of characteristics[Yang(1970);Ting(1969,1972,1973);Wu and Lin(1974);Lin and Ballmann(1993);Karagiozova(2004)];however,it is not easy to determine the speeds of characteristic lines which are dependent on the loading history for elastic-plastic material.Unfortunately,these analyses cannot be applied to the stress wave propagation when wave reflection occurs at the boundary[Tanimoto(2007)].

It is known that the elastoplastic equations are a set of differential algebraic equations(DAEs)with discontinuities,which take place at the points of transition from an elastic state to a plastic state and vice versa[Büttner and Simeon(2002);Eckert,Baaser,Gross and Scherf(2004)].Most practical problems of wave propagation in such a complicated situation cannot be solved analytically due to the shape of the boundaries,the hardening effect and the complicated loading functions.Owing to these facts,the elastic-plastic wave equation is not a pure conservation law.Although there exist many numerical schemes for hyperbolic conservation laws,which cannot be directly used to solve the elastic-plastic wave equations and some modifications are required[Giese and Fey(2003,2005)].

In the present theoretical computation of the elastic-plastic wave propagations in a solid material,we place a greater emphasis on a unifi ed point-of-view under the framework of a Lie-group differential algebraic equations(LGDAE)method,which is drastically different from the above mentioned methods.It is important to note that unlike most other numerical schemes,the present method can be easily extended to the approximations of any order in space and time and any complicated elasto-plastic model with complex hardening/softening effect[Liu(2005,2006,2007);Liu and Chang(2004)],which is extremely advantageous for a highaccuracy solution of the elastic-plastic wave propagation problems.

The remainder of this paper is structured as follows.In Section 2 we emphasize the complementary trio appearing in a system of six postulations for isotropic ally hardening plasticity,in tensor form.Then,in Section 3 we introduce the Fischer-Burmeister NCP-function,and propose a novel system of index-one DAEs for the isotropically-hardening plasticity equations,where we describe the governing equations,nonlinear complementarity problem,and momentum balance equation for a thin-walled tube problem.These equations are coined as a system of DAEs for all material points.For solving the resultant DAEs we develop a Lie-group method base donGL(n,R)in Section4,while the numerical algorithm for the elastic-plastic wave propagation problem of the isotropically-hardening plasticity is given in Section 5.In Section 6 we study the wave propagation problems of an exponentially saturated hardening material,an isotropically linear hardening material,as well as an isotropically softening material under different loading conditions,by using the corresponding LGDAE numerical algorithm.Finally,some conclusions are drawn in Section 7.

2 The isotropically-hardening model

In a small-deformation theory,the elastoplastic model for solid materials proposed by Prandtl(1924)and Reuss(1930)is re-formulated as follows:

whereIis the third-order identity tensor,the symbol tr denotes the trace of the tensor,and a dot between two tensors of the same order denotes their Euclidean inner product.The model has only two experimentally determined material constants and a material property,namely the bulk modulusK,the shear modulusG,and the shear yield strength τy,which are postulated to be

τyis supposed to be a function of Λ where˙Λ=λ,Λ(0)=0.The boldfaced ε,εe,εpand σ are,respectively,the strain,elastic strain,plastic strain,and stress tensors,all symmetric,whereas λ is a scalar.All the ε,εe,εp,σ and λ are functions of one and the same independent variable,which in most cases is taken as the ordinary time.

3 Axial-torsional plastic wave

3.1 Governing equations

For a two-dimensional axial-torsional deformation problem the following constitutive equations can be obtained from Eqs.(2)-(4)and(1),

and the yield condition reduces to

We may use Eqs.(9)and(10)to solve for σ11and σ12,respectively.However,they are coupled through Eq.(11),where τyis not a constant value,but a function of Λ with˙Λ=λ.

3.2 Nonlinear complementarity problem

A general complementarity problem is to fi nd a solutionx∈Rnof the following complementary trio system:

whereP,Q∈Rndenote vector functions.Many applications from engineering sciences,economics,game theory,etc.lead to problems of this kind;see Ferris and Pang(1997)for a survey.Most algorithms for the solution of nonlinear complementarity problem(NCP)are based on a suitable reformulation of Eq.(12)either as a system of algebraic equations,as an optimization problem,or as a fixed-point problem,etc.Refer the survey paper by Harker and Pang(1990)for some basic algorithms.

Letxbe a solution of an NCP,that is,x≥0,F(x)≥0,andxF(x)=0.Obviously,it is equivalent to the requirement thatxis a solution of the minimum problem:min(x,F(x))=0.The function φ is said to be an NCP-function:if φ :and φ(a,b)=0 if and only ifa≥0,b≥0,ab=0.

In addition to the above minimum function,there are many other NCP-functions,for example,the Fischer-Burmeister NCP-function:

The most interesting property of this merit function is that,as is easily veri fi ed:

Thus,for a general NCP of Eq.(12)in a component form we write it to be

wherePiandQiare respectively the components ofPandQ.Accordingly,we can reformulate the complementary trio in plasticity as an NCP by

It is easy to check that Eqs.(16)-(19)are index-one DAEs,which are better than those of the index-two DAEs formulation of plasticity[Büttner and Simeon(2003a,2003b);Liu(2013c)].On the other hand,in Eqs.(16)-(19)we do not need to treat the on-off switch,while for the original index-two DAEs we need to treat a two phase system and the on-off switching criteria of elasticity and plasticity.

3.3 Kinematic and momentum balance equation

Our method is easily extended to the three-dimensional elastic-plastic wave propagation problem by adjoining the momentum balance equation together:

where ρ is the material density whilewis the velocity field.However,we only consider the process of the propagation of an elastic-plastic wave in a finite length thin-walled tube with length ℓ consisting of an isotropically hardening material,subjected to dynamical boundary conditions and initial conditions.The problem will be discussed in a Lagrangian system where the axisxcoincides with the thin walled tube axis,and the originx=0 is assumed to be the left-end point of the thin-walled tube.It is also assumed that the thin-walled tube does not buckle in the course of deformation.We consider small strains and assume that the thin-walled tube density ρ does not change.The only stress components are σxx= σ and σxθ= τ,and the non-zero strain components are εxx= ε,εxθ= γ/2,and εrr= εθθ.The velocity fi eld has only the fi rst two components being nonzero and the third component is zero,i.e.,w=(u,v,0)T.Hence,the governing equations at a material pointxand a genric timetcan be written as

whereuandvare the axial and circumferential velocities.

If we want to know the field variables of strains ε(x,t)and γ(x,t)we can supplement the following kinematic equations:

4 A Lie-group DAE method

After a suitable discretization of the above equations at each material point,E-qs.(21)-(26)constitute a system of nonlinear differential algebraic equations(DAEs).Hereby,we give a general setting to treat the DAEs which govern the evolution ofn+qvariablesxi,i=1,...,nandyj,j=1,...,q,withnnonlinear ordinary differential equations(ODEs)andqnonlinear algebraic equations(NAEs):

4.1 Endowing the ODEs with a Lie-group GL(n,R)

The general linear group is a Lie group,whose manifold is an open subsetGL(n,R):={G∈ Rn×n|detG/=0}of the linear space of alln×nnon-singular matrices.Thus,GL(n,R)is ann×n-dimensional manifold.The group composition is given by the matrix multiplication.

Here we give a new form of Eq.(29)from theGL(n,R)Lie-group structure.The vector fi eldfon the right-hand side of Eq.(29)can be written as

is the coefficient matrix.The symbol⊗inu⊗ydenotes the dyadic operation ofuandy,i.e.,(u⊗y)z=y·zu.

Because the coefficient matrixAis well-defined,the Lie-group elementGgenerated from the above dynamical system(31)with˙G=AGsatisfies detG(t)/=0,such thatG∈GL(n,R).

Liu(2013a)has found the essential form in Eq.(31)for nonlinear ODEs,and developed a very effective Lie-groupGL(n,R)scheme to solve ODEs by only assuming thatf/‖x‖ is a constant vector within a small time step.Then,Liu(2013b)developed a Lie-groupGL(n,R)scheme to solve ODEs by assuming that bothf/‖x‖andx/‖x‖ are constant vectors in a small time incremental step.Liu(2013c)has developed a powerful numerical method to solve the nonlinear differential algebraic equations,based on the above Lie-groupGL(n,R)scheme.Liu(2014a)has used the Lie-group differential algebraic equations(LGDAE)method to solve the sliding mode control problem,and Liu(2014b)has solved the heat source identification problem by using the LGDAE.

4.2 An implicit GL(n,R)Lie-group scheme

Eq.(31)is a new starting point for the development of the Lie-groupGL(n,R)algorithm.In order to develop a numerical scheme from Eqs.(31)and(32),we suppose that the coefficient matrixAis constant with

being two constant vectors,which can be obtained by taking the values offandxat a suitable mid-point of¯t∈[t0=0,t],wheret≤t0+handhis a small time stepsize.

The variableyis suppose to be a constant vector in this small time interval.Thus from Eqs.(31)and(32)we have

Eq.(34)becomes

Then,from the above two equations we can derive

is a constant scalar in a small time step.From Eq.(37)it follows that

wherew0=b·x0.

Inserting Eq.(39)forw(t)into Eq.(36)and integrating the resultant equation we can obtain

where the superscript T denotes the transpose,x0is the initial value ofxat an initial timet=t0=0,and

LetGbe the coefficient matrix beforex0in Eq.(40):

which is one sort of elementary matrices.According to Liu(2013a,2013c)we can prove

such thatGis a Lie-group element ofGL(n,R).

Within a small time step we can suppose that the variablesyj,j=1,...,mare constant in the interval oftk＜t＜tk+1.Consequently,we can develop the following implicit scheme for solving the ODEs(29)whereyat thekth time step,denoted byyk,is viewed as parameters:

(i)Give 0≤θ≤1.

(ii)Give an initialx0at an initial timet=t0and a time stepsizeh.

(iii)Fork=0,1,...,we repeat the following computations to a terminal time:

wherefk:=f(xk,yk,tk).With the abovexk+1generated from an Euler step as an initial guess we can iteratively solve the newxk+1by

Ifzk+1converges according to a given stopping criterion,such that,

then go to(iii)to the next time step;otherwise,letxk+1=zk+1and go to the computations in Eqs.(45)-(46)again.In all the computations given below we will fix θ=1/2.

4.3 Newton iterative scheme for DAEs

Now,we turn our attention to the DAEs defined in Eqs.(29)and(30).Within a small time step we can suppose that the variablesyj,j=1,...,mare constant in the interval oftk＜t＜tk+1.We give an initial guess ofyj,j=1,...,m,and insert them into Eq.(29).Then we apply the above implicit scheme to find the nextxk+1,supposing thatxkis already obtained in the previous time step.Whenxk+1are available we insert them into Eq.(30),and then apply the Newton iterative scheme to solveyk+1by

till the following convergence criterion is satisfied:

Otherwise,go to Eq.(45).In the above the componentBijof the Jacobian matrixBis given by∂Fi/∂yj.

The numerical process as a combination of the Lie-group method based onGL(n,R)and the Newton method to solve the DAEs in Eqs.(29)and(30)is called the Liegroup DAE(LGDAE)method.

5 Numerical algorithm for elastic-plastic wave

Now we suppose that all the variables are discretized to be σi(t)= σ(xi,t),τi(t)=τ(xi,t),ui(t)=u(xi,t),vi(t)=v(xi,t),λi(t)= λ(xi,t)and Λi(t)= Λ(xi,t),wherexi=(i-1)∆x=(i-1)ℓ/(m-1),andmis the number of grid points.Then,we have totallyn=5mODEs andq=mconstraints:

We apply the implicitGL(n,R)scheme to solve σi, τi,ui,vi,and Λithrough E-q.(50)and then iteratively solve the unknown function λithrough Eq.(51)by the Newton iterative method.The numerical processes of this implicit LGDAE are given below,where we useQ=(σ1,...,σm,τ1,...,τm,u1,...,um,v1,...,vm,Λ1,...,Λm)Tandfis used to denote the right-hand sides of Eq.(50):

(i)Give an initial guess of λi0,for example,λi0=0.

(ii)Give an initial conditionQ0at an initial timet=0 and a time stepsize δt.

(iii)Fork=0,1,...,we repeat the following computations to a speci fi ed terminal timet=tf:

With the aboveQk+1generated from an Euler step as an initial guess we then iteratively solve the newQk+1by

then go to(iv);otherwise,letQk+1=zk+1and go to Eq.(53).

(iv)Forj=0,1,...,we repeat the following computations:

where the prime denotes the differential with respect to λi,and

where bothi(super or sub scripted)denote theith components.If λijconverges according to

then go to(iii)with λijas an initial guess of λifor the next time step;otherwise,let λij=λij+1and go to Eq.(53).

6 Examples of plastic waves

In order to assess the performance of the above numerical method of LGDAE,we consider elastic-plastic wave propagation problems of a thin-walled tube with length ℓ=0.1 m=100 mm,and the material constants areE=70,000 MPa,G=70,000/[2(1+ ν)]=26923.08 MPa with ν =0.3,and ρ =2700 kg/m3.We consider three hardening cases and one softening case,and different loading conditions and initial conditions.

6.1 Step-loading of a linearly hardening material,comparison with Clifton(1966)

First,we consider an isotropically-linearly hardening material with τybeing a linear function of Λ:

The thin-walled tube is initially at rest and pre-stressed over the initial yield stress,being subjected simultaneously at the endx=0 to constant stresses σ0and τ0.The initial conditions are,respectively,

The thin-walled tube is pre-stressed to a plastic state and the intial value of Λ can be solved from Eq.(58).For the pre-tension case we consider σ0=200 and τ0=300,while for the pre-shear case we consider σ0=350 and τ0=300.

Figure 1:For an isotropically linear hardening material showing the stress paths for(a)pre-tension stress,and(b)pre-shear stress.

We fi x ∆x=0.01/100 and solve this problem under a time stepsize δt=5×10-7sec in a time intervalt∈[0,0.0001].The stress paths for the above two cases are plotted in Figs.1(a)and 1(b)for a point atx=∆x.For the pre-tension case,from Fig.1(a)it can be seen that the material quickly tends from point a to point b with a fast plastic wave where normal stress falls to a much smaller normal stress and shear stress quickly increases from zero shear stress to a much larger shear stress,and then from point b to point c along a slow plastic wave.For the pre-shear case,from Fig.1(b)it can be seen that the material quickly tends from point a to point b,which propagates with fast plastic wave,and then from point b to point c along a slow plastic wave.The above behaviors were described by Clifton(1966)for a semi-in finite tube by using the characteristic theory,for which our results are similar.

6.2 Exponentially saturated hardening material

We let τybe a function of Λ,which is given by

For the isotropically hardening plastic wave problem we need to calculate Λiby

Ting(1973)has analyzed the plastic wave speeds of an isotropically work hardening material.

The boundary conditions are supposed to beu(0,t)=1+t+t2,v(0,t)=2(1+t+t2),and the initial conditions are given byu(x,0)=v(x,0)=0.01cos[(πx)/(2ℓ)],and σ(x,0)= τ(x,0)=50.The right-end of the tube is fixed.We consider three different terminal times:case 1 withtf=0.002 sec,case 2 withtf=0.01 sec,and case 3 withtf=0.1 sec.

We fix ∆x=0.002 for all cases.We solve case 1 under a time stepsize δt=10-6sec.In Fig.2(a)we plot the elastic and plastic zones in the plane(x,t),where the point corresponding to a plastic state is marked by a black square point.In Fig.3 we plot the waves at some sampled timest1=500δt=0.0005 sec,t2=1000δt=0.001 sec,t3=1500δt=0.0015 sec,andt4=2000δt=0.002 sec.Then we plot the time histories of wave at two pointsx1=0.002 andx2=0.08 in Fig.4.It can be seen that at the first point the material quickly tends to yielding and enters into the plastic state.At pointx2the material is elastically unloading.

We solve case 2 under a time stepsize δt=5×10-6sec.In Fig.2(b)we plot the elastic and plastic zones in the plane(x,t).In Fig.5 we plot the waves at some sampled timest1=10δt=0.00005 sec,t2=500δt=0.0025 sec,t3=1500δt=0.0075 sec,andt4=2000δt=0.01 sec.We plot the time histories of wave at two pointsx1=0.008 andx2=0.08 in Fig.6.It can be seen that at the first pointx1=0.008 the material quickly tends to plastic state and intervening by elastic unloading state.At pointx2the material is elastically unloading to compressive stress and negative shear stress.

Figure 2:The elastic-plastic zones for three different terminal times for the elastic plastic wave propagations of isotropically hardening material.

Figure 3:For case 1 the elastic-plastic wave propagation for isotropically hardening material:(a)axial stress,(b)shear stress,(c)axial velocity and(d)torsional velocity.

Figure 4:For case 1 of isotropically hardening material the elastic-plastic wave propagation with respect to time:(a)axial stress,(b)shear stress,(c)axial velocity,and(d)torsional velocity.

Figure 5:For case 2 the elastic-plastic wave propagation for isotropically hardening material:(a)axial stress,(b)shear stress,(c)axial velocity,and(d)torsional velocity.

Figure 6:For case 2 of isotropically hardening material the elastic-plastic wave propagation with respect to time:(a)axial stress,(b)shear stress,(c)axial velocity,and(d)torsional velocity.

Figure 7:For case 3 the elastic-plastic wave propagation for isotropically hardening material:(a)axial stress,(b)shear stress,(c)axial velocity,and(d)torsional velocity.

Figure 8:For case 3 of isotropically hardening material the elastic-plastic wave propagation with respect to time:(a)axial stress,(b)shear stress,(c)axial velocity,and(d)torsional velocity.

Figure 9:For case 2 of isotropically hardening material the stress wave pro fi les:(a)axial stress,and(b)shear stress.

Figure 10:For case 3 of isotropically hardening material the stress wave pro files:(a)axial stress,and(b)shear stress.

We solve case 3 under a time stepsize δt=4×10-5sec.In Fig.2(c)we plot the elastic and plastic zones in the plane(x,t).In Fig.7 we plot the waves at some sampled timest1=10δt=0.0004 sec,t2=1500δt=0.06 sec,t3=2000δt=0.085 sec,andt4=2500δt=0.1 sec.We plot the time histories of wave at two pointsx1=0.002 andx2=0.08 in Fig.8.It can be seen that at both two points the material beavior is very complex

In Fig.9 we plot the axial stress and shear stress over the plane(x,t)for case 2,while in Fig.10 we plot the axial stress and shear stress over the plane(x,t)for case 3.It can be seen that in case 3 the stress distributions are much complicated than that in case 2.

6.3 Linearly hardening material

We let τybe a linear function of Λ:

When the right-end of the tube is fixed,the left-end is loaded by impulse stresses.They simultaneously increase parabolic ally,after that,they keep constant stresses after the interval of rising timetr=n1δtwhere we fix δt=10-5sec:

The thin-walled tube is supposed to be pre-stressed to a plastic state with initial conditionsu(x,0)=0.3cos[(πx)/(2ℓ)],v(x,0)=0.25cos[(πx)/(2ℓ)],σ(x,0)= σ0=0,and τ(x,0)= τ0=150 MPa which is equal to the initial shear yield stress 150 MPa.The parameters used aren1=5,σf=300 MPa and τf=250 MPa.

We fix ∆x=0.002 and solve this problem under a time stepsize δt=10-6sec in a time intervalt∈[0,0.002].In Fig.11(a)we plot the elastic and plastic zones in the plane(x,t),where the point corresponding to a plastic state is marked by a black point.In Figs.11(b)and 11(c)we plot the stress waves at some sampled timest1=10δt=10-5sec,t2=500δt=0.0005 sec,t3=1500δt=0.0015 sec,andt4=2000δt=0.002 sec.Then we plot the time histories of stress waves at two pointsx1=0.002 andx2=0.006 in Figs.12(a)and 12(b).The stresses path is plotted in Fig.12(c)and the strains path at the second point is plotted in Fig.12(d).From Fig.12(c)it can be seen that at the first point the material quickly tends from point a to point b,which propagates with fast plastic wave,and then from point b to point c along a slow plastic wave.At the second point the material propagate along another fast plastic wave from point a to point d.In Fig.13 we plot the axial stress and shear stress over the plane(x,t),from which we can observe the stress distributions.

Figure 11:For an isotropically linear hardening material the stress wave propagation:(a)plastic zone,(a)axial stress,and(c)shear stress.

Figure 12:For an isotropically linear hardening material investigating stress wave propagation at two different points:(a)axial stress,(b)shear stress,(c)stress paths,and(d)strain path.

6.4 Isotropically-softening material

Finally we consider an isotropically-softening material with τybeing a function of Λ:

The initial yield stress is 200 MPa,and then fast tends to 50 MPa.

When the right-end of the tube is fixed,the left-end is loaded by impulse velocities.

They simultaneously increase parabolic ally,after that,they keep constant velocities after the interval of rising timetr=n1δtwhere we fix δt=10-5sec:

Figure 13:For an isotropically linear hardening material the stress wave pro files:(a)axial stress,and(b)shear stress.

Figure 14:For an isotropically softening material the stress wave propagation:(a)plastic zone,(a)axial stress,and(c)shear stress.

Figure 15:For an isotropically softening material investigating stress wave propagation at two different points:(a)axial stress,(b)shear stress,(c)stress paths,and(d)strain path.

Figure 16:For an isotropically softening material the stress wave pro files:(a)axial stress,and(b)shear stress.

where we takeu0=0.01,v0=0.05 andn1=10.The thin-walled tube is supposed to be pre-stressed to a plastic state with σ (x,0)=100 MPa,and τ(x,0)=200 MPa,and other initial conditions are zeros.

We fix ∆x=0.002 and solve this problem under a time stepsize δt=10-5sec in a time intervalt∈[0,0.05].In Fig.14(a)we plot the elastic and plastic zones in the plane(x,t),where the point corresponding to a plastic state is marked by a black point.In Figs.14(a)and 14(b)we plot the stress waves at some sampled timest1=1000δt=0.01 sec,t2=2000δt=0.02 sec,t3=3000δt=0.03 sec,andt4=5000δt=0.05 sec.Then we plot the time histories of stress waves at two pointsx1=0.002 andx2=0.006 in Figs.15(a)and 15(b).The stresses path is plotted in Fig.15(c)and the velocity path is plotted in Fig.15(d).From Fig.15(c)it can be seen that at the first point the material quickly tends from point a to point b,which propagates with fast plastic wave,and then from point b to point c and then to d along a slow plastic wave.In Fig.16 we plot the axial stress and shear stress over the plane(x,t),from which we can observe the stress distributions.

7 Conclusions

The conventional implicit index-two DAEs formulation of plasticity is successfully transformed into a set of explicit index-one DAEs in this paper,for the constitutive equations of an isotropically-hardening/softening material,which provides some advantages such that one can solve the nonlinear elastic-plastic wave propagation problems without resorting on the switching criteria and the two phase equations at each material point.A novel Lie-group differential algebraic equations(LGDAE)method was developed for the solutions of different elastoplastic wave propagation problems,and we found that the LGDAE can faithfully reveal the complex wave behavior of a finite thin-walled tube under different loading conditions and initial conditions,to produce different axial-torsional stress waves.The elastic-plastic wave for the isotropic-hardening-softening material has a complicated feature with many elastic-plastic loadings and unloadings,together with complicated moving elastic-plastic boundaries.Because the LGDAE can provide a very effective time marching solution,without using the two-phase equations and the on-off switch,it can compute the solution straightforward for much saving computational time.

Acknowledgement:Taiwan’s National Science Council project NSC-102-2221-E-002-125-MY3 and the 2011 Outstanding Research Award,granted to the fi rst author,are highly appreciated.The work of the second author is supported by a UCI/ARL colloborative research grant.

Büttner,J.;Simeon,B.(2002):Runge-Kutta methods in elastoplasticity.Appl.Numer.Math.,vol.41,pp.443-458.

Büttner,J.;Simeon,B.(2003a):Büttner J,Simeon B.Differential-algebraic equations in elasto-viscoplasticity.In Deformation and Failure in Metallic Materials,Hutter K,Baaser H,Eds.,pp.31-49,Springer,New York.

Büttner,J.;Simeon,B.(2003b):Numerical treatment of material equations with yield surfaces.In Deformation and Failure in Metallic Materials,Hutter K,Baaser H,Eds.,pp.3-30,Springer,New York.

Clifton,R.J.(1966):An analysis of combined longitudinal and torsional plastic waves in a thin-walled tube.In Proc.5th U.S.National Congress of Appl.Mech.ASME,1966,pp.465-480.

Cristescu,N.(1967):Dynamical Plasticity.North-Holland,Amsterdam.

Eckert,S.;Baaser,H.;Gross,D.;Scherf,O.(2004):A BDF2 integration method with step size control for elasto-plasticity.Compu.Mech.,vol.34,pp.377-386.

Ferris,M.C.;Pang,J.S.(1997):Engineering and economic applications of complementarity problems.SIAM Review,vol.39,pp.669-713.

Giese,G.;Fey,M.(2003):High-order simulation of the elastic-plastic wave equation in two space dimensions.Wave Motion,vol.38,pp.327-347.

Giese,G.;Fey,M.(2005):A high-resolution residual splitting scheme for the elastic-plastic wave equation.Wave Motion,vol.41,pp.29-44.

Goel,R.P.;Malvern,L.E.(1970):Biaxial plastic simple waves with combined kinematic and isotropic hardening.ASME J.Appl.Mech.,vol.37,pp.1100-1106.

Harker,P.T.;Pang,J.S.(1990):Finite-dimensional variational inequality and nonlinear complementarity problems:A survey of theory,algorithms and applications.Math.Program,vol.48,pp.161-220.

Hill,R.(1950):The mathematical theory of plasticity,Oxford University,New York.

Karagiozova,D.(2004):Dynamical buckling of elastic-plastic square tubes under axial impact–I:stress wave propagation phenomenon.Int.J.Impact Eng.,vol.30,pp.143-166.

Klisinski,M.;Mroz,Z.;Runesson,K.(1992):Structure of constitutive equations in plasticity for different choices of state and control variables.Int.J.Plasticity,vol.8,pp.221-243.

Lin,X.;Ballmann,J.(1993):A Riemann solver and a second-order Godunov method for elastic-plastic wave propagation in solids.Int.J.Impact Eng.,vol.13,pp.463-478.

Lipkin,J.;Clifton,R.J.(1970a):Plastic waves of combined stresses due to longitudinal impact of a pretorqued tube(Part 1:experimental results).ASME J.Appl.Mech.,vol.37,pp.1107-1112.

Lipkin,J.;Clifton,R.J.(1970b):Plastic waves of combined stresses due to longitudinal impact of a pretorqued tube(Part 2:comparison of theory with experiment).ASME J.Appl.Mech.,vol.37,pp.1113-1120.

Liu,C.-S.(2005):Computational applications of the Poincar´e group on the elastoplasticity with kinematic hardening.CMES:Computer Modeling in Engineering&Sciences,vol.8,pp.231-258.

Liu,C.-S.(2006):Computing Prager’s kinematic hardening mixed-control equations in a pseudo-Riemann manifold.CMES:Computer Modeling in Engineering&Sciences,vol.12,pp.161-179.

Liu,C.-S.(2007):Five different formulations of the fi nite strain perfectly plastic equations.CMES:Computer Modeling in Engineering&Sciences,vol.17,pp.73-93.

Liu,C.-S.(2013a):A method of Lie-symmetryGL(n,R)for solving non-linear dynamical systems.Int.J.Non-Linear Mech.,vol.52,pp.85-95.

Liu,C.-S.(2013b):A state feedback controller used to solve an ill-posed linear system by aGL(n,R)iterative algorithm.Commu.Numer.Anal.,vol.2013,Article ID cna-00181,22 pages.

Liu,C.-S.(2013c):Solving nonlinear differential algebraic equations by an implicitGL(n,R)Lie-group method.J.Appl.Math.,vol.2013,ID 987905,8 pages.

Liu,C.-S.(2014a):A new sliding control strategy for nonlinear system solved by the Lie-group differential algebraic equation method.Commun.Nonlinear Sci.Numer.Simulat.,vol.19,pp.2012-2038.

Liu,C.-S.(2014b):Lie-group differential algebraic equations method to recover heat source in a Cauchy problem with analytic continuation data.Int.J.Heat Mass Transfer,vol.78,pp.538-547.

Liu,C.-S.;Chang,C.W.(2004):Lie group symmetry applied to the computation of convex plasticity constitutive equation.CMES:Computer Modeling in Engineering&Sciences,vol.6,pp.277-294.

Nadai,A.(1950):Theory of Flow and Fracture of Solids.McGraw-Hill,New York.

Nowacki,W.K.(1978):Stress Waves in Non-elastic Solids.Pergamon Press,New York.

Tanimoto,N.(2007):An analysis of combined longitudinal and torsional elasticplastic-viscoplastic waves in a thin-walled tube.J.Solid Mech Mater.Eng.,vol.1,pp.1112-1127.

Ting,T.C.T.(1969):Interaction of shock waves due to combined two shear loadings.Int.J.Solids Struct.,vol.5,pp.415-435.

Ting,T.C.T.(1972):The initiation of combined stress waves in a thin-walled tube due to impact loadings.Int.J.Solids Struct.,vol.8,pp.269-293.

Ting,T.C.T.(1973):Plastic wave speeds in isotropically work-hardening materials.ASME J.Appl.Mech.,vol.40,pp.1045-1049.

Wu,H.C.;Lin,H.C.(1974):Combined plastic waves in a thin-walled tube.Int.J.Solids Struct.,vol.10,pp.903-917.

Yang,C.Y.(1970):Strain hardening effect on unloading spherical waves in an elastic plastic medium.Int.J.Solids Struct.,vol.6,pp.757-772.