An Accurate Limit Load Solution for an Anisotropic Highly Undermatched Tension Specimen with a Crack

: Plastic anisotropy signiﬁcantly inﬂuences the behavior of structures subjected to various loading conditions. The extremum principles in the theory of rigid plastic solids are convenient and reliable tools for plastic design. The present paper combines the upper bound theorem and Hill’s quadratic yield criterion for orthotropic materials to evaluate the plastic collapse load of a highly undermatched welded tensile panel with a crack in the weld. The base material is supposed to be rigid. The shape of the crack is quite arbitrary. The orientation of the principal axes of anisotropy varies through the thickness of the weld. The upper bound solution is based on an exact solution for a layer of an anisotropic material. This feature of the upper bound solution is advantageous for increasing its accuracy. A numerical treatment is only necessary to ﬁnd the solution for the uncracked specimen. This specimen has two axes of symmetry, which simpliﬁes the solution. Simple analytic formulae transform this solution into a solution for the cracked specimens with one axis of symmetry and no symmetry. It is shown that the through-thickness distribution of anisotropic properties signiﬁcantly affects the limit load.


Introduction
Flaw assessment procedures, such as ETM, R6, SINTAP, FITNET, BS 7910, and API 579, are widely used in engineering fracture mechanics (for example, [1][2][3][4][5]). The plastic limit load is one of the key parameters affecting the accuracy of predictions based on these procedures [6]. Compendiums of limit load solutions for structures containing cracks are available in the literature [7][8][9]. Most of these solutions are for isotropic materials. Several solutions for anisotropic materials are reviewed in [9]. These solutions and the solutions found in [10][11][12] reveal a great effect of plastic anisotropy on the limit load. An overview of the state-of-the-art application of engineering fracture mechanics to weldments are provided in [13]. This work emphasizes the inhomogeneous microstructure of the weld. The microstructure of materials, among other factors, affects the mechanical anisotropy. Some of these factors are reviewed in [14] for steels used for pipelines. Only the yield criterion affects the plastic limit load. In particular, elastic properties are immaterial for calculating the plastic limit load [15]. In the present paper, it is assumed that the yield criterion represents the current state of anisotropy, independently of the factors that have caused it. The solutions [9][10][11][12] are based on the assumption that the parameters involved in the yield criterion are constant in the weld and the base material. The review [13] strongly suggests that the distribution of anisotropic properties is not uniform within the weld. Taking into account that the effect of plastic anisotropy on various parameters resulting from elastic/plastic solutions is large [16], it is natural to expect that those non-uniform distributions of anisotropic properties may significantly affect the plastic limit load. A plane strain solution for a highly undermatched middle cracked tensile panel is proposed in [17]. This solution assumes that the weld is orthotropic and obeys Hill's quadratic yield criterion [18]. The thickness of the weld is constant. The present paper proposes a two-fold extension of the solution [17]. Firstly, quite an arbitrary through-thickness distribution of the parameters involved in Hill's yield criterion is assumed. Secondly, the shape of the crack is also quite arbitrary, though some restrictions apply. The method proposed in [19] is used to extend the solution for symmetric specimens to a class of the corresponding asymmetric specimens.
Other upper bound solutions for cracked specimens that account for plastic anisotropy are proposed in [10][11][12]. The solutions [11,12] are for over-matched specimens. The general method developed in the present paper is not applicable for such specimens. However, the solution [10] is for highly undermatched scarf joints. Therefore, the method developed is immediately applicable to generalize the solution [10] to non-uniform distributions of anisotropic properties and the crack of arbitrary shape.

Statement of the Problem
A schematic diagram of the weld joint under consideration is shown in Figure 1. The width of the specimen is 2B, and its thickness is W. The shape of a through-crack is quite arbitrary. The specific restrictions on this shape will be formulated later. The specimen is subject to tension by two forces Q. It is required to evaluate the value of Q at plastic collapse under plane strain deformation. It is assumed that the weld is much softer than the base material. Therefore, plastic yielding occurs in the weld, whereas the base material is rigid. For this reason, the plastic properties of the base material are immaterial. The elastic properties of both the weld and base material are immaterial at plastic collapse, as follows from the theorem proven in [15]. The weld is plastically anisotropic. Let σ αα , σ ββ , and σ αβ be the stress components referred to as the principal axes of anisotropy in planes of flow. Under plane strain conditions, Hill's quadratic yield criterion reads as follows [18]: Here, s αβ = σ αα − σ ββ 2, T is the yield stress in shear with respect to the (α, β) axes, and c is expressed in terms of the yield stresses with respect to the principal axes of anisotropy. The expression for c is provided in [18]. For the present paper, it is sufficient to know that the following holds: In the case of isotropic materials, c = 0. The upper bound theorem involves the plastic work rate. Let ξ αα , ξ ββ , and ξ αβ be the strain rate components referred to as the principal axes of anisotropy (α, β) . The general expression for the plastic work rate per unit volume is derived in [18]. This expression reduces to the following: under plane strain deformation. It is taken into account that the material is incompressible and, therefore, ξ αα = −ξ ββ . We introduce a Cartesian coordinate system (x, y) in the planes of flow. Let θ be the inclination of the α-axis to the x-axis, measured counterclockwise. Then, the following holds: Here, σ xx , σ yy and σ xy are the stress components referred to as the Cartesian coordinates. Equation (4) can be rewritten as follows: Analogously, Here, ξ xx , ξ yy , and ξ xy are the strain rate components referred to as the Cartesian coordinates.

