# Solution of a Problem Linear Plane Elasticity with Mixed Boundary Conditions by the Method of Boundary Integrals

Authors: Nahed S. Hussein
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

A numerical boundary integral scheme is proposed for the solution to the system ofeld equations of plane. The stresses are prescribed on one-half of the circle, while the displacements are given. The considered problem with mixed boundary conditions in the circle is replaced by two problems with homogeneous boundary conditions, one of each type, having a common solution. The equations are reduced to a system of boundary integral equations, which is then discretized in the usual way, and the problem at this stage is reduced to the solution to a rectangular linear system of algebraic equations. The unknowns in this system of equations are the boundary values of four harmonic functions which define the full elastic solution and the unknown boundary values of stresses or displacements on proper parts of the boundary. On the basis of the obtained results, it is inferred that a stress component has a singularity at each of the two separation points, thought to be of logarithmic type. The results are discussed and boundary plots are given. We have also calculated the unknown functions in the bulk directly from the given boundary conditions using the boundary collocation method. The obtained results in the bulk are discussed and three-dimensional plots are given. A tentative form for the singular solution is proposed and the corresponding singular stresses and displacements are plotted in the bulk. The form of the singular tangential stress is seen to be compatible with the boundary values obtained earlier. The efficiency of the used numerical schemes is discussed.


---

## Body

## 1. Introduction

The plane problem of the linear theory of elasticity has received considerable attention long ago as being a simplified alternative to the more realistic three-dimensional problems of practical interest. A large class of two-dimensional problems has been tackled using various analytical techniques. Due to the increasing mathematical difficulties encountered in the theoretical studies of problems involving arbitrary boundary shapes or complicated boundary conditions, many purely numerical techniques have been developed in the past few decades, which rely on finite difference or finite element techniques. In both methods, the natural boundary of the body is usually replaced by an outer polygonal shape which involves a multitude of corner points and necessarily adds or deletes parts to the region occupied by the body. This, in turn, necessitates the application of boundary conditions on artificial boundaries, a fact that introduces additional inaccuracies into the solution. Minimizing the error requires large computing times.Many of the disadvantages of the numerical techniques are overcome by the use of alternative, semianalytical treatments based on boundary integral formulations of the problem. Such approaches are usually classified under the general title of boundary integral methods. They have the advantage of reducing the volume of calculations by considering, at one stage, only the boundary values of the unknown functions and then using them to find the complete solution in the bulk. In addition, the procedure deals exclusively with the real boundary of the medium (restricted, though, to certain regularity conditions) and need not introduce artificial boundaries. An extensive account of integral equation methods in potential theory and in elastostatics may be found in [1, 2]. Natroshvili et al. [3] give a brief review of boundary integral methods as applied to the theory of micropolar elasticity. Constanda [4] investigates the use of integral equations of the first kind in plane elasticity. Atluri and Zhu [5] present a meshless local Petrov-Galerkin approach for solving problems of elastostatics. Sladek et al. [6] and Rui et al. [7] present meshless boundary integral methods for 2D elastodynamic problems. Elliotis et al. [8] present a boundary integral method for solving problems involving the biharmonic equation with crack singularities. Li et al. [9] present a numerical solution for models of linear elastostatics involving crack singularities.The solution to plane problems of elasticity for isotropic media with mixed boundary conditions is a difficult task. Boundary methods may be useful in providing such solutions, especially when the geometry of the domain boundary is not simple. Several papers deal with such problems, either for Laplace’s equation or for the biharmonic equation. Shmegera [10] finds exact solutions of nonstationary contact problems of elastodynamics for a half-plane with friction condition in the contact zone in a closed form. A new method of solution based on the use of Radon transform is used. Schiavone [11] presents integral solutions to mixed problems in plane strain elasticity with microstructure. Haller-Dintelmann et al. [12] consider three-dimensional elliptic model problems for heterogeneous media, including mixed boundary conditions. Helsing [13] studies Laplace’s equation under mixed boundary conditions and their solution by an integral equation method. Problems of elasticity are also considered. Lee et al. [14, 15] study singular solutions at corners and cracks in linear elastostatics under mixed boundary conditions. Explicit solutions are obtained. Khuri [16] outlines a general method for finding well-posed boundary-value problems for linear equations of mixed elliptic and hyperbolic type, which extends previous techniques. This method is then used to study a particular class of fully nonlinear mixed type equations.Abou-Dina and Ghaleb [17, 18] proposed a method to deal with the static, plane problems of elasticity in stresses for homogeneous isotropic media occupying simply connected regions. The method relies on the representation of the biharmonic stress function in terms of two harmonic functions and on the well-known integral representation of harmonic functions expressed in real variables. These authors applied their method to a number of examples with boundary conditions of the first or of the second type only but did not consider mixed conditions. Constanda [19] discusses Kupradze’s method of approximate solution in linear elasticity. The same author [20] explains the advantages and convenience of the use of real variables due to its generality in dealing with the different forms of the boundary, unlike the approach based on the use of complex variables “where the essential ingredients of the solution must be constructed in full for every individual situation.” Abd-Alla et al. [21] investigated the effect of the nonhomogeneity on the composite infinite cylinder of isotropic material. Abd-Alla and Farhan [22] investigated the effect of the nonhomogeneity on the composite infinite cylinder of orthotropic materials. Abd-Alla et al. [23] investigated the effect of magnetic field and nonhomogeneity on the radial vibrations in hollow rotating elastic cylinder. Abd-Alla et al. [24] studied the propagation of Rayleigh waves in a rotating orthotropic material elastic half-space under initial stress and gravity. The extensive literature on the topic is now available and we can only mention a few recent interesting investigations in [25–28].In the present paper, we propose a numerical scheme for the solution to a mixed boundary-value problem of plane, linear elasticity for homogeneous, isotropic elastic bodies occupying a circular domain. Part of the boundary is subjected to a given pressure, and the remaining part of the boundary is fixed. The initial problem with mixed boundary conditions is replaced by two subproblems with homogeneous boundary conditions, one of each type, having a common solution. Following the scheme presented in [17], the equations for each of these two subproblems are reduced to a system of boundary integral equations which are then discretized in the usual way, and the problem at this stage is reduced to the solution to a linear system of algebraic equations. The obtained results are thoroughly discussed and graphs are given. In particular, we put in evidence the singular behavior of the tangential stress component at the two separation boundary points. Three-dimensional plots for the stress function, the stress components, and the displacement components in the whole domain, obtained by the boundary collocation method, are also provided. A singular solution is proposed; the corresponding singular stresses and displacements are plotted. The ensuing form of the singular tangential stress is seen to be compatible with the boundary values obtained earlier. The efficiency of the used numerical schemes is discussed. All figures were produced using Mathematica 7.0 software.


