Application of the Generalized Method of Moving Coordinates to Calculating Stress Fields near an Elliptical Hole

The distribution of stresses near holes is of great importance in fracture mechanics and material modeling. The present paper provides a general stress solution near a traction-free surface for an arbitrary piecewise linear yield criterion, assuming plane-strain conditions. The generalized method of moving coordinates is proven efficient in this case. In particular, the solution reduces to evaluating one ordinary integral. The boundary value problem solved is a Cauchy problem for a hyperbolic system of equations. Therefore, the stress solution in the plastic region is independent of other boundary conditions, though the occurrence of plastic yielding at a specific point is path-dependent. The general solution applies to calculating the stress field near an elliptic hole. It is shown that the parameter that controls the pressure-dependency of the yield criterion affects the stress field significantly. The aspect ratio is less significant as compared to that parameter. However, for a given material, the aspect ratio should also be considered to predict the stress field accurately, especially in the near vicinity of the hole. The solution reduces to an available solution for the pressure-independent yield criterion, which is a particular yield criterion of the considered class of yield criteria.


Introduction
The distribution of stresses near holes and cracks is of great importance in fracture mechanics and material modeling. A vast amount of literature is devoted to developing methods for finding such distributions using the theories of elasticity and plasticity. The linearly elastic stress field at the base of a sharp stationary crack has been found in [1]. This analysis has been extended to the plastic range in [2][3][4].
The importance of stress solutions near holes is twofold. Firstly, such solutions are required to understand the fracture process. Paper [5] has presented a stress solution in the plastic region around an elliptic hole. This solution can be used in conjunction with the ductile fracture model proposed in [6]. Exact solutions for an elliptic hole embedded in a thermoelectric material have been presented in [7,8]. Paper [9] has found the stress distribution around an elliptic hole for non-local elastic material models. An elastic solution for a plate subject to biaxial loading and containing an inclined elliptic defect has been provided in [10]. The stress concentration generated by an elliptic hole in a functionally graded panel subjected to uniform tension has been analyzed in [11]. The effect of elastic anisotropy on the stress distribution around an elliptical crack has been revealed in [12], assuming that the medium is subjected to general loading at infinity. The influence of multiple pre-existing holes on the fracture process in granite subject to compressive loading has been studied in [13]. The stress concentration due to a circular hole in a plate subject to the uniaxial extension has been investigated in [14], using new constitutive equations with density-dependent elasticity moduli. This result has been extended to bi-axial deformation in [15]. A perturbation technique has been adopted in [16] for deriving analytical solutions for an elastic/plastic plate with a hole, including an elliptic hole. In particular, the elastic/plastic boundary has been found. Another analytical method for determining the elastic/plastic boundary in the vicinity of holes has been proposed in [17].
Another area of applications of the stress solutions near elliptic holes is that structural members contain such holes. Examples of such members have been considered in [18,19].
The present paper extends the solution [5] to an arbitrary piecewise linear yield criterion. Such yield criteria are widely used for describing metallic and nonmetallic materials [20][21][22]. The generalized method of moving coordinates is developed and used. The solution is semi-analytical. One needs to evaluate one ordinary integral to calculate the principal stresses at any point of the plastic region.