Solution for the Specimen with No Crack
It will be seen later that the solution for the specimen with no crack can be immediately adopted to evaluate the limit load of the specimen with quite an arbitrary crack. A schematic diagram of the specimen considered in this section, together with the Cartesian coordinate system, is shown in Figure 2. The axes of this coordinate system coincide with the symmetry axes of the specimen.

An Exact Solution for a Layer of Anisotropic Material
Even though the upper bound theorem does not require that the equilibrium equations are satisfied, it is advantageous to choose a kinematically admissible velocity field so that the associated stress field satisfies the equilibrium equations and associated flow rule. An exact solution for the compression of an anisotropic material layer between two parallel plates is provided in [20]. Since the material is incompressible, this solution can be immediately adopted for the tension of the layer. Because of symmetry, it is sufficient to consider the solution in the first quadrant, x ≥ 0 and y ≥ 0 . The rigid base material moves with velocity V along the y-axis. The magnitude of this velocity is immaterial.
Let u x and u y be the velocity components referred to the Cartesian coordinate system. The velocity field found in [20] is as follows: Here, u 0 is constant, χ is a dummy variable of integration, and cos 2ϕ is determined as follows: Here, s xy = σ xx − σ yy 2. The pair s xy and σ xy involved in (8) should satisfy the yield criterion. The solution [20] is valid if ϕ is independent of x. In what follows, it is convenient to use the dimensionless quantities: Then, Equation (7) becomes the following: Using (7), one can find the strain rate components as follows: Equation (1) is satisfied by the following substitution: Here, γ is an arbitrary function of ζ. Equations (5) and (12) combine to give the following: It is seen from (8) and (13) that the condition that ϕ is a function only of ζ is satisfied if θ is a function only of ζ. The latter is assumed in the present paper. It follows from (8) and (13) that the following holds: If σ xy has the same value at y = h and y = −h (Figure 2), the solution for the shear stress given in [20] reads as follows: Here, µ is a constant whose value is determined by the boundary condition at ζ = 1. The second equation in (13) and (15) combine to provide the following: Using the trigonometric identities sin 2γ = 2 tan γ 1 + tan 2 γ and cos 2γ = one transforms (16) to the following: It is a quadratic equation for tan γ. Its solution is the following: Four values of γ are determined from this equation. The unique solution satisfies the following inequalities: and The latter follows from the assumption that s xy ≤ 0 and the first equation in (13). Substituting γ found from (19) into (14), one finds ϕ as a function of ζ. Then, Equation (11) supplies the strain rate components as functions of ζ.