## 2. Problem Formulation and Basic Equations

Let us considerD as the circular region occupied by the isotropic elastic medium. Its contour C has parametric representation
(1)x=acos⁡⁡(θ),y=asin⁡(θ),0≤θ<2π,
where (x,y) denote a pair of orthogonal Cartesian coordinates in the plane of  D with origin O at the center of the circle and θ is the polar angle measured from the x-axis. For application, we will take a=1. Let n be the unit outwards normal to C and τ the unit vector tangent to C at any arbitrary point, in the positive sense associated with C. One has
(2)τ=x˙ωi+y˙ωj,n=y˙ωi-x˙ωj,where the dot over a symbol denotes differentiation with respect toθ, and
(3)ω=x.2+y.2=a=1.The general equations of the linear theory of elasticity for a homogeneous and isotropic material are well established and may found in standard references. In what follows, we will quote these equations as presented in [17] without proof, to be used throughout the paper. In the absence of body forces, the equations of equilibrium are automatically satisfied if the identically nonvanishing stress components are defined through the stress function U by the relations
(4)σxx=∂2U∂y2,σyy=∂2U∂x2,σxy=∂2U∂x∂y.With respect to polar coordinates, the stress components are(5)σrr=1r∂U∂r+1r2∂2U∂θ2,σθθ=∂2U∂r2,σrθ=1r2∂U∂θ-1r∂2U∂r∂θ.It is assumed that the stress functionU∈C4(D) and that its second order partial derivatives are univalued functions in the whole region D. Hooke’s law reads
(6)σxx=vE(1+v)(1-2v)[∂u∂x+∂v∂y]+E(1+v)∂u∂x,σyy=vE(1+v)(1-2v)[∂u∂x+∂v∂y]+E(1+v)∂u∂y,σxy=E2(1+v)[∂u∂y+∂v∂x],
where E and v are Young’s modulus and Poisson’s ratio, respectively, for the considered elastic medium.The compatibility condition for the solution to (6) for the displacement components leads to the following homogeneous biharmonic equation for the stress function U:
(7)Δ2U=0.The stress functionU solving (7) is
(8)U=xΦ+yΦc+Ψ,
where Φ and Ψ are two harmonic functions, the superscript c′ denotes the harmonic conjugate, and D¯ is the closure of D. Since the boundary integral representation is to be used, it seems adequate to suppose from the outset that the functions Φ and Ψ and their conjugates belong to the class of functions C2(D). The following representation for the mechanical displacement components may be easily deduced:
(9)E1+vu=-∂U∂x+4(1-v)Φ,E1+vυ=-∂U∂y+4(1-v)Φc.In terms of the harmonic functionsΦ, Φc, and Ψ, the stress and the displacement components are expressed as follows:
(10)vvσxx=x∂2Φ∂y2+2∂Φc∂y+y∂2Φc∂y2+∂2Ψ∂y2,(11)vvvσyy=x∂2Φ∂x2+2∂Φ∂x+y∂2Φc∂x2+∂2Ψ∂x2,(12)σxy=-x∂2Φ∂x∂y-y∂2Φc∂x∂y-∂2Ψ∂x∂y,(13)vvvvvE1+vu=(3-4v)Φ-x∂Φ∂x-y∂Φc∂x-∂Ψ∂x,(14)vvvvvvE1+vυ=(3-4v)Φc-x∂Φ∂y-y∂Φc∂y-∂Ψ∂y.


