Thermal stress singularity analysis for V-notches by natural boundary element method

Thermal stress singularity analysis for V-notches by natural boundary element method

Accepted Manuscript Thermal stress singularity analysis for V-notches by natural boundary element method Changzheng Cheng , Shenyu Ge , Shanlong Yao ...

770KB Sizes 0 Downloads 26 Views

Accepted Manuscript

Thermal stress singularity analysis for V-notches by natural boundary element method Changzheng Cheng , Shenyu Ge , Shanlong Yao , Zhongrong Niu PII: DOI: Reference:

S0307-904X(16)30278-5 10.1016/j.apm.2016.05.028 APM 11181

To appear in:

Applied Mathematical Modelling

Received date: Revised date: Accepted date:

3 November 2015 1 April 2016 19 May 2016

Please cite this article as: Changzheng Cheng , Shenyu Ge , Shanlong Yao , Zhongrong Niu , Thermal stress singularity analysis for V-notches by natural boundary element method, Applied Mathematical Modelling (2016), doi: 10.1016/j.apm.2016.05.028

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT

Highlights  The nearly singular integral is reduced by one order in the newly established natural BEM.  The natural BEM can be used to calculate the thermal

CR IP T

stress much closer to the boundary.  The thermal stress singularity is investigated basing on the near V-notch tip stress field.

AN US

 The mathematic and physic double singularity are dealt

AC

CE

PT

ED

M

with in the present paper.

1

ACCEPTED MANUSCRIPT

Thermal stress singularity analysis for V-notches by natural boundary element method Changzheng Cheng*, Shenyu Ge, Shanlong Yao, Zhongrong Niu (Department of Engineering Mechanics, Hefei University of Technology, Hefei 230009, China)

CR IP T

Abstract: The accuracy of thermal stresses at internal points by the conventional boundary element method becomes deteriorate when the points are approaching to the boundary due to the inaccuracy of the calculation of nearly singular integrals by the Gaussian integration. Herein, a thermal stress natural boundary integral equation is proposed, in which the nearly hyper-strongly

AN US

singular integral is reduced to a nearly strongly singular one and then dealt by the regularization method. Thus, it can be applied to model the high stress gradient in the vicinity very close to the V-notch vertex. The stress method is subsequently introduced to calculate the stress singularity

M

orders and thermal stress intensity factors once the thermal stresses along the bisector and very close to the notch tip are yielded. The mathematical and physical singularity difficulties, i.e., the

ED

evaluation of nearly singular integrals and singular stress fields, are both overcome in this paper. After a benchmark model being given out to verify the efficiency of the present method, the

PT

thermal stress singularities for a symmetrical V-notch and an inclined one are respectively

CE

analyzed. The benchmark example manifests that the thermal stress nature boundary integral equation can be successfully used to calculate the thermal stresses much closer to the boundary by

AC

the comparison with the conventional thermal stress boundary integral equation. The accuracy of the stress singularity orders and thermal stress intensity factors by the present method is confirmed and the computational effort is dramatically decreased by comparing with the finite element method. Key words: V-notch, stress singularity order, thermal stress intensity factor, nearly singular integral, boundary element method

*

Corresponding author. E-mail address: [email protected] 2

ACCEPTED MANUSCRIPT

1. Introduction The thermal stress concentration caused by the structural or/and material discontinuities usually arises in the electronic devices, welding structures and construction of composite materials [1-3]. When a plate is subjected to thermal loading which results in varying temperature T (see

CR IP T

Fig.1a), the stress concentration near the V-notch tip is not a priori equal to the one inducing by the equivalent mechanical loading  0  E  T (see Fig.1b), where 

is the thermal

expansion coefficient, E is the material elasticity modulus. The calculated stress intensity factors of the V-notches respectively shown in Fig.1a and Fig.1b are definitely different when the

AN US

displacement method is used, and they gradually deviate from each other seriously with the increase of notch depth when the stress method is utilized. In addition, two components weakened by blunt notches and scaled in geometrical proportion have the same values of the theoretical

M

stress concentration factor. However, two components weakened by sharp V-notches and scaled in

ED

geometrical proportion have different notch thermal stress intensity factors [4]. Although the research focused on the thermal stress singularity of the V-notch could be dated back to 1960s

PT

[5,6], it is still deserved to be studied nowadays.

CE

There are some numerical methods, such as the finite element method and Hamiltonian

AC

approach, developed to calculate the thermal stress intensity factors for cracks [7-10] and V-notches [11-17]. The boundary element method is a powerful tool for the numerical evaluation, which can also be used to analyze the singular thermal stress near the tip of V-notches [18-26]. When the source point is close to but not on the boundary element, the distance r between the source point and the integrated node (field point) is very small, namely, r  0 . In this case, the nearly hyper-strongly, strongly and weakly singular integrals will occur in the thermal stress

3

ACCEPTED MANUSCRIPT

boundary integral equation because the integral kernels contain 1 r 6 , 1 r 4 and 1 r 2 terms, respectively. The accuracy of the stresses becomes worse and worse with the inner point approaching to the boundary in the conventional boundary element method, due to the calculating difficulty of the nearly singular integrals. However, the accuracy of the stress filed very close to

CR IP T

the notch tip is a key issue for the determination of the stress singularity orders and stress intensity factors of the V-notch under thermal loading. The first task of using the boundary element method to the singularity analysis of the V-notch is to deal with the nearly singularity integrals in the

AN US

boundary integral equations.

The recently developed methods are separately dealing with the nearly hyper-strongly singular integral and nearly strongly singular one [27-38]. In the author’s opinion, the nearly

M

singular integrals are singular in the sense of mathematics, but not in the sense of physics. These two kinds of nearly singular integrals could be handled together [39]. Herein, the nearly

ED

hyper-strongly singular integral is transformed into a nearly strongly singular one by introducing

PT

two new natural boundary variables into the thermo-elasticity displacement derivative boundary integral equation. Then, the newly created nearly strongly singular integral is added together with

