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 j2 ( x ) Fikj* j1 ( 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[ m2 ( x ) nm1 ( 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 n1
1 exp N
N
lg(
jn
n1i )
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 105 C1 .
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 ) 4r
(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 mm21 )
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 mm21 ) 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 mm21 ) 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 mm22 )
K1/(N mm21 )
/ ( )
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 mm21 )
/ ( )
K 2 /(N mm22 )
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