## 3. Boundary Integral Representation of the Basic Equations

In what follows, we present the boundary integral representation of the basic equations and boundary conditions to be used in the sequel. We closely follow the guidelines of [17].


### 3.1. Boundary Integral Representation of Harmonic Functions

Let us consider thatf∈C2(D-) is harmonic in D. We use the well-known integral representation for f at an arbitrary field point (x,y) in D in terms of the boundary values of the function f and its complex conjugate fc in the form
(15)f(x,y)=12π∮s[f(s´)∂∂n´ln⁡⁡R+fc(s´)∂∂s´ln⁡⁡R]ds´,
where R is the distance between the point (x,y) in D and the current integration point (x(s´),y(s´)) on S. The representation of the conjugate function is given by
(16)fc(x,y)=12π∮s‍[fc(s´)∂∂n´ln⁡⁡R-f(s´)∂∂s´ln⁡⁡R]ds´.The integral representations (15) and (16) for the harmonic functions f and fc replace the usual Cauchy-Riemann conditions
(17)∂f∂x=∂fc∂y,∂f∂y=-∂fc∂x.
When the point (x,y) tends to a boundary point (x(s),y(s)), relation (15) yields
(18)f(s)=1π∮s‍[f(s´)∂∂n´ln⁡R+fc(s´)∂∂s´ln⁡R]ds´.Replacing(∂/∂n´)ln⁡⁡R by (∂/∂s´)Θ in (15), (16), and their boundary version (18), where Θ is the complex conjugate of ln⁡⁡R, it is readily seen that these integral relations are invariant under the transformation of parameter from the arc length s to any other suitable parameter. This property makes the method more flexible.


### 3.2. Conditions for the Uniqueness of the Solution

Before dealing with each of the two above-mentioned fundamental problems, we first turn to the conditions to be satisfied in order to determine the unknown harmonic functions in an unambiguous manner. This is of primordial importance for any numerical treatment of the problem, for a proper use of the solving algorithm. We will require the following supplementary conditions to be satisfied at the pointQ0(s=0) of the boundary, in order to determine the totality of the arbitrary integration constants appearing throughout the solution process. These additional conditions have no physical implications on throughout the problem.(i)
The vanishing of the functionU and its first order partial derivatives at Q0(19)U=∂U∂x=∂U∂y=0.Or, equivalently,
(20)U=∂U∂s=∂U∂n=0,
which, in terms of the boundary values of the unknown harmonic functions, give
(21)vvx(0)Φ(0)+y(0)Φc(0)+Ψ(0)=0,x(0)Φ.(0)+y(0)Φc.(0)+Ψ.(0)+x.(0)Φ(0)+y.(0)Φc(0)=0,x(0)Φc.(0)-y(0)Φ.(0)+Ψc.(0)+y.(0)Φ(0)-x.(0)Φc(0)=0,(ii)
The vanishing of the expression(22)x(0)Φc(0)-y(0)Φ(0)+Ψc(0)=0.This last additional condition amounts to determining the value ofΨc at Q0 and this is chosen to simplify the formule.


### 3.3. Boundary Condition for the First Fundamental Problem of Elasticity

In the first fundamental problem, we are given the force distribution on the boundaryS of the domain D.Let(23)f=fxi+fyj=fττ+fnn
denote the external force per unit length of the boundary. Then, at a general boundary point Q, the stress vector
(24)σn=for, in components,(25)σxxnx+σxyny=fx,σxynx+σyyny=fy.The stress functionU at the boundary point Q(26)∂U∂s(s)=-x.(s)Y(s)+y.(s)X(s),∂U∂n(s)=-y.(s)Y(s)-x.(s)X(s),
or, in terms of the unknown harmonic functions,
(27)x(s)Φ.(s)+y(s)Φc.(s)+Ψ.(s)+x.(s)Φ(s)+y.(s)Φc(s)=-x.(s)Y(s)+y.(s)X(s),x(s)Φc.(s)-y(s)Φ.(s)+Ψ.c(s)+y.(s)Φ(s)-x.(s)Φc(s)=-y.(s)Y(s)-x.(s)X(s).


### 3.4. Boundary Condition for the Second Fundamental Problem of Elasticity