CE

the original one and handled by the regularization method proposed by the authors before [40].

AC

Since the nearly singularity is decreased by one order, the newly built thermal stress natural boundary integral equation can be successfully used to calculate the thermal stress field much closer to the notch tip. The stress method is then applied to analyze the thermal stress singularity of the V-notch, basing on the obtained thermal stress field in the near tip region. After a benchmark model being proposed to investigate the efficiency of the present method, a symmetrical and an inclined V-notch models are built to calculate the stress singularity orders

4

ACCEPTED MANUSCRIPT

and thermal stress intensity factors, which are compared with the ones from the finite element method. The accuracy of the present method is verified and the computational effort is obviously decreased by the comparison with the finite element method. 2. Evaluation of natural boundary variables

CR IP T

For 2-D thermo-elasticity with a uniform temperature variation, the displacement boundary integral equation without considering the body force can be written as Cij ( y)ui ( y) 

 U t ( x)d   T u ( x)d   [R T ( x)  Q 

* ij j



* ij

j

i



i

T ( x ) ]d n

(1)

AN US

where i, j  1,2 , x is a field point and y is a source point; Cij ( y) is the singular coefficient which depends on the material constant and local geometry of the domain boundary 

at point

y ; u j (x ) and t j (x ) are respectively the displacement and traction component on the boundary

 ; T (x ) and T (x) n are respectively the temperature and normal temperature gradient on

M

the boundary  ; the expressions of U ij , Tij , Ri and Qi in the integral kernels are listed in

ED

Appendix A. The customary standard Euler notation for summation over repeated subscripts is

PT

used.

The singular coefficient Cij ( y) is equal to one when y is an inner point. In this case, Eq.(1)

CE

becomes a displacement boundary integral equation of an inner point and is written as below

 U t ( x)d   T u ( x)d   [R T ( x)  Q 

* ij j

* ij



j

i



i

T ( x ) ]d n

(2)

AC

ui ( y) 

By differentiating Eq.(2) with respect to the coordinate xk , the displacement derivative boundary integral equation of an interior point can be yielded as follows ui , k ( y) 

U 

t ( x )d 

* ij, k j

T 

* ij, k

u j ( x )d 

 [R 

i,k

T ( x )  Qi , k

T ( x ) ]d n

(3)

Basing on Eq.(3) and using 2-D thermo-elasticity Hooke’s law, the conventional thermal stress boundary integral equation for an interior point can be derived and written as

5

ACCEPTED MANUSCRIPT

 ik ( y) 

W

t ( x)d 

* ikj j



S 

u j ( x)d 

* ikj

 [ R T ( x)  Q ik



ik

T ( x)]d    ik0 ( y) n

(4)

The expressions of the integral kernels in Eqs.(3,4) and initial thermal stress  ik0 ( y) in Eq.(4) are listed in Appendix A, from which it can be observed that Qik has a weak singularity of order * has a hyper-strong 1 r 2 , Wikj* and Rik has a strong singularity of order 1 r 4 and S ikj

CR IP T

singularity of order 1 r 6 , and r is the distance from the source point to field point. The integrals containing these three kinds of singularity orders are called the nearly singular integrals when the source point is approaching to, but not on, the integral element, i.e., r  0 , because the

AN US

magnitude of the integrand may be quite large when the Gaussian integration is used.

After respectively multiplying the Kronecker operator  ik and permutation symbol eik (i.e., eik

i k

 0 , e12  1 , e21  1 ) to the two sides of Eq.(3), one can get r, j 1  2 1 1 [ t j ( x )d  (2r, j r,n  n j )u j ( x )d ]  2  r 2 (1   ) 2G  r r  (1   ) T ( x ) [ ,n T ( x )  (ln r  1) ]d 1  2  r n





1

1  2G [



ED

 eik ui ,k ( y) 



(5a)

M

 ik ui ,k ( y ) 



r,k r

ekjt j ( x )d 





1 (2r, j r,   j )u j ( x )d ] r2

(5b)

PT

where  is the Poisson’s ratio and G is the shear modulus. When the source point is approaching towards the boundary integral element, in other words, r  0 , the 1st and 2nd

CE

integrals at the right side in Eq.(5) respectively present the nearly weak singularity and nearly strong one. For distinguishing the boundary on which the nearly singular integrals occur or not,

AC

the boundary on which the integrals in Eq.(5) are nearly singular is denoted as 1 , and the other part is called  2 as shown in Fig.2, where x1 and x 2 are respectively the starting and terminal points of boundary 1 , n and  are respectively the normal and tangential direction along the boundary. Furthermore, the integral

 ()d 1

 ()d 

and the none nearly singular part

in Eq.(5) is departed to the nearly singular part

 ()d . 2

For replacing the partial derivative u j (x)  in Eq.(5), two new boundary variables 6

ACCEPTED MANUSCRIPT named natural boundary variables 1 ( x ) and 2 ( x ) are respectively defined as

1 ( x ) 

tn ( x ) u ( x ) t ( x ) un ( x )   , 2 ( x )   2G  2G 

(6)

After 1 ( x ) and 2 ( x ) being introduced in, Eq.(5) can be rewritten as follows r,n

r,k

1

1



2 ( y) 

r,

 [ r  ( x)  r  ( x)]d  r e u ( x) 2

 (1   ) 1  2 r,

kj

r,n

 [ r T ( x)  (ln r  1) 

j

r, j

r,n

1

2



 [

x  x1

2

r, j t j ( x ) r 2G



T ( x ) ]d n

 [ r  ( x)  r  ( x)]d  r u ( x) 1

x2

j

x2 x x1



 [

r,k

2

r

ekj

1 (2r, j r,n  n j )u j ( x )]d r2 (7a)

CR IP T

1 ( y ) 

t j ( x) 2G



1 (2r, j r,   j )u j ( x )]d (7b) r2

The natural boundary variables i ( x)(i  1,2) in the thermo-elasticity problem can be

AN US

