An Upper Bound Solution for the Compression of an Orthotropic Cylinder

The upper bound theorem is used in conjunction with Hill’s quadratic yield criterion for determining the force required to upset a solid cylinder. The kinematically admissible velocity field accounts for the singular behavior of the real velocity field in the vicinity of the friction surface if the maximum friction law is adopted. The regime of sticking is also taken into consideration. The effect of this regime on the upper bound limit load is revealed. In particular, the kinematically admissible velocity field that includes the regime of sticking may result in a lower upper bound than that with no sticking. The boundary value problem is classified by a great number of geometric and material parameters. Therefore, a systematic parametric analysis of the effect of these parameters on the compression force is practically impossible. An advantage of the solution found is that it provides a quick estimate of this force for any given set of parameters.


Introduction
The upsetting of solid cylinders is an important metal-forming operation [1]. Moreover, this process is used as a test for evaluating flow stress and friction [2]. An experimental technique for accurately determining the strains on the cylindrical surface and the flat ends of a cylindrical compression specimen was developed and applied in [3]. Another experimental method for evaluating strain inhomogeneity was proposed in [4]. Hot upsetting tests on steel cylindrical specimens were carried out in [5,6]. The development of damage in the upsetting of cylinders was studied in [7]. A measure of barreling appearing in the cylinders' upsetting was used to study frictional conditions and lubricant properties in [8][9][10][11][12][13].
The experimental studies above should be complemented with theoretical solutions. For this purpose, the finite element method is widely used (for example, [2,5,6,13]. However, the upper bound theorem also provides an efficient method for finding an approximate solution to these boundary value problems. It is worthy to note that the upper bound method applies to micro-forming [14] and novel metal forming processes (for example, the process for extruding curved profiles [15]). Sometimes, the upper bound method is even more efficient than the finite element method. In particular, solutions for several materials models are singular in the vicinity of maximum friction surfaces [16,17]. Finite element solutions based on standard shape functions do not converge in this case [18,19]. On the other hand, upper bound solutions incorporate this singularity with no difficulty (for example, [20,21]). Other advantages of the upper bound method over the finite element method were summarized in the recent paper [22]. A typical disadvantage of upper bound solutions is that the regime of sticking friction is ignored (for example, [23][24][25][26][27]). On the other hand, such a regime inevitably occurs over a certain portion of the friction surface in many metal forming processes.
Cylindrical orthotropy of material properties is frequently generated in the course of such processes as axisymmetric extrusion and drawing. This material property is often taken into account in structural analysis and design (for example, [28][29][30][31][32][33]). However, as follows from the review above, the analysis of cylinder's upsetting is usually based on isotropic models. An exception is the solution provided in [34] for the upsetting of a hollow cylinder. This solution is based on Hill's quadratic yield criterion [35]. However, it ignores the regime of sticking. Therefore, its practical value is questionable. The present paper also adopts Hill's quadratic yield criterion. The effect of plastic anisotropy on the force required to upset a solid cylinder is demonstrated. Moreover, the effect of the regime of sticking on the upper bound limit load is revealed. In particular, the kinematically admissible velocity field that includes this regime may result in a lower upper bound than that with no sticking.
The overall motivation of this research is to demonstrate, using a simple example, that plastic anisotropy, which is a very common property of metallic materials, should not be ignored in upper bound solutions for metal forming processes.

Statement of the Problem
A circular solid cylinder is upset between two parallel rough plates ( Figure 1). The radius of the cylinder is R 0 , and its height is 2H 0 . Each plate moves with velocity U. The force P applies to each plate. This force should be found from the solution. The material of the cylinder is plastically orthotropic. The principal axes of anisotropy coincide with the radial, circumferential, and axial directions. Hill's quadratic yield criterion is adopted [35]. Under the conditions of axial symmetry, this criterion reads as follows.
where σ rr , σ θθ , σ zz , and σ rz are the components of the stress tensor referred to a cylindrical coordinate system (r, θ, z). The z-axis coincides with the axis of symmetry of the cylinder, and the plane z = 0 coincides with the plane of symmetry of the cylinder. Since z = 0 is the plane of symmetry for the flow, it is sufficient to find the solution in the region z ≥ 0. The coefficients involved in (1) are material constants. Let Υ, Θ, and Z be the tensile yield stresses in the radial, circumferential, and axial directions, respectively, and S be the shear yield stress in the rz-plane. Then, we have the following: The equivalent strain rate is determined as follows [35]: where ξ rr , ξ θθ , ξ zz , and ξ rz are the components of the strain rate tensor referred to the cylindrical coordinate system. The plastic work rate per unit volume is as follows [35]: The velocity boundary conditions are as follows: for z = 0 and u z = −U (6) for z = H 0 . Moreover, we have the following: for r = 0. Here, u r and u z are the radial and axial velocities, respectively. The lateral surface of the cylinder is traction-free. Friction occurs at the surface z = H 0 . It is assumed that the friction stress is equal to a constant fraction of the local shear yield stress. Taking into account the orientation of the principal axes of anisotropy, one obtains the friction stress as τ f = mS where m is constant. Then, we have the following for z = H 0 : The friction law (8) is valid if the regime of sliding occurs.

Kinematically Admissible Velocity Field
A general kinematically admissible velocity field for a class of axisymmetric problems was proposed in [21]. This paper deals with isotropic materials. However, the same velocity field is kinematically admissible for anisotropic materials. The general kinematically admissible velocity field can be reduced to a form appropriate for the upsetting of cylinders. The general structure of this velocity field is shown in Figure 2. The rigid region moves together with the plate. Therefore, its velocity is U. The rigid and plastic regions are separated by a velocity discontinuity line ac. This line must pass through the origin of the coordinate system. The shape of this line should be found from the solution. The solution below is valid if the radial coordinate of point c is less or equal to the radial coordinate of point b. Since the axis of symmetry belongs to the rigid region, the boundary condition (7) is satisfied. The kinematically admissible velocity field in the plastic region is as follows: where ρ = r/R 0 , ζ = z/H 0 , h = H 0 /R 0 and f (ζ) is an arbitrary function of its argument. One can readily verify that the velocity field (9) satisfies the incompressibility equation ∂u r ∂r + u r /r + ∂u z ∂z = 0. In addition, the axial velocity satisfies the boundary conditions (5) and (6). The strain rate components involved in (3) are determined from (9) as follows: In what follows, it is convenient to use the following dimensionless strain rate components: Equations (10) and (11) combine to give the following: Let i and j be the unit base vectors of the axes r and z, respectively. Then, the velocity vector in the rigid region is represented as U r = −Uj, and in the plastic region as U p = u r i + u z j. The velocity normal to the velocity discontinuity line must be continuous. Then, U r · n = U p · n, or the following holds: The velocity components u r and u z are understood to be calculated at the velocity discontinuity line using (9). From the geometry of Figure 2, we have the following: where ϕ is the orientation of the tangent to the velocity discontinuity line relative to the r-axis. Substituting (9) and (14) into (13) and using the dimensionless coordinates yields the following: Using the substitution η = ρ 2 , one transforms (15) into the following linear ordinary differential equation: The solution of this equation supplies the shape of the velocity discontinuity line. Since this line must pass through the origin of the coordinate system, the boundary condition to Equation (16) is as follows: for ζ = 0. The general solution of (16) can be represented as follows: where C is constant and The denominator in (18) approaches zero as ζ → 1. Therefore, the velocity discontinuity line can reach the friction surface ζ = 1 only if C = 0. Then, Equation (18) becomes the following: This equation determines the velocity discontinuity line. Applying l'Hospital's rule and taking into account (19), one obtains the following: where η c is the value of η at point c ( Figure 2). Having found the velocity discontinuity line, one can readily determine an infinitesimal length element of this line from (19) and (20) as follows: The amount of velocity jump across the velocity discontinuity line is [u] = U r − U p . The velocity vectors are understood to be calculated at this line. Then, using (9) and (20), one arrives at the following:

Upper Bound Solution
It follows from the upper bound theorem that [35] where W V is the plastic work rate in the plastic region, W d is the plastic work rate at the velocity discontinuity line, and W f is the plastic work rate at the friction surface. The infinitesimal volume element in the cylindrical coordinate system is dV = rdrdθdz. Using the dimensionless coordinates, one arrives at the following: The plastic work rate in the plastic region is determined as follows: Substituting (3) into (4) and using (11), one obtains the following: where (2) is used. In (12), one can express ρ in terms of η. As a result, we have the following: Eliminatingξ rr ,ξ θθ ,ξ zz , andξ rz in (28) using (29), one arrives atω as a function of η and ζ. Equation (26) becomes the following: The integral here can be evaluated numerically. The general expression for the plastic work rate at the velocity discontinuity line is as follows: where τ s is the shear stress on the velocity discontinuity line. It is known that [35] where Using some trigonometry, one transforms (32) to the following: where the angle ϕ is understood to be calculated at the velocity discontinuity line. Therefore, from (14) and (16), we have the following: Equations (22) and (23) yield the following: where Substituting (34) and (36) into (31) gives the following: Taking into account (20), (36), and (37), one can evaluate the integral in (38) numerically.
Using (8), one can represent the plastic work rate at the friction surface as follows: The radial velocity is understood to be calculated at the friction surface. Therefore, using (9), one transforms (39) to the following: where Equation (21) is used. The integral in (40) can be immediately evaluated to give the following: Substituting (30), (38), and (41) into (24) results in the following: where P u is the upper bound of the force required to deform the cylinder, and p u is its dimensionless representation. The right-hand side of this equation can be found, using the procedure above for the arbitrary function f (ζ) and any set of parameters. The function f (ζ) may involve additional parameters. In this case, the right-hand side of (42) should be minimized with respect to these parameters to find the best upper bound based on the kinematically admissible velocity field chosen. The solution above is valid if the following holds: The solution with no rigid region can be obtained from the solution above, with the following: In this case, there is no velocity discontinuity line, and Equation (24) becomes as follows: The regime of sliding occurs over the entire friction surface. Therefore, Equations (9), (39) and (44) combine to give the following: Moreover, using (27), one transforms Equation (26) to the following: Substituting (44) into (29) leads to the following: Then, Equation (28) becomes the following: Substituting (49) into (47) and integrating gives the following: Equations (45), (46) and (50) supply the following dimensionless upper bound limit load: (51)

Numerical Examples
One should choose the function f (ζ) to evaluate the right-hand side of (42). This function specifies the kinematically admissible velocity field. In general, it is advantageous that kinematically admissible velocity fields exhibit some mathematical properties of the real velocity field. The process under consideration is symmetric relative to the plane ζ = 0. Therefore, the real velocity field is described by an even function of ζ. Moreover, the real velocity field is singular near the friction surface if m = 1 in (8). In particular, we have the following [17]: as ζ → 1. Here, A is independent of ζ. The general asymptotic analysis carried out in [17] is for plane strain deformation. However, a particular solution presented in [36] shows that (44) is also valid for axisymmetric problems. One of the simplest functions satisfying the symmetry condition and (44) is the following: Taking into account (21), one rewrites this equation as the following: Equations (19) and (54) combine to give the following: It follows from (17) and (20) that Then, Equations (55) and (56) combine to give β 1 = −2η c (πh). Eliminating β 1 in Equations (54) and (55), one obtains the following: Moreover, Substituting (49) and (50) into (42), one arrives at the right-hand side of this equation as a function of ζ. Moreover, this function involves one free parameter, η c . One should minimize the right-hand side of (42), taking into account (43).
The boundary value problem is classified by five dimensionless parameters: Υ/Z, Θ/Z, S/Z, h, and m. Therefore, its systematic parametric analysis is practically impossible. The numerical example below includes the isotropic case, the anisotropic properties determined experimentally and reported in [35], and three sets of arbitrarily chosen yield stresses. The isotropic solution is used for showing the effect of plastic anisotropy on the limit load. The solutions for the chosen sets of yield stresses allow one to gain insight into the possible effect of plastic anisotropy on the interpretation of the corresponding friction test.
The experimental results from [37] are summarized in Table 1. Table 2 shows the chosen sets of yield stresses. The numerical integration has used 20 Gaussian integration points and weights for the integrals involved in (30) and (38). With a given set of the parameters that classify the boundary value problem, the minimization process to find the optimal value of η c is conducted by using hundreds to thousands of η c -data points in the range between 0 and 1 to calculate p u . The number of data points determines the resolution and accuracy of the optimal value of η c . Once the η c , p u data list is calculated, it is straightforward to find the minimum value of p u . This value is denoted as p (1) u . The limit load found from (51) is denoted as p (2) u . Then, the solution to the boundary value problem is as follows: The optimal value of η c is also of interest because it may affect the interpretation of the friction test.
It is natural to expect that p (1) u = p (2) u for certain combinations of the parameters that classify the boundary value problem. Figure 3 illustrates this feature of the solution for the isotropic material at m = 0.7. In this case, the only parameter that varies is h. It is seen from Figure 3 that p (1) u = p (2) u at h = h c . It follows from (59) that the solution (42) is valid in the range h ≤ h c and the solution (51) in the range h ≥ h c . It is found that Case 2 results in p (1) u < p (2) u if 0 ≤ m ≤ 1. Therefore, in this case p u = p (1) u . Figure 4 depicts the dependence of h c on m for the isotropic material and three anisotropic materials whose properties are given in Tables 1 and 2. It is seen from this figure that the range of validity of the solution with no rigid region decreases as the friction factor increases for all the cases considered.   (1) u < p (2) u in the m range. Figure 5 illustrates the limit load solution at m = 1. The curves show the dependence of p u on h for the isotropic material and four anisotropic materials whose properties are given in Tables 1 and 2. As expected, p u increases with h. The limit load for anisotropic materials may be lower or higher than that for the isotropic material. The difference between the limit load for the isotropic material and that for the anisotropic material whose properties are determined experimentally (Table 1) is between 15% and 20% in the range of h considered.

Conclusions
A new upper bound solution for the upsetting of a circular cylinder is proposed. A distinguishing feature of the solution is that plastic anisotropy and the existence of a rigid region are taken into account. The existence of the rigid region automatically means that there is a region of sticking friction. From this work, the following conclusions can be drawn: 1.
Plastic anisotropy affects the limit load required to deform the specimen. It may increase or decrease the limit load as compared to the isotropic case ( Figure 5).

2.
The upsetting of a cylinder is often used as a friction test. Plastic anisotropy significantly affects the region of sticking friction ( Figure 6). Since Equation (8) is not valid in this region, this effect of plastic anisotropy should be taken into account in the interpretation of the friction test results.

3.
Five dimensionless parameters classify the boundary value problem. Therefore, its systematic parametric analysis is invisible. An advantage of the proposed solution is that it quickly estimates the upper bound limit load for a given set of parameters. 4.
The real velocity field is singular near the friction surface if m = 1 in Equation (8).
The solution proposed accounts for this singularity, which is impossible when using ordinary finite element solutions. Institutional Review Board Statement: Not applicable.
Informed Consent Statement: Not applicable.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author upon reasonable request.
Acknowledgments: This research was supported by the RUDN University Strategic Academic Leadership Program. We thank the National Center for High-Performance Computing (NCHC) for providing computational and storage resources.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Nomenclature
i and j unit base vectors m friction factor n unit vector normal to the velocity discontinuity line p u dimensionless upper bound limit load (t, θ, z) cylindrical coordinate system u r and u z radial and axial velocities F, G, H and M material parameters introduced in Equation (2) H 0 half-height of the cylinder P force R 0 radius of the cylinder S shear yield stress in the rz-plane U velocity of the plate W d plastic work rate at the velocity discontinuity line W f plastic work rate at the friction surface W v plastic work rate in the plastic region Υ, Θ and Z tensile yield stresses in the radial, circumferential, and axial directions ξ eq equivalent strain rate ξ rr , ξ θθ , ξ zz and ξ rz strain rate components referred to the cylindrical coordinate system ξ rr , ξ θθ , ξ zz and ξ rz dimensionless strain rate components introduced in Equation (11) ρ, ζ and h dimensionless quantities introduced after Equation (9) σ rr , σ θθ , σ zz and σ rz stress components referred to the cylindrical coordinate system τ f friction stress ϕ orientation of the tangent to the velocity discontinuity line ω plastic work rate per unit volume