In this problem, we are given the displacement vector on the boundaryS of the domain D. Let this vector be denoted
(28)d=dxi+dyj=dττ+dnn.Multiplying the restriction of expression (13) to the boundary S by x˙(s) and that of expression (14) by y˙(s) and adding, one gets
(29)(3-4v)(x.(s)ϕ(s)+y.(s)ϕc(s))-x(s)ϕ.(s)-y(s)ϕ˙c(s)-ψ.(s)=E1+v(x.(s)dx(s)+y.(s)dy(s))ω.Similarly, if one multiplies the restriction of expression (13) to the boundary S by y˙(s) and that of expression (14) by x˙(s) and subtracts, one obtains
(30)(3-4v)(y.(s)ϕ(s)-x.(s)ϕc(s))-x(s)ϕ.c(s)+y(s)ϕ.(s)-ψ.c(s)=E1+v(y.(s)dx(s)-x.(s)dy(s))ω.These last two relations may be conveniently rewritten as(31)(3-4v)(x.(s)ϕ(s)+y.(s)ϕc(s))-x(s)ϕ.(s)-y(s)ϕ.c(s)-ψ.(s)=E1+vdτ(s)ω,(3-4v)(y.(s)ϕ(s)-x.(s)ϕc(s))-x(s)ϕ.c(s)+y(s)ϕ.(s)-ψ.c(s)=E1+vdnω.


### 3.5. Boundary Conditions for the Third Fundamental Problem of Elasticity

This is a problem with mixed boundary conditions. For definiteness, we will restrict further considerations to the case where one-half of the boundary has a prescribed pressure on it, while the other half of the boundary is fixed. This problem will be replaced by two subproblems, each with homogeneous boundary condition. The first subproblem is of the first kind. It involves the given known pressure on the same half of the boundary as the initial problem and an unknown stress on the other half. This stress is expressed through its normal and tangential components, respectively, denotedP~n, P~τ. The second subproblem is of the second type. It involves zero displacement on the same half of the boundary as the initial problem and an unknown displacement on the other half. This displacement is expressed through its normal and tangential components, respectively, denoted u~n, u~τ.In what follows, we will apply this idea to solve three problems for the ellipse, the nearly circular boundary, and the Cassini oval. In choosing these boundaries, we have tried to keep away from boundaries involving singular points, as our main task is to deal with the mixed boundary conditions, which, on its own, includes two boundary points of separation which need special attention.


### 3.6. Calculation of the Harmonic Functions at Internal Points

Having determined the boundary values of the harmonic functions, formulae (15) and (16) may now be used to calculate the values of these functions at any point (x,y) inside the domain. For this we write
(32)R=(x-x(s′))2+(y-y(s′))2,∂ln⁡⁡(R)∂n=n·∇(ln⁡⁡R),∂ln⁡⁡(R)∂s=τ·∇(ln⁡⁡R).
We can also proceed otherwise. In fact, if we write down expansions of the four harmonic functions involved in the calculations in terms of some adequately chosen basis, we can then determine the expansion coefficients using the well-known boundary collocation method (BCM). This is in fact the method we have used to calculate the unknown functions in the circular domain. The expansions of the four basic harmonic functions in terms of polar harmonics are as follows:
(33)Φ=∑n=1Nrn(Ancos⁡nθ+Bnsinnθ),Φc=∑n=1Nrn(Ansinnθ-Bncos⁡nθ),Ψ=∑n=1Nrn(Encos⁡nθ+Dnsinnθ),Ψc=∑n=1Nrn(Ensinnθ-Dncos⁡nθ),while the stress function is(34)U=xΦ+yΦc+Ψ
and the quantities of practical interest are
(35)2μu=∑n=1Nrn(Ancos⁡nθ+Bnsinnθ)-∑n=1Nnrn(Ancos⁡⁡(n-2)θ+Bnsin⁡(n-2)θ)-∑n=1Nnrn-1(Encos⁡(n-1)θ+Dnsin(n-1)θ),2μυ=∑‍n=1Nrn(Ansin⁡nθ-Bncos⁡⁡nθ)-∑‍Nn=1nrn(-Ansin⁡(n-2)θ+Bncos⁡⁡(n-2)θ)-∑‍n=1Nnrn-1(-Ensin⁡(n-1)θ+Dncos⁡⁡(n-1)θ).The equations for the normal and the tangential stresses on any given element of area with unit normal(nx,ny) inside the body or on its boundary are given by the following formulae:
(36)σnn=∑‍n=1Nrn-1An(cos⁡⁡(n-1)θ((3n-n2)nr2+(n2+n)nθ2)vvvvvvvv∑+2(n2-n)nrnθsin⁡(n-1)θ)+∑‍n=1Nrn-1Bn(sin⁡(n-1)θ((3n-n2)nr2+(n2+n)nθ2)vvvvvvvvvv∑+2(n2-n)nrnθcos⁡⁡(n-1)θ)+∑‍n=1Nrn-2Cn(cos⁡⁡nθ((n-n2)nr2+(n2-n)n-θ2)vvvvvvvvvv∑+2(n2-n)nrnθsin⁡nθ)+∑‍n=1Nrn-2Dn(sin⁡nθ((n-n2)nr2+(n2-n)nθ2)vvvvvvvv∑+2(n-n2)nrnθcos⁡⁡nθ),σnτ=∑‍n=1Nrn-1An(cos⁡⁡(n-1)θ((-2n+2n2)nrnθ)vvvvvvvvvvvv+(n2-n)(nr2-nθ2)sin⁡(n-1)θ)+∑‍n=1Nrn-1Bn(sin⁡(n-1)θ((-2n+2n2)nrnθ)vvvvvvvvvvvv+(n-n2)(nr2-nθ2)cos⁡⁡(n-1)θ)+∑‍n=1Nrn-2Cn(cos⁡⁡nθ(-2n+2n2)nrnθ+(n2-n)nθ2vvvvvvvvvvvv+2(n2-n)(nr2-nθ2)sinnθ)+∑‍n=1Nrn-2Dn(sinnθ(-2n-2n2)nrnθvvvvvvvvvv+(n-n2)(nr2-nθ2)cos⁡nθ).The relevant boundary relations are discretized in the usual way by considering a partition of the boundary. As a result, the actual boundary is replaced by a contour formed by broken lines. The differential and integral equations thus reduce to a rectangular system of linear algebraic equations which are solved by the least squares method. The convergence of the solution to the discretized system of equations to the solution to the initial problem was discussed elsewhere [20]. Here, we only notice the existence of removable singularities in the formulae of integral representation of  harmonic functions. These are dealt with in the manner explained in [18]. Also, the tangential derivatives of the unknown harmonic functions have to be evaluated carefully. We have calculated these derivatives using 15 points.