yielded from Eq.(7) when it is discretized on the boundary.

3. Establishment of thermal stress natural boundary integral equation In this section, the thermal stress boundary integral equation expressed with the natural boundary variables 1 ( x ) and 2 ( x ) will be given out.

where

Tij*,k  

Eij*,k

(8)



PT

ED

with respect to  as follows,

M

The integral kernel Tij*,k in Eq.(3) can be written as the negative partial differential of Eij*,k

CE

Eij*,k  

1 [2(1   )emk r,m ij  2e jl r,i r,k r,l  (1  2 )eij r,k  e jl (r,l ik  r,i lk )] 4 (1   )r

(9)

AC

According to the integration by parts, there is



Eij*,k

1



u j ( x )d  Eij,k u j ( x )

x2 x x1





Eij*,k





Eij*,k

1

u j ( x ) 

d

(10)

d

(11)

Introducing Eq.(8) into Eq.(10), one can get



1

Tij*,k u j ( x )d   Eij,k u j ( x )

x2 1

x x

When Eq.(11) is substituted in, Eq.(3) can be rewritten as

7

1

u j ( x ) 

ACCEPTED MANUSCRIPT

ui , k ( y ) 

U 

t ( x )d 

 ij ,k j



2

Tij,k u j ( x )d 



1

u j ( x )

Eij,k



d  Eij,k u j ( x )

x2 x  x1

T ( x )  [ Ri ,kT ( x )  Qi ,k ]d  n



(12)

Furthermore, after Eq.(12) being introduced into the Hooke’s law, one can get the thermal stress boundary integral equation at an interior point which is expressed as



Wikj t j ( x )d 



[ RikT ( x )  Qik









2

 Sikj u j ( x )d 



1

Fikj

u j ( x ) 

d  Fikj u j ( x )

x2 x x1

CR IP T

 ik ( y) 

T ( x ) ]d   ik0 ( y) n

(13)

The expression of Fikj* in Eq.(13) is the function of 1 / r 4 which can be found in Appendix

AN US

* A. The nearly hyper-strongly singular integral kernel S ikj in Eq.(4) has been transformed into a

nearly strongly singular one Fikj* in Eq.(13), which means that the singularity has been reduced by one order. However, there appears a new variable u j (x)  in Eq.(13), which is never

M

solved before and should be replaced with the known boundary variables. By noting the relationship below

u j ( x )

ED



un ( x ) u ( x )  nl e jl   

  l e jl

(14)

PT

two nearly strongly singular integrals at the right side in Eq.(13) can be added together and expressed by 1 ( x ) and 2 ( x ) as follows

 [W

CE 1

AC

t ( x )  Fikj*

* ikj j

u j ( x ) 

]d 

 [F 1

* ikj

n j2 ( x )  Fikj*  j1 ( x )]d

t ( x) t ( x)  [F  j n  Fikj n j   Wikj* n j tn ( x )  Wikj*  j t ( x )]d 1 2G 2G



(15)

* ikj

Inserting Eq.(15) back into Eq.(13) yields  ik ( y ) 



2

Wikj t j ( x )d 



2

 Sikj u j ( x )d  Fijk u j ( x )

x2 x  x1





1

Fikj* e jm[ m2 ( x )  nm1 ( x )]d

T ( x )  H t ( x )d  [ RikT ( x )  Qik ]d   ik0 ( y ) 1  n



* ikj j



(16)

Eq.(16) is called the thermal stress natural boundary integral equation. The expression of * H ikj can be found in Appendix A, which is the function of 1 / r 2 and the corresponding integral

8

ACCEPTED MANUSCRIPT

only presents the nearly weakly singularity. Since the maximum singularity is the nearly strong singularity from the integral kernel Fikj* , the singularity in Eq.(16) is reduced by one order in comparison with the nearly hyper-strong singularity in the conventional thermal stress boundary integral equation Eq.(4).

CR IP T

Nevertheless, if the boundary integrals in Eq.(16) are directly calculated by the 8-point Gaussian integration, the results of thermal stresses would become deteriorate when the interior point is approaching to the boundary due to the existence of nearly strongly singular integrals. We

AN US

proposed a regularization algorithm, in which the nearly singular terms are moved outside from the sign of integration by repeating the integration by parts, to deal with the nearly singular integrals [40]. This regularization algorithm is introduced here to evaluate the nearly strongly

M

singular integrals occurred in the thermal stress natural boundary integral equation.

ED

4. Determination of thermal stress intensity factors In the polar coordinate system o centered at a V-notch tip (see Fig.1a), the singular

PT

thermal stress fields of   and   near the vertex of a V-notch can be respectively expressed

AC

CE

with the series asymptotic expansions as follows   ( , ) 

M

K  i

i 1 ~ i( i,



 )

i 1, 3, 5

  (  ,  ) 

M 1

K  i

(17) i 1 ~  i ( i ,



 )

i 2, 4, 6

where  is the radial distance to the notch tip, i is the stress singularity order, M is the number of truncated series item and here is set odd, ~ i (i ,  ) and ~ i (i ,  ) are stress characteristic angular functions, K i denotes the thermal stress intensity factor. For the symmetrical deformation mode ( 1 ),

9

ACCEPTED MANUSCRIPT

K1  lim 2  1 1  (  , ) |  0

(18)

 0

For the anti-symmetrical deformation mode ( 2 ), K 2  lim 2  1 2   (  , ) |  0  0

(19)

Usually, only the first two dominant stress singular orders are located in (0,1) which make the

CR IP T

thermal stress field singular. Taking the logarithm to the dominant terms at two sides of Eq. (17) respectively yields

(20a)

lg    (2  1) lg   lg K 2  lg ~ 2 (2 , )

(20b)

AN US

lg    (1  1) lg   lg K1  lg ~ 1 (1 , )

It can be found from Eq.(20) that lg   ~ lg  and lg   ~ lg  are linear in the stress singularity domain. The slopes of the straight lines expressed by Eq.(20a) and Eq.(20b) are 1  1 and 2  1 , respectively. The intersections of the straight lines lg   ~ lg  and lg   ~ lg 