Basic Equations
Considered is plane-strain deformation of a material obeying an arbitrary piecewise linear yield criterion. Such yield criteria are represented by a linear function of the two principal stresses in a typical plane of flow. These stresses are denoted as σ 1 and σ 2 . The yield criterion can be represented as: Here, λ is a constitutive parameter and Note that Equation (1) represents the totality of piecewise linear yield criteria if the stresses are non-dimensionalized using an appropriate reference stress. For example, a standard form of the Mohr-Coulomb yield criterion is [23]: where the cohesion k and angle of internal friction χ are constants. Equation (3) reduces to Equation (1) if λ = − sin χ and the reference stress is k cos χ. It is also possible to choose without loss of generality that σ 1 > σ 2 . Then, q > 0. Let ψ be the inclination of the direction of the stress σ 1 relative to the x 1 −axis of a Cartesian coordinate system (x 1 , y 1 ), measured anticlockwise. Then, the stress components referred to the Cartesian coordinates are: Equations (1), (2) and (4) combine to give: Equation (1) and the equilibrium equations comprise a statically determinate system. This system, without the assumption that the yield criterion is piecewise linear, has been analyzed in [24]. In the case of Equation (1), this analysis provides the following characteristics: which are termed the α− and β−lines, respectively, where The relations along the characteristics are: on the α− and β−lines, respectively. It is seen from Equations (6) and (7) that the angle between the direction of the stress σ 1 and each characteristic line is constant. Therefore, Equation (8) can be rewritten as: is the angle between the x−axis and the α−characteristic line measured from the axis anticlockwise. If both families of characteristics are curved, the equations in (9) can be immediately integrated to give: where The constants ϕ 0 and p 0 have been introduced for further convenience. Note that Equation (12) reduces to P = (p − p 0 )/2 as λ → 0 . Solving the equations in (11) for P and ϕ yields: P = (β − α) sin 2φ and ϕ = (β + α) sin 2φ + ϕ 0 .

Generalized Moving Coordinates
The moving coordinates (x, y) are used for solving boundary value problems in pressure-independent plasticity. This method has been described in many monographs, for example, [24]. In short, x and y are the coordinates of the point C under consideration, referred to axes passing through the origin O and parallel to the characteristic directions at C. The geometric interpretation of the moving coordinates is shown in Figure 1a. In the case of pressure-independent plasticity, the characteristic directions are orthogonal. Therefore, the moving coordinates are orthogonal. It is seen from Equations (6) and (7)  consideration, referred to axes passing through the origin O and parallel to the characteristic directions at C. The geometric interpretation of the moving coordinates is shown in Figure 1a. In the case of pressure-independent plasticity, the characteristic directions are orthogonal. Therefore, the moving coordinates are orthogonal. It is seen from Equations (6) and (7)  The geometric interpretation of the generalized moving coordinates, which are also denoted as ( ) , , x y is shown in Figure 1b. In this figure, ( ) , x y is a skew rectilinear coordinate system with the origin at O. The x − axis is tangent to the base α − characteristic line passing through O, and the y − axis is tangent to the base β − characteristic Let e x and e y be the unit base vectors of the (x, y)−coordinate system. The unit vectors e α and e β are directed along the α− and β−characteristic lines at the point C under consideration. Therefore, the angle between the vectors e x and e α is ϕ. Then (Figure 1b): e α · e x = cos ϕ, e α · e y = cos(ϕ − 2φ), e β · e x = cos(ϕ + 2φ), e β · e y = cos ϕ, e y · e x = cos 2φ.
The position vector of any point can be represented as Using Equations (14) and (15), one can find the scalar products, R · e x and R · e y , as: Solving these equations for x and y yields: Alternatively, the equations in (16) can be solved for x and y. As a result, Differentiating the first equation in (17) with respect to β and the second with respect to α leads to: It is seen from Equation (13) that ∂ϕ/∂α = ∂ϕ/∂β = sin 2φ. Therefore, Equation (19) Using Equation (18), one can get after some trigonometry: Equation (6) can be rewritten using Equation (10) as: One can situate the origin of the (x 1 , y 1 )−coordinate system at the origin of the (x, y)−coordinate system and direct the x 1 −axis along the x−axis ( Figure 2). Then, Consider the terms of the equations in (20) containing the derivatives x α ∂ ∂ , y α ∂ ∂ , x β ∂ ∂ , and y β ∂ ∂ . Using Equations (22) and (24), one can represent these terms as: Using trigonometric identities, one can show that the coefficient of ishes in each of these equations. Therefore, equations (20) and (22) combine to give: Consider the terms of the equations in (20) containing the derivatives ∂x/∂α, ∂y/∂α, ∂x/∂β, and ∂y/∂β. Using Equations (22) and (24) Using trigonometric identities, one can show that the coefficient of ∂y 1 /∂α vanishes in each of these equations. Therefore, Equations (20) and (22) combine to give: ∂x ∂β − (y + x cos 2φ) = 0 and ∂y ∂α + y cos 2φ + x = 0.
Introduce new variables X and Y as Substituting Equation (27) into Equation (26) yields: It is seen from this equation that the quantities X and Y separately satisfy the equation of telegraphy. The latter is integrated by the method of Riemann.