## 3.1. Boundary Integral Representation of Harmonic Functions

Let us consider thatf∈C2(D-) is harmonic in D. We use the well-known integral representation for f at an arbitrary field point (x,y) in D in terms of the boundary values of the function f and its complex conjugate fc in the form
(15)f(x,y)=12π∮s[f(s´)∂∂n´ln⁡⁡R+fc(s´)∂∂s´ln⁡⁡R]ds´,
where R is the distance between the point (x,y) in D and the current integration point (x(s´),y(s´)) on S. The representation of the conjugate function is given by
(16)fc(x,y)=12π∮s‍[fc(s´)∂∂n´ln⁡⁡R-f(s´)∂∂s´ln⁡⁡R]ds´.The integral representations (15) and (16) for the harmonic functions f and fc replace the usual Cauchy-Riemann conditions
(17)∂f∂x=∂fc∂y,∂f∂y=-∂fc∂x.
When the point (x,y) tends to a boundary point (x(s),y(s)), relation (15) yields
(18)f(s)=1π∮s‍[f(s´)∂∂n´ln⁡R+fc(s´)∂∂s´ln⁡R]ds´.Replacing(∂/∂n´)ln⁡⁡R by (∂/∂s´)Θ in (15), (16), and their boundary version (18), where Θ is the complex conjugate of ln⁡⁡R, it is readily seen that these integral relations are invariant under the transformation of parameter from the arc length s to any other suitable parameter. This property makes the method more flexible.


## 3.2. Conditions for the Uniqueness of the Solution

Before dealing with each of the two above-mentioned fundamental problems, we first turn to the conditions to be satisfied in order to determine the unknown harmonic functions in an unambiguous manner. This is of primordial importance for any numerical treatment of the problem, for a proper use of the solving algorithm. We will require the following supplementary conditions to be satisfied at the pointQ0(s=0) of the boundary, in order to determine the totality of the arbitrary integration constants appearing throughout the solution process. These additional conditions have no physical implications on throughout the problem.(i)
The vanishing of the functionU and its first order partial derivatives at Q0(19)U=∂U∂x=∂U∂y=0.Or, equivalently,
(20)U=∂U∂s=∂U∂n=0,
which, in terms of the boundary values of the unknown harmonic functions, give
(21)vvx(0)Φ(0)+y(0)Φc(0)+Ψ(0)=0,x(0)Φ.(0)+y(0)Φc.(0)+Ψ.(0)+x.(0)Φ(0)+y.(0)Φc(0)=0,x(0)Φc.(0)-y(0)Φ.(0)+Ψc.(0)+y.(0)Φ(0)-x.(0)Φc(0)=0,(ii)
The vanishing of the expression(22)x(0)Φc(0)-y(0)Φ(0)+Ψc(0)=0.This last additional condition amounts to determining the value ofΨc at Q0 and this is chosen to simplify the formule.


## 3.3. Boundary Condition for the First Fundamental Problem of Elasticity

