Introduction

The evaluation of the stress state in the proximity of the crack tip represents an important topic in structural engineering. In this regard, the Stress Intensity Factor (SIF) plays a significant role in the investigation of the behaviour of cracked brittle components and structures. As an example, the SIF together with the crack growth rate provides important information on the fatigue behaviour of the material.

SIF can be assessed by adopting different experimental techniques, such as DIC (Digital Image Correlation) [1, 2], photoelasticity [3], and TSA (Thermoelastic Stress Analysis) [4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19]. This latter presents some advantages that make it suitable for in-service applications on real components when cycling loading conditions are applied. In particular, TSA requires a suitable surface preparation, unlike the other experimental techniques (an opaque paint is sufficient to ensure high and uniform emissivity).

In recent years, many researchers have adopted the TSA for obtaining the SIF value and locating the crack tip. In his work [4], Stanley described a procedure to obtain the Paris law parameters. This procedure is based on Westergaard’s equations for elastic plane stress and elastic plane strain conditions. In particular, the SIF can be obtained using a simple procedure that relates the first stress invariant and the distance from the crack line, thus avoiding the need for the exact identification of the crack tip. The main advantage is that only a series of signal profiles (horizontal or vertical), scanning around the crack-tip stress field and away from the plastic zone, are required to obtain the SIF measurement. The main drawback of the approach regards the omission of the T-stress term in Westergaard’s equations, which can lead to an error in SIF evaluation [19, 20]. Finally, a difference of less than 5% was estimated with the SIF obtained by an independent numerical technique [9].

In Lesniak et al. [5], Williams’s stress solution was employed for measuring the mixed-mode stress intensity factor of isotropic materials. Errors of up to 20% were obtained for mixed-mode stress intensity measurements. To reduce these errors, Tomlinson et al. [6] presented an alternative methodology to determine the SIF for cracks under mixed-mode displacements. In this case, the stress field was described using a Fourier series according to Muskhelishvili’s complex potentials approach.

Dulieu-Barton et al. [10] proposed a new curve-fitting routine based on a genetic algorithm to obtain information about SIFs and the crack tip in a mixed-mode configuration. In this way, a more accurate evaluation of the crack tip was obtained. In the work of Diaz et al. [12], a fitting algorithm based on the Downhill Simplex Method was employed to improve the SIF determination in combination with a method to monitor fatigue crack growth. In particular, the influence of residual stresses and crack closure phenomena were also experimentally investigated on as-welded and stress-relieved steel specimens. Finally, Pitarresi et al. [19], adopted Williams’s approach to investigate the influence of the series terms number, the selection of the fitting area, and the crack tip location on the SIF evaluation. In particular, the effect of the T-stress term was considered and a quantitative comparison with the Stanley approach was also provided.

All the cited works are based on the classical theory of TSA, in which temperature changes are related to the stress changes through the thermoelastic constant for an isotropic material, under linear elastic and local adiabatic conditions [21, 22]. This classical approach is only valid if the changes in the mechanical properties of the material with temperature can be considered negligible [21, 22]. Experimental work [23,24,25,26] and a review of the thermoelastic theory presented by Wong [27,28,29] demonstrated a relationship between the thermoelastic signal and the mean stress. These important results pushed researchers towards two new topics: the possibility of evaluating the influence of residual stresses [30, 31] and the development of a novel formulation to account for the mean stress effect [32, 33].

In the works of Palumbo and Galietti et al. [32, 33], a new procedure to process thermoelastic data and assess the correct value of the measured uniaxial stress was proposed. In particular, in the work of Palumbo et al. [33], an error analysis of titanium was carried out to investigate the mean stress effect on thermoelastic data, and a procedure for residual stress evaluation was proposed. This work demonstrated that significant errors can occur during stress evaluation if the mean stress effect is neglected, especially with the presence of high-stress gradients.

More recently, Di Carolo et al. [30, 31] proposed a general TSA formulation to account for the influence of biaxial residual stress on aluminium and titanium alloys. Based on a statistical analysis, the minimum residual stress value that leads to significant and measurable variations in the TSA signal was estimated. Moreover, the error in stress amplitude evaluation was assessed in cases where residual stresses were neglected. This error depended on the modulus and angle of the principal residual stresses with respect to the applied stresses.

Jones et al. [34], in an appendix to their work, developed thermoelastic equations considering the second-order effects in the case of a mode I fracture. They found that the thermoelastic temperature variation also had a harmonic component at twice the loading frequency. Moreover, the thermoelastic signal presents a change in the order of singularity, from 0.5 to 1, due to the mean stress effect. However, they did not investigate the effects on the SIF evaluation.

The aim of the present work is to show a new TSA formulation that considers the variations of the material properties (i.e., Young’s modulus) with temperature in the proximity of a crack tip loaded in mode I on titanium. The proposed formulation accounts for the second-order effects (mean stress effect and component at twice the mechanical frequency) that become significant for some materials like titanium and aluminium, affecting the TSA results in fracture mechanics applications, such as SIF determination.

Starting from the revised TSA theory, a TSA equation has been developed by describing the stress state, in terms of principal stresses, at the crack tip using Westergaard’s and Williams’s solutions, including the T-stress term. Then, the new model has been compared with the classical TSA equation to investigate the main parameters affecting the SIF measurement. In addition, a new approach has been employed to quantitatively and qualitatively evaluate its experimental implications for titanium through the difference in the TSA signal and the SIF inferred from Westergaard’s and Williams’s models, with and without second-order effects considered.

Theory of Thermoelastic Stress Analysis

The thermoelastic effect [35] relates the temperature change to the change in the sum of principal stresses for an isotropic material in linear elastic and adiabatic conditions. In particular, temperatures and stresses are related by the thermoelastic constant, which is generally assumed to be constant and independent of the applied stress.

The pioneering works of Belgen [36], Machin et al. [23], and Dunn et al. [24] showed, for some materials, a linear dependence of the thermoelastic signal on the mean stress value. The physical interpretation of this effect and a review of the thermoelastic theory was presented by Wong et al. [27,28,29]. Following this approach, a thermoelastic equation was derived, imposing the conservative laws that govern the mechanics of small quasi-static deformations [21] and the second law of thermodynamics for a reversible process. For an isotropic material, without any internal heat source besides the thermoelastic heat source, the temperature variations can be described by the following equation:

$${\rho }_{0}{C}_{\epsilon }\dot{T}=T\frac{\partial {\sigma }_{ij}}{\partial T}{\dot{\varepsilon }}_{ij}-{Q}_{i,i},$$
(1)

where ρ0 is the density of the material, Cε is the specific heat under constant strain, T is the temperature, σij is the stress tensor, εij is the strain tensor (summing over i, j with i, j = 1–3), and Qi,i is the heat flux through a unit surface area whose outward direct normal is ni. The dotted symbols represent derivatives with respect to time [27, 28]. The constitutive law is:

$${\sigma }_{ij}=2\mu {\varepsilon }_{ij}+\left(\lambda {\varepsilon }_{kk}-\beta\Delta T\right){\delta }_{ij},$$
(2)

where

$$\beta =\left(3\lambda +2\mu \right)\alpha ,\Delta T=T-{T}_{0}$$
$$\mu =\frac{E}{2(1+\nu )},\lambda =\frac{\nu E}{(1+\nu )(1-2\nu )},$$
(3)

where α is the coefficient of linear thermal expansion, λ and μ are the Lamé constants, T0 is the reference temperature, εkk is the first strain invariant, and δij is the Kronecker delta.

Deriving the constitutive law with respect to the temperature and substituting it into equation (1) yields [27, 28]:

$$T\left[\left(-\beta -\frac{\partial \beta }{\partial T}\Delta T+\frac{\partial \lambda }{\partial T}{\varepsilon }_{kk}\right){\delta }_{ij}+2\frac{\partial \mu }{\partial T}{\varepsilon }_{ij}\right]\dot{{\varepsilon }_{ij}}-{\rho }_{0}{C}_{\epsilon }\dot{T}={Q}_{i,i}.$$
(4)

Under adiabatic conditions and expressing equation (4) in terms of principal strain εi [27, 28]:

$${\rho }_{0}{C}_{\epsilon }\frac{\dot{T}}{T}=-\left(\beta +\frac{\partial \beta }{\partial T}\Delta T-\frac{\partial \lambda }{\partial T}{\varepsilon }_{kk}\right)\dot{{\varepsilon }_{kk}}+2\frac{\partial \mu }{\partial T}{\varepsilon }_{i}\dot{{\varepsilon }_{i}}.$$
(5)

Substituting equation (3) into equation (5), neglecting the high-order terms and the term (∂β/∂T)ΔT, and expressing it in terms of the principal stresses [27, 28], we obtain the following equation:

$$\rho {C}_{\varepsilon }\frac{\dot{T}}{T}=-\left[\alpha +\left(\frac{\upsilon }{{E}^{2}}\frac{\partial E}{\partial T}-\frac{1}{E}\frac{\partial \upsilon }{\partial T}\right)s\right]\dot{s}+\left[\frac{\left(1+\upsilon \right)}{{E}^{2}}\frac{\partial E}{\partial T}-\frac{1 }{E}\frac{\partial \upsilon }{\partial T}\right]{\sigma }_{i}{\dot{\sigma }}_{i},$$
(6)

where s represents the first stress invariant.

It is interesting to note that, for a one-dimensional stress state where the stress changes over time with a sinusoidal waveform, by integrating equation (6), we obtain:

$$\rho {C}_{\varepsilon }\frac{\Delta T}{{T}_{0}}=-\left(\alpha -\frac{1}{{E}^{2}}\frac{\partial E}{\partial T}{\sigma }_{m}\right){\sigma }_{a}\mathrm{sin} \; \omega t-\frac{1}{4{E}^{2}}\frac{\partial E}{\partial T}{\left({\sigma }_{a}\right)}^{2}\mathrm{cos} \; 2\omega t,$$
(7)

where σm and σa represent the mean and amplitude stresses, respectively. Equation (7) shows that the mean stress effect is strictly related to the Young’s modulus variation with the temperature. Moreover, it is interesting to observe that the thermoelastic signal also occurs at twice the loading frequency.

In the next section, the proposed approach will be presented by starting from equation (6) and considering a cracked material.

The Proposed Approach: A New Formulation for Describing the Thermoelastic Effect in Proximity of a Crack Tip

In this section, a new equation is obtained to describe the thermoelastic behaviour of materials in the presence of a crack.

In a similar way to the works of Patterson et al. [37], Palumbo et al. [33], and Di Carolo et al. [30], two material constants can be defined as reported in equation (8):

$$a=\frac{\alpha }{\rho {C}_{\varepsilon }}, b=\frac{1}{{E}^{2}}\frac{\partial E}{\partial T}\frac{1}{\rho {C}_{\varepsilon }}.$$
(8)

By considering the plane stress conditions, substituting equation (8) into equation (6) and neglecting the variations of the Poisson’s ratio with temperature [38, 39], we obtain:

$$\frac{\dot{T}}{T}=-\left[a+\upsilon b\left({\sigma }_{1}+{\sigma }_{2}\right)\right]\dot{\left({\sigma }_{1}+{\sigma }_{2}\right)}+\left(b+\upsilon b\right)\left({\sigma }_{1}\dot{{\sigma }_{1}}+{\sigma }_{2}\dot{{\sigma }_{2}}\right).$$
(9)

Under a generic sinusoidal loading in which the load changes between its maximum and minimum values, Pmin and Pmax, the loading ratio R can be defined as follows:

$$R=\frac{{P}_{min}}{{P}_{max}}=\frac{{\sigma }_{i,min}}{{\sigma }_{i,max}}, {\sigma }_{mi}={\sigma }_{ai}\frac{1+R}{1-R}={\sigma }_{ai}{R}_{f}, \qquad \text{where } i=\mathrm{1,2}$$
(10)

Thus, under a generic sinusoidal loading, the principal stresses and their rate of change are:

$${\sigma }_{1}={\sigma }_{a1}({R}_{f}+\mathrm{sin} \; \omega t), {\sigma }_{2}={\sigma }_{a2}({R}_{f}+\mathrm{sin} \; \omega t)$$
$$\dot{{\sigma }_{1}}={\sigma }_{a1}\omega \; \mathrm{cos} \; \omega t, \dot{{\sigma }_{2}}={\sigma }_{a2}\omega \; \mathrm{cos} \; \omega t.$$
(11)

Substituting equations (10) and (11) into equation (9) gives:

$$\begin{aligned}\frac{\dot{T}}{T}=&-\left[a+\upsilon b\left({\sigma }_{a1}+{\sigma }_{a2}\right)({R}_{f}+\mathrm{sin} \; \omega t)\right]\left({\sigma }_{a1}+{\sigma }_{a2}\right)\omega \; \mathrm{cos} \; \omega t \\ &+\left(b+\upsilon b\right)\left({{\sigma }_{a1}}^{2}+{{\sigma }_{a2}}^{2}\right)({R}_{f}+\mathrm{sin} \; \omega t)\omega \; \mathrm{cos} \; \omega t. \end{aligned}$$
(12)

After a few simple mathematical steps, equation (12) can be written as follows:

$$\begin{aligned}\frac{\dot{T}}{T}&=\left[-a\left({\sigma }_{a1}+{\sigma }_{a2}\right)+b{R}_{f}\left(-2\nu {\sigma }_{a1}{\sigma }_{a2}+{{\sigma }_{a1}}^{2}+{{\sigma }_{a2}}^{2}\right)\right]\omega \; \mathrm{cos} \; \omega t \\ &+\left[b\left(1+2\nu \right)\left({{\sigma }_{a1}}^{2}+{{\sigma }_{a2}}^{2}\right)++2\nu b{\sigma }_{a1}{\sigma }_{a2}\right]\omega \; \mathrm{sin} \; \omega t \; \mathrm{cos} \; \omega t. \end{aligned}$$
(13)

Integrating equation (13), considering ΔT < < T0, and neglecting the static components, the general expression of the thermoelastic signal in proximity of the crack is:

$$\frac{\Delta T (t)}{{T}_{0}}={g}_{1}\mathrm{sin} \; \omega t+{g}_{2} \; \mathrm{cos} \; 2\omega t.$$
(14)

where:

$${g}_{1}=\left[-a\left({\sigma }_{a1}+{\sigma }_{a2}\right)+b{R}_{f}\left(-2\nu {\sigma }_{a1}{\sigma }_{a2}+{{\sigma }_{a1}}^{2}+{{\sigma }_{a2}}^{2}\right)\right]$$
(15)
$${g}_{2}=-\frac{1}{2}\left[b\left(1+2\nu \right)\left({{\sigma }_{a1}}^{2}+{{\sigma }_{a2}}^{2}\right)+2\nu b{\sigma }_{a1}{\sigma }_{a2}\right].$$
(16)

It is important to highlight that the mean stresses in equation (9) could include the presence of residual stresses. As shown in [30], residual stresses affect the thermoelastic signal and can significantly affect the SIF evaluation. In this work, residual stress-free samples were employed, so residual stress was not included in the model.

Another important issue is represented by the plastic zone at the crack tip. It is well known that the stress values are limited by the yield stress of the material and then by the plastic behaviour, which generates a stress redistribution around the crack tip. However, it is important to highlight that the aim of this work is to investigate the effect of the mean stress and how second-order terms change the thermoelastic equation in the proximity of the crack in the new formulation. The effect of the plastic zone in TSA application has been extensively treated in the literature [14, 16, 18, 19, 40] and, in this regard, methods and procedures based on the classical TSA equation (and its validity hypothesis), used for describing the stress state at the crack tip in presence of the plastic area, can be extended to the proposed formulation. In the present research, it is not the aim of the authors to make any speculation on the plastic zone and its extension.

In the presence of a crack in a flat plate under plane loading, the state of the stress is characterized by two SIF values, mode I and mode II, respectively, KI and KII. In this regard, several theoretical models have been developed and can be employed to describe the state of stress around the crack. In this work, two models were considered: Westergaard’s equations [41] and Williams’s series expansion [42], truncated at the third term.

Derivation of TSA Temperature Variation Using Westergaard’s Solution

Westergaard’s equations [41] can be used to describe the state of stress around the crack in polar coordinates, as is shown in Fig. 1 and equation (17).

Fig. 1
figure 1

Polar coordinates used to describe the stress state around the crack tip

Full size image

It is important to highlight that equation (17) includes the T-stress [19] for the stress σx. Its omission can lead to an error in SIF evaluation that can be difficult to quantify since the T-stress may vary significantly with varying structure geometries and loadings [15, 19].

$$\begin{aligned}\left(\begin{array}{c}{\sigma }_{x}\\ {\sigma }_{y}\\ {\tau }_{xy}\end{array}\right)=&\frac{{K}_{I}}{\sqrt{2\pi r}}\mathrm{cos}\frac{\theta }{2}\left[\begin{array}{c}1-\mathrm{sin}\frac{\theta }{2}\mathrm{sin}\frac{3\theta }{2}\\ 1+\mathrm{sin}\frac{\theta }{2}\mathrm{sin}\frac{3\theta }{2}\\ \mathrm{sin}\frac{\theta }{2}\mathrm{cos}\frac{3\theta }{2}\end{array}\right]\\ &+\frac{{K}_{II}}{\sqrt{2\pi r}}\left[\begin{array}{c}-\mathrm{sin}\frac{\theta }{2}\left(2+\mathrm{cos}\frac{\theta }{2}\mathrm{cos}\frac{3\theta }{2}\right)\\ \mathrm{sin}\frac{\theta }{2}\mathrm{cos}\frac{\theta }{2}\mathrm{cos}\frac{3\theta }{2}\\ \mathrm{cos}\frac{\theta }{2}\left(1-\mathrm{sin}\frac{\theta }{2}\mathrm{sin}\frac{3\theta }{2}\right)\end{array}\right]+\left[\begin{array}{c}{\sigma }_{0x}\\ 0\\ 0\end{array}\right]\end{aligned}$$
(17)

The principal stresses σ1 and σ2, can be obtained by applying equation (18):

$${\sigma }_{1}, {\sigma }_{2}=\frac{{\sigma }_{x}+{\sigma }_{y}}{2}\pm \sqrt{{\left(\frac{{\sigma }_{x}-{\sigma }_{y}}{2}\right)}^{2}+{\tau }_{xy}^{2}}.$$
(18)

Considering only mode I of loading, and expressing the stress component in equations (15) and (16) by Westergaard’s equation (17), a few mathematical steps give:

$$\begin{aligned}{g}_{1}=&-a\left(\frac{2{K}_{Ia}}{\sqrt{2\pi r}}\mathrm{cos}\frac{\theta }{2}+{\sigma }_{0x}\right)+b{R}_{f}\left\{\frac{{K}_{Ia}^{2}}{4\pi r}\left[2\left(1-\upsilon \right)\right.\right. \\ &\left.\left.\left(\mathrm{cos}\theta +1\right) +(1+\upsilon ){\mathrm{sin}}^{2}\theta \right] +\frac{{\sigma }_{0x}{K}_{Ia}}{2\sqrt{2\pi r}}\left[\left(3-5\upsilon \right)\mathrm{cos}\frac{\theta }{2}\right.\right. \\ &\left.\left.+(1+\upsilon )\mathrm{cos}\frac{5\theta }{2}\right]+{{\sigma }_{0x}}^{2}\right\}\end{aligned}$$
(19)

and

$$\begin{aligned}{g}_{2}=&-\frac{1}{2}b\left\{\frac{{{K}_{Ia}}^{2}}{2\pi r}\left[\left(1+3\nu \right)\mathrm{cos}\theta -\frac{\left(1+\nu \right)}{4}\mathrm{cos}\left(2\theta \right)+\frac{1}{4}\left(5+13\nu \right)\right] \right. \\ &\left. +\frac{{K}_{Ia}\sqrt{2}}{\sqrt{\pi }}\frac{{T}_{s}}{4\sqrt{r}}\left(\left(\left(3+11\nu \right)\mathrm{cos}\left(\frac{\theta }{2}\right)+\left(1+\nu \right)\mathrm{cos}\left(\frac{5}{2}\cdot \theta \right)\right)\right)\right. \\ &\left.+{{\sigma }_{0x}}^{2}\left(2\nu +1\right)\right\}.\end{aligned}$$
(20)

Equations (19) and (20) show the thermoelastic signal and terms in which the order of singularity is 1, induced by the presence of the mean stress, as already noted by Jones et al. [34].

By setting b to 0, with equations (14) and (19) we can obtain the classical solution used for relating the thermoelastic signal and KIa:

$${\Delta T}_{nc}\left(t\right)=-{T}_{0}\left(\frac{2a{K}_{Ia}}{\sqrt{2\pi r}}\mathrm{cos}\frac{\theta }{2}+{\sigma }_{0x}\right)\mathrm{sin} \; \omega t,$$
(21)

where ΔTnc stands for the noncorrected value of ΔT.

Comparing equations (14) and (19) with equation (21), it is worth noting that the thermoelastic temperature variation ΔT also depends on the stress ratio R and on material constants b and υ. This means that an error in ΔT evaluation and then in SIF evaluation can be made by using equation (21) instead of equation (19).

Derivation of TSA Temperature Variation by Using Williams’s Series Expansion

The Williams’s series expansion, truncated to the third term, was combined with the TSA second-order equation to obtain a new formulation.

Considering a polar coordinate system with its centre at the crack tip, the first three terms of the Williams’s series expansion describing the elastic stress field surrounding the crack for mode I are:

$$\begin{aligned}{\sigma }_{x}=&\frac{{K}_{I}}{\sqrt{2\pi r}}\mathrm{cos}\frac{\theta }{2}\left(1-\mathrm{sin}\frac{\theta }{2}\mathrm{sin}\frac{3\theta }{2}\right)+{T}_{s}\\ &+{A}_{I3}\sqrt{r} \; \mathrm{cos}\frac{\theta }{2}\left(1+{\mathrm{sin}}^{2}\frac{\theta }{2}\right)+O\left(r\right).\end{aligned}$$
(22)
$$\begin{aligned}{\sigma }_{y}=&\frac{{K}_{I}}{\sqrt{2\pi r}}\mathrm{cos}\frac{\theta }{2}\left(1+\mathrm{sin}\frac{\theta }{2}\mathrm{sin}\frac{3\theta }{2}\right)\\ &+{A}_{I3}\sqrt{r} \; \mathrm{cos}\frac{\theta }{2}\left(1-{\mathrm{sin}}^{2}\frac{\theta }{2}\right)+O\left({r}^{3/2}\right)\end{aligned}$$
(23)
$${\tau }_{xy}=\frac{{K}_{I}}{\sqrt{2\pi r}}\mathrm{cos}\frac{\theta }{2}\mathrm{sin}\frac{\theta }{2}\mathrm{cos}\frac{3\theta }{2}-{A}_{I3}\sqrt{r} \; \mathrm{sin}\frac{\theta }{2}{\mathrm{cos}}^{2}\frac{\theta }{2}+O\left(r\right),$$
(24)

where Ts is the T-stress and AI3 is the third term coefficient.

Expressing the stress component in equations (15) and (16) by Williams’s series, truncated to the third term (equations (22)–(24)), after a few mathematical steps, gives:

$$\begin{aligned}{g}_{1}=a&\left[-\mathrm{cos}\left(\frac{\theta }{2}\right)\left(2{A}_{I3}\sqrt{r}+\frac{{K}_{Ia}\sqrt{2}}{\sqrt{\pi }}\right)-{T}_{s}\right] \\ &+bRf\Bigg\{\left(\frac{{{K}_{Ia}}^{2}}{2\pi r}+{{A}_{I3}}^{2}r\right)\left[\left(1-\nu \right)\mathrm{cos}\theta -\frac{\left(1+\nu \right)}{4}\mathrm{cos}\left(2\theta \right)+\frac{1}{4}\left(5-3\nu \right)\right] \\ &+\frac{{K}_{Ia}\sqrt{2}}{\sqrt{\pi }}\bigg[{A}_{I3}\left(\frac{1}{8}\left(\left(7-9\nu \right)\mathrm{cos}\left(\theta \right)+\left(1+\nu \right)\mathrm{cos}\left(3\theta \right)\right)+\left(1-\nu \right)\right) \\ &+\frac{{T}_{s}}{4\sqrt{r}}\left(\left(\left(3+5\nu \right)\mathrm{cos}\left(\frac{\theta }{2}\right)+\left(1+\nu \right)\mathrm{cos}\left(\frac{5}{2}\theta \right)\right)+\left(1-\nu \right)\right)\bigg] \\ &+\frac{{A}_{I3}{T}_{s}\sqrt{r}}{2}\left[\left(5-3\nu \right)\mathrm{cos}\left(\frac{\theta }{2}\right)-\left(1+\nu \right)\mathrm{cos}\left(\frac{3}{2}\theta \right)\right]+{T}_{s}^{2}\Bigg\}. \end{aligned}$$
(25)

and