Solution near a Traction-Free Surface
One of the principal stresses vanishes on any traction-free surface. The other principal stress and p are determined from Equations (1) and (2). It is necessary to consider two cases. One of these cases demands: on a traction-free surface. The other case demands: on a traction-free surface. Consider Equation (29). The general structure of the solution is shown in Figure 3. It is required to find X and Y at a point C. BC is an α−characteristic line, and AC is a β−characteristic line. Riemann's method for the telegraph equation results in    Here f is X or Y and.
Also, J 0 (z) is the Bessel function of the first kind of zero order, α C is the value of α at the point C, and β C is the value of β at the point C. Note that J 0 (0) = 1. and J 1 (0) = 0 (33) where J 1 (z) = −dJ 0 (z)/dz is the Bessel function of the first kind of first order. Since dβ = 0 on BC and dα = 0 on AC, Equation (31) can be rewritten as: Moreover, β = β C on BC and α = α C on AC. Therefore, z = 0 on these lines. Then, it follows from Equations (33) and (34) that: Here f A , f B , and f C are the value of f at A, B, and C, respectively. Put p 0 = 1/(1 + λ) in (12). Then, it follows from Equations (12), (13), and (29) that: , and dα = dϕ 2 sin 2φ (36) on AB.
Using Equation (36), one can rewrite the integral in Equation (35) as It follows from Equation (13) Equations (37) and (38) combine to give: The derivative ∂ f /∂P can be eliminated using Equation (28). In particular, it follows from Equations (28) and (38) that: Moreover, by replacing α and β in Equation (32) with ϕ and P using Equation (13), one gets: where ϕ C is the value of ϕ at C, P C is the value of P at C, and It follows from Equations (41) and (42) that: Substituting Equation (40) into Equation (39) and integrating by parts, one arrives at: Eliminating here the derivatives ∂G/∂P and ∂G/∂ϕ using Equation (43) yields: Substituting Equation (45) into Equation (35) leads to: The coefficients of X and Y in the integrands are known functions of ϕ due to Equations (36) and (42).
Consider Equation (30). In this case, AC is an α−characteristic line, and BC is a β−characteristic line (Figure 3). Putting p 0 = −1/(1 + λ) in Equation (12) and repeating the line of reasoning that has led to Equation (46) in the previous case, one gets: It follows from Equation (11) that: if Equation (29) is valid and Equation (30) is valid. Using Equations (48) and (49), one can eliminate ϕ C and P C in Equations (46) and (47). The free surface geometry determines the quantities X and Y involved in Equations (46) and (47). Since P = 0 on the free surface, it is seen from Equations (13) and (27) that: on AB. One of the principal stress directions is orthogonal to the free surface. Therefore, replacing ϕ with ψ in Equations (46) and (47) is convenient. It is seen from Equation (10) that dψ = dϕ. Then, using Equations (48)-(50), one rewrites Equations (46) and (47) as: and respectively.

