1 Introduction

Titanium and titanium alloys are materials widely used for the biomedical applications. For the more demanding implants like those of hip or knee joints, Ti6Al4V alloy is preferred over pure titanium due to its superior strength properties. However, aluminium and vanadium ions released into the human organism have a toxic effect [1]. Commercially pure (CP) titanium does not have this flaw, but it has lower strength than Ti6Al4V alloy. Therefore, it is imperative to increase its properties. According to Hall–Petch law, this can be achieved by reduction of the grain size [2, 3].

The methods of attaining micro- or nano-grains in metals can be divided into two groups: thermo-mechanical processing (TMP) and severe plastic deformation (SPD) methods [4].

During TMP processes, material behaviour is associated with complex mechanisms such as thermal activation, dynamic recovery, strain hardening and dynamic recrystallization, leading to grain fragmentation and creation of the micro- and subgrains [5, 6]. The main advantage of those techniques is that they make it possible to produce usable ingots suitable for mill operations.

SPD methods like equal channel angular pressing (ECAP), high-pressure torsion (HPT), accumulative roll bonding (ARB) and others [7,8,9,10,11,12] allow to device ultrafine and nano-grained structure, greatly increasing material strength, hardness and ductility [13]. However, severe plastic deformation methods utilized to develop submicron or even nanoscale structures can be used only on laboratory scale. The major problems encountered with these techniques are the limitations placed on size, dimensions and geometry of the manufactured components. These constraints have restricted the industrial application of SPD methods, including the production of bioprostheses, which generally are of complex geometries [14].

Complex nature of phenomena occurring during thermoplastic deformation can be described using constitutive equations. Furthermore, application of such mathematical models into finite element simulation allows one to control microstructure and mechanical properties [15,16,17]. Multiple models and constitutive equations have been proposed over time. The most commonly used were phenomenological models like Arrhenius equation, Johnson–Cook (JC) and Khan–Huang–Liang (KHL) models, or physical-based models, such as Zerrili–Armstrong model and others [13, 18, 19]. In our previous work, Arrhenius equation was chosen as the potentially most promising. However, our experimental and analytical studies showed that this equation was accurate in the case of the pure titanium only in certain circumstances [20, 21]. Therefore, other models were taken into consideration. We have chosen JC, KHL and Anand models to mathematically describe the phenomena taking place in Ti during the TMP process. Not only do the models describe deformation of metals under high strain rates and temperatures but also they take into consideration the material strain level [22]. KHL model, however, handles material constants differently. Anand equation is an elastic–viscoplastic model mostly used in works regarding deformation of various solders. According to various authors [23,24,25,26,27], it is possible to apply the equation for different purposes/materials as the material behaviour during deformation is similar to that observed during titanium deformation.

The goal of this paper was to apply the aforementioned models to mathematically describe the processes occurring in the conducted thermo-mechanical experiments. The practical aim of the study was to adopt a proper constitutive model for the TMP process conducted on pure Ti to enhance its mechanical properties. The required material constants were calibrated on the basis of experimental tests with the future intention to implement the model into a finite element method (FEM) package, in order to simulate and optimize the titanium thermo-mechanical deformation.

2 Methods

2.1 Experimental procedure

The experiments of thermoplastic deformation were conducted on samples of the annealed Grade 2 titanium of chemical composition presented in Table 1 with the use of the Gleeble 3500 thermoplastic simulator—the equipment for dynamic thermal–mechanical testing of materials and physical simulation of processes.

Cylindrical specimens 12 mm in height and 10 mm in diameter were quickly heated to the planned temperatures and then compressed until strain equal to 0.6 was reached. Tests were conducted with the constant strain rates varied from 10\(^{\mathrm {-2}}\) s\(^{\mathrm {-1}}\) to 10\(^{\mathrm {1}}\) s\(^{\mathrm {-1 }}\) (four different strain rates values per one temperature set). The detailed plan of the experiment is given in Table 2.

Table 1 Maximum content of impurities in Grade 2 titanium in mass percentages
Full size table
Table 2 Parameters of the experiments
Full size table

Stress–strain graphs for each sample were created to utilize them to curve-fitting procedure of the constant identification. A simple correlation was be observed, i.e. the maximal stress for 775 K temperature was lower than that for 875 K.