Upper Bound Solution
The general structure of the solution is illustrated in Figure 3. The velocity field (10) is kinematically admissible in the plastic region. In particular, it is seen from (11) that the incompressibility equation, ξ xx + ξ yy = 0, is satisfied. It follows from (10) that u y = 0 at ζ = 0 and u y = V at ζ = 1. The velocity component u x does not satisfy the condition u x = 0 at η = 0. Therefore, there must exist a rigid region that contains the y-axis (Figure 3). This region moves along the y-axis, together with the rigid base material. Therefore, its velocity is represented as follows: where j is one of the unit base vectors of the Cartesian coordinate system. The other unit base vector is denoted as i. The velocity vector in the plastic region is represented as follows: Here, Equation (10) supplies u x and u y . Let φ be the inclination of the tangent to velocity discontinuity line OA to the x-axis, measured counterclockwise ( Figure 3). Then, the unit vector normal to OA is represented as follows: The velocity component normal to the velocity discontinuity line OA must be continuous. Therefore, U r · n = U p · n. Substituting (22), (23), and (24) into this equation and using (10), one obtains the following: It follows from the geometry of Figure 3 and (9) that the following holds: Equations (25) and (26) combine to give the following: It is a linear differential equation for determining the shape of the velocity discontinuity line. This line should pass through the origin of the Cartesian coordinate system. Therefore, the boundary condition to Equation (27) is η = 0 for ζ = 0. The solution of Equation (27) satisfying this boundary condition is as follows: Integrating by parts, one obtains the following: Equations (28) and (29) combine to provide the following: The velocity discontinuity line may have a common point with the line ζ = 1 only if the numerator in (30) vanishes at ζ = 1. Therefore, the following holds: Substituting (31) into (30) and applying l'Hospital's rule, one determines the value of η at point A (Figure 3) as the following: The present solution is valid if the following holds: The final expression for the velocity discontinuity line follows from (30) and (31) as follows: Using (14) and (19), one can evaluate the integrals numerically. The infinitesimal length element of the velocity discontinuity line is dL = (dx) 2 + (dy) 2 = dx dy 2 + 1dy.
Here, the derivative dx/dy is determined from (26) and (27). Then, using (9) and (30), one can represent the infinitesimal length element as follows: where Ω 1 (y) is a known function of y. Alternatively, using (9), one transforms the expression dL = (dx) 2 + (dy) 2 into the following: The derivative dη OA dζ is determined from (27) and (31) as follows: Using (34) and (37), one obtains the coefficient of dζ in (36) as a function of ζ. Then, the infinitesimal length element is represented as follows: where Ω 2 (ζ) is a known function of ζ. The amount of velocity jump across the velocity discontinuity line is represented as [u] = U r − U p . Substituting (22) and (23) into this equation, one arrives at the following: The velocity components u x and u y are understood to be calculated at the velocity discontinuity line. Therefore, Equations (10), (31) and (39) combine to give the following: The strain rate components involved in (3) are determined from (6) and (11) as follows: Equations (3) and (41) combine to give the following: It follows from the upper bound theorem and Equation (35) that the following holds: The velocity component u x is understood to be calculated at ζ = 1. The first integral in (43) represents the plastic work rate in the plastic region, the second at the velocity discontinuity line OA, and the third at the velocity discontinuity line AC (Figure 3). Then, τ is the actual shear stress over the velocity discontinuity line OA and τ s is the actual shear stress over the velocity discontinuity line AC.
Using (9) and (38), one can rewrite (43) as follows: Using (42), one transforms the first integral in (44) to the following: The integrand is a known function of ζ. Therefore, the integral can be evaluated numerically.
Using (36) and (40), one transforms the second integral in (44) to the following: Since τ is the actual shear stress over the velocity discontinuity line, its value is unknown. Therefore, the integral in (46) cannot be evaluated. In the case of isotropic materials, this difficulty is resolved by replacing τ with the shear yield stress, which is a material constant. Since τ is less or equal to the shear yield stress, this replacement cannot reduce the value of the integral in (46). Based on this argument, one concludes that the resulting right-hand side of (44) provides an upper bound on Q. However, the shear yield stress in anisotropic materials depends on the plane's orientation over which it applies. The second equation in (13) gives the shear stress over any plane when the yield criterion is satisfied. It follows from this equation that we have the following: One can replace τ with τ m in (44) and (46) to arrive at the following: Here, Q u is the upper bound limit load and q u is its dimensionless representation. The essential difference between the OA and AC velocity discontinuity lines is that the bi-material interface is an envelope of characteristics in the exact solution. Therefore, τ s involved in the third integral in (48) is as follows [18]: Here, θ s is the value of θ at ζ = 1. It follows from (10), (31), and (32) that the following holds: Substituting (49) and (50) into (48), one obtains the following: (51) Using (45), (47) and the second equation in (48), one can calculate q u from (51). It follows from (15) and (49) that the following holds:

Numerical Examples
It is seen from the previous section that a numerical technique is only necessary for evaluating ordinary integrals. The Wolfram Language function NIntegrate is used for this purpose. It is assumed that b = 5 in all calculations presented in this section. Figures 4 and 5 depict the variation of the dimensionless limit load with c if θ is a constant function and a linear function of ζ, respectively. The linear function is defined as θ = θ 0 ζ, where θ 0 is constant. Figure 6 illustrates the effect of the through-thickness variation of θ on the dimensionless limit load. In this case, θ is represented by second-order polynomials as follows: where δ is constant. The particular case δ = 0 corresponds to one of the solutions shown in Figure 4. The range of δ is restricted by the condition that sin 2θ is a monotonic function of ζ in the range 0 ≤ ζ ≤ 1. All the curves in Figures 4 and 5 intersect at c = 0 (isotropic material). The effect of θ ( Figure 4) and θ 0 ( Figure 5) on the limit load increases as |c| increases. The qualitative behavior of the solution is different in the ranges of c > 0 and c < 0. In particular, q u is an increasing function of |c| for all θ (Figure 4) and all θ 0 Figure 5) in the range c < 0. However, it may be a decreasing or increasing function of c in the range c > 0. Therefore, predictions of flaw assessment procedures based on the assumption of isotropic material may be non-conservative if c > 0.  Figure 6 shows that the limit load is sensitive to the through-thickness distribution of anisotropic properties. It is particularly sensitive if δ is small. Small δ-values mean that θ is close to π/2 everywhere, as follows from (53).