Numerical Example
The distributions of x and y along AB are required to evaluate the integrals in Equations (51) and (52). The free surface shape determines these distributions. The numerical solution below is for an elliptic hole. Its equation can be written as: where x 1 = ρ cos θ and y 1 = ρ sin θ. Equations (17) and (18) have been derived upon the assumption that the x−axis coincides with the x 1 −axis. Therefore, the origin of the (x, y)−coordinate system should be situated at a point of the ellipse where the tangent to the α−is parallel to the x 1 −axis. This condition is equivalent to: at the origin of the (x, y)−coordinate system. Let θ 0 be the value of θ at which this condition is satisfied. The components of the unit vector orthogonal to the ellipse are readily determined from Equation (53). As a result, Then, the condition Equation (54) yields: The ellipse, (x, y)−coordinate system, and (x 1 , y 1 )−coordinate system are shown in Figure 4. The unit base vectors of the (x, y)−coordinate system are denoted as e x and e y , and the unit base vectors of the (x 1 , y 1 )−coordinate system as i and j. The position vector of a generic point is: Using Equation (56), one can represent the vector R 0 as: It follows from the geometry of Figure 4 that: e x · i = 1, e x · j = 0, e y · i = cos 2φ, e y · j = sin 2φ.
shown in Figure 4. The unit base vectors of the ( )   Using Equations (58) and (59), one can find the scalar products R · i and R · j from Equation (57) as: x 1 = x 0 + x + y cos 2φ and y 1 = y 0 + y sin 2φ.
One can eliminate r in these equations using Equation ( (51) to calculate X C and Y C . Having X C and Y C , one can calculate the x 1 − and y 1 −coordinates of point C using Equations (18), (27), and (60). The principal stresses at this point are found from Equations (1), (2), (12), and (48) (or (49)). The numerical results below are for Equation (29). Figure 5 depicts the variation of the stress σ 1 along the x 1 −axis for several values of λ and a. Figure 6 shows the variation of the stress σ 2 . In these figures, s is the distance from the hole surface, s = x 1 − a. Figures 5 and 6 reveal a significant effect of λ on the stresses. In particular, the effect of a is invisible in these figures. However, the hole shape does affect the stresses near the surface. Figure 7 depicts the variation of the principal stresses for several a-values at λ = −1/4. It is seen from this figure that the value of |σ 1 | near a circular hole (a = 1) is larger than its value near an elliptic hole. However, elliptical holes produce higher values of |σ 1 | than the corresponding circular hole at some distance from the hole's surface. The same tendency occurs with the stress |σ 2 |.

Conclusions
From this work, the following conclusions can be drawn: (1) The generalized method of moving coordinates is efficient for determining stress fields near traction-free surfaces. In particular, determining the stress components at any point of the plastic region requires evaluating one ordinary integral in Equations (51) or (52). (2) The solutions in Equations (51) and (52) are independent of other boundary conditions, though the shape of the plastic region is. (3) The solution found is identical to the solution in [5] if λ = 0 in Equation (1). (4) The parameter λ involved in Equation (1) has a much greater effect on the stress field than the aspect ratio ( Figures 5 and 6). However, it is a consequence of the yield criterion. For a given material, the aspect ratio affects the stress distribution in a small region near the hole (Figure 7). (5) The solution found is for perfectly plastic material. However, slip-line solutions for the class of problems considered agree well with experimental results for strain-hardening materials [25].

Conclusions
From this work, the following conclusions can be drawn: (1) The generalized method of moving coordinates is efficient for determining stress fields near traction-free surfaces. In particular, determining the stress components at any point of the plastic region requires evaluating one ordinary integral in Equations (51) or (52). (2) The solutions in Equations (51) and (52) are independent of other boundary conditions, though the shape of the plastic region is. (3) The solution found is identical to the solution in [5] if λ = 0 in Equation (1). (4) The parameter λ involved in Equation (1) has a much greater effect on the stress field than the aspect ratio ( Figures 5 and 6). However, it is a consequence of the yield criterion. For a given material, the aspect ratio affects the stress distribution in a small region near the hole (Figure 7). (5) The solution found is for perfectly plastic material. However, slip-line solutions for the class of problems considered agree well with experimental results for strain-hardening materials [25].