In order to confirm the grain fragmentation and measure the grain size, we performed microstructure observations using scanning electron microscope (SEM) Hitachi Su70 and Hitachi SU-8000 and scanning transmission electron microscope (STEM) Hitachi HD2700. The electron backscatter diffraction (EBSD) tests were performed on the HKL Nordlys II detector coupled to the Hitachi Su70 microscope.

2.2 Constitutive equations of TMP

2.2.1 Arrhenius-type equation

The proposed hyperbolic sine Arrhenius equation is the most commonly used constitutive model due to its simplicity and low number of required constitutive constants [28,29,30]. As a form of modification, the temperature compensated strain rate, known as Zener–Hollomon parameter, was incorporated into the Arrhenius equation expressed by formula (1):

$$\begin{aligned} \dot{\varepsilon }=AF\left( \sigma \right) \exp \left( -\frac{Q}{RT} \right) , \end{aligned}$$
(1)

where:

$$\begin{aligned} F\left( \sigma \right) =\left\{ {\begin{array}{lll} \sigma ^{n_{1}} ,&{} if &{} \quad \alpha \sigma <0.8, \\ \exp \left( \beta \sigma \right) ,&{} if &{} \quad \alpha \sigma >1.2, \\ \left[ \sinh \left( \alpha \sigma \right) \right] ^{n}, &{} &{}\quad \text { for all } \sigma , \\ \end{array}} \right. \end{aligned}$$
(2)

and \(\dot{\varepsilon }\) is strain rate [\(s^{-1}\)], T stands for the absolute temperature [K], \(\sigma \) for flow stress [MPa] and Q for activation energy [kJ/mol\(^{\mathrm {-1}}\)]. The symbols \(A, \beta , n_{1}, n,\alpha =\beta /n_{1}\) denote material constants. Effects of temperature and strain rate on deformation can be described by Zener–Hollomon parameter as follows:

$$\begin{aligned} Z=\dot{\varepsilon }\exp \left( \frac{Q}{RT} \right) . \end{aligned}$$
(3)

Substituting Eq. (1) into (3) and using the definition of the hyperbolic law, flow stress can be also presented as:

$$\begin{aligned} \sigma =\frac{1}{\alpha }\ln \left\{ \left( \frac{Z}{A} \right) ^{1/n}+\left[ \left( \frac{Z}{A} \right) ^{1/n}+1 \right] ^{1/2} \right\} . \end{aligned}$$
(4)

2.2.2 Anand equation

Anand model was proposed by Anand and Brown [25] to describe the large, isotropic elastic–viscoplastic deformations of metals in high temperatures. Over time, the proposed constitutive equation was successfully used for various materials, mostly solders [25, 31, 32]. Analyses of the multiple research on the subject indicate that this model might be used to predict the material flow during thermo-mechanical compression of the titanium, especially in higher temperatures, where saturation stress, which is vital parameter in this model, is easy to determine. The Anand constitutive equations can reflect the physical phenomena of strain rate and temperature sensitivity, taking into account also strain hardening and dynamic recovery, which are critical in titanium plastic deformation.

The parameters of the Anand constitutive model are typically determined from uniaxial stress–strain tests at several strain rates and temperatures using a standard multistep model parameter determination procedure. Strengthening mechanisms such as dislocation density, solid solution strengthening, subgrain and grain size effects are characterized by the variable s, which is related to the equivalent stress    \(\sigma \):

$$\begin{aligned} \sigma =cs;c<1, \end{aligned}$$
(5)

where c is a function of strain rate and temperature which can be expressed as:

$$\begin{aligned} c=\frac{1}{\xi }\sinh ^{-1}\left\{ \left[ \frac{\dot{\varepsilon _{p}}}{A}\exp \left( \frac{Q}{RT} \right) \right] ^{m} \right\} . \end{aligned}$$
(6)

In Eq. (6), \(\xi \) is the multiplier of stress,\( \dot{\varepsilon _{p}}\) is the inelastic strain rate, A is the material constant (physical interpretation of the constant is given below), Q is the activation energy, R is the universal gas constant, T is the absolute temperature and m is the strain rate sensitivity.

The following Anand flow equation defines strain rate as

$$\begin{aligned} \dot{\varepsilon _{p}}=A\exp \left( -\frac{Q}{RT} \right) \left[ \sinh \left( \xi \frac{\sigma }{s} \right) \right] ^{\frac{1}{m}} . \end{aligned}$$
(7)

The rate of variable s changes is described in the following evolution equation:

$$\begin{aligned} \dot{s}=h(\sigma , s, T)\dot{\varepsilon _{p}}. \end{aligned}$$
(8)

The function h is associated with dynamic recovery and strain hardening. The specific form of this equation is given:

$$\begin{aligned} \dot{s}=\left[ h_{0}\left( A-\frac{s}{s^{*}} \right) ^{a}\mathrm {sign}\left( 1-\frac{s}{s^{*}} \right) \right] \dot{\varepsilon _{p}}; a>1, \end{aligned}$$
(9)

where

$$\begin{aligned} s^{*}=\hat{s}\left[ \frac{\dot{\varepsilon _{p}}}{A}\exp \left( \frac{Q}{RT} \right) \right] ^{n}. \end{aligned}$$
(10)

The variable \(s^{*}\) represents a saturation value of s for the given strain rate and temperature. The \(h_{0}\) is a hardening/softening constant, A is the strain rate sensitivity of hardening/softening, n is the strain rate sensitivity for the saturation value, and \(\hat{s}\) is coefficient of deformation resistance.

By analogy with Eq. (5), we formulate the relation between \(\sigma ^{*}\) and \(s^{*}\) in the form \(\sigma ^{*}=cs^{*}\). Taking into account Eqs. (6) and (10), we obtain the following formula:

$$\begin{aligned} \sigma ^{*}=\frac{\hat{s}}{\xi }\left[ \frac{\dot{\varepsilon _{p}}}{A}\exp \frac{Q}{RT} \right] ^{n}\sinh ^{-1}\left\{ \left[ \frac{\dot{\varepsilon _{p}}}{A}\exp \frac{Q}{RT} \right] ^{m} \right\} . \end{aligned}$$
(11)

Finally, Eqs. (5) and (9) lead to:

$$\begin{aligned} \frac{d\sigma }{d\varepsilon _{p}}=ch_{0}\left| 1-\frac{\sigma }{\sigma ^{*}} \right| \mathrm {sign}\left( 1-\frac{\sigma }{\sigma ^{*}} \right) ; a>1 \end{aligned}$$
(12)

the integration form of Eq. (12):

$$\begin{aligned} \sigma =\sigma ^{*}-\left[ \left( \sigma ^{*}-cs_{0} \right) ^{1-a}+\left( a+1 \right) \left\{ \left( ch_{0} \right) \left( \sigma ^{*} \right) ^{-a} \right\} \varepsilon _{p} \right] ^{1/(1-a)}, \end{aligned}$$
(13)

where \(s_{0}\) is the initial value of s.

For the Anand model, 9 parameters must be defined, i.e. \(A, Q, \xi , m \hat{s}, s_{0}, n, h_{0}, a\). The following procedure is applied to calibrate the constants:

  1. (1)

    firstly, the saturation stress must be obtained from the stress–strain graphs obtained from the experiments under constant temperatures and strain rates;

  2. (2)

    then the variables \(A, Q, \xi , m \hat{s}, n\) from Eq. (11) must be calculated by means of nonlinear least square method using the data acquired in step 1. The fraction \(\hat{s}/\xi \) is treated as one parameter. The lacking values of \(\dot{\varepsilon _{p}}\) were acquired from the stress–strain rate graph (Fig. 1);

  3. (3)

    variables \(\hat{s} \mathrm {and} \xi \) can be attained from the step 2. The \(\xi \) parameter is selected in such a way that c variable in Eq. (6) is less than unity. Then \(\hat{s}\) is determined from the \(\hat{s}/\xi \) form;

  4. (4)

    finally, constants \({ch}_{0}, {cs}_{0}, a\) from Eq. (13) are determined by the nonlinear square fit from the constant strain rate data. With c previously calculated \(h_{0} \mathrm {and} s_{0}\) parameters are acquired.

Fig. 1
figure 1

Stress–strain rate relationship. Additional strain rate values (beside experimental ones) can be determined from this graph

Full size image

2.2.3 Johnson–Cook equation

Original Johnson–Cook model [33] was proposed to handle large deformations of metals under high strains, strain rates and temperatures. The original form of the equation is as follows:

$$\begin{aligned} \sigma =\left( A+B\varepsilon ^{n} \right) \left( 1+C\ln \dot{\varepsilon }^{*} \right) \left( 1-T^{*m} \right) , \end{aligned}$$
(14)

where \(\sigma \) is the flow stress, \(\varepsilon \) is plastic strain, A is the yield stress at reference temperature and reference strain rate, B is the coefficient of strain hardening, n is the strain hardening exponent, C is coefficient of strain hardening and m is thermal softening exponent. The \(\dot{\varepsilon }^{*}=\dot{\varepsilon /}\dot{\varepsilon _{0}}\) is the dimensionless strain rate, where \(\dot{\varepsilon _{0}}\) is the arbitrarily chosen reference strain rate. \(T^{*}{=}\left( T-T_{r} \right) \!{/}(T_{m}-T_{r})\), where T stands for current temperature, \(T_{m}\) is meeting temperature and \(T_{r}\) is the arbitrarily chosen reference temperature.

Over time multiple modifications of the original Johnson–Cook model were made in order to tailor the equation for specific needs [22, 33, 34]. Modified JC model shown in Eq. (15) proposed by Lin et al. [35] is characterized by the fact that it takes into consideration coupled effects of temperature and strain rate:

$$\begin{aligned} \sigma =\left( A_{1}+B_{1}\varepsilon +B_{2}\varepsilon ^{2} \right) \left( 1+C_{1}\ln \dot{\varepsilon }^{*} \right) \exp \left[ \left( \lambda _{1}+\lambda _{2}\ln \dot{\varepsilon }^{*} \right) \left( T-T_{r} \right) \right] , \end{aligned}$$
(15)

where \(A_{1},B_{1},B_{2}, C_{1},\lambda _{1},\lambda _{2}\) are material constants and the rest parameters have the same meaning as in Eq. (14).

JC model as opposed to Arrhenius and Anand models has direct relationship between stress and strain, and that allows us to use nonlinear least square method to find desired constants.

2.2.4 Khan–Huang–Liang equation

Khan and co-workers developed a model with main goal of predicting the work hardening behaviour at large strain rates [35]. It also takes into consideration coupled effects of temperature and strain rate. The equation is as follows:

$$\begin{aligned} \sigma =\left[ A+B\left( 1-\frac{\ln \dot{\varepsilon }}{\ln D_{0}^{p}} \right) ^{n_{1}}\left( \varepsilon ^{p} \right) ^{n_{0}} \right] \left( \frac{\dot{\varepsilon }}{\dot{\varepsilon }^{*}} \right) ^{C}\left( \frac{T_{m}-T}{T_{m}-T_{r}} \right) ^{m}, \end{aligned}$$
(16)

where \(\sigma \) is the stress, \(\varepsilon ^{p}\) is plastic strain, \(T_{m}, T \)and \(T_{r}\) are the same as in JC equation, \(A,B, n_{1},n_{0},C,m\) are material constants and \(D_{0}^{p}={10}^{6} s^{-1}\) is arbitrarily chosen upper bound strain rate. Like in the case of JC model, direct relationship between stress and strain allowed us to use nonlinear least square fit to determine the model constants.

3 Results

3.1 Experiment results and Arrhenius model calculations

The experimental results of the thermo-mechanical tests as well as the theoretical values of stress calculated by means of the Arrhenius equation are presented in Fig. 2. The required constants in the Arrhenius equation and the procedure of obtaining them was presented in the previous work [20]. Data obtained from the experiment allowed for calculation parameters necessary to attain theoretical stresses. For specific values of the strain rates, stresses were calculated and compared to experimental results (Fig. 2).

Fig. 2
figure 2

Comparison of experimental data and calculated values of stresses according to applied Arrhenius model for 775 K temperature (left) and 875 K temperature (right)

Full size image
Fig. 3
figure 3

Comparison of experimental data and calculated values of stresses according to applied Anand model for 775 K temperature (left) and 875 K temperature (right)

Full size image

We can observe that the Arrhenius model is not adequate to describe mathematically the processes occurring in the material during the tests. This might be caused by the fact that the value of activation energy Q was based on the available literature or the model is too simple to be utilized to describe correctly the titanium behaviour in the thermo-mechanical tests.

3.2 Anand model fitting

We obtained much better results of the curve fitting for the Anand model (Fig. 3). The fitting is relatively precise with the mean adjusted R-square equal 0.935. Confidence bounds for the coefficients are equal 0.95. The model describes much better processes in higher temperatures (875 K), where the saturation stress is easier to determine like in samples 3b and 4b (Fig. 3 right). Unfortunately, Anand model does not describe well the change in the material during passing from elastic to plastic phase. In samples 2a, 3a and 4a, saturation stress was not reached, therefore, for the calculations, the highest value of stress was chosen. However, model may be used in the FEM analysis for prediction of the final material shape after thermoplastic processing.

Two sets of the parameters were calculated (Table 3) with accordance to experiment temperature. Those values are in good agreement with literature [36,37,38].

Table 3 Anand model parameters
Full size table

3.3 Johnson–Cook model fitting

Nonlinear least square method was used to determine the parameters of modified JC model. Calculations were performed with use of Eq. (15). The chosen reference temperature and reference strain rate were, respectively, 760 K and 0.009. Titanium melting temperature is equal \(T_{m }=\) 1938.16 K.

However, this model is usable only for the inelastic part of the material deformation. Fit within that region is very precise with adjusted R-square on the level of 0.98. Confidence bounds for the fit are equal 0.95. The results of this model fitting are presented on the Fig. 4.

Fig. 4
figure 4

Comparison of experimental data and calculated values of stresses according to applied modified Johnson–Cook model for 775 K temperature (left) and 875 K temperature (right)

Full size image

Two sets of the parameters were calculated (Table 4) with accordance to experiment temperature.

Table 4 Mean values of modified Johnson–Cook model parameters
Full size table

3.4 Khan–Huang–Liang fitting

The procedure of KHL model fitting is similar to JC fitting. The chosen reference temperature is equal \(T_{r} \quad =\) 760 K, and arbitrary chosen upper bound strain rate is \(D_{0}^{p}={10}^{6} s^{-1}\).

Similarly to JC model, fitting procedure works best for the inelastic deformation. In this case, however, accuracy of fit is lower with mean adjusted R-square equal 0.91. Calculated material constants (Table 5) are in good agreement with literature [39]. Figure 5 presents KHL model fitting results.

Fig. 5
figure 5

Comparison of experimental data and calculated values of stresses according to applied Khan–Huang–Liang model for 775 K temperature (left) and 875 K temperature (right)

Full size image
Table 5 Mean values of Khan–Huang–Liang model parameters
Full size table

3.4.1 Material observations

The sample deformed at the temperature of 875 K, at \(\dot{\varepsilon }={10}^{1}\), was selected for a scanning electron microscope (SEM) and scanning transmission electron microscope (STEM) observations. Such a choice was dictated by the fact that with such parameters of plastic forming, the grain disintegration should be the highest [20].

The observation revealed changes in the material. Initial structure consisted of mainly the grains size in the range of 2–8 \(\upmu \)m, some large grains of the size in the range 10–20 \(\upmu {\hbox {m}}\), and comparable to the large grains, a small amount of subgrains (Fig. 6).

Fig. 6
figure 6

SEM image of the initial state of the structure: a BSE image, b mapping of Euler angles with grain boundaries [20]

Full size image

Plastic forming at the temperature of 875 K has caused the disintegration of grains to the subgrains of the size of approx. 0.5 \(\upmu \)m (Fig. 7).

Fig. 7
figure 7

SEM image. The micrometre grains surrounded by the subgrains are twinned initially and secondary. The state of the structure: a Backscatter electron (BSE) image, b mapping of Euler angles with grain boundaries [20]

Full size image
Fig. 8
figure 8

Histogram of the frequency of equivalent grains diameter of initial state with the determined size of 2–8 \(\upmu \)m and sample deformed at the temperature of 875 K with the maximum frequency size in the range below 2 \(\upmu \)m [20]

Full size image

The electron backscatter diffraction (EBSD) method allowed us to generate a more accurate image of the distribution of grain sizes (Fig. 8), indicating a major reduce in the grain size and therefore an increase in structure strength.

4 Conclusions and discussion

In the paper, the TMP applied to CP titanium was mathematically modelled. We adopted four models to describe the process. The research showed that the Anand constitutive equation is able to model the process in the whole range of deformation.

The TMP was applied to enhance the mechanical properties of pure titanium. As a result of the process, grain fragmentation to the desired size was achieved. We managed to decrease the grain size down to even 0.5 \(\upmu \)m. Refined microstructure of titanium greatly influences mechanical properties of the material and its biomedical performance. The authors of [40] proved that smaller grain size improves bioactivity of titanium. This was verified by measuring hydroxyapatite amount formation on the surface. Similar studies were described in [41] where the authors noted higher mechanical properties and better biological response of ultrafine grained titanium. The titanium grain refinement that we achieved significantly increased mechanical properties of the material. Similar results have been reported in [41,42,43,44] where the authors improved microhardness, tensile strength, yield strength, fracture performance and fatigue strength. The grain refinement processes adopted in the aforementioned papers were different than that utilized in our studies. The proposed method of plastic deformation of titanium proved to be effective and allowed us to achieve the desired grain fragmentation. This is confirmed by the SEM and EBSD analyses.

However, there is still needed to optimize the main parameters of the process, i.e. temperature, strain rate and sample deformation. The best way to accomplish the optimization of the parameters is to perform finite element analyses of TMP. Therefore, we need a proper constitutive model describing the process. In the paper, we proposed and verified several constitutive equations.

The Arrhenius, Johnson–Cook and Khan–Huang–Liang equations precisely described the process but only in the plastic deformation zone [39, 45,46,47,48,49]. The Arrhenius model, in spite of the fact that it is widely used [50, 51], turned out to be less precise than JC and Anand models in our studies.

Johnson–Cook and Khan–Huang–Liang models are similar to each other in terms of fitting the experimental data. However, JC equation gives overall better results, especially for titanium thermo-mechanical processing under higher strain rates. According to [47], there is one of the KHL models’ variation which might be interesting to test in the future. This one incorporates the material grain size as one of the parameters, thus allowing for grain fragmentation prediction.

The Anand model gave good prediction of the material behaviour during thermo-mechanical deformation, especially in later phases of the process, i.e. during material flow. One of the main advantages over previously used models is that it does not require more experimental data for correct predictions. The best fits were obtained for the tests done in higher temperatures (where saturation stress is easier to determine). Moreover, this model is simple to implement and use, which is advantageous for further optimization and finite element model simulations [52]. The model, in opposition to the other ones presented in this paper, allowed us to predict the transition between elastic and plastic deformation, however, not in a very precise way. The values of Anand model parameters are in good agreement with literature (taking into consideration differences between the studied materials) [36,37,38].

The further work should include the finite element modelling of the pure titanium processing and optimization of work parameters in order to find the best input parameters leading to desired material structure and properties. Expansion of current experimental data is also required, as it allows to build more precise models. The future studies in the area will focus on implementation of the constitutive model into FE analyses and optimization of the TMP parameters. Some other studies [52,53,54] show that different approaches in the mathematical description of the continuum deformation might be worth of considering.

In recent years, another kind of models mimicking materials of grain structure has been proposed. Those models can be adopted to simulate and predict the creep, creep recovery and rate-dependent monotonic response [55]. Another model [56] can be utilized to simulate reduction in grain size resulting in a stiffer stress–strain behaviour or a transition from brittle to ductile type of failure.

It is also worth verifying whether the models we adopted to describe the titanium processing are able to model asymmetric behaviour of the material. The studies [57] show that materials of granular structure can induce various behaviour under compressive and tensile loads. It might be interesting whether or not combination of tensile and compression tests during the thermoplastic processing of titanium influences grain size in the material.

Another aspect is the presence of impurities in titanium. As shown in [58], they can affect mechanical properties of the material. Modelling of such effects will require, however, to consider higher-order or higher-gradient models [59,60,61,62,63,64].