Solution for the Specimen with a Crack
This section demonstrates that the solution derived in Section 3 can be immediately adopted to find an upper bound limit load for a class of the specimens with a crack.

Crack on the x-Axis
As shown in Figure 2, the boundary value problem is symmetric relative to the x-axis. Therefore, it is sufficient to consider the domain y ≥ 0. A schematic diagram of the weld and the general structure of a kinematically admissible velocity field are shown in Figure 7. It is seen from this figure that the solution from Section 3 is valid on the right to line l 1 if B is replaced with B 1 . Analogously, the solution from Section 3 is valid on the left to line l 2 if B is replaced with B 2 . Then, we have the following: Here, q u and q (2) u are determined from (48) at b = b 1 = B 1 h and b = b 2 = B 2 h, respectively. The factor 1/2 appears because (51) gives the limit load for the specimen whose width is 2B.
Both the velocity discontinuity lines (Figure 7) should satisfy the inequality (33). Of practical importance is the inequality η A ≤ 1. Putting η A = 1 in (32), one obtains the following: The solution is valid if b 1 ≥ b m and b 2 ≥ b m . Figure 8 depicts the variation of b m with c at several values of θ, assuming that θ is constant.

Arbitrary Crack
The boundary value problem is not symmetric relative to any axis. Therefore, it is necessary to consider the entire domain occupied by the weld. A schematic diagram of the weld and the general structure of a kinematically admissible velocity field are shown in Figure 9. The kinematically admissible velocity fields on the right to l 1 and left to l 2 are taken the same as in the specimen shown in Figure 7. For the entire velocity field to be kinematically admissible, it is necessary to introduce two additional velocity discontinuity lines, a 1 c 1 and a 2 c 2 . Then, the only restriction on the shape of the crack is that it is compatible with the motion of the two rigid regions along the y-axis in the opposite directions. The contribution of the additional velocity discontinuity lines to the right-hand side of (43) involves τ. As in the case of the velocity discontinuity line OA (Figure 3), one can replace τ with τ m given by (47). Then, the contribution of the additional velocity discontinuity lines to the right-hand side of (43) can be represented as follows: Figure 9. Schematic diagram of the weld and the general structure of a kinematically admissible velocity field in the case of an arbitrary crack.
Here, L 1 and L 2 are the lengths of lines a 1 c 1 and a 2 c 2 , respectively. Using (47), (54), and (56), one arrives the following expression for the dimensionless limit load:

Numerical Examples
Given the mechanical properties, Equations (54) and (57) allow the limit load to be easily calculated for the specimen having any crack if the dependence of q  Figure 10 depicts the variation of q (i) u with b i at several values of θ, assuming that θ is constant. Figure 11 shows the variation of q   The limit load for isotropic materials is an increasing function of b. An important feature of the solution illustrated in Figures 10 and 11 is that the limit load considered a function of b may have a minimum if the material is anisotropic. This feature should affect the prediction of flaw assessment procedures.

Conclusions
A new upper bound limit load solution is proposed for a highly undermatched tension specimen containing a crack. The weld is supposed to be plastically anisotropic. The distinguishing features of this solution are (i) through-thickness variation in the plastic anisotropy, and (ii) crack of an arbitrary shape. Some restrictions of the variation in the plastic anisotropy and the shape of the crack apply, but they are not restrictive for practical applications.
An advantage of the solution is that the stress field associated with the kinematically admissible velocity field chosen satisfies the equilibrium equations. Therefore, this solution should be more accurate than other upper bound solutions of the same level of complexity.
The solution is practically analytic. A numerical technique is only required for evaluating ordinary integrals. It is particularly important because the material properties are classified by an arbitrary function rather than a set of parameters. Therefore, no parametric study is possible. However, the calculation of the limit load for any given function is fast and straightforward.
The numerical example has revealed the following interesting features of the solution: (i) predictions of flow assessment procedures based on the assumption of isotropic materials may be non-conservative for certain anisotropic materials, and (ii) the limit load is represented by a non-monotonic function of the relative thickness of the weld for certain anisotropic properties.
The solution is ready for use in conjunction with flaw assessment procedures. Such procedures are widely accepted in practical applications [1].

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author upon reasonable request.