In the first fundamental problem, we are given the force distribution on the boundaryS of the domain D.Let(23)f=fxi+fyj=fττ+fnn
denote the external force per unit length of the boundary. Then, at a general boundary point Q, the stress vector
(24)σn=for, in components,(25)σxxnx+σxyny=fx,σxynx+σyyny=fy.The stress functionU at the boundary point Q(26)∂U∂s(s)=-x.(s)Y(s)+y.(s)X(s),∂U∂n(s)=-y.(s)Y(s)-x.(s)X(s),
or, in terms of the unknown harmonic functions,
(27)x(s)Φ.(s)+y(s)Φc.(s)+Ψ.(s)+x.(s)Φ(s)+y.(s)Φc(s)=-x.(s)Y(s)+y.(s)X(s),x(s)Φc.(s)-y(s)Φ.(s)+Ψ.c(s)+y.(s)Φ(s)-x.(s)Φc(s)=-y.(s)Y(s)-x.(s)X(s).


## 3.4. Boundary Condition for the Second Fundamental Problem of Elasticity

In this problem, we are given the displacement vector on the boundaryS of the domain D. Let this vector be denoted
(28)d=dxi+dyj=dττ+dnn.Multiplying the restriction of expression (13) to the boundary S by x˙(s) and that of expression (14) by y˙(s) and adding, one gets
(29)(3-4v)(x.(s)ϕ(s)+y.(s)ϕc(s))-x(s)ϕ.(s)-y(s)ϕ˙c(s)-ψ.(s)=E1+v(x.(s)dx(s)+y.(s)dy(s))ω.Similarly, if one multiplies the restriction of expression (13) to the boundary S by y˙(s) and that of expression (14) by x˙(s) and subtracts, one obtains
(30)(3-4v)(y.(s)ϕ(s)-x.(s)ϕc(s))-x(s)ϕ.c(s)+y(s)ϕ.(s)-ψ.c(s)=E1+v(y.(s)dx(s)-x.(s)dy(s))ω.These last two relations may be conveniently rewritten as(31)(3-4v)(x.(s)ϕ(s)+y.(s)ϕc(s))-x(s)ϕ.(s)-y(s)ϕ.c(s)-ψ.(s)=E1+vdτ(s)ω,(3-4v)(y.(s)ϕ(s)-x.(s)ϕc(s))-x(s)ϕ.c(s)+y(s)ϕ.(s)-ψ.c(s)=E1+vdnω.


## 3.5. Boundary Conditions for the Third Fundamental Problem of Elasticity

This is a problem with mixed boundary conditions. For definiteness, we will restrict further considerations to the case where one-half of the boundary has a prescribed pressure on it, while the other half of the boundary is fixed. This problem will be replaced by two subproblems, each with homogeneous boundary condition. The first subproblem is of the first kind. It involves the given known pressure on the same half of the boundary as the initial problem and an unknown stress on the other half. This stress is expressed through its normal and tangential components, respectively, denotedP~n, P~τ. The second subproblem is of the second type. It involves zero displacement on the same half of the boundary as the initial problem and an unknown displacement on the other half. This displacement is expressed through its normal and tangential components, respectively, denoted u~n, u~τ.In what follows, we will apply this idea to solve three problems for the ellipse, the nearly circular boundary, and the Cassini oval. In choosing these boundaries, we have tried to keep away from boundaries involving singular points, as our main task is to deal with the mixed boundary conditions, which, on its own, includes two boundary points of separation which need special attention.


## 3.6. Calculation of the Harmonic Functions at Internal Points