M

with the vertical axis are lg K1  lg ~ 1 (1 , ) and lg K 2  lg ~ 2 (2 , ) , respectively. By

ED

coupling the singularity asymptotic expansion technique with the interpolating matrix method [41],

PT

we have got the characteristic angular functions ~ 1 (1, ) and ~ 2 (2 , ) for a V-notch. The stress intensity factors K1 and K 2 can be therefore determined from the intersections of the

CE

aforementioned straight lines with the vertical axis.

AC

Once the thermal stresses   and   in the stress singular region are calculated by the

method proposed in Section 2, the stress singularity orders and thermal stress intensity factors can be interpolated from Eq.(20). In order that the nodes selected for calculating the accurate thermal stress intensity factors are within the dominant region of the singular stress field, we firstly give out the figure of lg   ~ lg  and lg   ~ lg  near the notch tip, and then select a number of successive nodes in the interval where lg   ~ lg  and lg   ~ lg  are linear. After that, the

10

ACCEPTED MANUSCRIPT

least square method [42] is used to determine the stress singularity order and related thermal stress intensity factors, where

i  1 

 (lg 

jn

lg  n ) 

n 1

n

n 1

 (lg  )  ( lg  ) 2

n

(i  1, j   ; i  2, j   )

(21a)

2

n

n 1

Ki 

N

jn

n 1 N

N

N

N

 lg  lg  n1

1 exp   N



N

 lg(

jn

 n1i )

n 1

~ ji (i ,  0 )



CR IP T

N

N

(i  1, j   ; i  2, j   )

(21b)

where N is the number of selected inner points;  n and  jn ( j   ,  ) are respectively the

stresses of the inner points are evaluated. 5. Numerical examples

AN US

 and  j ( j   ,  ) of the nth inner point;  0 is the selected direction on which the thermal

The procedure for evaluating the thermal stress intensity factors of the V-notch can be

M

concluded as follows. Firstly, the unknown displacement component ui ( x)(i  1,2) and traction

ED

component ti ( x)(i  1,2) on the structure boundary can be calculated according to the

PT

conventional displacement boundary integral equation Eq.(1). Secondly, the natural boundary variable i ( x)(i  1,2) is solved from the natural boundary integral equation Eq.(7). Then, the

CE

thermal stress components of interior points are calculated by introducing the displacement

AC

component, traction component and natural boundary variable into the thermal stress natural boundary integral equation Eq.(16). Lastly, the stress singularity orders and thermal stress intensity factors of a V-notch are interpolated from Eq.(20). Example 1. A plate subjected to uniform temperature increase To check the efficiency of the proposed method in calculating the thermal stresses of interior points close to the boundary, a vertical-edge-constrained square plate subjected to the uniform

11

ACCEPTED MANUSCRIPT

temperature increase shown in Fig.3 is tested as a benchmark model. The temperature increase

T is set as 1 C ; the elasticity modulus E  210GPa and Poisson's ratio   0.3 ; the thermal expansion coefficient   1.2 10 5  C 1 . There 32 uniform linear boundary elements are discretized on the boundary. The thermal stresses of inner points approaching to the corner point

CR IP T

A(2m,2m) along the diagonal AB are tested. Three different projects are applied to calculate the thermal stress field of interior points. The first one is named CBEM without regularization, in which the Gaussian integration is used to

AN US

calculate the nearly hyper-strongly singular integrals in the conventional stress boundary integral equation Eq.(4). The second one is called CBEM with regularization, in which the regularization algorithm [40] is used to cope with the nearly hyper-strongly singular integrals in the conventional

M

stress boundary integral equation Eq.(4). The last one is denoted as NBEM with regularization, in which the regularization algorithm [40] is used to calculate the nearly strongly singular integrals in

ED

the thermal stress natural boundary integral equation Eq.(16).

PT

The calculated thermal stresses close to the corner point A are illustrated in Fig.4, from which it can be seen that the results from the CBEM without regularization begin to decompose

CE

when x1  x2  1.89m , and the results obtained by the CBEM with regularization become invalid

AC

when x1  x2  1.99999m . On the contrary, the results evaluated by the NBEM with regularization method are still accurate when x1  x2  1.9999999m . It can be concluded that the thermal stress natural boundary integral equation can be used to calculate the thermal stresses of inner points much closer to the boundary by comparison with the conventional boundary element method.

12

ACCEPTED MANUSCRIPT

Example 2. Thermal stress singularity analysis for a symmetrical V-notch Let’s consider a plate weakened by a single edge V-notch, which is constrained at its two ends and subjected to uniform temperature variation ΔT  1 C (see Fig.1a). The thickness of the plate is set at 1mm, the height h  200mm , width w  40mm , the material constants

CR IP T

E  210GPa ,   0.3 and   1.2  105  C1 .

There 180 linear elements are discretized in the boundary element method (BEM) model and 19 interior points along the bisector of the V-notch are chosen for calculating the thermal stresses.

AN US

The finite element method (FEM) is also applied to analyze this model for providing the results for reference. FE analyses were carried out with the commercial FE code Ansys, release 12.1. There 7606 finite elements are totally used in the FEM model, in which the triangle elements are

M

chosen near the notch tip and four node plane elements are adopted far from the notch tip. When the Intel Core i7-4790 CPU is used, the cost time of BEM is 0.852 second, however, FEM cost

ED

3.230 seconds. As it can be seen, the computational effort of the BEM is much smaller than the one of the FEM basing on the accuracy declared hereinafter.

PT

The double logarithmic distributions of lg   versus lg  along the bisector (  0  0 ) and

CE

closing to the V-notch tip are plotted in Fig.5, where three different opening angles are chosen for

AC

examples and the thermal stresses of 19 inner points are calculated by the present boundary element method. It can be seen from Fig.5 that the relationships between lg   and lg  are linear. As have been known, the slopes of these straight lines are the stress singularity orders