$$\begin{aligned}{g}_{2}=-\frac{1}{2}b&\Bigg\{\left(\frac{{{K}_{Ia}}^{2}}{2\pi r}+{{A}_{I3}}^{2}r\right)\left[\left(1+3\nu \right)\mathrm{cos}\theta -\frac{\left(1+\nu \right)}{4}\mathrm{cos}\left(2\theta \right)+\frac{1}{4}\left(5+13\nu \right)\right] \\ & +\frac{{K}_{Ia}\sqrt{2}}{\sqrt{\pi }}\bigg[{A}_{I3}\left(\frac{1}{8}\left(\left(7++23\nu \right)\mathrm{cos}\left(\theta \right)+\left(1+\nu \right)\mathrm{cos}\left(3\theta \right)\right)+\left(1+3\nu \right)\right) \\ &+\frac{{T}_{s}}{4\sqrt{r}}\left(\left(\left(3+11\nu \right)\mathrm{cos}\left(\frac{\theta }{2}\right)+\left(1+\nu \right)\mathrm{cos}\left(\frac{5}{2} \cdot \theta \right)\right)\right)\bigg] \\ &+\frac{{A}_{I3}{T}_{s}\sqrt{r}}{2}\left[\left(5+13\nu \right)\mathrm{cos}\left(\frac{\theta }{2}\right)-\left(1+\nu \right)\mathrm{cos}\left(\frac{3}{2}\theta \right)\right]+{T}_{s}^{2}\left(2\nu +1\right)\Bigg\}. \end{aligned}$$
(26)

For b = 0, we can find the classical solution for the thermoelastic expression combined with Williams’s equations:

$$\frac{{\Delta T}_{nc,}(t)}{{T}_{0}}=a\left[-\mathrm{cos}\left(\frac{\theta }{2}\right)\left(2{A}_{I3}\sqrt{r}+\frac{{KI}_{a}\sqrt{2}}{\sqrt{\pi }}\right)-{T}_{s}\right]\sin\omega t,$$
(27)

where ΔTnc stands for the noncorrected value of ΔT.

Again, comparing equations (14) and (25) with equation (27), we find that the thermoelastic temperature variation at the angular velocity ω also depends on the stress ratio R and on material constants b and υ, while the thermoelastic variation at twice the angular velocity depends only on the material constants.

Methods: Experimental Implications for SIF Evaluation

The proposed formulations will be verified and compared with the classical theory for both Williams’s and Westergaard’s approaches. Moreover, experimental tests are carried out titanium grade 2 at two different stress ratios and the SIF value is assessed via the overly deterministic method [15, 19, 43]. The results will be compared with the synthetic data to show the statistical significance of including the second-order effects in the TSA formulation.

TSA Analysis with Synthetic Data

Equation (14), combined with equations (19), (20), (25), and (26), can be employed to obtain the synthetic thermoelastic data derived from the two considered models, if the mechanical and thermophysical constants of the material are known.

In the proposed equations, the thermoelastic temperature variation now consists of two harmonic components at ω and 2ω. Generally, experimental thermographic data are processed via hardware or software [4, 18, 21] to extract, separately, the amplitude and phase images related to the first and second harmonic of the thermoelastic signal. In this regard, the temperature variation, ΔTc, obtained by the new formulation can be represented as follows:

$${\Delta T}_{c}=\left({g}_{1}{T}_{0}\right)\mathrm{sin} \; \omega t+{(g}_{2}{T}_{0})\mathrm{cos} \; 2\omega t={\Delta T}_{c1} \; \mathrm{sin} \; \omega t+{\Delta T}_{c2} \; \mathrm{cos} \; 2\omega t.$$
(28)

Equations (19), (20), (25), and (26) can be rewritten as a function of KImax and the stress ratio R. In this case, we can write:

$${K}_{Ia}={K}_{Imax}\left(\frac{1-R}{2}\right)=\frac{{\Delta K}_{I}}{2}$$
(29)

and then, by using Westergaard’s equations, the components are:

$$\begin{aligned}{\Delta T}_{c1,west}(r,\theta )={T}_{0} &\Bigg\{-a\left(\frac{{K}_{Imax}}{\sqrt{2\pi r}}\left(1-R\right)\mathrm{cos}\frac{\theta }{2}+{\sigma }_{0x}\right)\\&+b{R}_{f} \Bigg[\frac{{K}_{Imax}^{2}}{16\pi r}{\left(1-R\right)}^{2}\left(2\left(1-\upsilon \right)\left(\mathrm{cos}\theta +1\right)\right.\\&\left.+(1+\upsilon ){\mathrm{sin}}^{2}\theta \right) +\frac{{\sigma }_{0x}{K}_{Imax}}{4\sqrt{2\pi r}}\left(1-R\right)\\ &\left(\left(3-5\upsilon \right)\cos\frac{\theta }{2}+(1+\upsilon )\cos\frac{5\theta }{2}\right)+{{\sigma }_{0x}}^{2}\Bigg]\Bigg\}\end{aligned}$$
(30)
$$\begin{aligned}{\Delta T}_{c2.west}\left(r,\theta \right)=&-\frac{1}{2}{T}_{0}b\Bigg\{\frac{{{{K}_{Imax}}^{2}\left(1-R\right)}^{2}}{8\pi r}\left[\left(1+3\nu \right)\cos\theta \right.\\&\left.-\frac{\left(1+\nu \right)}{4}\cos\left(2\theta \right)+\frac{1}{4}\left(5+13\nu \right)\right] \\ &+\frac{{K}_{Imax}\left(1-R\right)}{\sqrt{2}\sqrt{\pi }}\frac{{T}_{S}}{4\sqrt{r}}\left(\left(\left(3+11\nu \right)\cos\left(\frac{\theta }{2}\right)\right.\right.\\&\left.\left.+\left(1+\nu \right)\cos\left(\frac{5}{2}\theta \right)\right)\right)+{{T}_{S}}^{2}\left(2\nu +1\right)\Bigg\}.\end{aligned}$$
(31)

In the same way, the thermoelastic components for Williams’s series are:

$$\begin{aligned}{\Delta T}_{c1, will}(r,\theta )={T}_{0}&\Bigg\{a\left[-\mathrm{cos}\left(\frac{\theta }{2}\right)\left(2{A}_{I3}\sqrt{r}+\frac{{K}_{Imax}\sqrt{2}}{2\sqrt{\pi }}\left(1-R\right)\right)-{T}_{S}\right] \\ & +bRf \Bigg\{\left(\frac{{{K}_{Imax}}^{2}}{8\pi r}{\left(1-R\right)}^{2}+{{A}_{I3}}^{2}r\right)\left[\left(1-\nu \right)\mathrm{cos}\theta -\frac{\left(1+\nu \right)}{4}\mathrm{cos}\left(2\theta \right)+\frac{1}{4}\left(5-3\nu \right)\right] \\ &+\frac{{K}_{Imax}}{\sqrt{2}\sqrt{\pi }}\left(1-R\right)\Bigg[{A}_{I3}\left(\frac{1}{8}\left(\left(7-9\nu \right)\mathrm{cos}\left(\theta \right)+\left(1+\nu \right)\mathrm{cos}\left(3\theta \right)\right)+\left(1-\nu \right)\right) \\ &+\frac{{T}_{S}}{4\sqrt{r}}\left(\left(\left(3-5\nu \right)\mathrm{cos}\left(\frac{\theta }{2}\right)+\left(1+\nu \right)\mathrm{cos}\left(\frac{5}{2}\theta \right)\right)+\left(1-\nu \right)\right)\Bigg] \\ &+\frac{{A}_{I3}{T}_{S}\sqrt{r}}{2}\left[\left(5-3\nu \right)\mathrm{cos}\left(\frac{\theta }{2}\right)-\left(1+\nu \right)\mathrm{cos}\left(\frac{3}{2}\theta \right)\right]+{{T}_{S}}^{2}\Bigg\}\Bigg\} \end{aligned}$$
(32)
$$\begin{aligned}{\Delta T}_{c2, will}(r,\theta )=&{-\frac{1}{2}T}_{0}b\Bigg\{\left(\frac{{{{K}_{Imax}}^{2}\left(1-R\right)}^{2}}{8\pi r}+{{A}_{I3}}^{2}r\right)\left[\left(1+3\nu \right)\cos\theta -\frac{\left(1+\nu \right)}{4}\cos\left(2\theta \right)+\frac{1}{4}\left(5+13\nu \right)\right] \\ & +\frac{{K}_{Imax}\left(1-R\right)}{\sqrt{2}\sqrt{\pi }}\bigg[{A}_{I3}\left(\frac{1}{8}\left(\left(7+23\nu \right)\cos\left(\theta \right)+\left(1+\nu \right)\cos\left(3\theta \right)\right)+\left(1+3\nu \right)\right) \\ &+\frac{{T}_{S}}{4\sqrt{r}}\left(\left(\left(3+11\nu \right)\cos\left(\frac{\theta }{2}\right)+\left(1+\nu \right)\cos\left(\frac{5}{2}\theta \right)\right)\right)\bigg] \\ &+\frac{{A}_{I3}{T}_{S}\sqrt{r}}{2}\left[\left(5+13\nu \right)\cos\left(\frac{\theta }{2}\right)-\left(1+\nu \right)\cos\left(\frac{3}{2}\theta \right)\right]+{{T}_{S}}^{2}\left(2\nu +1\right)\Bigg\}. \end{aligned}$$
(33)

Equations (30)–(33) can be used for obtaining the synthetic corrected and noncorrected thermoelastic temperature variations, which will be compared with the measured thermoelastic data, as will be explained in the following sections, in order to study the second-order effects for each model.

TSA Experimental Calibration

In this work, the thermoelastic parameters a and b were evaluated by applying the experimental approach proposed by Palumbo et al. [33].

The procedure involved the implementation of dynamic tests on samples made of the same material as the fracture mechanics test samples, with a known stress distribution. Thus, dog bone samples with a monoaxial and uniform stress distribution in the useful section were employed.

In this case, σa2 is null and the first harmonic amplitude of the thermoelastic signal became (equation (15)):

$$\frac{{\Delta T}_{1}}{{T}_{0}}=-a{\sigma }_{a1}+b{R}_{f}{{\sigma }_{a1}}^{2}.$$
(34)

Substituting this into equation (10) and dividing by σa1 gives:

$$\frac{{\Delta T}_{1}}{{T}_{0}{\sigma }_{a1}}=-a+b{\sigma }_{m1}.$$
(35)

The parameters a and b represent the intercept and slope, respectively, of the linear relation expressed by equation (35) and can be evaluated by fitting the experimental data.

The calibration involved the following steps:

  1. 1.

    Acquisition with an IR camera of thermoelastic signal during sinusoidal loading tests on dog bone samples at different levels of mean load.

  2. 2.

    Lock-in analysis and calculation of the first harmonic amplitude ΔT1 and reference temperature T0.

  3. 3.

    Evaluation of the thermoelastic parameters a and b by fitting the experimental data (equation (35)).

Material and Experimental Setup

A calibration procedure was performed, following the procedure suggested by Palumbo et al. [33] in order to find the thermoelastic parameters a and b.

In this regard, TSA dynamical uniaxial tensile test with different levels of mean load were performed on a dog-bone sample in pure titanium grade 2 material, with a useful section of 14.55 mm2. The dynamical tensile tests were performed using a loading frame MTS model 370 with a 25 kN capacity (Fig. 2).

Fig. 2
figure 2

Setup and equipment: loading frame and infrared camera

Full size image

Table 1 shows the test plan performed for calibration: the tests were carried out with a loading frequency of 15 Hz, applying three different levels of loading amplitudes Pa = 700, 1000, or 1200 N, and three levels of mean stress Pm = 1400, 1600, or 1800 N. Three sequences of 1000 frames each were recorded for each load condition.

Table 1 Calibration load conditions
Full size table

TSA was performed on pure titanium grade 2 compact tensile (CT) specimens 1 mm thick, whose geometry is reported in Fig. 3. In Table 2, the mechanical characteristics of the material are reported.

Fig. 3
figure 3

Sample geometry with dimensions in mm. The thickness is 1 mm

Full size image
Table 2 Pure titanium grade 2 mechanical characteristics
Full size table

All the samples were prepared by painting the surface with a black matte spray in order to enhance the surface emissivity and make it uniform (Fig. 3).

Table 3 shows the loading conditions applied to the titanium CT samples during the dynamical test. In particular, two loading ratios were investigated, R = 0.1 and R = 0.5, keeping the maximum load (750 N) and the loading frequency (17 Hz) constant.

Table 3 Test parameters used for CT specimens
Full size table

Thermal sequences of 300 frames each were acquired with different numbers of cycles following a crack with a crack growth of about 0.5 mm between one acquisition and the next. The thermal acquisition was performed with a FLIR X6581sc IR camera, with the acquisition parameters set as reported in Table 4.