Having determined the boundary values of the harmonic functions, formulae (15) and (16) may now be used to calculate the values of these functions at any point (x,y) inside the domain. For this we write
(32)R=(x-x(s′))2+(y-y(s′))2,∂ln⁡⁡(R)∂n=n·∇(ln⁡⁡R),∂ln⁡⁡(R)∂s=τ·∇(ln⁡⁡R).
We can also proceed otherwise. In fact, if we write down expansions of the four harmonic functions involved in the calculations in terms of some adequately chosen basis, we can then determine the expansion coefficients using the well-known boundary collocation method (BCM). This is in fact the method we have used to calculate the unknown functions in the circular domain. The expansions of the four basic harmonic functions in terms of polar harmonics are as follows:
(33)Φ=∑n=1Nrn(Ancos⁡nθ+Bnsinnθ),Φc=∑n=1Nrn(Ansinnθ-Bncos⁡nθ),Ψ=∑n=1Nrn(Encos⁡nθ+Dnsinnθ),Ψc=∑n=1Nrn(Ensinnθ-Dncos⁡nθ),while the stress function is(34)U=xΦ+yΦc+Ψ
and the quantities of practical interest are
(35)2μu=∑n=1Nrn(Ancos⁡nθ+Bnsinnθ)-∑n=1Nnrn(Ancos⁡⁡(n-2)θ+Bnsin⁡(n-2)θ)-∑n=1Nnrn-1(Encos⁡(n-1)θ+Dnsin(n-1)θ),2μυ=∑‍n=1Nrn(Ansin⁡nθ-Bncos⁡⁡nθ)-∑‍Nn=1nrn(-Ansin⁡(n-2)θ+Bncos⁡⁡(n-2)θ)-∑‍n=1Nnrn-1(-Ensin⁡(n-1)θ+Dncos⁡⁡(n-1)θ).The equations for the normal and the tangential stresses on any given element of area with unit normal(nx,ny) inside the body or on its boundary are given by the following formulae:
(36)σnn=∑‍n=1Nrn-1An(cos⁡⁡(n-1)θ((3n-n2)nr2+(n2+n)nθ2)vvvvvvvv∑+2(n2-n)nrnθsin⁡(n-1)θ)+∑‍n=1Nrn-1Bn(sin⁡(n-1)θ((3n-n2)nr2+(n2+n)nθ2)vvvvvvvvvv∑+2(n2-n)nrnθcos⁡⁡(n-1)θ)+∑‍n=1Nrn-2Cn(cos⁡⁡nθ((n-n2)nr2+(n2-n)n-θ2)vvvvvvvvvv∑+2(n2-n)nrnθsin⁡nθ)+∑‍n=1Nrn-2Dn(sin⁡nθ((n-n2)nr2+(n2-n)nθ2)vvvvvvvv∑+2(n-n2)nrnθcos⁡⁡nθ),σnτ=∑‍n=1Nrn-1An(cos⁡⁡(n-1)θ((-2n+2n2)nrnθ)vvvvvvvvvvvv+(n2-n)(nr2-nθ2)sin⁡(n-1)θ)+∑‍n=1Nrn-1Bn(sin⁡(n-1)θ((-2n+2n2)nrnθ)vvvvvvvvvvvv+(n-n2)(nr2-nθ2)cos⁡⁡(n-1)θ)+∑‍n=1Nrn-2Cn(cos⁡⁡nθ(-2n+2n2)nrnθ+(n2-n)nθ2vvvvvvvvvvvv+2(n2-n)(nr2-nθ2)sinnθ)+∑‍n=1Nrn-2Dn(sinnθ(-2n-2n2)nrnθvvvvvvvvvv+(n-n2)(nr2-nθ2)cos⁡nθ).The relevant boundary relations are discretized in the usual way by considering a partition of the boundary. As a result, the actual boundary is replaced by a contour formed by broken lines. The differential and integral equations thus reduce to a rectangular system of linear algebraic equations which are solved by the least squares method. The convergence of the solution to the discretized system of equations to the solution to the initial problem was discussed elsewhere [20]. Here, we only notice the existence of removable singularities in the formulae of integral representation of  harmonic functions. These are dealt with in the manner explained in [18]. Also, the tangential derivatives of the unknown harmonic functions have to be evaluated carefully. We have calculated these derivatives using 15 points.


## 4. Numerical Results and Discussion

The force acting on one-half of the boundary is a pressure of intensityf given by
(37)f=-p0(sin⁡θ)6π<θ≤2π.The other half of the boundary is completely fixed:(38)u=υ=0,0<θ≤π.For definiteness, we have takenp0=1. The motivation for the above choice of the pressure on half the boundary is to make the pressure distribution tend to zero smoothly enough at both ends of its interval of definition.The above boundary integral equations are solved numerically, from which we have obtained the boundary values of the harmonic functionsϕ, ϕc, ψ, ψc, σ~n, σ~τ, u~n, and u~τ and, accordingly, of the stress function U. The boundary was discretized by placing a number of nodal points on it as explained; 240 boundary nodes were needed in order to get the present results. The nodal points were distributed uniformly on the boundary. The results are shown below.Figure1 gives the boundary values of the four basic harmonic functions.Boundary values of the basic harmonic functions.
(a)
(b)
(c)
(d)All the four figures show a weak discontinuity at the pointθ=π. It goes without saying that the same takes place for θ=0 from symmetry considerations. The emplacements of these discontinuities are referred to by arrows on the figures.Figure2 shows the normal and the tangential components of the unknown displacement and stress on the relevant parts of the boundary. As may be noticed three of figures, some difficulties were encountered when performing the calculations in the vicinity of the singular points at θ=0 and at θ=π. It was not possible to increase the number of nodal points beyond 240 for stability reasons in order to improve the results.Boundary values of the normal and tangential stress and displacement components.
(a)
(b)
(c)
(d)The normal displacement component reached zero value at both separation points as should be, while the tangential displacement component failed to do so, but the results improved as the number of nodes was augmented, up to a certain limit. Unwanted oscillations appeared on the curves for stress near these two points. Curve fitting techniques by polynomial functions were used to improve the curves. Additionally, two logarithmic functions based on the singularities were used for fitting, only the functionσ~τ. These are the smooth curves on the figures, based on these observations. One suspects the presence of a logarithmic behavior of the function σ~τ at the singular points.Figure3 shows the boundary values of the stress function. The curve is skewed towards the second half of the boundary, but thus asymmetry should not raise any concerns as it depends on the additional conditions imposed on this function.Figure 3
Boundary values of the stress function.We have used the well-known boundary collocation method to directly compute the unknown functionsσnn, σnτ, u, and υ on concentric circles centered at the origin inside the domain using (33)–(36), together with the given boundary conditions. Thus, there is no ambiguity in the meaning of the normal unit vector appearing in the equations. A maximum number of 120 nodes were used. The results are shown on the following three-dimensional plots, where we have also shown on each of them the circular region in which the unknown functions are plotted (Figure 4).Displacement and stress components inside the circle by BCM.
(a)
u
(b)
v
(c)
σ
n
n
(d)
σ
n
τThe surfaces for the displacement components are in conformity with our expectations. The normal stress component is regular, while the tangential stress component shows the singular behavior mentioned above (Figure5).Stress components inside the circle by BCM.
(a)
σ
x
x
(b)
σ
y
y
(c)
σ
x
yHere again, we notice a singular behavior of the stress componentsσyy and σxy at the two singular boundary points. In the other cases, the comparison is made with the correction of the analytical results obtained in [29] as a special case. The numerically obtained results are compared with those obtained analytically in [30].