1  1 and the intersections of these lines with the vertical axis reflect the thermal stress intensity factors K1 . The calculated stress singularity orders and thermal stress intensity factors for different notch opening angles are listed in Table 1, where the relative error defined by

13

ACCEPTED MANUSCRIPT



absolute (BEM result - FEM result)  100 FEM result

(21)

By using the interpolating matrix method to solve the singularity characteristic differential equations, Ref.[41] also provided the stress singularity orders of a V-notch, which are listed in Table 1 as well for reference. It can be observed that the results by the BEM and FEM are accurate

CR IP T

to at least two digits by contrast with the referenced ones. It can also be found from Table 1 that the maximum relative error between the thermal stress intensity factors from the BEM and the ones from the FEM is smaller than 0.8%.

AN US

Table 2 gives out the calculated thermal stress intensity factors for the V-notches in different notch depth, where two kinds of notch opening angle   30 and   60 are chosen for examples. It can be observed that the maximum relative error between the thermal stress intensity

than 1.0% when   60 .

M

factors by the BEM and the ones by the FEM is smaller than 1.3% when   30 and not bigger

ED

The stress intensity factors of the V-notch under the thermal load (see Fig.1a) are compared

PT

with the ones under the equivalent mechanical load (see Fig.1b) in Table 3, from which it can be found that the difference between the results from thermal load and the ones from equivalent

CE

mechanical load are slight when the notch is not deep, however, this difference is growing with the

AC

increase of the notch depth. When the notch depth is half width of the plate, the stress intensity factors from the thermal load is only half of the ones under the equivalent mechanical load. Example 3. Thermal stress singularity analysis for an inclined V-notch The geometry of an inclined V-notch plate constrained at two ends and subjected to uniform temperature variation ΔT  1 C is shown in Fig.6, where the angle between the notch bisector and horizontal line is denoted as  . The material constants, thickness, width and height of the

14

ACCEPTED MANUSCRIPT

plate are the same as the ones in Example 2. Because the notch is inclined with the load direction, it is a mix-mode model. Table 4 lists the calculated thermal stress intensity factors K1 and K 2 when l w  0.2 , where two kinds of notch opening angles   30 and   60 are respectively chosen for examples. There 180 linear

CR IP T

elements are discretized in the boundary element model and 7852 elements are used in the finite element model. It can be observed from Table 4 that K1 decreases with the inclined angle  , while K 2 increases with  . The relative errors between the results from the BEM and FEM are

AN US

smaller than 1.55% and 0.92% for K1 and K 2 , respectively.

The calculated stress intensity factors for the inclined V-notch under the thermal load and equivalent mechanical load when l / w  0.2 are shown in Table 5, from which it can be seen that the stress intensity factors under the mechanical load are bigger than the corresponding ones

M

under the thermal load. It’s not right to replace the thermal load with the equivalent mechanical

PT

6. Conclusions

ED

load when the thermal stress intensity factors of a sharp V-notch is considered.

The natural boundary integral equations, from which the natural boundary variables can be

CE

yielded, is firstly established basing on the displacement derivative boundary integral equations.

AC

The conventional thermal stress boundary integral equation is then transformed into the thermal stress natural one, in which only the nearly strongly singular integrals exist. After the regularization method is introduced to deal with the remained nearly strongly singular integrals, the thermal stress natural boundary integral equation is successfully applied to calculate the thermal stress field very close to the boundary. At last, the stress singularity orders and thermal stress intensity factors of the V-notch are effectively interpolated basing on the singular thermal

15

ACCEPTED MANUSCRIPT

stress field calculated by the thermal stress natural boundary integral equation. Some conclusions are drawn as follows, 1. The thermal stress natural boundary integral equation can be used to calculate the thermal stress field much closer to the boundary by comparison with the conventional boundary

CR IP T

integral equation. The accurate stress field very close to the notch tip, which is the basis of the stress intensity factor calculation, can be determined by the thermal stress natural boundary integral equation.

AN US

2. As to the calculation of the stress singularity orders and thermal stress intensity factors, the accuracy of the present method is confirmed, while the computational amount is sharply reduced by comparing with the finite element method.

M

3. For a symmetrical V-notch, the thermal stress intensity factors increase with the notch opening angle slightly, but increase with the notch depth sharply.

ED

4. For an inclined V-notch, the thermal stress intensity factor K1 decreases with the inclined

PT

angle between the notch bisector and horizontal line, while K 2 increases with this angle. 5. The stress intensity factor of the V-notch under the thermal load is approximately equal to the

CE

one under the equivalent mechanical load when the notch is not deep, while the difference

AC

between them is growing with the notch depth. Acknowledgments This work was supported by the National Natural Science Foundation of China, China

(No.11372094, No.11272111).

16

ACCEPTED MANUSCRIPT

Appendix A. Expressions of integral kernels The expressions of U ij* , Tij* , Ri and Qi in Eq.(1) are respectively written as U ij*  

(A.1)

1 {(1  2 )( r,i n j  r, j ni )  r,n [(1  2 ) ij  2r,i r, j ]} 4π(1   )r Ri 

 (1   ) 1 1 [(ln  )ni  r,i r, n ] 4π(1   ) r 2

Qi 

 (1   ) 1 1 (ln  )ri 4π(1   ) r 2

where r and its derivatives can be respectively presented as

(A.2) (A.3)

CR IP T

Tij* 

1 [(3  4 ) ln r ij  r,i r, j ] 8π(1   )G

(A.4)

AN US

r  ri ri , ri  xi  yi , r,i  r xi  ri r , r, n  r n  r,i ni

(A.5)

The derivative of U ij* and Tij* with respective to xk can be respectively expressed as 1 [(3  4 )r,k  ij  r,i jk  r, j ki  2r,i r, j r,k ] 8πG(1   )r

(A.6)