Table 4 IR camera acquisition parameters
Full size table

Experimental TSA Data Analysis

In this work, experimental data were acquired to compare the ability of the two considered thermoelastic models (with and without second-order correction) to describe the elastic stress field around the crack tip, and thus to characterize the fracture behaviour of the material.

Experimental data were collected by recording with a cooled IR camera the cracked samples’ thermal response to a dynamic load. The thermal data were then processed using the IRTA [44] software, which performed the lock-in analysis to obtain the first harmonic amplitude ΔT1 and phase ϕ1 maps.

The data processing involved the application of the Least Squares Overdeterministic method [12] to obtain the value of the unknowns of the two systems, which are the ΔKI and T-stress, respectively, and, only for Williams’s series expansion, the third term coefficient AI3.

The crack tip coordinates are also optimized. The algorithm is based on a crack tip search in an area close to the initial guess, choosing the solution that minimizes the sum of the deviations.

Following the approach proposed by several authors [14, 16, 18], the phase data were employed to identify the crack tip initial value for the optimization and to estimate the plastic zone extension in order to select the minimum radius for the data point selection; the crack tip was predicted as being at the point of sign inversion of a horizontal profile passing for the minimum (Fig. 4). The plastic zone can be approximated to the area with negative phase, which can overestimate the size due to dissipative phenomena, depending on the load frequency. Thus, the minimum radius for the selected data points can be selected as the point where the phase returns to zero (Fig. 4).

Fig. 4
figure 4

Phase (ϕ1) map (a) and data along the profile passing through the minimum (b)

Full size image

The maximum radius was selected by moving 2 mm from the minimum radius, so as to remain close enough to the crack tip to obtain a signal affected by the singularity.

The data processing workflow is shown in Fig. 5 and can be summarized by the following steps:

  1. 1.

    Experimental data were acquired with a cooled IR camera during dynamical tests on CT titanium samples.

  2. 2.

    A lock-in analysis was performed: the acquired signal was decomposed into harmonics for each pixel, and the reference temperature T0, the first amplitude ΔT1, and its phase ϕ1 maps were obtained (Figs. 4 and 6).

  3. 3.

    The map of the normalized first amplitude ΔT1/T0 was evaluated.

  4. 4.

    The ϕ1 map was employed to identify the crack tip initial value and the elastic stress field boundaries.

  5. 5.

    A semiannular path of selected data points was created in the elastic zone, having as its centre the coordinates of the crack tip and radius determined as explained above. The annular path angle ranges from –4π/9 and 4π/9 with step π/36 and 10 points for each angle (Fig. 6).

  6. 6.

    The experimental calibration was performed on a titanium dog bone sample and the thermoelastic parameters a and b were evaluated.

  7. 7.

    For each considered thermoelastic model, an overdeterministic system of equation was written by substituting the ΔT1/T0 values of the selected data points, together with the calibration parameters a and b, the loading conditions, and the polar coordinates in equations (30) and (32), respectively. The unknowns of this system are the parameters that describe the two models (ΔKI and Ts in Westergaard’s equations and ΔKI, Ts, and AI3 in Williams’s series expansion).

  8. 8.

    The crack tip was located with a trial-and-error procedure, by assuming the initial guess as coordinates for the crack tip and searching the optimum in a close area. The overdeterministic system was solved for each point and the solution was obtained by minimizing deviations using a direct research method based on the Nelder–Mead method [45].

Fig. 5
figure 5

First amplitude ΔT1 (a), maps of reference temperature T0 (b), and normalized first amplitude ΔT1/T0 (c) with selected data points

Full size image
Fig. 6
figure 6

Workflow of the data processing

Full size image

Results and Discussion

In the present section, the assessment of the thermoelastic constants from thermoelastic data is presented. The effect of the correction of Westergaard’s equation (30) and Williams’s models (equation (32)) is shown in terms of thermoelastic maps, comparing the synthetic and experimental data, and in terms of evaluated SIFs (ΔKI_c,west and ΔKI_c,will) for the tests carried out at R = 0.1 and R = 0.5.

Assessing the Thermoelastic Constants

The parameters a and b were found by calculating the intercept and slope of the experimental curve ΔT1/T0 versus the mean load. By focusing on the semi-amplitude of the temperature running at the same loading frequency ΔT1 and for a uniaxial loading condition, one can describe the temperature variations through equation (35).

The results obtained by applying the linear regression of data are as reported in Fig. 7:

  • a = (3.59e-6 ± 1.47e-8) [1/MPa]

  • b = (–2.825e-9 ± 1.36e-10) [1/MPa2]

Fig. 7
figure 7

Thermoelastic data as a function of mean stress and linear fitting. The temperature amplitude is here reported with a positive sign

Full size image

Effect of Correction on Thermoelastic Data and SIFs Values

A first qualitative demonstration of the importance of considering the second-order effects of the thermoelastic signal in proximity to the crack tip in two loading conditions at R = 0.1 and R = 0.5 can be observed in the TSA maps in Fig. 8(a-d). The maps compare the contour lines of experimental thermoelastic temperature data ΔT1 divided by T0 and the synthetic ΔT1 /T0 data obtained by considering equations (30) and (21).

Fig. 8
figure 8

Thermoelastic temperature ΔT1 divided by T0: experimentally assessed (solid contours) and synthetic data with (dotted contour) and without (dashed contour) second-order correction for Westergaard’s model with (a) R = 0.1 or (b) R = 0.5, and Williams’s model with (c) R = 0.1 or (d) R = 0.5

Full size image

