Analysis of Two- and Three-dimensional Steady-state Thermo-mechanical Problems Including Curved Line/Surface Heat Sources Using the Method of Fundamental Solutions

M. Mohammadi1,∗, M. R. Hematiyan2 and A. Khosravifard2

1Department of Mechanical Engineering, College of Engineering, Shiraz Branch,
Islamic Azad University, Shiraz, Iran
2Department of Mechanical Engineering, Shiraz University, Shiraz 71936, Iran
E-mail: mehr4457@gmail.com; mohammadi mehrdad@iaushiraz.ac.ir
*Corresponding Author

Received 05 March 2019; Accepted 28 June 2019;
Publication 05 August 2019

Abstract

In this work, the method of fundamental solutions (MFS) and the method of particular solutions (MPS) are used to solve two and three dimensional steady-state thermoelastic problems involving curved-shape heat sources. The geometrical shape of the heat sources can be very complicated. Each curved heat source is modelled by assembling several simple sources with quadratic shapes. The particular solutions for temperature and stress are presented in simple forms and they are used without considering any internal points or internal cells. Several examples are analysed to demonstrate the efficiency of the presented formulation. Numerical results show that the presented MFS-MPS formulation is very efficient and useful. Unlike the finite element method, only a small number of collocation and source points are sufficient to achieve accurate results in the proposed MFS formulation.

Keywords: Method of fundamental solutions, Curved line heat source, Curved surface heat source, Thermo-mechanical, Meshfree

1 Introduction

In practical thermo-mechanical problems such as infrared heating, electrical heating and laser beam welding, it is necessary to consider concentrated heat sources for accurate modelling of such systems. Only a few analytical solutions are presented yet for these types of problems. For example, [1] presented a general analytical solution for an annular involving a point heat source. [2] presented Green’s functions for some boundary value problems in an infinite plane with a point heat source, and [3] presented Green’s functions for a multi-material domain including a point heat source. Despite these analytical solutions, for complicated practical problems involving concentrated heat sources, numerical methods should be used.

The finite element method (FEM) is a popular numerical method for solving large majority of thermo-mechanical problems. But in the FEM, a concentrated heat source should be modeled as a separated region with a small area. The size of this area affects the accuracy of the analysis and a large number of internal elements and nodes should be used, especially near the heat source region. As a powerful alternative, the boundary element method (BEM) only requires boundary discretization. The BEM, which is based on fundamental solutions, has been used effectively to solve direct and inverse thermo-mechanical problems containing concentrated heat sources. [4] proposed a BEM approach for identification of intensity of point heat sources in time-dependent heat conduction. He tested his proposed method by conducting some 2D experiments [5]. [6] formulated line and distributed heat sources in the BEM analysis of nonlinear heat conduction problems in two and three dimensions. They also used the proposed method for inverse non-linear heat conduction analysis [7]. There is no need for domain discretization in their proposed method. [8] used the BEM for thermoelastic analysis of two-dimensional anisotropic media involving concentrated heat sources. They converted domain integrals to boundary integrals in their formulation. In another work, [9] presented a BEM approach for heat conduction in anisotropic non-homogeneous media involving internal point heat sources. They used the direct domain mapping (DDM) technique to enable the exact volume integral transformation. [10] presented a BEM formulation with no need to any internal points or cells for 2D and 3D thermo-elastic problems including point, straight line and flat surface heat sources. They considered a linear variation for the strength of heat sources in their formulation. [11] implemented curved line heat sources with arbitrary variations in 2D and 3D isotropic domains by considering several quadratic heat sources.

Another efficient numerical method is the method of fundamental solutions (MFS), which is used in the present work. Similar to the BEM, the MFS is a boundary-type technique and it is applicable when a fundamental solution of the problem is known. However, the MFS does not require any mesh and it is an integration-free method. Therefore, developing computer codes for the MFS is much simpler than the BEM and it provides very good results with low computational cost [12, 13]. The MFS is a simple method in which the solution of the problem is expressed in terms of fundamental solutions with some unknown coefficients. The singular points of the fundamental solutions are positioned on a boundary outside the body. The intensities of the MFS sources are unknown and are determined by satisfying boundary conditions on some collocation points on the boundary of the body. After the first numerical implementation of the MFS by [14], various applications of this method have been presented in the literature. An overview can be found in the works of [15–18].

The governing equations of the steady-state thermo-elasticity involving heat sources are the heat conduction equation including the heat source term (a standard Poisson equation) and the Navier equations with thermal effects. To solve these equations using the MFS, particular solutions corresponding to heat source terms in addition to homogeneous solutions of the equations should be found. The MFS in conjunction with the dual reciprocity method (DRM) [19] was used by [20] for solving Poisson’s equation, by [21] to solve two-dimensional thermo-elasticity with general body forces, by [22] for thermo-mechanical analysis of functionally graded materials and by [23] to solve three-dimensional thermo-elasticity with arbitrary body forces. Another similar approach is the MFS in combination with the so called method of particular solutions (MPS). The MFS-MPS was used by [24] to solve Poisson’s equation, by [25] for thermo-elastic problems in axisymmetric domains, by [26] to solve two-dimensional thermoelasticity, by [27] for three-dimensional inverse thermoelasticity, by [28] to solve two-dimensional thermoelasticity with the non-singular MFS (NMFS), and by [13] for steady state heat conduction in 2D and 3D isotropic domains involving curved heat sources.