1 {2r, n [r,i jk  r, j ki  (1  2 )r, k  ij  4 r, i r, j r, k ]  4π(1   )r 2 (1  2 )( jk ni   ij nk   kin j )  2 r, i r, j nk  2(1  2 )( r, i r, k n j  r, j r, k ni )}

(A.7)

U ij*,k 

M

Tij*, k 

ED

* The expressions of Wikj* , Sikj , Fikj* , Rik , Qik and the initial thermal stress  ik0 ( y) in

Eq.(13) are respectively written as

PT

* Wikj 



1 (1  2 )( r,k  ij  r,i kj  r, j ki )  2r,i r, j r,k 4π(1   )r



G {2 r, n [(1  2 )r, j ki   (r,i jk  r, k ij )  4r, i r, j r, k ]  (1  2 )( 2r,i r, k n j   jk ni 2π(1   )r 2   ij nk )  2 (r, i r, j nk  r, j r, k ni )  (1  4 ) kin j }

CE

* Sikj 

2G [(1   )(emk r,m ij  emi r,m kj )  2e jmr,i r,k r,m   (eij r,k  ekj r,i ) 4 (1   )r  (1  2 )emj  ik r,m ]

AC

Fikj*  

(A.8)

(A.9)

(A.10)

Rik 

G (1   ) r 1 [ (  ik  2r,i r,k )  ni r,k  nk r,i ] 2π(1   )r n 1  2

(A.11)

Qik 

G (1   ) 1 1 1  2 [r,i r,k  (ln  ) ik ] 2π(1   ) 1  2 r 2

(A.12)

 ik0 ( y)  2G

1 T ik 1  2

* The expression of integral kernel H ikj in Eq.(16) is written as

17

(A.13)

ACCEPTED MANUSCRIPT

1 (eli r,l ekm  elk r,l eim  r,k  im  r,i km ) 4r

(A.14)

AC

CE

PT

ED

M

AN US

CR IP T

* H ikj 

18

ACCEPTED MANUSCRIPT

References [1] N.S. Lu, Z. Zhang, J. Yoon, Z.G. Suo, Singular stress fields at corners in flip-chip packages, Eng. Fract. Mech. 86 (2012) 38-47. [2] M. Heinzelmann, Singular thermal stress fields in brazed three layer joints, Eng. Fract. Mech. 55(4) (1996) 647-655. [3] A. Barroso, D. Vicentini, F. Paris, V. Mantic, Representativity of thermal stresses in designing composite joints

CR IP T

based on singular stress states at multimaterial corners, Composites: Part A 42 (2011) 1084-1092.

[4] P. Ferro, F. Berto, P. Lazzarin, Generalized stress intensity factors due to steady and transient thermal loads with applications to welded joints, Fatigue Fract. Engng. Mater. Struct .29 (2006) 440-453.

[5] A.F. Emery, J.A. Williams, J. Avery, Thermal stress concentration caused by structural discontinuities, Exp.

AN US

Mech. 2 (1969) 558-564.

[6] W.K. Wilson, I.W. Yu, The use of the J-integral in thermal stress crack problems, Int. J. Fracture 15(4) (1979) 377-387.

[7] Z.H. Zhou, A.Y.T. Leung, X.S. Xu, X.W. Luo, Mixed-mode thermal stress intensity factors from the finite

M

element discretized symplectic method, Int. J. Solids Struct. 51 (2014) 3798-3806.

[8] M. Nagai, T. Ikeda, N. Miyazaki, Stress intensity factor analysis of a three-dimensional interface crack

ED

between dissimilar anisotropic materials under thermal stress, Eng. Fract. Mech. 91 (2012) 14-36. [9] A.Y.T. Leung, X.S. Xu, Z.H. Zhou, Hamiltonian approach to analytical thermal stress intensity factors-part 2

PT

thermal stress intensity factor, J. Therm. Stresses 33 (2010) 279-301. [10] M. Nagai, T. Ikeda, N. Miyazaki, Stress intensity factor analysis of a three-dimensional interface crack

CE

between dissimilar anisotropic materials under thermal stress, Eng. Fract. Mech. 91(2012) 14-36. [11] X.C. Ping, M.C. Chen, B.B. Zheng, B. Xu, An effective numerical analysis of singular stress fields in

AC

dissimilar material wedges under thermo-mechanical loads, Eng. Fract. Mech. 106(2013) 22-37. [12] Y. Nomura, T. Ikeda, N. Miyazaki, Stress intensity factor analysis at an interfacial corner between anisotropic bimaterials under thermal stress, Eng. Fract. Mech. 76(2009) 221-235.

[13] H.J. Ding, N.L. Peng, The analysis of thermal residual stresses near the apex in bonded dissimilar materials, Int. J. Solids Struct. 36 (1999) 5611-5637. [14] Z. Yosibash, Thermal generalized stress intensity factors in 2-D domains, Comput. Methods Appl. Mech. Engrg. 157 (1998) 365-385. [15] E.V. Glushkov, N.V. Glushkova, D. Munz, Y.Y. Yang, Analytical solution for bonded wedges under thermal 19

ACCEPTED MANUSCRIPT

stresses, Int. J. Fracture 106 (2000) 321-339. [16] Y.Y. Yang, D. Munz, Stress singularities in a dissimilar materials joint with edge tractions under mechanical and thermal loadings, Int. J. Solids Struct. 34(10) (1997) 1199-1216. [17] Y. Nomura, T. Ikeda, N. Miyazaki, Stress intensity factor analysis of a three-dimensional interfacial corner between anisotropic bimaterials under thermal stress, Int. J. Solids Struct. 47 (2010) 1775-1784.

boundary element method, Eng. Anal. Bound. Elem. 37 (2013) 116-127.

CR IP T

[18] Y. Ochiai, V. Sladek, J. Sladek, Three-dimensional unsteady thermal stress analysis by triple-reciprocity

[19] J. Sladek, V. Sladek, Evaluation of T-stress and stress intensity factors in stationary thermoelasticity by the conservation integral method, Int. J. Fracture 86 (1997) 199-219.