Westergaard’s data (the contour lines reported in Fig. 8(a), show that the synthetic corrected data are closer to the experimental data than noncorrected ones, and the difference (corrected vs. non-corrected) is marked. This consideration is also valid for the data of the test at R = 0.5 (Fig. 7(b)). It is worth noting that there is a difference in the shape of contour lines between experimental and synthetic data, especially far away from the plastic area (ΔT1 /T0 is less than 6 10–4 at R = 0.1 and 5 10–4 at R = 0.5). This effect depends on the model adopted for representing the stress distribution around the crack tip; as also demonstrated in [19], the use of more terms of the Williams’s series expansion leads to a better matching of the isopach lines farthest from the crack tip.

The synthetic maps obtained by using corrected Williams’s solution (equation (32)), shown in Fig. 8(c,d), are close to the experimental data at both R = 0.1 and R = 0.5 near the plastic area. The difference in the shape of the contour lines between corrected (the solid black line) and noncorrected data (the dotted red line) is more marked in the case of R = 0.5.

Far away from the plastic area, especially at R = 0.5, for ΔT1 /T0 less than 5 10–4, both the corrected model and noncorrected contours fail to correctly approximate the experimental data.

The effect of the correction on the SIF evaluation can be easily obtained by observing Table 4(a,b), in which ΔKI_c,west and ΔKI_c,will are compared to the corresponding value without correction, respectively: ΔKI_nc,west and ΔKI_nc,will. These values represent the best fitting values obtained from the minimization of the sum of the squares of the residuals of the overdetermined system.

In Table 4(a), the difference between ΔKI_c,west and ΔKI_nc,west appears qualitatively significant for both stress ratios. Clearly, as the ΔKI value depends on the stress amplitude, the SIF values are higher for R = 0.1 than for R = 0.5.

The statistical significance of the difference between ΔKI_c,west and ΔKI_nc,west was quantitatively assessed by calculating the percent difference Δ% (Table 4), according to the following formula:

$${\Delta }_{\%}=\frac{{\Delta K}_{I c}-{\Delta K}_{I nc}}{{\Delta K}_{I c}}*100.$$
(36)

The estimated difference is roughly 20% at R = 0.1, while it increases from 21 to 29% at R = 0.5. ΔKI_nc,west is always higher than ΔKI_c,west, at both R = 0.1 and R = 0.5.

Table 5(a) also reports the residuals (Resc_west, Resnc_west), a quantitative measurement of the goodness of fit of the adopted model (lower values mean a better fit). In either case, R = 0.1 or R = 0.5, the residuals of the corrected Westergaard’s model are lower than the related noncorrected values. This means that the correction is necessary as it improves the data analysis.

Table 5 Evaluation of corrected and noncorrected SIFs and related statistics parameters using Westergaard’s (a) and Williams’s solutions (b)
Full size table

As a further confirmation of the necessity to consider the second-order effects, the latter statistical parameter is reported in Table 5(a): the coefficient of determination, R2. The R2c,west values for the two stress ratios are higher than R2nc,west.

As for Westergaard’s model, there is a marked difference between the ΔKI_c,will and ΔKI_nc,will of both test conditions (R = 0.1 and R = 0.5). In general, second-order effects produce lower SIF values with respect to those obtained with the noncorrected Williams’s model.

In this case, the percent difference Δ% changes significantly between corrected and corrected SIFs. Between corrected and noncorrected data at R = 0.1, Δ % varies in the range 5–26%, while for the test at R = 0.5 it varies between 17 and 32%.

The Resc_will values are smaller than the corresponding Resnc_will, meaning that the corrected model is capable of better ΔKI estimation. As for the coefficients of determination, both R2c,will and R2nc,will are close to 1; in particular, R2c,will is always higher than R2nc,will.

To summarize, the quantitative analyses focused on investigating the effect of the second-order effects correction on Westergaard’s and Williams’s models with T-stress led us to draw the following conclusions:

  • ΔKI_c,will values are similar to ΔKI_c,west; specifically, the effect of the correction involves lower values of SIF;

  • The evaluated residuals and R2 coefficients demonstrate the improvement in the analysis when adopting a corrected model fitting using both Westergaard’s and Williams’s models.

In order to make more consistent considerations, in the following graphs of Fig. 9(a–d), the SIF values at different crack lengths are compared to the values obtained by using the formula proposed by the Standard ASTM [46], considered as a reference curve.

Fig. 9
figure 9

ΔKI was assessed by using Westergaard’s solution for the test at R = 0.1 (a) and R = 0.5 (b), and Williams’s solution for the test at R = 0.1 (c) and R = 0.5 (d) compared to ΔKI assessed by using the Standard ASTM [46]

Full size image

Figure 8 reports the ΔKI measured by ASTM (the solid black line) and the ΔKI, corrected and noncorrected, at R = 0.1 (Fig. 9(a)) and R = 0.5 (Fig. 9(b)). The small difference in crack length between the corrected (ΔKI_c,west, ΔKI_c,will) and noncorrected (ΔKI_nc,west, ΔKI_nc,will) data points (see Fig. 8) can be explained by considering that the crack length depends on the crack tip, found by implementing the fitting procedure described in the previous section, which, in turn, depends on the model adopted.

As previously shown in Table 5, for Westergaard’s data, the difference between corrected and noncorrected data (Fig. 9(a,b)) is almost constant; in effect, ΔKI_nc,west seems to have shifted upwards to the reference curve. The ΔKI_c,west data fit the ASTM reference curve very well.

Figure 9(c,d) reports the data of ΔKI_c,will and ΔKI_nc,will for the two tests at different stress ratios (R = 0.1 and R = 0.5), each compared to the ΔKI provided by the Standard. In this case, the proposed correction is in total agreement with the reference curve for the tests at both R = 0.1 and 0.5, while, as seen from the results in Table 5(b), the ΔKI_nc,will values are higher than the reference curve and very scattered.

Conclusions

In this work, a new formulation has been proposed for describing the TSA signal in proximity to the crack tip on materials characterized by a non-negligible second-order effect on the thermoelastic signal. The proposed approach starts from the revised theory of the TSA, where the effect of the mean stress on the thermoelastic signal is considered.

The thermoelastic equation has been rewritten for describing the stress distribution around a crack by using Westergaard’s and Williams’s solutions with T-stress. The main results can be summarized as follows:

  • A part of the thermoelastic equation occurs at twice the loading frequency. This component depends on the square of the SIF and the stress ratio for both Westergaard’s and Williams’s solutions.

  • The component of the thermoelastic signal that occurs at the load frequency depends on the material properties and the stress ratio.

The effects of the proposed formulation have been investigated experimentally by performing tests on small CT samples of pure titanium in two loading conditions: R = 0.1 and R = 0.5. The experimental results have been compared with the synthetic thermoelastic data. The major outcome from this analysis is a better ability to describe the stress distribution in proximity to the crack tip by the proposed formulation. The effect of the correction is pronounced in both Williams’s and Westergaard’s models.

Furthermore, the possible implications of SIF evaluation using the classical thermoelastic equation have been investigated. The SIFs obtained as the output of an algorithm based on an over-deterministic approach were compared to those found by following the ASTM Standard procedure.

The proposed equations can be a useful tool for understanding the limits of the applicability of the classical theory and well-known solutions and the error that arises from neglecting the second-order effects.