To the authors’ best knowledge, the MFS formulation for thermo-elastic problems involving internal curved heat sources has not been presented yet. The method proposed by [13] for heat conduction analysis is developed here for thermo-mechanical problems. In the present work, the MFS-MPS is used for 2D and 3D analysis of steady-state thermo-mechanical problems involving curved line/surface heat sources.The particular solutions associated with the heat source effect are presented in simple forms. They are calculated without considering any internal points or internal cells and therefore the attractiveness of the MFS as a boundary-type mesh-free method is preserved. Several numerical example problems are presented to demonstrate the advantages of the proposed MFS in comparison with the BEM and the finite element method (FEM).

2 Basic Equations and Formulation

images

Figure 1: Domain, main boundary, pseudo boundary, and a curved line heat source.

Consider a medium Ω with its boundary Γ (Figure 1). In the presence of heat sources, the governing equation of steady-state heat conduction can be expressed as

2τ(x) = -s(x) k inΩ (1)

The boundary conditions can be written as

f1(xb)τ(xb) + f2(xb)τ(xb) n = f3(xb)on Γ (2)

where f1, f2 and f3 are given functions on the boundary, 2 is the Laplace operator, n is the direction normal to the boundary, x is a point in the domain, xb is a point on the boundary, τ is the temperature, s(x) is a known function describing the heat source distribution in the domain and k is the thermal conductivity.

Navier equations for two and three dimensional problems with consideration of the thermal effects are given as:

G2u i(x) + G1 - 2ν̄ j=1d 2u j(x)xixj -γ̄τ(x)xi = 0 (3)

where d is the dimension of the problem, i.e. d = 2 for a 2D problem and d = 3 for a 3D problem, and i = 1, 2 for 2D problems and i = 1, 2, 3 for 3D problems.

Mechanical boundary conditions can be expressed as follows:

u i ( x b ) = g i ( x b ) o n   Γ 1                             ( 4 ) t i ( x b ) = h i ( x b ) o n   Γ 2

where Γ 1 Γ 2 = Γ . u i       a n d       t i are the components of displacement vector and traction vector, respectively, G is the modulus of rigidity and g i and hi are prescribed functions on the boundary. In Equation (3) γ ¯ = 2 G α ¯ ( 1 + ν ¯ ) / 2 G α ¯ ( 1 + ν ¯ ) ( 1 2 ν ¯ ) ( 1 2 ν ¯ ) , where

