- Email: [email protected]
Steady state creep analysis of a cracked body using the superposition method
Nuclear Engineering and Design 54 (1979) 79-89 © North-Holland Publishing Company
STEADY STATE CREEP ANALYSIS OF A CRACKED BODY USING THE SUPERPOSITION METHOD
G. YAGAWA Department of Nuclear Engineering, University of Tokyo, Hongo, Bunkyoku, Tokyo, Japan and N. MIYAZAKI Japan Atomic Energy Research Institute, Tokai Research Establishment, Tokaimura, Ibaraki, Japan Received 16 February 1979
The superposition method of the finite element and analytical displacement models is formulated to be applicable to the steady state creep analysis in the plane stress condition. The effectiveness of the method is demonstrated by solving several numerical examples which include the steady state creep analysis of a cracked plate under tension force. 1. Introduction
The finite element method has been widely used in various kinds of stress analyses, including the singularity analysis of cracked structure. It is well known, however, that extremely fine mesh discretization is necessary in order to obtain the accurate stress intensity factor (K-value) of the structure with crack if the conventional finite element method is used. Therefore, various techniques have been developed in the K-value analysis, namely the singular element method [ 1 - 6 ] , the superposition method [7,8] and so on [9-11]. In nonlinear analysis, on the other hand, the computer running time is considerable because of unavoidable iterative calculation, so some kind of effective finite element analysis is desirable in this case. Recently, the superposition method has been applied to the elastic plastic problem [12], where the conventional finite element method is combined with the analytical displacement function by using the incremental form of the principle of virtual work and several numerical examples are solved. A generalization of the superposition method has been also attempted with the use of the Lagrange multiplier technique [13], In the present paper, the superposition method is formulated to be applicable to the plane stress steady state creep analysis and several numerical examples are sohted in order to show the effectiveness of the superposition method. These numerical examples include the steady-:~tate creep problem of a cracked plate with the r -1/(n+1) singularity of stress around the crack tip. This problem is important with regard to fracture mechanics under creep condition and to some kinds of elastoplastic fracture mechanics as well.
2. Method of analysis In this paper *, we employ the Norton type equation as the uniaxial expression of steady state creep strain rate dc , i.e., de = A o
~
(1)
.
* Indical notation and usual summation convention with repeated indices are used to describe equations. 79
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked body
80
Here, A and n are material constants and o represents the uniaxial stress. Under the assumption of incompressibility of creep strain, eq. (1) can be extended into multiaxial state-as follows
~. = ~--~'~ 3 ~ _ . - I , ~ m/.
(2)
Here, ~. is the creep strain rate tensor, oe and S o are equivalent and deviatoric stresses, respectively, which are defined as 3 1/2 a~-- - (~oo%.) ,
(3)
1
(4)
Si/ = oi] - ~ % ~ i / ,
where cri/is the stress tensor and 6q the Kronecker delta. The governing equations of the steady state creep problem can be written in small displacement theory as follows
(5)
in V,
oi],/=0
4=½(a~,;+h;,;)
(6)
in V,
~. = %,-1 S 0. in V,
(7)
hi = ai
on Su,
(8)
on S o ,
(9)
and Ti = 7"i
with
(lo)
on S .
Ti = Gqni
Here, ai is the velocity, [-] implies a prescribed quantity and [ ] ,i represents the differentiation with respect to orthogonal Cartesian coordinates x i. V is the volume occupied by the body, Su and $a are the parts of surface S of body V, in which the velocity ui and the surface traction Ti are prescribed, respectively, ni is the unit normal drawn outward on S. Body and inertia forces have been neglected for simplicity. By using the Hoff analogy, the governing equations of steady state creep problem, eqs. (5)-(10), are converted into the following nonlinear elasticity problem uii,] = 0
in V,
1
ei/=~(ui, i + u i , i) -
3 A ~n--l~,
el/- ~e ui =ui
(11) in V,
(12)
in V
"i/
(13)
onSu,
(14)
onSo,
(15)
and
Ti=T; with
(16)
Ti = oi/ni .
As for a variational formulation equivalent to the above nonlinear elasticity problem, we consider the modified principle of virtual work as follows
faijSeiidV- f >isu,dS- f s[x~(ui-~,)l d S = 0 , V
Sa
Su
(17)
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked b o d y
81
where ;~i is the Lagrange multiplier defined on Su. The vanishing of eq. (17) for arbitrary parameters will provide eqs. (11), (14) and (15) as the Euler equations. Introducing eq. (4) and incompressibility condition of creep strain eli = 0 into eq. (17) leads to the following equation
fsij6e q d V - f Ti6ui dS - f 8 [~/(ui- ~i)] dS = 0 . V
SO
(18)
Su _
/1
By using eq. (13) together with the relation e e - A %, we obtain the following equation
f2A-'/"e('-")/"e,jSe,j dV- f V
~,Su,
So
dS-
_f 6 [;~,(ui - ~/)1 (is = 0 .
(19)
Su
Here, ee is equivalent strain defined by t-2e..e..~l/2 q tll
6e = k3
(20)
•
Next, we assume that the displacement ui in eq. (19) is expressed as follows
(21)
ui = ai + uT.
Here, fii denotes the displacement of the conventional finite element model. On the other hand, u~ is the analytical solution which represents singularity of stress or special deformation of higher mode. Both functions t~i and u[ can be written as follows
ai = Niradm ,
(22)
ui =Nimdm,
(23)
J~rim and N/*m are shape functions, and d m and dm are the nodal displacement and the generalized coordinates, respectively. The Lagrange multiplier ?~i can be written by using the generalized coordinates qm as follows
where
Xi = M i m q m
.
(24)
Introducing eqs. (22) and (23) into eq. (12), strain components can be expressed as follows
eii= l~milbrn + Bmi/dm ,
(25)
where B m i ] = ~(Nim,/ 1 ^
+ N/re,i),
B* 1 (•r* + mr* "t mij = ~VVimd ~Vjm,i)
(26)
,
(27)
Using eqs. (20)--(27) in eq. (19) with the consideration of arbitrariness of the independent parameters, we obtain the following matrix equation
G 11 a= sym.
a2~lld* t = F2 . ass] ~,q ) Fs
(28)
where
GI~t ) = fCBkt/B,q dV,
(29a)
V 13 a(kl) = -- fNi~'Iit dS, 8
Gt20 = fCJO~ilB~/dV,
(29b)
V
(29c1
a~2,) = f C B ~*j B u*/ d V , g
(29d1
82
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked body
23 = G(kt)
f N~kM i t d S ,
(29e)
a~al) = 0 ,
(290
Su F~)= =
f ~,~ dS,
(29g)
=
So
f -
T ' ~ i k dS
Scr
"
,
(29h) .
(29i)
F~k) = - f -uiMik dS , Su
with
c : ~)("+')/~" A -'/" [(~,pC~r + Br'~e~;)(gp.a, + Bbqa3] ('-")/~"
(30)
As ~,/in eq. (19) is physically equivalent to the reaction force Ti on S u, we can also use the following equation
fi.-'/"-('-")/"-~ ~e ~ij v~o ~i1
d V -
V
f ~,~u,dS f s[r,(u, -
So
-
Ui)] dS=O.
(31)
Su
In this equation, Ti is no longer independent parameter, but calculated by the use of eq. (16). It should be noted here that the displacement ui given in eq. (21) must be compatible on the geometrical boundary Su or on the boundary between the superposed part and the remaining finite element part in the whole structure if the superposed domain is partial in the whole structure [ 13 ]. In this regard, the modified principle of virtual work as shown in eq. (19) or eq. (31) is particularly convenient. In the case of the plane stress condition, constitutive eq. (13) can be written as * (32)
oil = gijklekl
or explicitly
ol(,,q [sym,
to12)
lJ~e12)
Using eqs. (32) and (12) in eq. (16), T i in eq. (31) can be calculated from displacement as follows (34)
T i = EiiklUk, l n i .
Introducing eq. (34) into eq. (31) and taking the similar procedure as in deriving eq. (28), we have the following matrix equation
sym.
G-22-1(d*~= (E-~J '
(35/
where
-11 G oa)
= fcj~i~u i dV - f Eiipqnl(Bkp¢~ri, +,SlpqfiCik)dS , v
(36a)
su
* In the case of three dimensional and plane strain creep problems, eq. (13) can not be written as eq. (32) because of the incompressibility of creep strain. Therefore, the formulation presented in this paper is restricted to the plane stress condition. The special variational principle for the incompressible solid may be necessary to solve the three dimensional and the plane strain problems [14].
83
G. Yagawa, N, Miyazaki / Steady state creep analysis o f a cracked body
f
= f ca ,j ,tj dr'V
G{~,) = /CB;iiB,~ i d V -
/ E#pqni(B~pctNi*l+ BI~qNi*k) dS
v
Ftk) =
f Tdg,k dS So
us,
+
(36b)
Su
(36c)
,
Su -
f ~EiipqniBk~ dS ,
(36d)
Su
(36e)
sa
Su
It is noted that eqs. (28) and (35) are nonlinear, since the coefficient matrices in these equations are the functions o f d and d °. Direct iteration method is here employed to solve these equations.
3. Numerical examples The computer program for plane stress steady state creep analysis was developed using the theory described above. The finite element utilized here is of conventional type, i.e., constant strain triangle. 3.1. Analysis o f can tilever b earn
The first numerical example is the steady state creep bending analysis of a cantilever beam with one edge clamped and the other subjected to shear force as shown in fig. 1. It is well known that the conventional finite element model with constant strain provides very poor results even in the case of simple elasticity problem of this
I
I
I
I
I
b LIJ ~Q::
0.3
Z o_ I---~ O+ WE: ._1 -J EL = LLI 0..
0
o
0.2 -
P Z
/
o
- - L
z
~.
~ ~ ~
h=THICKNES~ ;
L
O Z
I-~ 13_
0.1
o --: 0
~ 0.0
ELEMENTARY BEAM THEO(~ : PRESENT METHOD
13 : CONVENTIONAL FINITE ELEMENT METHOD I I I 2 3 4
I t
z,u,
CREEP EXPONENT
Fig. l. Element mesh of cantilever beam with shear force P. Fig. 2. Variation of nondimensional deflection rate with creep exponent n.
I 5
n
84
G. Yagawa,N. Miyazaki / Steady state creep analysis of a cracked body
structure. In this relation, we apply the present method to the creep analysis of the structure. The finite element mesh discretization is given in the same figure. The analytical solution introduced into eq. (21) is assumed as follows
+
{ . } | [ - - x 2 , - X a . . . . , -xn+2 ~uz)
n+IZ]lql ]ii
(37)
,
I Iqn+l
where Ux and Uz represent the analytical displacements in the x and z directions, respectively, and q 1, q2, ..., qn+l are the generalized coordinates. Provided the Bernoulli assumption is valid, eq. (37) yields the exact solution for the present problem [ 15 ]. In other words, the role of the finite element model ai in this problem is relatively small in eq. (21). Fig. 2 shows the nondimensional deflection rate at point A in fig. 1 versus creep exponent n, where B =A/(ln) n with In = 2hnb2+(1/n)/(2n + 1). In this figure, the elementary beam theory refers to the book by Odqvist and Hult [15] and the same mesh as shown in fig. 1 is used for the solution of the conventional finite element method. The present method is seen from the figure to coincide fairly well with the elementary beam theory for small values of n. Some discrepancy between these two theories for larger values of n might be explained from the fact that the nonlinearity of deformation along the z-direction becomes significant as the creep exponent n increases. 3.2. Analysis o f square plate with circular hole
The second numerical example is the steady state creep analysis of a square plate with a central hole under uniform tension ao as shown in fig. 3a. The analytical solution is introduced in the superposition method as follows ' . - ( rrx x a - 3 r y 2] x r ~ ( x a - 3 x y 2) ' ux + r4 ] , r2 , r6 q1
'= ao:o ~ * Uy
q2 , (r-S 4
Ya-3_x2y~ r4 ]
y ,r 2 ,
r~(ya-3x2y) r6
(38)
qa
where tc = (3 - v)/(3 + v), G =E/2(1 + v) and r = ~ ) . Poisson's ratio v is taken equal to 0.5 because of incompressibility of creep strain. Eq. (38) is obtained using the theory of two dimensional elasticity with the complex stress function of infinite plate with a circular hole under uniform tension force. The detail of the derivation of the equation is given in ref. [12]. In this calculation, the diameter of the circular hole 2to is taken equal to one third of the plate width, and the finite element mesh is given as depicted in fig. 3b. In fig. 4, the variation of stress concentration factor oy/oo at point A in fig. 3a is shown with the reciprocal of creep exponent n. From the figure, the present solution is found to agree well with Hata's analytical one [16], which is obtained using the stress function method. Fig. 5 shows the variation of the ratio of the analytical displacement in the x direction Ux to the total one Ux + Ux with creep exponent n at point A in fig. 3a. The trend obtained as shown in this figure is reasonable in the sense that the role of the analytical part of the solution gradually decreases with n because the analytical solution utilized here merely represents the displacement field of the simple elasticity problem around a circular hole in an infinite plate. It is interesting, however, to note that the ratio is more than 60% even in n = 5 as shown in the figure. This may imply that, even in the creep exponent n as large as 5, the mixed use of the analytical solution of the simple elasticity problem is still effective to express such a high strain rate around a circular hole as shown in fig. 6.
G. Yagawa,N. Miyazaki / Steady state creep analysis o[ a cracked body
85
y, Uy
Y L
) 2L-
- -
--X,Ux
a o
b
;X
0o - ro hA Fig. 3. Square plate with a circular hole and its element mesh.
3.3. Analysis o f square plate with crack
The numerical example presented here is the steady state creep analysis of a square plate with a center crack whose length is half of the plate width and which is subjected to a uniform tension stress uo as shown in fig. 7. First of all, we consider an elastic-plastic problem whose material model is given as follows e = ale
(o< ay),
(39a)
e=Ao n
(a>~ov),
(39b)
where E and o v are Young's modulus and yield stress, respectively. Rice and Rosengren [ 17], and Hutchinson
i
n--
4.0
o I(_)
i
i
i
i
PRESENT METHOD / o SOLUTION BY HA~A//
Z
_o ~
3.0
Z U.I
g (.)
o
I--
/oJ
2.0
/
I.C *'~ 08 ÷
0
0
0 0 0
o
o Oa O2 PRESENTMETHOD
1"8.0
i
I
I
I
I
02
0.4
0.6
0.8
1.0
INVERSE OF CREEP EXPONENT 1/n
OO
~
1
I
2
I
3
I
'4
5
CREEP EXPONENT n
Fig. 4. Variation of stress concentration factor with inverse of creep exponent 1/m Fig. 5. Variation of ratio Ux/(U x + u x) with creep exponent n.
86
G. Yagawa, N. Miyazaki / Steady state creep analysis of a cracked body
, m 0 t-
,
,
]
o
r
i i
Z
o
Bi
~ o m o
J
20.C
o
§
i
L
o
~
10.(3
o
7
(to PRESENT
Pr
~
I I i
30£
oo
METHOD
i
i
t
i
2
3
4
5
CREEP
EXPONENT
n
// 0
CRACK
Fig. 6. Variation of strain rate concentration factor with creep exponent n. Fig. 7. Square plate with a center crack. [18,19] have determined the character of the plastic crack tip singularity for the material model given in eq. (39) under the conditions of the plane stress and the plane strain. From their studies, the asymptotic stress ai/, the equivalent stress ae, the strain el~, and the displacement ui can be written as follows ai/ = Kor-l/(n+l)~q(O) ,
(40a)
Oe = Kor-1/(n+t)~e(O) ,
(40b)
eq = Ker-n/(n+l)~q(O) ,
(40c)
ui = Kerl/(n+l)ffi(O) ,
(40d)
where Ke = K~ and (r, 0) are the polar coordinates centered at the crack tip. Plastic stress intensity factor Ka is defined as the amplitude of ~e(0) which is normalized by the maximum value of ~e(0). The function ~i/(O), ae(O), "di/(O) and fi'i(0) are determined through the numerical solution of nonlinear ordinary differential equation [18,19] However, the explicit form of these functions has not been obtained. According to the Hoff analogy, the elastic-plastic problem with eq. (39) is equivalent to the following steady state creep problem: when the creep strain rate in a steady state is given by ~c =Aon, the distributions of stress, equivalent stress, strain rate and displacement rate around the crack tip are expressed as follows
oii=Kol'-l/(n+l)oil(O
),
(41a)
oe=Kor-t/(n+l)~e(O ) .
(41b)
~ =Rer-nl(n+l)~q(O) ,
(41c)
h i =I(erll(n+l)ffi(O) ,
(41d)
Therefore, instead of solving the steady state creep problem directly, we may solve the elastic-plastic problem with the stress-strain relation as given in eq. (39). In view of the symmetry of the structure as shown in fig. 7 and the character of eq. (40d), the analytical solution superposed on the finite element solution in this problem can be written as follows m
u x* = ~
m
amr '1("+1) cos ½mO ,
(42a)
u~, = m=, ~ bmr 1/(n+I) sin lmO
(42b)
m=l
where am and bin, m = 1, 2, ..., M, are the generalized coordinates, and M is the number of truncation which is set
G. Yagawa,N. Miyazaki / Steady state creep analysis of a cracked body
87
a
CRACK
E
E 0 F
CRACK
£)
F
Fig. 8. Element mesh of square plate with a center crack.
15ia co co UJ
I
'
I
,
0
1.5
, co
~
o')
t.-(.~
z iJ.J
zU.I
F, z
i
I
0.5
5O3
ub
i
I.C
I--co d <
._1 <[ z
b
OC
Z
-0~
0
PRESENT METHOD n=3 I
-1.0
[3=7
I
0
-OE
L
~
ANGLE
I
-1(
~/2
0
e
I
I
:~/2
.~
ANGLE
e
Fig. 9. O-variation of nondimensional stresses around crack tip (present method).
U3
co
1.5
13/ L~
1.0
...I
< z co z
a
Z 0
- 10
rr t-o3
Z 1.1.1
OC -05
CO CO UJ
,< Z 0 co
05
:~
z
'
HUTCHINSON n=3 I
i
i
-,~/2
ANGLE
e
I
1.0 05 O0
HUTCHINSON ~ -0.5 SLIP LINE THEORY d3 Z n=oo i -t,9.~ 0 z - I.C 0 ~/2" ANGLE
O
Fig. 10.0-variation of nondimensional stresses around crack tip (Hutchinson).
. I
, .K
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked body
88 1 5
(/3 co LIJ r'r I-o9 _j z 0 (,o z
,
i
100
i
0 Z
PRESENT METHOD
10 03 03 LIJ n," I-O0
05
-05
I0
~/
-'SLOPE=-I/2
J
O.C
Z
re 123 Z
•
50
o_ 09 Z U.I
- CONVENTIONAL FINITE ELEMENT METHOD n=3 i
-1C 0
i
ANGLE
0 Z
i
ix/2
J[ e
1
0.001
0'005 0101
0105 011
NONDIMENSIONAL DISTANCE FROM CRACK TIP r / o
Fig. 11.0-variation of nondimensional stresses around crack tip (conventional finite element method). Fig. 12. Logarithmic representation of nondimensional stress distributions near crack tip along x-axis.
to be 10 in this calculation. Fig. 8 is the finite element mesh discretization, in which figs. 8a and 8b show, respectively, a quarter of the plate and the enlargement of the crack tip region. The numbers of the elements and the nodes in this calculation are 242 and 139, respectively. The distributions of ~i/and °e calculated for n = 3 and n = 7 are, respectively, given in figs. 9a and 9b. These figures are depicted by using the stress at the center of gravity of the elements around the crack tip. Fig. 10 shows the distributions of oij and ue for n = 3 and n = o. obtained analytically by Hutchinson. The results shown in figs. 9a and 10a for n = 3 coincide well with each other except near the crack surface 0 = 7r, and the results given in Fig. 9b for n = 7 compare favourably at least qualitatively with those given in fig. 10b for n = ~o based on slip line theory. Compared with the results of figs. 9a and 10a, the results using the conventional finite element method in fig. 11 seems to be cosiderably poor. The reason for the discrepancy for the solutions near the crack surface between the superposition method and the analytical one may attribute rather coarse mesh in the former near 0 = zr as shown in fig. 8a in view of very complex stress distributions around there. The same reason can apply to the numerical oscillation of stresses in the conventional finite element method as shown in fig. 11. Fig. 12 shows the logarithmic representation of the nondimensional stress ay/ao versus nondimensional distance from the crack tip r/a, where the straight lines show the inclinations o f - 1 / ( n + 1) originating from the stress estimate points nearest the crack tip. It can be seen from the figure that the stress singularities near the crack tip are well furnished in the superposition method.
4. Conclusion A superposition method of the analytical and the finite element solutions has been presented for analysis of the plane stress steady state creep problems. Effectiveness of the proposed method has been shown by solving several numerical examples. It may be concluded that the superposition method proposed in this paper is a powerful numerical technique to the steady state creep problems, particularly to the problem of a cracked plate with the singularity around the crack tip.
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked body
89
References [1] [2] [3] [4] [5] [6] [7] [8]
E. Byskov, Int. J. Fracture Mech. 6 (1970) 159. A.K. Rao, I.S. Raju and A.V. Krishna Murthy, Int. J. Numerical Methods Eng. 3 (1971) 389. S.N. Atluri, A.S. K0bayashi and M. Nakagaki, Int. J. Fracture 11 (1975) 257. T.H. Pian, P. Tong and C.H, Luk, in: Proc. 3rd Conf. Matrix Methods in Structural Mechanics (1973) p. 661. D.M. Tracey, Eng. Fracture Mech. 3 (1971) 255. R.D. Henshell and K.G. Shaw, Int. J. Numerical Methods Eng. 9 (1975) 495. Y. Yamamoto and N. Tokuda, Int. J. Numerical Methods Eng. 6 (1973) 427. G. Yagawa, T. Nishioka, Y. Ando and N. Ogura, in! Computational Fracture Mechanics, eds. E.F. Rybicki and S.E. Benzley (ASME, New York, 1975) p. 15. [9] R.H. Gallagher, in: Numerical Methods in Fracture Mechanics, eds. A.R. Luxmoore and D.RJ. Owen (University College Swansea, Swansea, 1978) p. 1. [10] S.E. Benzley and D.M. Parks, in: Structural Mechanics Computer Programs, eds. W. Pilkey, K. Saczalski and H. Schaeffer, (University Press of Virginia, 1974) p. 81. [11] M.C. Apostal, S. Jordan and P.V. Marcal, EPRI SR-22 Special Report (1975). [12] G. Yagawa, T. Nishioka and Y. Ando, Nucl. Eng. Des. 34 (1975) 247. [131 G. Yagawa, N. Miyazaki and Y. Ando, Int. J. Numerical Methods Eng. 11 (1977) 1107. [14] L.R. Herrmann, AIAA J. 3 (1965) 1896. [ 15] F.K.G. Odqvist and J. Hult, Kriechfestigkeit Metalliseher Werkstoffe (Springer-Verlag, Berlin, 1962). [16] T. Hata, Nucl. Eng. Des. 32 (1975) 325. [17] J.R. Rice and G.F. Rosengren, J. Mech. Phys. Solids 16 (1968) 1. [18] J.W. Hutchinson, J. Mech. Phys. Solids 16 (1968) 13. [19] J.W. Hutchinson, J. Mech. Phys. Solids 16 (1968) 337.
STEADY STATE CREEP ANALYSIS OF A CRACKED BODY USING THE SUPERPOSITION METHOD
G. YAGAWA Department of Nuclear Engineering, University of Tokyo, Hongo, Bunkyoku, Tokyo, Japan and N. MIYAZAKI Japan Atomic Energy Research Institute, Tokai Research Establishment, Tokaimura, Ibaraki, Japan Received 16 February 1979
The superposition method of the finite element and analytical displacement models is formulated to be applicable to the steady state creep analysis in the plane stress condition. The effectiveness of the method is demonstrated by solving several numerical examples which include the steady state creep analysis of a cracked plate under tension force. 1. Introduction
The finite element method has been widely used in various kinds of stress analyses, including the singularity analysis of cracked structure. It is well known, however, that extremely fine mesh discretization is necessary in order to obtain the accurate stress intensity factor (K-value) of the structure with crack if the conventional finite element method is used. Therefore, various techniques have been developed in the K-value analysis, namely the singular element method [ 1 - 6 ] , the superposition method [7,8] and so on [9-11]. In nonlinear analysis, on the other hand, the computer running time is considerable because of unavoidable iterative calculation, so some kind of effective finite element analysis is desirable in this case. Recently, the superposition method has been applied to the elastic plastic problem [12], where the conventional finite element method is combined with the analytical displacement function by using the incremental form of the principle of virtual work and several numerical examples are solved. A generalization of the superposition method has been also attempted with the use of the Lagrange multiplier technique [13], In the present paper, the superposition method is formulated to be applicable to the plane stress steady state creep analysis and several numerical examples are sohted in order to show the effectiveness of the superposition method. These numerical examples include the steady-:~tate creep problem of a cracked plate with the r -1/(n+1) singularity of stress around the crack tip. This problem is important with regard to fracture mechanics under creep condition and to some kinds of elastoplastic fracture mechanics as well.
2. Method of analysis In this paper *, we employ the Norton type equation as the uniaxial expression of steady state creep strain rate dc , i.e., de = A o
~
(1)
.
* Indical notation and usual summation convention with repeated indices are used to describe equations. 79
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked body
80
Here, A and n are material constants and o represents the uniaxial stress. Under the assumption of incompressibility of creep strain, eq. (1) can be extended into multiaxial state-as follows
~. = ~--~'~ 3 ~ _ . - I , ~ m/.
(2)
Here, ~. is the creep strain rate tensor, oe and S o are equivalent and deviatoric stresses, respectively, which are defined as 3 1/2 a~-- - (~oo%.) ,
(3)
1
(4)
Si/ = oi] - ~ % ~ i / ,
where cri/is the stress tensor and 6q the Kronecker delta. The governing equations of the steady state creep problem can be written in small displacement theory as follows
(5)
in V,
oi],/=0
4=½(a~,;+h;,;)
(6)
in V,
~. = %,-1 S 0. in V,
(7)
hi = ai
on Su,
(8)
on S o ,
(9)
and Ti = 7"i
with
(lo)
on S .
Ti = Gqni
Here, ai is the velocity, [-] implies a prescribed quantity and [ ] ,i represents the differentiation with respect to orthogonal Cartesian coordinates x i. V is the volume occupied by the body, Su and $a are the parts of surface S of body V, in which the velocity ui and the surface traction Ti are prescribed, respectively, ni is the unit normal drawn outward on S. Body and inertia forces have been neglected for simplicity. By using the Hoff analogy, the governing equations of steady state creep problem, eqs. (5)-(10), are converted into the following nonlinear elasticity problem uii,] = 0
in V,
1
ei/=~(ui, i + u i , i) -
3 A ~n--l~,
el/- ~e ui =ui
(11) in V,
(12)
in V
"i/
(13)
onSu,
(14)
onSo,
(15)
and
Ti=T; with
(16)
Ti = oi/ni .
As for a variational formulation equivalent to the above nonlinear elasticity problem, we consider the modified principle of virtual work as follows
faijSeiidV- f >isu,dS- f s[x~(ui-~,)l d S = 0 , V
Sa
Su
(17)
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked b o d y
81
where ;~i is the Lagrange multiplier defined on Su. The vanishing of eq. (17) for arbitrary parameters will provide eqs. (11), (14) and (15) as the Euler equations. Introducing eq. (4) and incompressibility condition of creep strain eli = 0 into eq. (17) leads to the following equation
fsij6e q d V - f Ti6ui dS - f 8 [~/(ui- ~i)] dS = 0 . V
SO
(18)
Su _
/1
By using eq. (13) together with the relation e e - A %, we obtain the following equation
f2A-'/"e('-")/"e,jSe,j dV- f V
~,Su,
So
dS-
_f 6 [;~,(ui - ~/)1 (is = 0 .
(19)
Su
Here, ee is equivalent strain defined by t-2e..e..~l/2 q tll
6e = k3
(20)
•
Next, we assume that the displacement ui in eq. (19) is expressed as follows
(21)
ui = ai + uT.
Here, fii denotes the displacement of the conventional finite element model. On the other hand, u~ is the analytical solution which represents singularity of stress or special deformation of higher mode. Both functions t~i and u[ can be written as follows
ai = Niradm ,
(22)
ui =Nimdm,
(23)
J~rim and N/*m are shape functions, and d m and dm are the nodal displacement and the generalized coordinates, respectively. The Lagrange multiplier ?~i can be written by using the generalized coordinates qm as follows
where
Xi = M i m q m
.
(24)
Introducing eqs. (22) and (23) into eq. (12), strain components can be expressed as follows
eii= l~milbrn + Bmi/dm ,
(25)
where B m i ] = ~(Nim,/ 1 ^
+ N/re,i),
B* 1 (•r* + mr* "t mij = ~VVimd ~Vjm,i)
(26)
,
(27)
Using eqs. (20)--(27) in eq. (19) with the consideration of arbitrariness of the independent parameters, we obtain the following matrix equation
G 11 a= sym.
a2~lld* t = F2 . ass] ~,q ) Fs
(28)
where
GI~t ) = fCBkt/B,q dV,
(29a)
V 13 a(kl) = -- fNi~'Iit dS, 8
Gt20 = fCJO~ilB~/dV,
(29b)
V
(29c1
a~2,) = f C B ~*j B u*/ d V , g
(29d1
82
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked body
23 = G(kt)
f N~kM i t d S ,
(29e)
a~al) = 0 ,
(290
Su F~)= =
f ~,~ dS,
(29g)
=
So
f -
T ' ~ i k dS
Scr
"
,
(29h) .
(29i)
F~k) = - f -uiMik dS , Su
with
c : ~)("+')/~" A -'/" [(~,pC~r + Br'~e~;)(gp.a, + Bbqa3] ('-")/~"
(30)
As ~,/in eq. (19) is physically equivalent to the reaction force Ti on S u, we can also use the following equation
fi.-'/"-('-")/"-~ ~e ~ij v~o ~i1
d V -
V
f ~,~u,dS f s[r,(u, -
So
-
Ui)] dS=O.
(31)
Su
In this equation, Ti is no longer independent parameter, but calculated by the use of eq. (16). It should be noted here that the displacement ui given in eq. (21) must be compatible on the geometrical boundary Su or on the boundary between the superposed part and the remaining finite element part in the whole structure if the superposed domain is partial in the whole structure [ 13 ]. In this regard, the modified principle of virtual work as shown in eq. (19) or eq. (31) is particularly convenient. In the case of the plane stress condition, constitutive eq. (13) can be written as * (32)
oil = gijklekl
or explicitly
ol(,,q [sym,
to12)
lJ~e12)
Using eqs. (32) and (12) in eq. (16), T i in eq. (31) can be calculated from displacement as follows (34)
T i = EiiklUk, l n i .
Introducing eq. (34) into eq. (31) and taking the similar procedure as in deriving eq. (28), we have the following matrix equation
sym.
G-22-1(d*~= (E-~J '
(35/
where
-11 G oa)
= fcj~i~u i dV - f Eiipqnl(Bkp¢~ri, +,SlpqfiCik)dS , v
(36a)
su
* In the case of three dimensional and plane strain creep problems, eq. (13) can not be written as eq. (32) because of the incompressibility of creep strain. Therefore, the formulation presented in this paper is restricted to the plane stress condition. The special variational principle for the incompressible solid may be necessary to solve the three dimensional and the plane strain problems [14].
83
G. Yagawa, N, Miyazaki / Steady state creep analysis o f a cracked body
f
= f ca ,j ,tj dr'V
G{~,) = /CB;iiB,~ i d V -
/ E#pqni(B~pctNi*l+ BI~qNi*k) dS
v
Ftk) =
f Tdg,k dS So
us,
+
(36b)
Su
(36c)
,
Su -
f ~EiipqniBk~ dS ,
(36d)
Su
(36e)
sa
Su
It is noted that eqs. (28) and (35) are nonlinear, since the coefficient matrices in these equations are the functions o f d and d °. Direct iteration method is here employed to solve these equations.
3. Numerical examples The computer program for plane stress steady state creep analysis was developed using the theory described above. The finite element utilized here is of conventional type, i.e., constant strain triangle. 3.1. Analysis o f can tilever b earn
The first numerical example is the steady state creep bending analysis of a cantilever beam with one edge clamped and the other subjected to shear force as shown in fig. 1. It is well known that the conventional finite element model with constant strain provides very poor results even in the case of simple elasticity problem of this
I
I
I
I
I
b LIJ ~Q::
0.3
Z o_ I---~ O+ WE: ._1 -J EL = LLI 0..
0
o
0.2 -
P Z
/
o
- - L
z
~.
~ ~ ~
h=THICKNES~ ;
L
O Z
I-~ 13_
0.1
o --: 0
~ 0.0
ELEMENTARY BEAM THEO(~ : PRESENT METHOD
13 : CONVENTIONAL FINITE ELEMENT METHOD I I I 2 3 4
I t
z,u,
CREEP EXPONENT
Fig. l. Element mesh of cantilever beam with shear force P. Fig. 2. Variation of nondimensional deflection rate with creep exponent n.
I 5
n
84
G. Yagawa,N. Miyazaki / Steady state creep analysis of a cracked body
structure. In this relation, we apply the present method to the creep analysis of the structure. The finite element mesh discretization is given in the same figure. The analytical solution introduced into eq. (21) is assumed as follows
+
{ . } | [ - - x 2 , - X a . . . . , -xn+2 ~uz)
n+IZ]lql ]ii
(37)
,
I Iqn+l
where Ux and Uz represent the analytical displacements in the x and z directions, respectively, and q 1, q2, ..., qn+l are the generalized coordinates. Provided the Bernoulli assumption is valid, eq. (37) yields the exact solution for the present problem [ 15 ]. In other words, the role of the finite element model ai in this problem is relatively small in eq. (21). Fig. 2 shows the nondimensional deflection rate at point A in fig. 1 versus creep exponent n, where B =A/(ln) n with In = 2hnb2+(1/n)/(2n + 1). In this figure, the elementary beam theory refers to the book by Odqvist and Hult [15] and the same mesh as shown in fig. 1 is used for the solution of the conventional finite element method. The present method is seen from the figure to coincide fairly well with the elementary beam theory for small values of n. Some discrepancy between these two theories for larger values of n might be explained from the fact that the nonlinearity of deformation along the z-direction becomes significant as the creep exponent n increases. 3.2. Analysis o f square plate with circular hole
The second numerical example is the steady state creep analysis of a square plate with a central hole under uniform tension ao as shown in fig. 3a. The analytical solution is introduced in the superposition method as follows ' . - ( rrx x a - 3 r y 2] x r ~ ( x a - 3 x y 2) ' ux + r4 ] , r2 , r6 q1
'= ao:o ~ * Uy
q2 , (r-S 4
Ya-3_x2y~ r4 ]
y ,r 2 ,
r~(ya-3x2y) r6
(38)
qa
where tc = (3 - v)/(3 + v), G =E/2(1 + v) and r = ~ ) . Poisson's ratio v is taken equal to 0.5 because of incompressibility of creep strain. Eq. (38) is obtained using the theory of two dimensional elasticity with the complex stress function of infinite plate with a circular hole under uniform tension force. The detail of the derivation of the equation is given in ref. [12]. In this calculation, the diameter of the circular hole 2to is taken equal to one third of the plate width, and the finite element mesh is given as depicted in fig. 3b. In fig. 4, the variation of stress concentration factor oy/oo at point A in fig. 3a is shown with the reciprocal of creep exponent n. From the figure, the present solution is found to agree well with Hata's analytical one [16], which is obtained using the stress function method. Fig. 5 shows the variation of the ratio of the analytical displacement in the x direction Ux to the total one Ux + Ux with creep exponent n at point A in fig. 3a. The trend obtained as shown in this figure is reasonable in the sense that the role of the analytical part of the solution gradually decreases with n because the analytical solution utilized here merely represents the displacement field of the simple elasticity problem around a circular hole in an infinite plate. It is interesting, however, to note that the ratio is more than 60% even in n = 5 as shown in the figure. This may imply that, even in the creep exponent n as large as 5, the mixed use of the analytical solution of the simple elasticity problem is still effective to express such a high strain rate around a circular hole as shown in fig. 6.
G. Yagawa,N. Miyazaki / Steady state creep analysis o[ a cracked body
85
y, Uy
Y L
) 2L-
- -
--X,Ux
a o
b
;X
0o - ro hA Fig. 3. Square plate with a circular hole and its element mesh.
3.3. Analysis o f square plate with crack
The numerical example presented here is the steady state creep analysis of a square plate with a center crack whose length is half of the plate width and which is subjected to a uniform tension stress uo as shown in fig. 7. First of all, we consider an elastic-plastic problem whose material model is given as follows e = ale
(o< ay),
(39a)
e=Ao n
(a>~ov),
(39b)
where E and o v are Young's modulus and yield stress, respectively. Rice and Rosengren [ 17], and Hutchinson
i
n--
4.0
o I(_)
i
i
i
i
PRESENT METHOD / o SOLUTION BY HA~A//
Z
_o ~
3.0
Z U.I
g (.)
o
I--
/oJ
2.0
/
I.C *'~ 08 ÷
0
0
0 0 0
o
o Oa O2 PRESENTMETHOD
1"8.0
i
I
I
I
I
02
0.4
0.6
0.8
1.0
INVERSE OF CREEP EXPONENT 1/n
OO
~
1
I
2
I
3
I
'4
5
CREEP EXPONENT n
Fig. 4. Variation of stress concentration factor with inverse of creep exponent 1/m Fig. 5. Variation of ratio Ux/(U x + u x) with creep exponent n.
86
G. Yagawa, N. Miyazaki / Steady state creep analysis of a cracked body
, m 0 t-
,
,
]
o
r
i i
Z
o
Bi
~ o m o
J
20.C
o
§
i
L
o
~
10.(3
o
7
(to PRESENT
Pr
~
I I i
30£
oo
METHOD
i
i
t
i
2
3
4
5
CREEP
EXPONENT
n
// 0
CRACK
Fig. 6. Variation of strain rate concentration factor with creep exponent n. Fig. 7. Square plate with a center crack. [18,19] have determined the character of the plastic crack tip singularity for the material model given in eq. (39) under the conditions of the plane stress and the plane strain. From their studies, the asymptotic stress ai/, the equivalent stress ae, the strain el~, and the displacement ui can be written as follows ai/ = Kor-l/(n+l)~q(O) ,
(40a)
Oe = Kor-1/(n+t)~e(O) ,
(40b)
eq = Ker-n/(n+l)~q(O) ,
(40c)
ui = Kerl/(n+l)ffi(O) ,
(40d)
where Ke = K~ and (r, 0) are the polar coordinates centered at the crack tip. Plastic stress intensity factor Ka is defined as the amplitude of ~e(0) which is normalized by the maximum value of ~e(0). The function ~i/(O), ae(O), "di/(O) and fi'i(0) are determined through the numerical solution of nonlinear ordinary differential equation [18,19] However, the explicit form of these functions has not been obtained. According to the Hoff analogy, the elastic-plastic problem with eq. (39) is equivalent to the following steady state creep problem: when the creep strain rate in a steady state is given by ~c =Aon, the distributions of stress, equivalent stress, strain rate and displacement rate around the crack tip are expressed as follows
oii=Kol'-l/(n+l)oil(O
),
(41a)
oe=Kor-t/(n+l)~e(O ) .
(41b)
~ =Rer-nl(n+l)~q(O) ,
(41c)
h i =I(erll(n+l)ffi(O) ,
(41d)
Therefore, instead of solving the steady state creep problem directly, we may solve the elastic-plastic problem with the stress-strain relation as given in eq. (39). In view of the symmetry of the structure as shown in fig. 7 and the character of eq. (40d), the analytical solution superposed on the finite element solution in this problem can be written as follows m
u x* = ~
m
amr '1("+1) cos ½mO ,
(42a)
u~, = m=, ~ bmr 1/(n+I) sin lmO
(42b)
m=l
where am and bin, m = 1, 2, ..., M, are the generalized coordinates, and M is the number of truncation which is set
G. Yagawa,N. Miyazaki / Steady state creep analysis of a cracked body
87
a
CRACK
E
E 0 F
CRACK
£)
F
Fig. 8. Element mesh of square plate with a center crack.
15ia co co UJ
I
'
I
,
0
1.5
, co
~
o')
t.-(.~
z iJ.J
zU.I
F, z
i
I
0.5
5O3
ub
i
I.C
I--co d <
._1 <[ z
b
OC
Z
-0~
0
PRESENT METHOD n=3 I
-1.0
[3=7
I
0
-OE
L
~
ANGLE
I
-1(
~/2
0
e
I
I
:~/2
.~
ANGLE
e
Fig. 9. O-variation of nondimensional stresses around crack tip (present method).
U3
co
1.5
13/ L~
1.0
...I
< z co z
a
Z 0
- 10
rr t-o3
Z 1.1.1
OC -05
CO CO UJ
,< Z 0 co
05
:~
z
'
HUTCHINSON n=3 I
i
i
-,~/2
ANGLE
e
I
1.0 05 O0
HUTCHINSON ~ -0.5 SLIP LINE THEORY d3 Z n=oo i -t,9.~ 0 z - I.C 0 ~/2" ANGLE
O
Fig. 10.0-variation of nondimensional stresses around crack tip (Hutchinson).
. I
, .K
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked body
88 1 5
(/3 co LIJ r'r I-o9 _j z 0 (,o z
,
i
100
i
0 Z
PRESENT METHOD
10 03 03 LIJ n," I-O0
05
-05
I0
~/
-'SLOPE=-I/2
J
O.C
Z
re 123 Z
•
50
o_ 09 Z U.I
- CONVENTIONAL FINITE ELEMENT METHOD n=3 i
-1C 0
i
ANGLE
0 Z
i
ix/2
J[ e
1
0.001
0'005 0101
0105 011
NONDIMENSIONAL DISTANCE FROM CRACK TIP r / o
Fig. 11.0-variation of nondimensional stresses around crack tip (conventional finite element method). Fig. 12. Logarithmic representation of nondimensional stress distributions near crack tip along x-axis.
to be 10 in this calculation. Fig. 8 is the finite element mesh discretization, in which figs. 8a and 8b show, respectively, a quarter of the plate and the enlargement of the crack tip region. The numbers of the elements and the nodes in this calculation are 242 and 139, respectively. The distributions of ~i/and °e calculated for n = 3 and n = 7 are, respectively, given in figs. 9a and 9b. These figures are depicted by using the stress at the center of gravity of the elements around the crack tip. Fig. 10 shows the distributions of oij and ue for n = 3 and n = o. obtained analytically by Hutchinson. The results shown in figs. 9a and 10a for n = 3 coincide well with each other except near the crack surface 0 = 7r, and the results given in Fig. 9b for n = 7 compare favourably at least qualitatively with those given in fig. 10b for n = ~o based on slip line theory. Compared with the results of figs. 9a and 10a, the results using the conventional finite element method in fig. 11 seems to be cosiderably poor. The reason for the discrepancy for the solutions near the crack surface between the superposition method and the analytical one may attribute rather coarse mesh in the former near 0 = zr as shown in fig. 8a in view of very complex stress distributions around there. The same reason can apply to the numerical oscillation of stresses in the conventional finite element method as shown in fig. 11. Fig. 12 shows the logarithmic representation of the nondimensional stress ay/ao versus nondimensional distance from the crack tip r/a, where the straight lines show the inclinations o f - 1 / ( n + 1) originating from the stress estimate points nearest the crack tip. It can be seen from the figure that the stress singularities near the crack tip are well furnished in the superposition method.
4. Conclusion A superposition method of the analytical and the finite element solutions has been presented for analysis of the plane stress steady state creep problems. Effectiveness of the proposed method has been shown by solving several numerical examples. It may be concluded that the superposition method proposed in this paper is a powerful numerical technique to the steady state creep problems, particularly to the problem of a cracked plate with the singularity around the crack tip.
G. Yagawa, N. Miyazaki / Steady state creep analysis o f a cracked body
89
References [1] [2] [3] [4] [5] [6] [7] [8]
E. Byskov, Int. J. Fracture Mech. 6 (1970) 159. A.K. Rao, I.S. Raju and A.V. Krishna Murthy, Int. J. Numerical Methods Eng. 3 (1971) 389. S.N. Atluri, A.S. K0bayashi and M. Nakagaki, Int. J. Fracture 11 (1975) 257. T.H. Pian, P. Tong and C.H, Luk, in: Proc. 3rd Conf. Matrix Methods in Structural Mechanics (1973) p. 661. D.M. Tracey, Eng. Fracture Mech. 3 (1971) 255. R.D. Henshell and K.G. Shaw, Int. J. Numerical Methods Eng. 9 (1975) 495. Y. Yamamoto and N. Tokuda, Int. J. Numerical Methods Eng. 6 (1973) 427. G. Yagawa, T. Nishioka, Y. Ando and N. Ogura, in! Computational Fracture Mechanics, eds. E.F. Rybicki and S.E. Benzley (ASME, New York, 1975) p. 15. [9] R.H. Gallagher, in: Numerical Methods in Fracture Mechanics, eds. A.R. Luxmoore and D.RJ. Owen (University College Swansea, Swansea, 1978) p. 1. [10] S.E. Benzley and D.M. Parks, in: Structural Mechanics Computer Programs, eds. W. Pilkey, K. Saczalski and H. Schaeffer, (University Press of Virginia, 1974) p. 81. [11] M.C. Apostal, S. Jordan and P.V. Marcal, EPRI SR-22 Special Report (1975). [12] G. Yagawa, T. Nishioka and Y. Ando, Nucl. Eng. Des. 34 (1975) 247. [131 G. Yagawa, N. Miyazaki and Y. Ando, Int. J. Numerical Methods Eng. 11 (1977) 1107. [14] L.R. Herrmann, AIAA J. 3 (1965) 1896. [ 15] F.K.G. Odqvist and J. Hult, Kriechfestigkeit Metalliseher Werkstoffe (Springer-Verlag, Berlin, 1962). [16] T. Hata, Nucl. Eng. Des. 32 (1975) 325. [17] J.R. Rice and G.F. Rosengren, J. Mech. Phys. Solids 16 (1968) 1. [18] J.W. Hutchinson, J. Mech. Phys. Solids 16 (1968) 13. [19] J.W. Hutchinson, J. Mech. Phys. Solids 16 (1968) 337.