Computational Materials Science 36 (2006) 303–309 www.elsevier.com/locate/commatsci
Creep sagging analysis of pressure pipes D.G. Pavlou a
a,*
, N.V. Vlachakis a, Ch. Tsitouras a, M.G. Pavlou b, A. Fatsis
c
Technological Institute of Halkidas TEI, School of Applied Technology, GR34400 Psahna, Halkida, Evoia, Greece b ERGOSE SA, Works of Greek Federal Railway Organization, Karolou 27, 10437 Athens, Greece c RINA Hellas Ltd., Gruppo REGISTRO ITALIANO NAVALE, 47-49 Akti Miaouli str., 18536 Piraeus, Greece Received 25 November 2004; received in revised form 18 February 2005; accepted 8 April 2005
Abstract A new concept is developed to predict the mechanical behavior of pipes under the combination of internal pressure and elevated temperature. A non-linear differential equation is formulated to approximate the stress/displacement distribution in these pipes, which are reinforced by rigid rings uniformly distributed along their axis. The proposed model incorporates a non-linear time-dependent stress–strain material behavior rule derived by uniaxial tests. Taking into account equilibrium equations and boundary conditions, a general analytical expression for the case of biaxial creep loading is formulated. Numerical solution and results of a typical example are included. 2005 Elsevier B.V. All rights reserved. Keywords: Creep; Pipe; Sagging; Non-linear differential equation; Finite differences method
1. Introduction Pipe-lines used in power plants are mainly operating at high temperature and pressure conditions resulting in permanent deformations as they age [1]. Often, the prediction of the mechanical behavior of the material under these conditions is restricted for the case of uniaxial loading. Due to the elevated values of temperature, the relation of the resulted displacements by the stress is given by a time-dependent non-linear functional form because of the creep mechanism [2,3]. To perform stress/ displacement analysis for the case of axisymmetric biaxial loading of creeping pressurized pipes, a non-linear bending theory must be developed. This work focuses on the derivation of the differential equation in terms of displacement for a pipe, which is simulated as an axisymmetrically loaded cylindrical shell *
Corresponding author. Tel.: +30 22280 99647/525; fax: +30 2280 99664. E-mail addresses:
[email protected],
[email protected] (D.G. Pavlou). 0927-0256/$ - see front matter 2005 Elsevier B.V. All rights reserved. doi:10.1016/j.commatsci.2005.04.003
with the two ends fixed. Due to the boundary effects, geometric discontinuities and non-linear material behavior, the displacement analysis results to a complex non-linear sixth-order differential equation. Using mathematical transformations and phenomenological assumptions, the above sixth-order non-linear differential equation is reduced to a fourth-order one. It is approved that the existing linear fourth-order differential equation [4] describing the displacement distribution of a pipe at room temperature can be derived as a particular case of the general non-linear differential equation developed in present work. The derived non-linear equation is solved numerically using the finite differences method (FDM).
2. Formulation of the problem A long axisymmetric pipe of diameter D and wall thickness s is considered to be reinforced by rigid rings. The rings are uniformly distributed along the pipe-line and the distance between two rings is L. The axes of
304
D.G. Pavlou et al. / Computational Materials Science 36 (2006) 303–309
Fig. 1. Coordinate system of the problem.
the coordinates (x, w) as well as the geometrical parameters are shown in Fig. 1. The pipe is loaded by internal pressure Po at temperature T. Under these conditions, at every point x along the wall, a radial displacement w is developed progressively. Therefore, it can be written as w ¼ wðx; P o ; T ; tÞ
ð1Þ
Due to the relation between stress–strain and strain– displacement, the stress and strain fields within the material of the pipe can be expressed as a function of w. The boundary conditions of the problem are
Fig. 2. Uniaxial mechanical behavior of material under creep conditions.
wð0Þ ¼ 0 wðLÞ ¼ 0
ð2Þ ð3Þ
where n is the creep exponent, which usually lies between 3 and 8. The expression (9) is called ‘‘power-law’’ or Norton–Bailey rule [2]. Plotting the experimental values of the natural logarithm of de/dt against the reciprocal of the absolute temperature (1/T) the following relation can be obtained:
w0 ð0Þ ¼ 0 w0 ðLÞ ¼ 0
ð4Þ ð5Þ
de ¼ CeQc =Rg T dt
w0 ðL=2Þ ¼ 0
ð6Þ
wðaÞ ¼ wðL aÞ; 0 6 a 6 L w0 ðaÞ ¼ w0 ðL aÞ; 0 6 a 6 L
ð7Þ ð8Þ
where Rg is the universal gas constant (8.31 J/mole K) and Qc is the activation energy for creep (it has units of J/mole). Combining the dependences of de/dt by the stress (Eq. (9)) and temperature (Eq. (10)) it can be written as de ¼ AeQc =Rg T rn dt
3. Mechanical behavior of the material Materials operating at elevated temperatures have different mechanical behavior compared to those operating at room temperatures due to the creep mechanisms. From macroscopic point of view, creep is slow continuous permanent deformation activated when T > 0.3 to 0.4 Tm for metals (e.g. [2]), where Tm is the melting point. At these conditions, the strain e measured by uniaxial loading tests is related by the time t for several constant stress values r showing the general shape of Fig. 2. Although the primary (stage I) and the tertiary (stage III) creep periods cannot be neglected, they occur quickly. The longer creep period is the steady-state period when the strain increases steadily versus time. In designing against creep, it is usually this steady accumulation of strain that is considered as the most important. Plotting the experimental values of the logarithm of the steady creep-rate log(de/dt) against logarithm of the uniaxial stress log(r), at constant temperature T, it can be established that de ¼ Brn dt
ð9Þ
ð10Þ
ð11Þ
where A is the creep constant. The values of the constants A, n, Qc vary from material to material and have to be found experimentally by uniaxial creep tests.
4. Bending theory In order to develop the fundamental equation describing the radial displacement distribution w(x) of the pipeÕs wall, the equilibrium of a longitudinal strip of length L and width RDh = 1 of the pipe is considered firstly (Fig. 3). Due to bending, a cross section BC of above strip (Fig. 4) is moved to the position B 0 C 0 . After the bending it is assumed that the cross section BC remains plane [4]. The displacement u(x, y) of a point A(x, y) (Fig. 4) can be determined by the simple geometric formula uðx; yÞ ¼ yh
ð12Þ
where h is the slope of the cross section BC after bending. Taking into account Fig. 4 it can be written as h ¼ w0
ð13Þ
D.G. Pavlou et al. / Computational Materials Science 36 (2006) 303–309
305
where J is called herein ‘‘non-linear moment of inertia’’ for a section with width equal unity and it is given by the equation Z 1 y nþ1 dF ð21Þ J¼ F
Taking into account Eq. (20) and the following equilibrium equations (22) and (23) for the shear force Q and the distributed load q Q¼
dM dx
ð22Þ
q¼
d2 M dx2
ð23Þ
Fig. 3. Equilibrium of a longitudinal strip of pipeÕs wall.
the following differential equation is obtained: J =n 00 1n1 000 0 ðw Þ w ¼q ðBtÞ1=n
ð24Þ
The distributed load q is the superposition of the action q** due to the internal pressure Po and the reaction q* due to the wall resistance, so that q ¼ q q
ð25Þ
4.1. Action of the internal pressure Fig. 4. Slope of a cross section due to bending.
Therefore uðx; yÞ ¼ yw0
ð14Þ
The resulted strain at the point A(x, y) can be obtained [4] by the expression
q ¼ 2pRLP 0 =L
ð26Þ
or q ¼ 2pRP 0
ou e¼ ox
ð27Þ
ð15Þ 4.2. Reaction of the wall
According to Eq. (14), Eq. (15) can be written as e ¼ yw00
ð16Þ
Due to creep conditions, the strain e can be obtained approximately by the Norton-Bailey rule given by Eq. (9). Then, Eq. (16) can be written as Brn t ¼ yw00
ð17Þ
or ð18Þ
Combination of Eq. (18) with the following equilibrium condition for the bending moment M Z M¼ yr dF ð19Þ F
results to J ðBtÞ
1=n
ðw00 Þ
The increase of the circumference of the pipe due to the radial displacement DR can be calculated as Dð2pRÞ ¼ 2pDR
1=n
ð20Þ
ð28Þ
Then, the resulting circumferential strain can be expressed as et ¼
00 1=n w r¼ y 1=n Bt
M¼
The force acting in the wall due to the internal pressure Po can be expressed by the product 2pRLPo. Then, the distributed load q** can be written as
Dð2pRÞ 2pR
ð29Þ
DR w ¼ R R
ð30Þ
or et ¼
Taking into account the Norton–Bailey rule given by Eq. (9), Eq. (30) can be written as w ð31Þ Btrnt ¼ R Then, the resulting circumferential stress rt due to the radial displacement DR can be obtained as
306
rt ¼
D.G. Pavlou et al. / Computational Materials Science 36 (2006) 303–309
1 ðBRtÞ1=n
w1=n
ð32Þ
The radial resistance force F of the wall, developed due to the radial displacement w, can be determined by the equilibrium of the elementary section of the wall shown in Fig. 5. Therefore dh dF ¼ 2ðrt sLÞ sin ð33Þ 2 Because of the small values of the dh/2, the following approximation can be used: dh dh ¼ sin 2 2
ð34Þ
Then the total radial resistance F due to the radial displacement w can be obtained as Z 2p F ¼ dF ð35Þ 0
or ð36Þ
Therefore, the corresponding distributed force derived as F L
q*
can be
ð37Þ
1 0 1=n ðJ =nÞ w00n1 w000 þ 2psR1=n w1=n ¼ 2pRP o ðBtÞ ð42Þ Substituting the values n = 1 and 1/Bt = E, above equation leads to the well known linear differential equation [4] describing the sagging of pipe at room temperature w0000 þ cw ¼ d
ð43Þ
where c, d are constants containing geometrical and material parameters. Using the transformation w00 ¼ ðf 00 Þ
n
ð44Þ
Eq. (42) can be expressed as Jf 0000 þ 2psR1=n w1=n ¼ 2pRP o ðBtÞ
1=n
ð45Þ
f 0000 þ kw1=n ¼ m
ð46Þ
where 2ps k ¼ 1=n JR and
ð47Þ
1=n
or
2pRP o ðBtÞ ð48Þ J The solution of the differential equation (46) is very difficult due to the non-linearity. In order to simplify this equation, the non-linear reaction q* of the wall, given by Eq. (39), is assumed to be approximated by the substitute linear reaction m¼
q ¼ 2psrt
ð38Þ
With the aid of Eq. (32), Eq. (38) can be written as q ¼
or
or
F ¼ 2psrt L
q ¼
Then, the differential equation (24) can be written as J =n 001n1 000 0 2ps w w þ w1=n ¼ 2pRP o ð41Þ 1=n 1=n ðBtÞ ðBRtÞ
2ps ðBRtÞ
1=n
w1=n
ð39Þ
Taking into account Eqs. (27) and (39), Eq. (25) takes the form q ¼ 2pRP o
2ps ðBRtÞ
w1=n 1=n
ð40Þ
qm ¼ k m w
ð49Þ
of a simulated wall which accumulates the same energy as the true wall. To this scope, the behavior of the original and the simulated wall is demonstrated in Fig. 6. The equivalence between the accumulated energy of the original and the simulated wall can be expressed as Z wf 1 wf ðk m wf Þ ¼ q dw ð50Þ 2 0 where wf is the radial displacement corresponding to the material failure. Taking into account Eq. (39), Eq. (50) can be written as 1 2ps n ðnþ1Þ=n k m w2f ¼ wf 1=n 2 ðBRtÞ n þ 1
ð51Þ
or km ¼ Fig. 5. Superposition of the internal pressure and the wallÕs resistance.
n 4ps ð1nÞ=n w n þ 1 ðBRtÞ1=n f
ð52Þ
D.G. Pavlou et al. / Computational Materials Science 36 (2006) 303–309
307
4.3. Boundary conditions Combining Eqs. (2) and (56) it can be written as f 0000 ð0Þ ¼ m
ð63Þ
Taking into account Eq. (60), Eq. (63) results to g00 ð0Þ ¼ m
ð64Þ
Considering Eqs. (3) and (7), and following the same procedure, the boundary conditions g00 ðLÞ ¼ m
ð65Þ
and g00 ðaÞ ¼ g00 ðL aÞ;
can be derived. Differentiation of Eq. (56) leads to
Fig. 6. Resistance of the original and the simulated wall.
When the material failures, i.e. e = ef, Eq. (30) results to wf ð53Þ ef ¼ R or
f ðVÞ þ bw0 ¼ 0
ð67Þ
Combination of Eqs. (66), (4)–(6), (8), (60) results to the following boundary conditions: g000 ð0Þ ¼ 0
ð68Þ
000
wf ¼ Ref
ð54Þ
where ef is the tensile ductility of the material at the service temperature T. Then, Eq. (52) takes the following form: km ¼
n 4ps ð1nÞ=n e ðn þ 1ÞR ðBtÞ1=n f
ð55Þ
With the aid of Eqs. (25), (27), (44), (49) the fundamental differential equation (24) can be simplified to f 0000 þ bw ¼ m
ð56Þ
where the parameter b can be derived with the aid of Eq. (55) b¼
ð66Þ
06a6L
n 4ps 1n en n þ 1 RJ f
ð57Þ
Two times differentiation of Eq. (56) results to f ðVIÞ þ bw00 ¼ 0
ð58Þ
Considering Eq. (44), the above equation yields to n
f ðVIÞ þ bðf 00 Þ ¼ 0
ð59Þ
Using the transformation f 00 ¼ g
ð60Þ
above equation can be written as g0000 þ bgn ¼ 0
ð61Þ
where, according to Eq. (44), the function g can be expressed as g ¼ ðw00 Þ
1=n
ð62Þ
g ðLÞ ¼ 0 g000 ðL=2Þ ¼ 0 g000 ðaÞ ¼ g000 ðL aÞ;
ð69Þ ð70Þ ð71Þ
06a6L
5. Numerical procedure The non-linear differential equation (61) is solved numerically using the finite difference method (FDM) for a practical problem. Details for using this type of methods with MATLAB [5], are given in [6]. The considered pipe has length L = 10 m, radius R = 0.25 m and thickness t = 0.008 m. The ends of the pipe are fixed. Therefore at these locations (x = 0, x = L) the deflection and the slope of the wall are zero (w = 0, w 0 = 0). The used material is the steel X8CrNiMoNb1616 which is subjected at creep conditions p = 3 · 105 Pa and T = 973 K. For this material the Norton–Bailey parameters have been obtained by Ref. [7] and take the values n = 6.43 and B = 3.85 · 1033 (for A = 3.5 · 1019 h1, Qc = 260 kJ/mole, Rg = 8.314 J/mole.grad, T = 973 K). The geometrical parameter J for the section of width R 1 2pR and thickness t takes the value J ¼ 2pR F y nþ1 dF ¼ 1.01 105 . Then the value m = 297284 can be obtained by Eq. (48). Using sixth-order differences with a special operation at the boundaries we arrive to the non-linear FDM of the form Mz þ h2 ml þ h4 bzn ¼ 0;
ð72Þ T
where h = L/k, z [g(0)g(h)g(2h) g(L h)g(L)] is a vector containing k + 1 approximations of g(x) at the
308
D.G. Pavlou et al. / Computational Materials Science 36 (2006) 303–309
corresponding points and zn is a componentwise power of vector z. The vector l is given by the equation l ½8.2835 2.4386 0.5218 0 0 0 0 0 T
0.5218 2.4386 8.2835 ;
ð73Þ
while the (k + 1) · (k + 1) matrix M can be written as
pipe. However, it can be shown that the distribution w(x) has two additional local maximum values near the locations x = 0.5 m, x = 9.5 m and two local minimum values near the locations x = 1.5 m and x = 8.5 m. The physical interpretation of these deformations near the ends is not known. However, are in agreement with the results shown schematically in Fig. 9,
2
13.7427 35.6562 33.6425 16.2538 5.6115 1.2068 0.1201 0 6 4.8716 12.8268 11.8209 4.7934 1.0869 0.1730 0.0137 0 6 6 6 2.5582 9.1009 12.8440 9.0248 3.1680 0.4824 0.0378 0 6 1 13 28 13 1 6 6 2 2 0 6 6 2 3 2 6 6 6 .. 6 0 0 . M 6 6 .. .. 6 . 0 . 6 6 1 13 28 13 6 0 0 2 6 6 2 3 2 6 6 0 0 0.0378 0.4824 3.1680 9.0248 12.8440 6 6 4 0 0 0.0137 0.1730 1.0869 4.7934 11.8209 0 0 0.1201 1.2068 5.6115 16.2538 33.6425
..
3
0 0
.
2 9.1009
7 7 7 7 0 7 7 7 0 7 7 7 7 7 7 7 0 7 1 7 7 7 6 7 2.5582 7 7 7 4.8716 5
12.8268 35.6562 13.7427
ð74Þ We solved the FDM using MATLAB, and at a particular choice of k = 20, the resulted distribution of the parameter w00 = (Bt/Jn)Mn along the wall of the pipe is shown in Fig. 7. Numerical integration of this result leads to the wallÕs radial deformations w(x) along the pipe (Fig. 8). From the results shown in (Fig. 7) it is concluded that the maximum values of bending moment are located at the fixed ends. The derived safety factor N at these locations can be calculated by the equation Su rmax ¼ ð75Þ N where rmax is the applied maximum stress given by Eq. (18) for y = s/2, while Su is the ultimate stress at the temperature T = 973 K. Using the value Su = 10 Ksi = 68060 Pa obtained by Ref. [1] the value N 2 of the safety factor can be determined by Eq. (75). Fig. 8 indicates that the deflection distribution w(x) takes the maximum value at the middle x = 5 m of the Bt J^n M^n 15
obtained by Ref. [8], of existing solution concerning similar thin-walled pipes at room temperature. An attempt for mathematical interpretation of the pipeÕs deformation shown in Fig. 8 is possible when the temperature of the material is reduced to room temperature values. As mentioned earlier, for this case the developed non-linear differential equation given by Eq. (42) is simplified to the linear differential equation in Eq. (43). The solution of the latter equation with the boundary conditions of the bending of a long cylindrical shell subjected to a concentrated line-load, uniformly distributed along a circular section (due to a ring), is given in [9] wðxÞ ¼ AuðbxÞ
ð76Þ
w cm 2 1.5
10
1
5 0.5 2
4
6
8
10
x m
-5 Fig. 7. Normalized bending moment distribution along the pipe.
x m 2
4
6
8
10
Fig. 8. Radial deflection distribution of the wall along the pipe.
D.G. Pavlou et al. / Computational Materials Science 36 (2006) 303–309
309
not similar, because the solution expressed by Eq. (76) corresponds to infinitely long pipe.
6. Conclusions
Fig. 9. Schematic representation of the deformed pipe reinforced by rings at room temperature, obtained by Ref. [8].
0.4
ϕ (βx)
0.2 2
4
6
8
10
βx
-0.2
1. A new concept was developed to describe the stress and deformation behavior of power pipes under creep conditions. 2. The pipes under consideration are reinforced by rigid rings, uniformly distributed along the pipe-line. 3. The developed model is a non-linear differential equation, which can be solved by numerical methods only. 4. A practical problem of a pipe with fixed ends, under combination of internal pressure and elevated temperature was solved numerically by the finite differences method. 5. The results indicated that the maximum bending stresses were developed at the pipeÕs ends while the predicted deformed solid is in agreement with corresponding results at room temperature.
-0.4 -0.6
References
-0.8 -1 Fig. 10. Graphical representation of the function u(bx).
where A, b are constant coefficients related to the loading, as well as the material and geometrical properties of the structure, and u(bx) is the function uðbxÞ ¼ ebx ðcos bx þ sin bxÞ
ð77Þ
Near the location x = 0 (ringÕs location), the similarity between the graphical representation of the function u(bx) shown in Fig. 10 and the numerical results shown in Fig. 8 is obvious. However, far away from the value x = 0, the distributions shown in Figs. 8 and 10 are
[1] ASME B31.1, Power Piping 2001: Standard for the identification, design and welding of piping, ASME/ANSI, 2001. [2] F.H. Norton, Creep of Steel at High Temperatures, McGraw-Hill, NY, 1929. [3] D.G. Pavlou, Creep life prediction under stepwise constant uniaxial stress and temperature conditions, J. Eng. Struct. 23 (2001) 656– 662. [4] A.P. Boresi, O.M. Sidebottom, Advanced Mechanics of Materials, Wiley, 1985. [5] MATLAB 6, The MathWorks, Inc., Natick, MA, 2000. [6] L.F. Shampine, I. Gladwell, S. Thompson, Solving ODEs with MATLAB, Cambridge University Press, 2003. [7] S.G. Pantelakis, Kriechverhalten metalliscer Werkstoffe bei zeitveraenderlicher Spannung, Ph.D. Thesis, RWTH-Aachen, 1983. [8] H. Oery, Leichtbau, Rheinisch-Westfalischen Technis chen Hochschule-Aachen, 1990. [9] S.P. Timoshenko, S. Woinowsky-Krieger, Theory of Plates and Shells, McGraw-Hill, NY, 1959.