[20] C.Z. Cheng, Z.R. Niu, N. Recho, Analysis of the stress singularity for a bi-material V-notch by the boundary

AN US

element method, Appl. Math. Model. 37 (2013) 9398-9408.

[21] S.S. Lee, Boundary element analysis of singular hygrothermal stresses in a bonded viscoelastic thin film, Int. J. Solids Struct. 38 (2001) 401-412.

[22] P.H. Tehrani, A.R.H. Godarzi, M. Tavangar, Boundary element analysis of stress intensity factor KI in some

M

two-dimensional dynamic thermoelastic problems, Eng. Anal. Bound. Elem. 29 (2005) 232-240. [23] K. Yang, W.Z. Feng, H.F. Peng, J. Lv, A new analytical approach of functionally graded material structures

ED

for thermal stress BEM analysis, Int. Commun. Heat Mass Transfer 62 (2015) 26-32. [24] C.Z. Cheng, Z.L. Han, S.L. Yao, Z.R. Niu, N. Recho, Analysis of heat flux singularity at 2D notch tip by

1-9.

PT

singularity analysis method combined with boundary element technique, Eng. Anal. Bound. Elem. 46 (2014)

CE

[25] M. Prukvilailert, H. Koguchi, Boundary element analysis of the stress field at the singularity lines in three-dimensional bonded joints under thermal loading, J. Mech. Mater. Struct. 2(1) (2007) 149-166.

AC

[26] L.K. Keppas, N.K. Anifantis, Fatigue life prediction under cyclic thermal loads using the boundary elements method for two-dimensional problems, Comp. Struct. 89 (2011) 590-598.

[27] C.Z. Cheng, Z.L. Han, H.L. Zhou, Z.R. Niu, Analysis of the temperature field in anisotropic coating-structures by the boundary element method, Eng. Anal. Bound. Elem. 60 (2015) 115-122. [28] J.H. Lv, Y. Miao, H.P. Zhu, The distance sinh transformation for the numerical evaluation of nearly singular integrals over curved surface elements, Comput. Mech. 53(2014) 359-367. [29] Y. Gu, W. Chen, B. Zhang, W.Z. Qu, Two general algorithms for nearly singular integrals in two dimensional anisotropic boundary element method, Comput. Mech. 53 (2014) 1223-1234. 20

ACCEPTED MANUSCRIPT

[30] G.Z. Xie, J.M. Zhang, X.Y. Qin, G.Y. Li, New variable transformations for evaluating nearly singular integrals in 2D boundary element method, Eng. Anal. Bound. Elem. 35 (2011) 811-817. [31] Y.J. Liu, Analysis of shell-like structures by the boundary element method based on 3-D elasticity: formulation and verification, Int. J. Numer. Methods Eng. 41(3) (1998) 541-558. [32] M. Dehghan, H. Hosseinzadeh, Calculation of 2D singular and near singular integrals of boundary elements method based on the complex space C, Appl. Math. Model. 36 (2012) 545-560.

Meth. Engng. 60 (2004) 1075-1102.

CR IP T

[33] L.B. Sills, C. Ishbir, A conservative integral for bimaterial notches subjected to thermal stresses, Int. J. Numer.

[34] M. Dehghan, H. Hosseinzadeh, Improvement of the accuracy in boundary element method based on high-order discretization, Comput. Math. Appl. 62 (2011) 4461-4471.

AN US

[35] H. Hosseinzadeh, M. Dehghan, A new scheme based on boundary elements method to solve linear Helmholtz and semi-linear Poisson's equations, Eng. Anal. Bound. Elem. 43 (2014) 124-135.

[36] H. Hosseinzadeh, M. Dehghan, A simple and accurate scheme based on complex space C to calculate boundary integrals of 2D boundary elements method, Comput. Math. Appl. 68 (2014) 531-542.

M

[37] M. Dehghan, M. Shirzadi, The modified dual reciprocity boundary elements method and its application for solving stochastic partial differential equations, Eng. Anal. Bound. Elem. 58 (2015) 99-111.

ED

[38] M. Dehghan, A. Ghesmati, Solution of the second-order one-dimensional hyperbolic telegraph equation by using the dual reciprocity boundary integral equation (DRBIE) method, Eng. Anal. Bound. Elem. 34 (2010)

PT

51-59.

[39] C.Z. Cheng, Z.R. Niu, N. Recho, Z.Y. Yang, R.Y. Ge, A natural stress boundary integral equation for

CE

calculating the near boundary stress field, Comput. Struct. 89 (2011) 1449-1455. [40] Z.R. Niu, C.Z. Cheng, H.L. Zhou, Z.J. Hu, Analytic formulations for calculating nearly singular integrals in

AC

two-dimensional BEM, Eng. Anal. Bound. Elem. 31 (2007) 949-964. [41] Z.R. Niu, D.L. Ge, C.Z. Cheng, J.Q. Ye, Evaluation of the stress singularities of plane V-notches in bonded dissimilar materials, Appl. Math. Model. 33 (2009) 1776-1792.

[42] Y.H. Liu, Z.G. Wu, Y.C. Liang, X.M. Liu, Numerical methods for determination of stress intensity factors of singular stress field, Eng. Fract. Mech. 75 (2008) 4793-4803.

21

ACCEPTED MANUSCRIPT

Figures

 0  E  T

T h/2



l

 





O

h/2

x2 l

  x1 (bisector )



CR IP T

x2

O

h/2

x1 (bisector )

h/2

w

(a) Under thermal load

AN US

w

(b) Under equivalent mechanical load

AC

CE

PT

ED

M

Fig.1 Single edge V-notched plate

22

ACCEPTED MANUSCRIPT x2

 x2

n

y

2

1 x1

O

CR IP T

x1

AC

CE

PT

ED

M

AN US

Fig.2 Nearly singular integral occurred boundary 1

23

ACCEPTED MANUSCRIPT

x2

A(2m,2m)  22

 12

4m

O

 11

x1

B(2m,2m)

4m

CR IP T

T

AC