## 5. On the Singular Solution

We propose to add a functionψ(s) with boundary singularities to the function ψ, in order to get the required logarithmic behavior of the function σnτ at the singular points (±a,0). This function was proposed by Abou-Dina and Ghaleb [31] in connection with the solution to some boundary-value problems for Laplace’s equation in rectangular domains; here it is used in a special setting.Figure6 shows the emplacements of the singularities of function ψ(s), as well as the variables θ1, θ2, ρ1, and ρ2.Figure 6
Singular points.As to the functionψs, is of the form
(39)ψs=12π(∑ρ12(sin⁡2θ1ln⁡ρ1+θ1cos⁡⁡2θ1))+(ρ22(sin2θ2ln⁡ρ2+θ2cos⁡2θ2)),where(40)ρ1=r2-2arcos⁡⁡θ+a2,ρ2=r2+2arcos⁡⁡θ+a2,θ1=tan-1a-rcos⁡⁡θ-rsin⁡θ,θ2=tan-1rcos⁡⁡θ+a-rsin⁡θ.Figure7 shows the singular stress and displacement components obtained from the singular function ψs. These functions are also labeled “s.” The singular behavior of function σnτ is clear. It is recalled that the normal and the tangential components of stress are calculated on concentric circles centered at the origin inside the domain; thus the meaning of the unit normal vector is clear.Singular displacement and stress components inside the circle.
(a)
u
(
s
)
(b)
v
(
s
)
(c)
σ
n
n
(
s
)
(d)
σ
n
τ
(
s
)The details of the calculations and the final results for the stresses and displacements when the singular function is added will be considered in a separate publication for other types of boundaries.


## 6. Conclusions

The following conclusions are due.(1)
We have considered a boundary-value problem of the plane theory of elasticity with mixed boundary conditions in a circle. Half of the boundary is subjected to pressure, while the other half is completely fixed. The shape of the boundary was chosen as the simplest, in order to focus on problems raised by the boundary singularities. We have also chosen smooth boundary pressure that decreases smoothly enough to zero towards the points of separation.(2)
The correct calculation of this type of problems requires a large number of boundary points at which the unknowns are to be calculated. For the present case, 240 points could be reached without obtaining satisfactory results in the whole boundary. The reason for this is the presence of singular boundary points at the separation points of the boundary conditions.(3)
To get the solution on the boundary, we have replaced our problem by two subproblems, each with homogeneous boundary condition of one type, having a common solution.(4)
The calculations on the boundary were performed using a known boundary integral technique involving harmonic functions only, including regularization and a careful calculation of the tangential derivatives of functions using 15 points.(5)
The boundary calculations indicated a logarithmic behavior of the tangential stress component on the fixed part of the boundary.(6)
The solution inside the domain was obtained by the collocation method directly using the prescribed boundary conditions. A maximum number of 120 uniformly distributed nodes were used.(7)
In solving the arising systems of linear algebraic equations, we have used least squares and QR-factorization techniques; both yielded the same results. Each time a linear system of equations was solved, we verified that the obtained results satisfy the system. The errors did not exceed 1%.(8)
We have proposed a solution having a logarithmic boundary singularity to improve the solution. The absolute errors in satisfying the boundary conditions on an interval and including the separation points could thus be reduced from nearly 3.(9)
The obtained results show stress concentration and thus indicate the need to introduce domains of possible plastic behavior of the material around the two boundary separation points.(10)
The method extends to other geometries of the boundary. In the presence of corner points, a smoothing process must be applied. Numerical experiments have clearly indicated that the best results are obtained for boundaries of smoothness of the fourth order.(11)
Future work will involve more complicated shapes of the boundary and other types of boundary conditions.



---

*Source: 323178-2014-05-08.xml*