ν ¯ = {     ν p l a n e s t r a i n a n d 3 D ν / ( 1 + ν ) p l a n e s t r e s s α ¯ = {     α p l a n e s t r a i n a n d 3 D α ( 1 + ν ) / ( 1 + 2 ν ) p l a n e s t r e s s (5)

where α represents the coefficient of linear thermal expansion and ν is the Poisson’s ratio.

The solution of Equation (1) in the MFS is represented as follows:

τ(x) = l=1Na lτ*(x,ξ l)) + τp(x) (6)

where ξl and al are the location and intensity of the lth source point on the pseudo boundary Γ, respectively, and N is the number of source points. The fundamental solution of the Laplace equation is denoted by τ* and is expressed as

τ*(x,ξ l) = -1 2π ln(r(x,ξl)) for 2D (7) τ*(x,ξ l) = 1 4πr(x,ξl) for 3D

where r is the magnitude of the vector r, which represents the distance between the source point ξl to the point x.

τp(x) in Equation (6) is the particular solution and it can be obtained by constructing the associated Newton potential in the following domain integral form [15]:

τp(x) = 1 kΩs(ξ) τ*(x,ξ) dV (ξ) (8)

Considering N collocation points on the main boundary and satisfying the thermal boundary conditions at these points, the coefficients al can be determined.

Considering thermal effects in the Navier equation, the solution of this equation can be represented as follows:

ui(x) = uih(x) + u ip(x) (9)

The components of the homogeneous displacement vector (solution of the Navier equation without thermal terms) are expressed as follows [21]:

uih(x) = l=1N c lU1i*(x,ξ l) + cl+NU2i*(x,ξ l)(for 2D) (6) uih(x) = l=1N c lU1i*(x,ξ l) + cl+NU2i*(x,ξ l) + cl+2NU3i*(x,ξ l)(for 3D)

(10)

where:

Uij* = 1 8π(1-ν̄)G (3 - 4ν̄)ln 1 r δij + r,ir,j (for 2D) (8)U˙ijˆ* = 1 16π(1 -ν̄)Gr (3 - 4ν̄)δij + r,ir,j (for 3D)

(11)

in which

r,i = ri r (12)

where ri are the elements of the vector r.

The coefficients cl in Equation (10) are unknowns, which can be found by applying boundary conditions. The components of the particular displacement vector (solution of the Navier equations with thermal terms) can be represented by the following equation:

uip(x) = 1 Ā l=1N a lri(x,ξl)τ*(x,ξ l)+ 1 kĀΩs(ξ) ri(x,ξ)τ*(x,ξ) dV (ξ ) (13)

where:

Ā = 4G(1 -ν̄) γ̄(1 - 2ν̄) (14)

The components of the traction vector at boundary points can be calculated as follows:

ti(xb) = tih(x b) + tip(x b) (15)

where tih and tip are components of the homogeneous traction vector and the particular traction vector, respectively. The components of homogeneous traction vector can be derived from the following equations [21]:

tih(x b) = l=1N clP1i*(xb,ξl) + cl+NP2i*(xb,ξl)(for 2D) (14)t˙iˆh (x˙b) = l=1N c lP1i*(x b,ξl) + cl+NP2i*(x b,ξl) + cl+2NP3i*(x b,ξl) (15) (for 3D) (16)

where: P˙ijˆ*=-1 2d π(1-ν̄)r(d-1) ((1-2ν̄)δij+dr,jr,i) r n-(1-2ν̄) (r,inj-r,jni)

(17)

in which, ni are the elements of the unit vector n, normal to the boundary.

The components of the particular traction vector can be calculated as follows:

tip(x b) = σijp(x b)nj(xb) (18)

In order to calculate the components of the particular stress tensor σijp at a point, the following stress-strain relation can be used:

σijp(x) = 2G ε ijp(x) + ν̄ 1 - 2ν̄εkkp(x)δ ij -γ̄τ(x)δij (19)

The components of the particular strain tensor εijp at a point can be determined using the strain-displacement relation as

εijp(x) = 1 2 ui,jp(x) + u j,ip(x) (20) = 1 Ā l=1N a l τ*(x,ξ l)δij -ri(x,ξl)rj(x,ξl) 2(d-1)πrd(x,ξl) (21) + 1 kĀΩs(ξ) τ*(x,ξ)δ ij - ri(x,ξ)rj(x,ξ) 2(d-1)πrd(x,ξ) dV (ξ) (20)

For calculation of stress tensor components at a point of the domain, the following equation can be used:

σij(x) = σijh(x) + σ ijp(x) (21)

In which, the components of the homogeneous stress tensor σijh can be calculated as follows [21]:

σijh(x) = l=1N c lΨ1ij*(x,ξ l) + cl+NΨ2ij*(x,ξ l) (for 2D) (24)σij h (x ) = l=1N c lΨ1ij*(x,ξ l) + cl+NΨ2ij*(x,ξ l) + cl+2NΨ3ij*(x,ξ l) (25) (for 3D) (22)
where: Ψijk* = -1 2dπ(1-ν̄)r(d-1) (1 - 2ν̄) r,jδik + r,kδij - r,iδjk + dr,ir,jr,k (23)

As it can be seen in Equations (8), (13) and (20), the particular solutions associated with the heat source effect are represented by domain integrals. In the next section, these integrals are calculated for curved line/surface heat sources without any need for internal cells or internal points.

3 Formulations for Heat Sources Concentrated on a Curved Line/Surface

In this section, the particular solutions associated with curved line/surface heat sources in the MFS are presented. For the 2D case, curved line sources and for the 3D case, both curved line and curved surface sources are considered.

3.1 Curved Line Heat Source in Two-dimensional Problems

Any curved line heat source can be discretized by several quadratic segments. As shown in Figure 2, each segment has a quadratic shape and is represented by three points. The domain integrals in Equations (8), (13) and (20) for the 2D case, are converted to integrals over the quadratic segment as follows:

IT(x) = -1 2πklsources ξln r x,ξdl Iu(x) = -1 2πkĀlsources(ξ) ri(x,ξ)ln r x,ξ dl Iε(x) = -1 2πkĀlsources(ξ) ln r x,ξδij + ri(x,ξ)rj(x,ξ) r2(x,ξ) dl (24)

where dl is an infinitesimal element on the quadratic segment, r(x,ξ) = (x - xs )2 + (y - ys )2 is the distance between the field point x = (x,y) and the source points ξ = (xs,ys) on the quadratic line heat source, r1(x,ξ) = (x - xs), and r2(x,ξ) = (y - ys).

images

Figure 2: A curved line heat source with quadratic shape and its dimensionless local coordinate.

The coordinates of the point ξ = (xs,ys) on the quadratic segment can be expressed as

xs = N1x1 + N2x2 + N3x3 ys = N1y1 + N2y2 + N3y3 (25)

where (x1,y1), (x2,y2) and (x3,y3) are respectively, the coordinates of the starting point, the middle point and the end point of the quadratic segment and N1, N2 and N3 are the quadratic Lagrangian shape functions. These shape functions can be expressed as follows [29]:

N1 = 1 2η(η - 1) N2 = -(η + 1)(η - 1) N3 = 1 2η(η + 1) (26)

where - 1 η +1 is a dimensionless local coordinate shown in Figure 2. Using Equation (25), the integrals in Equation (24) are converted to the following integrals:

IT(x) = -1 2πk-11s(η)ln r(η)Jdη Iu(x) = -1 2πkĀ-11s(η)ri(η)ln r(η)Jdη Iε(x) = -1 2πkĀ-11s(η) ln r(η)δij + ri(η)rj r2(η) Jdη (27)

where J is the Jacobian of transformation.

The integrals in Equation (27) can be calculated using conventional numerical integration methods such as the Gaussian quadrature method (GQM). It should be noted that if the field point x = (x,y) is exactly on the line source, the integrals in Equation (27) will be weakly singular with a finite value. In other words, in the two-dimensional case, stress, displacement and temperature have finite values at points exactly on the curved line heat source. In this case, the integral in Equation (27) can be calculated by various methods such as the weighted Gaussian integration [30], transformation of variable [31] and subtraction of singularity [32] methods. In this research, weighted Gaussian integration method has been used.

3.2 Curved Line Heat Source in 3D Problems

Similar to the 2D case, any curved line heat source in 3D is discretized by some heat sources with quadratic shapes. The domain integrals in Equations (8), (13) and (20) for the 3D case, are converted to integrals on the quadratic segment as follows:

IT(x) = 1 4πklsource sξ rx,ξdl Iu(x) = 1 4πkĀlsources(ξ) ri(x,ξ) rx,ξ dl Iε(x) = 1 4πkĀlsources(ξ) δij rx,ξ -ri(x,ξ)rj(x,ξ) r3(x,ξ) dl (28)

where r(x,ξ) = (x - xs )2 + (y - ys )2 + (z - zs )2 is the distance between the field point x = (x,y,z) and the source points ξ = (xs,ys,zs) on the quadratic line heat source, r1(x,ξ) = (x - xs), r2(x,ξ) = (y - ys), and r3(x,ξ) = (z - zs).

Each quadratic segment is discretized by three points, i.e. (x1,y1,z1) at the starting point, (x2,y2,z2) at the middle point, and (x3,y3,z3) at the end point. The coordinates of the points ξ = (xs,ys,zs) on the quadratic segment can be expressed as follows:

xs = N1x1 + N2x2 + N3x3 ,ys = N1y1 + N2y2 + N3y3 , zs = N1z1 + N2z2 + N3z3 (29)

where N1, N2 and N3 are the quadratic Lagrangian shape functions as given in Equation (26).

Using Equation (29), the integrals in Equation (28) are converted to the following integrals:

IT(x) = 1 4πk-11s(η) r(η) JdηIu(x) = 1 4πkĀ-11s(η)ri(η) r(η) Jdη Iε(x) = 1 4πkĀ-11s(η) δij r(η) -ri(η)rj(η) r3(η) Jdη (30)

The integrals in Equation (30) can be calculated using conventional numerical integration methods such as the Gaussian quadrature method (GQM). It should be noted that if the field point x = (x,y,z) is exactly on the curved line heat source, the first and the third integrals in Equation (30) will be strongly singular without any finite values, but the second integral in Equation (30) will be weakly singular with a finite value. In other words, in the three-dimensional case, at a point on a curved line heat source, the displacement has a finite value but the temperature and the stress go to infinity.

3.3 Curved Surface Heat Source in 3D Problems

We consider an arbitrary heat source distributed over a curved surface in a 3D domain. Each surface heat source can be divided into several quadrilateral surface heat sources. Each quadrilateral surface heat source has a quadratic shape. A curved surface heat source and one of its quadratic segments are shown in Figure 3.

images

Figure 3: An arbitrary curved surface heat source and a quadratic segment with its local coordinates.

The domain integrals in Equations (8), (13) and (20) for the 3D case, are converted to integrals on the surface of the quadratic segment as follows:

IT(x) = 1 4πkAsource sξ rx,ξdA Iu(x) = 1 4πkĀAsources(ξ) ri(x,ξ) rx,ξ dA Iε(x) = 1 4πkĀAsources(ξ) δij rx,ξ -ri(x,ξ)rj(x,ξ) r3(x,ξ) dA (31)

where dA is an infinitesimal area element on the quadratic surface heat source and r(x,ξ) is the distance between the field point x = (x,y,z) and the source points ξ = (xs,ys,zs) on the quadratic surface heat source. As shown in Figure 3, each quadratic surface heat source is represented by 8 points and a set of local coordinates (ξ1,ξ2), which vary between -1 and 1. xs, ys, and zs can be expressed in terms of the 8 Lagrangian shape functions as follows:

xs(ξ1,ξ2) = i=18Ni(ξ1,ξ2)xi, ys(ξ1,ξ2) = i=18Ni(ξ1,ξ2)yi, zs(ξ1,ξ2) = i=18Ni(ξ1,ξ2)zi

(32)

where (xi,yi,zi) are the coordinates of the 8 points on the heat source. The shape functions N i in terms of ξ1 and ξ2 are expressed as follows [29]:

N1 = -1 4 (1 - ξ1)(1 - ξ2)(1 + ξ1 + ξ2), N2 = 1 2(1 - ξ12)(1 - ξ 2), (6) N˙3 = - 1 4 (1 + ξ1)(1 - ξ2)(1 - ξ1 + ξ2), N4 = 1 2(1 + ξ1)(1 - ξ22), (7) N˙5 = - 1 4 (1 + ξ1)(1 + ξ2)(1 - ξ1 - ξ2), N6 = 1 2(1 - ξ12)(1 + ξ 2), (8) N˙7 = - 1 4 (1 - ξ1)(1 + ξ2)(1 + ξ1 - ξ2), N8 = 1 2(1 - ξ1)(1 - ξ22). (33)

Using Equation (32), the integrals in Equation (31) are converted to the following integrals:

IT(x) = 1 4πk-11-11s(ξ1,ξ2) r(ξ1,ξ2)Jdξ1dξ2 Iu(x) = 1 4πkĀ-11-11s(ξ1,ξ2)ri(ξ1,ξ2) r(ξ1,ξ2) Jdξ1dξ2 Iε(x) = 1 4πkĀ-11-11s(ξ1,ξ2) δij rξ1,ξ2 -ri(ξ1,ξ2)rj(ξ1,ξ2) r3(ξ1,ξ2) Jdξ1dξ2 (34)

The integrals in Equation (34) can be calculated using standard 2D numerical integration methods such as the GQM. In the cases that the field point X is exactly on the surface of the heat source, the integrals in Equation (34) will be weakly singular with a finite value. In other words, in the three-dimensional case, the temperatures, displacements and stresses at points on a surface heat source have finite values. In this case, the integrals in Equation (34), which is weakly singular, can be calculated by various methods such as the transformation of variable and subtraction of singularity methods [32]. In this research, the method of transformation of variable has been used for such cases.

4 Numerical Examples

To verify the proposed MFS formulation, four thermo-mechanical examples are presented. In two-dimensional examples, curved line heat sources and in three-dimensional examples curved line and curved surface heat sources are considered. In the MFS modelling of all examples, the location of source points and collocation points are determined according to the procedure proposed by [33]. In this procedure, a “location parameter” for a source point is defined as the ratio of the distance from the source point to its corresponding collocation point to the distance from the same source point to the neighbouring collocation point [33]. It is recommended that the value of the location parameter be greater than a lower bound to obtain solutions without undesired oscillations. In this work, the locations of source points are determined in a way that the value of the location parameter of source points would be greater than 0.85.

In all examples, the results computed by the present MFS are compared with those of the FEM with a fine mesh as a reference solution. Also in the first example, the BEM results are included. The software package MATLAB is used for writing MFS and BEM source codes and ANSYS package is used for the FE analyses. In all examples, it is assumed that k = 60 WmC, α = 11.7 × 10-6 1C, E = 200 GPa and ν = 0.3.

images

Figure 4: The boundary conditions and the curved line heat sources in a circular domain.

4.1 Two-dimensional Examples

In the two-dimensional examples the state of plane strain is considered. The first 2D example deals with a circular domain with radius R = 0.5 m, which is centred at (0,0). According to Figure 4, boundary conditions are τB = 10C and ux = uy = 0. The domain involves two semi-circular heat sources with radius of r = 0.125m and centres at (0.125,0) and (-0.125,0). Their intensity functions are as follows:

s1 = 16000 x + 4000 Wm, s2 = -16000 x + 4000 Wm (35)

In the MFS modelling of the problem, only 8 boundary points and 8 source points are used. The pseudo boundary is considered as a circle with radius of R = 5 m and each semi-circular heat source is modelled with 2 quadratic segments. The MFS results are compared with those of the FEM with a fine mesh as a reference solution. The BEM results are obtained according to the procedure described in the authors’ previous work [11]. In the BEM modelling of the problem, 8 linear boundary elements are used and each semi-circular heat source is modelled with 2 quadratic segments. In other words, the number of unknowns in the MFS and BEM modelling of the problem is the same. In the finite element modelling, each curved line heat source should be modelled as an area heat source. The thickness of the area heat source should be selected small enough to reduce the error of the modelling. Also this thickness should be chosen so that there is a possibility of generating an appropriate mesh. In this example, each semi-circular segment is modelled as an area heat source over a semi-circular strip with the thickness 0.01 m. The whole domain is meshed by 3463 quadratic quadrilateral elements, which have 10502 nodes. The finite element mesh is shown in Figure 5.

images

Figure 05: The finite element discretization of the circular domain with curved line heat sources: (a) 3463 elements, (b) 10502 nodes.

images

Figure 06: Variations of field variables along the vertical diameter: (a) temperature, (b) displacement in y direction, (c) normal stress in y direction.

The obtained results for field variables along the vertical diameter of the circle are shown in Figure 6. As can be seen, with only 8 source points and without using any internal points, the MFS results are better than the BEM results with the same number of unknowns (8 nodes), and are in a close agreement with the FEM results with a large number of elements. Also, as presented in Section 3.1, it can be seen from the figures that the stress, displacement and temperature have finite values at points exactly on the curved line heat source.

In the second 2D example, a 1 × 1 m square domain with a large number of circular heat sources is considered. The domain involves 30 circular heat sources in different locations. The radius of each heat source is 0.05 m and the intensity of each heat source is 2500 W/m. The boundary conditions and the layout of the heat sources are shown in Figure 7.

images

Figure 7: A square with 30 circular heat sources.

In the proposed MFS, each circular heat source is modelled with 4 quadratic segments and 36 collocation points and 36 source points are considered on the main and pseudo boundaries, respectively. The configuration of these points is shown in Figure 8.

The results of the present MFS are compared with those of the FEM. In the FEM, 8524 quadratic elements with 25631 nodes are used and each circular heat source is considered as an area source over a ring with thickness of 0.005 m. The finite element mesh is shown in Figure 9.

images

Figure 8: The MFS modelling of the square with 30 circular heat sources.

images

Figure 9: The finite element discretization of the square domain with 30 circular heat sources: (a) 8524 elements, (b) 25631 nodes.

The obtained results for temperature, displacement and stress variation along the line AB of the square are shown in Figure 10. As can be seen, with only 36 source points and without using any internal points, the MFS results are in very good agreement with the FEM results with a large number of elements and nodes.

images

Figure 10: Variations of field variables along the line AB of the square: (a) temperature, (b) vertical displacement and (c) x-direction normal stress.

4.2 Three-dimensional Examples

In the following three-dimensional examples, a spherical domain with radius R = 1 m, and centred at (0, 0, 0) is considered. The surface of the sphere is kept at τ = 0C. In the first example, A circular heat source with the radius r = 0.4 m and the centre at ( -0.15, 0.2, 0.15) is considered. The plane of the circular heat source is parallel to the xz plane as shown in Figure 11. For better visualization, only a quarter of the spherical domain is shown in Figure 11. The intensity per unit length of the source is considered to be a function of x and the z coordinates with the following form:

s = 20000(1 + x2 + z) Wm-1 (36)

In the MFS modelling, the circular heat source is divided into 4 quadratic segments. Only 98 source points and 98 collocation points are used. Source points are distributed on surface of a sphere with radius R = 1.4 m.

images

Figure 11: A cut out of the spherical domain including a circular heat source.

The MFS solutions are compared with the FEM solutions. In the FE analysis, the heat source, which is concentrated over a circular ring, should be modelled as a torus with a small minor radius. The major and minor radii of this torus are considered as ro = 0.4 m and ri = 0.02 m, respectively. The torus and a quarter of the spherical domain are shown in Figure 12.

The whole domain is discretized with 40,094 3D quadratic elements and 54,040 nodes. As can be seen in Figure 13, a large number of nodes should be used, especially near the heat source region.

Variations of field variables along the y and z axes are depicted in Figures 14 and 15, respectively. As it can be observed, the MFS results are in a close agreement with the FEM solutions with a note that the modeling of the problem in the proposed MFS is much simpler than that in the FEM.

To verify the accuracy of the proposed formulation for curved surface heat sources, another 3D example is analysed. In this example, the previous sphere with a cylindrical surface heat source is considered. The base of the cylinder is centred at (0, 0, 0). The radius and the height of the cylindrical heat source are 0.2 m and 0.8 m, respectively. The heat source intensity is assumed a function of x, y and z coordinates as follows:

images

Figure 12: Finite element modelling of the circular heat source by a torus with a small thickness.

images

Figure 13: The nodes in the finite element discretization of the sphere with circular heat source.

images

Figure 14: The MFS and FEM results for (a) the temperature (b) horizontal displacement and (c) x-direction normal stress along the y-axis in the sphere with a circular heat source.

s = 20y 1 + exz kWm2 (37)

In the MFS modelling, 8 quadrilateral quadratic segments are used for discretization of the cylindrical surface heat source. Pseudo boundary, number of source points and collocation points are the same as the previous example. The results of the proposed MFS are compared with the results obtained by the FEM. In the finite element modelling, the surface heat source is considered as a cylindrical volume heat source with a small thickness of 0.04. The cylindrical volume and a quarter of the spherical domain are shown in the Figure 16.

images

Figure 15: The MFS and FEM results for (a) the temperature (b) horizontal displacement and (c) x-direction normal stress along the z-axis in the sphere with a circular heat source.

images

Figure 16: Finite element modelling of the curved surface heat source by a thin-walled cylinder.

images

Figure 17: The nodes in the finite element discretization of the sphere with cylindrical heat source.

images

Figure 18: The MFS and FEM results for (a) the temperature (b) vertical displacement and (c) y-direction normal stress along the y-axis in the sphere with a cylindrical heat source.

The total number of 3D quadratic elements is 33448 which include 46376 nodes. The nodal configuraion of the finite element mesh is shown in Figure 17.

images

Figure 19: The MFS and FEM results for (a) the temperature (b) vertical displacement and (c) y-direction normal stress along the z-axis in the sphere with a cylindrical heat source.

The variations of field variables along the y and z axes are shown in Figures 18 and 19, respectively. These figures show very good agreement between the proposed MFS results and the FEM solutions. Also, as presented in Section 3.3, it can be seen from Figure 19 that the stress, displacement and temperature have finite values at points exactly on the curved surface heat source.

5 Conclusions

A formulation based on the MFS-MPS was presented for 2D and 3D thermo-mechanical problems involving curved line/surface heat sources. The particular solutions associated with the arbitrary heat sources were presented in simple integral forms. These integrals were evaluated efficiently without considering any internal cells or internal points. Unlike the FEM, modelling of the concentrated heat sources is very simple in the proposed MFS. For reliable FEM modelling, a 2D curved line heat source should be modelled as an area heat source over a narrow region of the domain and 3D curved line/surface heat sources should be modelled as volume heat sources with small thicknesses. The accuracy of the analysis by the FEM depends on the mesh density especially near the concentrated heat source. On the other hand, in the presented MFS, a small number of source and collocation points are sufficient to achieve accurate results. By presenting four examples, the effectiveness and efficiency of the proposed formulation were demonstrated. The obtained results showed that the proposed MFS is very efficient.

References

[1] Chao, C. K., & Tan, C. J. (2000). On the general solutions for annular problems with a point heat source. Journal of Applied Mechanics, 67, 511–518.

[2] Han, J. J., & Hasebe, N. (2002). Green’s functions of point heat source in various thermoelastic boundary value problems. Journal of Thermal Stresses, 25, 153–167.

[3] Rogowski, B. (2016). Green’s function for a multifield material with a heat source. Journal of Theoretical and Applied Mechanics, 54, 743–755.

[4] Le Niliot, C. (1998). The boundary element method for the time-varying strength estimation of point heat sources: application to a two-dimensional diffusion system. Numerical Heat Transfer, part B, 33, 301–321.

[5] Le Niliot, C., Rigollet F., & Petit, D. (2000). An experimental identification of line heat sources in a diffusive system using the boundary element method. International Journal of Heat and Mass Transfer, 43, 2205–2220.

[6] Karami, G., & Hematiyan, M. R. (2000a). Accurate implementation of line and distributed sources in heat conduction problems by the boundary-element method. Numerical Heat Transfer, Part B, 38, 423–447.

[7] Karami, G., & Hematiyan, M. R. (2000b). A boundary element method of inverse non-linear heat conduction analysis with point and line heat sources. International Journal for Numerical Methods in Biomedical Engineering, 16, 191–203.

[8] Shiah, Y. C., Guao, T. L., & Tan, C. L. (2005). Two-dimensional BEM thermoelastic analysis of anisotropic media with concentrated heat sources. CMES Computer Modeling in Engineering and Sciences, 7, 321–338.

[9] Shiah, Y. C., Hwang P. W., & Yang, R. B. (2006). Heat conduction in multiply adjoined anisotropic media with embedded point heat sources. Journal of Heat Transfer, 128, 207–214.

[10] Hematiyan, M. R., Mohammadi, M., & Aliabadi, M. H. (2011). Boundary element analysis of two and three-dimensional thermo-elastic problems with various concentrated heat sources. Journal of Strain Analysis for Engineering Design, 46, 227–242.

[11] Mohammadi, M., Hematiyan, M. R., & Khosravifard, A. (2016). Boundary element analysis of 2D and 3D thermoelastic problems containing curved line heat sources. European Journal of Computational Mechanics, 25, 147–164.

[12] Burgess, G., & Mahajerin, E. (1984). A comparison of the boundary element and superposition methods. Computers & Structures, 19, 697–705.

[13] Mohammadi, M., Hematiyan, M. R., & Shiah, Y. C. (2018). An efficient analysis of steady-state heat conduction involving curved line/surface heat sources in two/three dimensional isotropic media. Journal of Theoretical and Applied Mechanics, 56, 1123–1137.

[14] Mathon, R., & Johnston, R. L. (1977). The approximate solution of elliptic boundary value problems by fundamental solutions. SIAM Journal of Numerical Analysis, 14, 638–650.

[15] Fairweather, G., & Karageorghis, A. (1998). The method of fundamental solutions for elliptic boundary value problems. Advances in Computational Mathematics, 9, 69–95.

[16] Fairweather, G., Karageorghis, A., & Martin, P. A. (2003). The method of fundamental solutions for scattering and radiation problems. Engineering Analysis with Boundary Elements, 27, 759–769.

[17] Karageorghis, A., Lesnic, D., & Marin, L. (2011). A survey of applications of the MFS to inverse problems. Inverse Problems in Science and Engineering, 19, 309–36.

[18] Liu, Q. G., & Šarler, B. (2013). Non-singular method of fundamental solutions for two-dimensional isotropic elasticity problems. Computer Modeling in Engineering and Sciences, 91, 235–266.

[19] Partridge, P. W., Brebbia, C. A., & Wrobel, L. C. (1992). The dual reciprocity boundary element method. Southampton, Computational Mechanics Publications.

[20] Golberg, M. A. (1995). The method of fundamental solutions for Poisson’s equation. Engineering Analysis with Boundary Elements, 16, 205–213.

[21] Medeiros, G. C., Partridge, P. W., & Brand?o, J. O. (2004). The method of fundamental solutions with dual reciprocity for some problems in elasticity. Engineering Analysis with Boundary Elements, 28, 453–461.

[22] Wang, H., & Qin, Q. H. (2008). Meshless approach for thermo-mechanical analysis of functionally graded materials. Engineering Analysis with Boundary Elements, 32, 704–712.

[23] Tsai, C. C. (2009). The method of fundamental solutions with dual reciprocity for three dimensional thermoelasticity under arbitrary forces. Engineering Computations, 26, 229–244.

[24] Poullikkas, A., Karageorghis, A., & Georgiou, G. (1998). The method of fundamental solutions for inhomogeneous elliptic problems. Computational Mechanics, 22, 100–107.

[25] Karageorghis, A., & Smyrlis, Y. S. (2007). Matrix decomposition MFS algorithms for elasticity and thermo-elasticity problems in axisymmetric domains. Journal of Computational and Applied Mathematics, 206, 774–795.

[26] Marin, L., & Karageorghis, A. (2013). The MFS–MPS for two-dimensional steady-state thermoelasticity problems. Engineering Analysis with Boundary Elements, 37, 1004–1020.

[27] Marin, L., Karageorghis, A., & Lesnic, D. (2016). Regularized MFS solution of inverse boundary value problems in three-dimensional steady-state linear thermoelasticity. International Journal of Solids and Structures, 91, 127–142.

[28] Liu, Q. G., & Šarler, B. (2017). A non-singular method of fundamental solutions for two-dimensional steady-state isotropic thermoelasticity problems. Engineering Analysis with Boundary Elements, 75, 89–102.

[29] Becker, A. A. (1992). The Boundary Element Method in Engineering: A Complete Course. McGraw-Hill Book Company.

[30] Stroud, A. H., & Secrest, D. (1966). Gaussian quadrature formulas. New York, Prentice-Hall.

[31] Telles, J. C. F. (1987). A self-adaptive coordinate transformation for efficient numerical evaluation of general boundary element integrals. International Journal for Numerical Methods in Engineering, 24, 959–973.

[32] Aliabadi, M. H. (2002). The boundary element method, applications in solids and structures (Vol.2). Chichester: John Wiley & Sons.

[33] Hematiyan, M. R., Haghighi, A., & Khosravifard, A. (2018). A two-constrained method for appropriate determination of the configuration of source and collocation points in the method of fundamental solutions for 2D Laplace equation. Advances in Applied Mathematics and Mechanics, 10, 554–580.

Biographies

image

M. Mohammadi received his Ph.D. in Mechanical Engineering from Shiraz University, Shiraz, Iran in 2010. His thesis was about the boundary element analysis of thermo-mechanical problems. He is currently an Assistant Professor at the Department of Mechanical Engineering of Islamic Azad University, Shiraz branch, where he has been a faculty member since 2010. His research interests include computational mechanics, boundary element method, method of fundamental solutions, elasticity and thermo-elasticity.

image

M. R. Hematiyan received his Ph.D. in Mechanical Engineering from Shiraz University, Iran in 2000. He is currently a full professor at the Department of Mechanical Engineering of Shiraz University, where he has been a faculty member since 2001. He has authored more than 75 journal papers and has been supervisor of more than 70 Ph.D. and M.Sc. thesis. His research interests include computational mechanics, inverse problems, finite element and mesh-free methods, boundary element method, elasticity, thermo-elasticity, hyper-elasticity, and visco-elasticity.

image

A. Khosravifard received his Ph.D. in Mechanical Engineering from Shiraz University, Shiraz, Iran, in 2012. He joined the Department of Solid Mechanics at School of Mechanical Engineering, Shiraz University as an Assistant Professor, in 2013. He has authored around 50 scientific publications. His areas of expertise include computational mechanics, meshfree methods, fracture mechanics, inverse methods, solidification and moving boundary problems, static and dynamic analysis of structures, and functionally graded materials.

Abstract

1 Introduction

2 Basic Equations and Formulation

images

3 Formulations for Heat Sources Concentrated on a Curved Line/Surface

images

3.1 Curved Line Heat Source in Two-dimensional Problems

3.2 Curved Line Heat Source in 3D Problems

3.3 Curved Surface Heat Source in 3D Problems

images

4 Numerical Examples

4.1 Two-dimensional Examples

images

images

images

images

images

images

images

4.2 Three-dimensional Examples

images

images

images

images

images

images

images

images

images

5 Conclusions

References

Biographies