CE

PT

ED

M

AN US

Fig.3 A vertical-edge-constrained square plate subjected to uniform temperature increase

24

ACCEPTED MANUSCRIPT

-1.8

CBEM without Regularization CBEM with Regularization NBEM with Regularization Exact solutions

-2.0

-2.2

(MPa)

-2.4

-2.6

x1=x2(m)

1.9999999

1.9999998

1.9999990

1.9999900

1.9999000

1.9990000

1.9900000

1.9800000

1.9000000

1.8900000

1.8800000

-3.2

1.8000000

-3.0

CR IP T

-2.8

AC

CE

PT

ED

M

AN US

Fig.4 Thermal stresses of interior points approaching to point A along diagonal AB

25

ACCEPTED MANUSCRIPT

2.8 0

30 60 90

2.7

0



2.6

lg

2.5 2.4 2.3

2.1 -4.0

-3.9

-3.8

-3.7

-3.6

-3.5

lg

-3.4

-3.3

-3.2

CR IP T

2.2

-3.1

-3.0

Fig.5 ln  vs. ln  along the bisector of symmetrical V-notches by the BEM

AC

CE

PT

ED

M

AN US

(Unit of   is MPa and  is mm )

26

ACCEPTED MANUSCRIPT

x2 bisector





T

w

h/2

x1

h/2

CR IP T

l

AC

CE

PT

ED

M

AN US

Fig.6 Inclined V-notch plate subjected to uniform temperature increase

27

ACCEPTED MANUSCRIPT

Tables Table 1 Stress singularity orders and thermal stress intensity factors of symmetrical V-notches when l w  0.2

1  1

( ) 

K1 (N  mm21 )

BEM

FEM

 (%)

BEM

FEM

 (%)

10

-0.4999

-0.4984

-0.4960

0.4839

15.8807

16.0031

0.7648

20

-0.4996

-0.4966

-0.4952

0.2827

15.9684

16.0175

0.3065

30

-0.4985

-0.4965

-0.4964

0.0201

16.0342

16.0755

0.2569

40

-0.4965

-0.4952

-0.4935

0.3445

16.1166

16.2014

0.5234

50

-0.4931

-0.4923

-0.4918

0.1017

16.2256

16.2178

0.0481

60

-0.4878

-0.4872

-0.4868

0.0822

16.3547

16.3571

0.0147

70

-0.4801

-0.4801

-0.4800

0.0208

16.5146

16.4816

0.2002

80

-0.4696

-0.4697

-0.4698

0.0213

16.6849

16.6401

0.2692

90

-0.4555

-0.4557

-0.4545

0.2640

16.8461

16.8980

0.3071

AC

CE

PT

ED

M

AN US

CR IP T

Ref.[41]

28

ACCEPTED MANUSCRIPT

Table 2 K1 (N  mm21 ) of symmetrical V-notches in different l w

  30

l w

  60

FEM

 (%)

BEM

FEM

 (%)

0.10

10.9333

10.8763

0.5243

11.1578

11.0968

0.5497

0.20

16.0342

16.0755

0.2569

16.3547

16.3571

0.0147

0.30

20.6557

20.7930

0.6606

20.9924

21.0999

0.5095

0.40

24.5816

24.8331

1.0123

24.9258

25.1226

0.7895

0.50

27.3247

27.6819

1.2903

27.6943

CR IP T

BEM

AC

CE

PT

ED

M

AN US

27.9649

29

0.9676

ACCEPTED MANUSCRIPT Table 3 K1 (N  mm21 ) of symmetrical V-notches under thermal load and equivalent mechanical load

  30 l w

  60

Equivalent mechanical load

Thermal load

Equivalent mechanical load

0.10

10.9333

10.7523

11.1578

11.1621

0.20

16.0342

17.3225

16.3547

17.7663

0.30

20.6557

25.4850

20.9924

26.1012

0.40

24.5816

37.5084

24.9258

38.2458

0.50

27.3247

56.0148

27.6943

56.8520

AC

CE

PT

ED

M

AN US

CR IP T

Thermal load

30

ACCEPTED MANUSCRIPT Table 4 Thermal stress intensity factors of inclined V-notches when l w  0.2

K 2 /(N  mm22 )

K1/(N  mm21 )

 /  ( )

FEM

 (%)

BEM

FEM

 (%)

30/ 0

16.0342

16.0853

0.3175

/

/

-

30/15

15.4860

15.6573

1.0941

2.8418

2.8499

0.2842

30/30

13.9573

14.1080

1.0682

5.1897

5.1854

0.0829

30/45

11.6124

11.7949

1.5473

6.6205

6.6086

0.1801

60/ 0

16.3547

16.4783

0.7501

/

60/15

15.7356

15.8415

0.6685

60/30

13.9551

14.0569

60/45

11.1470

11.0783

CR IP T

BEM

/

-

2.9189

2.9290

0.3448

0.7242

5.1842

5.1915

0.1406

0.6201

6.2435

6.1868

0.9165

AC

CE

PT

ED

M

AN US

Note: the slash means no stress intensity factors and the bar means no relative error.

31

ACCEPTED MANUSCRIPT

Table 5 Stress intensity factors of inclined V-notches under thermal load and equivalent mechanical load when

l w  0.2

K1/(N  mm21 )

 /  ( )

K 2 /(N  mm22 )

Equivalent mechanical load

Thermal load

Equivalent mechanical load

30/ 0

16.0342

17.2071

/

/

30/15

15.4860

16.7895

2.8418

3.0404

30/30

13.9573

15.1525

5.1897

5.5506

30/45

11.6124

12.7308

6.6205

7.1875

60/ 0

16.3547

17.7663

/

/

60/15

15.7356

17.3051

2.9189

3.1823

60/30

13.9551

15.4156

5.1842

5.6413

60/45

11.1470

12.6368

6.2435

6.7957

AC

CE

PT

ED

M

AN US

Note: the slash means no stress intensity factors.

CR IP T

Thermal load

32