Introduction

Martensitic transformation is a first-order solid-state displacive (diffusionless) ferroelastic phase transition that breaks the crystal symmetry of high-temperature austenite phase by inducing spontaneous anisotropic lattice strain upon cooling and producing multiple low-symmetry orientation variants of low-temperature martensite phase [1,2,3,4]. Martensitic transformations are widely observed in a large group of materials including metal alloys and ceramics [5]. They have been extensively studied [1, 2] for their fundamental importance in solid-state materials science and physics. They also provide the physical basis of shape memory effect, making some martensitic alloys important functional materials [6]. Despite successful applications of martensitic alloys as structural and functional materials, our fundamental understanding of martensitic transformations is still incomplete, and some important issues remain to be clarified [7,8,9,10]. One intriguing issue is pre-martensitic phenomena, also called martensite precursor effects. Prior to martensitic transformation, the high-symmetry cubic austenite phase usually undergoes various anomalies in a wide temperature range 10-100 K above the martensite start temperature, which are unexpected in cubic phase [4, 5, 10, 11]. In particular, strong diffuse scattering is observed in diffraction, which depends on the temperature and exhibits different characteristics around different Bragg reflection peaks [11]. As an example, Fig. 1 shows the experimentally measured three-dimensional (3D) diffuse scattering in Ni49.90Mn28.75Ga21.35 single crystal above the martensite start temperature 323 K, around (800) Bragg reflection peak, respectively, at 490 K in (a) and 327 K in (b), while around different Bragg reflection peaks \(\left( {HKL} \right)\) at 330 K in (c). Details of the 3D diffuse scattering experiment using in-situ high-energy synchrotron X-ray single-crystal diffraction are described in "Methodology" section.

Fig. 1
figure 1

High-energy synchrotron X-ray single-crystal diffraction measurement of 3D diffuse scattering in Ni49.90Mn28.75Ga21.35 single crystal above the martensite start temperature TC = 323 K: around (800) Bragg reflection peak at a 490 K and b 327 K, and c around different Bragg reflection peaks (HKL) at 330 K. The x, y and z axes are aligned along [100], [010] and [001] crystal axes, respectively

The diffuse scattering exemplified in Fig. 1 manifests structural deviations from the cubic austenite crystal structure above the martensite start temperature and also, according to reciprocal relation, implies an atomic lengthscale of the structural deviations with characteristic wavelengths of a few lattice parameters. Furthermore, repeated experiments reproduce the same diffuse scattering results at corresponding temperatures irrespective of prior thermal and mechanical histories, e.g., forward and reverse martensitic transformations, stress-induced martensitic transformation, annealing at temperatures up to 500 K, etc. These observations indicate that the diffuse scattering and underlying structural deviations are intrinsic properties of the cubic austenite phase in a thermodynamic equilibrium pre-martensitic state. In this paper, we focus on the intrinsic behaviors of a defect-free crystal that undergoes cubic ↔ tetragonal martensitic transformation. In particular, we carry out material modeling and computer simulation to correlate the diffuse scattering to the short-range ordering of atomic-scale heterogeneous displacement field developed in the pre-martensitic austenite crystal lattice prior to the development of long-range order during martensitic transformation.

It is worth noting recent progress in strain glass that provides valuable insights into the defect-driven behavior of disordered state in martensitic system, which complements our findings on the intrinsic behavior of pre-martensitic state in defect-free crystal. While the concept of strain glass was first introduced in 2005 [12] and has been intensively explored since then [13,14,15,16,17] (see Ref. [18] for review of recent progress in strain glass), this important idea was discussed earlier in 1990s [19,20,21]. In particular, the disordered state was called frustrated glassy phase rather than strain glass; nevertheless, the compositional disorder intrinsic to alloys and the resultant elastic strain frustration were discussed in direct correspondence with spin-glass phase [19, 20]. Despite the frustrated glassy phase concept (equivalent to strain glass) discussed by Krumhansl and coworkers in 1990s to explain pre-martensitic phenomenon and martensitic transformation, in 2000 Krumhansl listed martensitic precursor phenomenon as one of the unsolved fundamental problems at the frontier of multiscale materials science in the twenty-first century [10] (the last publication of his life). In this paper, we employ atomic-scale Monte Carlo simulations to cover the required time and length scales for study of diffuse scattering and displacement short-range ordering at pre-martensitic temperature.

Atomic-scale Monte Carlo method covers sufficient length and time scales with affordable computing requirements and thus is ideal for simulation study of diffuse scattering (time-averaging over seconds or longer in experimental measurements) and displacement short-range ordering at varying temperatures in the pre-martensitic austenite state. Various other methods previously used for martensitic transformations but unsuitable for this task include phase field microelasticity model [22, 23] limited by its mesoscale nature, Density Functional Theory [24] limited by supercell size and zero Kelvin computation, and Molecular Dynamics [25] limited by timescale (nanoseconds or shorter) and demanding computing requirements. In this work, a quasi-spin Ising model of ferroelastic phase transition is employed to perform atomic-scale Monte Carlo simulation of martensitic transformation, which are compared with complementary synchrotron X-ray single-crystal diffraction experiments. The computational and experimental methodologies are described in "Methodology" section, and the simulation results are presented and discussed in "Results and Discussions" section.

Methodology

Quasi-Spin Ising Model and Monte Carlo Simulation

In analogy to the spin variable in Ising model of ferromagnetism, the proposed Ising model of ferroelastic phase transition and thermoelastic martensitic transformation employs a multivalued quasi-spin variable \(s\left( {\mathbf{r}} \right)\) defined at each lattice site \({\mathbf{r}}\) to describe the system state. To be specific while without loss of generality, a martensitic system of tetragonal ground state is considered in the following. In such an exemplary system, the quasi-spin variable \(s\left( {\mathbf{r}} \right)\) assumes one of three values at each lattice site \({\mathbf{r}}\), namely, \(s = 1\), 2 or 3 characterizing the unit cell of tetragonal orientation variant 1, 2 or 3 of the ground-state martensite phase, respectively, each of which is associated with a stress-free lattice misfit strain \(\varepsilon_{ij}^{0} \left( s \right)\):

$$\varepsilon_{ij}^{0} \left( 1 \right) = \left( {\begin{array}{*{20}c} {\varepsilon_{c} } & 0 & 0 \\ 0 & {\varepsilon_{a} } & 0 \\ 0 & 0 & {\varepsilon_{a} } \\ \end{array} } \right),\;\varepsilon_{ij}^{0} \left( 2 \right) = \left( {\begin{array}{*{20}c} {\varepsilon_{a} } & 0 & 0 \\ 0 & {\varepsilon_{c} } & 0 \\ 0 & 0 & {\varepsilon_{a} } \\ \end{array} } \right),\;\varepsilon_{ij}^{0} \left( 3 \right) = \left( {\begin{array}{*{20}c} {\varepsilon_{a} } & 0 & 0 \\ 0 & {\varepsilon_{a} } & 0 \\ 0 & 0 & {\varepsilon_{c} } \\ \end{array} } \right)$$
(1)

The tetragonality of the ground-state martensite phase is characterized by \(\gamma_{0} = \varepsilon_{c} - \varepsilon_{a}\). The three ground-state tetragonal martensite variants possess the same volume, thus the martensitic transformation in this model is volume conserved, which is a reasonable approximation since martensitic transformation is dominated by shear strain while volume strain is relatively small. For a given quasi-spin configuration described by the field variable \(s\left( {\mathbf{r}} \right)\), the lattice misfit strain field is:

$$\begin{gathered} \varepsilon_{ij}^{0} \left( {\mathbf{r}} \right) = \frac{1}{2}\varepsilon_{ij}^{0} \left( 1 \right)\left[ {2 - s\left( {\mathbf{r}} \right)} \right]\left[ {3 - s\left( {\mathbf{r}} \right)} \right] + \varepsilon_{ij}^{0} \left( 2 \right)\left[ {s\left( {\mathbf{r}} \right) - 1} \right]\left[ {3 - s\left( {\mathbf{r}} \right)} \right] \\ + \frac{1}{2}\varepsilon_{ij}^{0} \left( 3 \right)\left[ {s\left( {\mathbf{r}} \right) - 1} \right]\left[ {s\left( {\mathbf{r}} \right) - 2} \right] \\ \end{gathered}$$
(2)

The elastic energy generated by this lattice misfit strain field \(\varepsilon_{ij}^{0} \left( {\mathbf{r}} \right)\) is [3]:

$$E_{{{\text{el}}}} = \frac{1}{2}\rlap{--} \smallint \frac{{d^{3} k}}{{\left( {2\pi } \right)^{3} }}B_{{ijkl}} \left( {\mathbf{n}} \right)\;\tilde{\varepsilon }_{{ij}}^{0} \left( {\mathbf{k}} \right)~\tilde{\varepsilon }_{{kl}}^{0} \left( {\mathbf{k}} \right)^{*}$$
(3)

where \(B_{ijkl} \left( {\mathbf{n}} \right) = C_{ijkl} - n_{p} \;C_{pqij} \;\Omega_{qr} \left( {\mathbf{n}} \right)C_{rskl} \;n_{s}\), \(C_{ijkl}\) is elastic modulus tensor, k is a vector in reciprocal space, \({\mathbf{n}} = {{\mathbf{k}} \mathord{\left/ {\vphantom {{\mathbf{k}} k}} \right. \kern-0pt} k}\) is a unit directional vector in k-space, \(\Omega_{ij} \left( {\mathbf{n}} \right)\) is the elastic Green function tensor inverse to \(\Omega_{ij}^{ - 1} ({\mathbf{n}}) = C_{ikjl}^{{}} n_{k} n_{l}\), \({\widetilde{\varepsilon }}_{ij}^{0}\left(\mathbf{k}\right)\) is the Fourier transform of the lattice misfit strain field \(\varepsilon_{ij}^{0} \left( {\mathbf{r}} \right)\), the superscript asterisk * indicates complex conjugate, summation convention over repeated indices is assumed, and \(\rlap{--} \smallint\) is a principal-value integral in k-space that is evaluated by excluding the point \({\mathbf{k}} = {\mathbf{0}}\). When elastic modulus tensor of cubic symmetry is considered, instead of the three cubic elastic constants \(C_{11}\), \(C_{12}\) and \(C_{44}\), it is more convenient to formulate in terms of shear modulus \(C^{\prime} = {{\left( {C_{11} - C_{12} } \right)} \mathord{\left/ {\vphantom {{\left( {C_{11} - C_{12} } \right)} 2}} \right. \kern-0pt} 2}\) on (110) plane, Poisson’s ratio \(\nu = {{C_{12} } \mathord{\left/ {\vphantom {{C_{12} } {\left( {C_{11} + C_{12} } \right)}}} \right. \kern-0pt} {\left( {C_{11} + C_{12} } \right)}}\) under [100] uniaxial stress, and elastic anisotropy parameter \(A = {{C_{44} } \mathord{\left/ {\vphantom {{C_{44} } {C^{\prime}}}} \right. \kern-0pt} {C^{\prime}}}\) defined as the ratio between shear modulus \(C = C_{44}\) on (100) plane and shear modulus \(C^{\prime}\) on (110) plane. In the case of isotropic elasticity, \(A = 1\). In the simulations, a typical value of Poisson’s ratio ν = 1/3 is assumed, various values of elastic anisotropy parameter A is considered, and the shear modulus \(C = C_{44}\) is used to define the reference energy unit for nondimensionalization as described in "Results and Discussions" section.

The elastic energy in Eq. (3) is the Hamiltonian function, which is used to evolve the quasi-spin configuration and explore equilibrium state at given temperature via Monte Carlo simulation using Metropolis algorithm [26]. In each sampling trial, one lattice site is randomly selected, its current orientation variant is switched randomly into one of the other two variants, and the acceptance ratio of this attempt is given by the probability ratio of the Boltzmann factors of the attempted and current configurations: \(\alpha = \exp \left( { - {{\Delta E_{{{\text{el}}}} } \mathord{\left/ {\vphantom {{\Delta E_{{{\text{el}}}} } {k_{{\text{B}}} T}}} \right. \kern-0pt} {k_{{\text{B}}} T}}} \right)\), where \(k_{{\text{B}}}\) is Boltzmann constant, \(T\) is absolute temperature, and \(\Delta E_{{{\text{el}}}}\) is the elastic energy change associated with this attempted variant switching on the selected lattice site which determines the acceptance ratio of the attempt in each Monte Carlo sampling trial.

Here it is worth making a brief conceptual comparison between the above formulated model with the Ising model of ferromagnetism. The role of exchange interaction between nearest-neighbor spin pairs in the Ising model of ferromagnetism is now replaced by the long-range dipole–dipole-like elastic interactions among the orientation variants over all lattice sites. At sufficiently low temperature, the ferromagnetic exchange interactions (with positive sign) bring the system into an ordered spin configuration of aligned spins, i.e., ferromagnetic phase; similarly, the elastic interactions (which always produce a nonnegative total elastic energy [27]) bring the system into an ordered quasi-spin configuration described by the field variable \(s\left( {\mathbf{r}} \right)\), which characterizes martensitic phase. At sufficiently high temperature, thermal fluctuations disturb the spins into disordered configurations and each site possesses a zero average magnetic dipole moment, i.e., paramagnetic phase; similarly, disordered quasi-spin configurations described by the time-dependent field variable \(s\left( {{\mathbf{r}};t} \right)\) are expected to result from thermal fluctuations, where each lattice site frequently switches among the three orientation variants and thus the unit cell possesses an average cubic symmetry, i.e., austenitic phase. It is worth noting that, since the high-temperature austenitic phase is characterized by an ergodic disordered state, where the average (both spatially and temporally) symmetry is cubic, a separate state is not needed to describe the austenitic state. In addition to the phase transitions occurring in both models, there is another aspect of fundamental importance in the analogy: a strong correlation develops among the spins at temperature close to but above the Curie point of ferromagnetic transition, which reflects the short-range ordering in the paramagnetic state; similarly, displacement short-range ordering is expected in the austenitic phase, which could shed light on the pre-martensitic phenomena. In particular, this proposed quasi-spin Ising model of ferroelastic phase transition and thermoelastic martensitic transformation will be employed to investigate the displacement short-range ordering phenomenon. To this end, diffuse scattering produced by the short-range ordered displacement field will be computed for the simulated quasi-spin configuration \(s\left( {{\mathbf{r}};t} \right)\), as described next.

The heterogeneous displacement field \(v_{i} \left( {\mathbf{r}} \right)\) produces local heterogeneous strain in the lattice, which describes the structural deviations from the homogeneous average lattice structure. The Fourier transform \({\widetilde{v}}_{i}\left(\mathbf{k}\right)\) of the heterogeneous displacement field \(v_{i} \left( {\mathbf{r}} \right)\) is given by [3]:

$${\widetilde{v}}_{i}\left(\mathbf{k}\right)= - \frac{i}{k}\Omega_{ij} \left( {\mathbf{n}} \right)C_{jklm} n_{k} {\widetilde{\varepsilon }}_{lm}^{0}\left(\mathbf{k}\right)$$
(4)

where \({\widetilde{v}}_{i}\left(\mathbf{k}=0\right)=0\). This heterogeneous displacement field \(v_{i} \left( {\mathbf{r}} \right)\) in the lattice produces diffuse scattering, where the scattering intensity is distributed in the reciprocal space at positions away from the Bragg reflection peaks. In particular, the diffuse scattering intensity distribution around Bragg peak centered at fundamental reciprocal lattice vector K is given by [28]:

$$I_{{{\text{diff}}}} \left( {{\mathbf{K}};{\mathbf{k}}} \right) = V_{{{\text{cell}}}}^{ - 2}\,{\left|{\widetilde{n}}_{0}\right|}^{2}{\left|(\mathbf{K}+\mathbf{k})\bullet \widetilde{\mathbf{v}}(\mathbf{k})\right|}^{2}\propto {\left|(\mathbf{K}+\mathbf{k})\bullet \widetilde{\mathbf{v}}(\mathbf{k})\right|}^{2}$$
(5)

where \(V_{{{\text{cell}}}}\) is the unit cell volume of the lattice, \({\widetilde{n}}_{0}\) is the structure factor, the reciprocal lattice vector K can be expressed via the Bragg Peak Index \(\left( {HKL} \right)\), and k is defined within the first Brillouin zone centered at K. For a cubic array of lattice sites r, \({\mathbf{K}} \propto \left( {H,K,L} \right)\). The diffuse scattering exhibits different intensity distribution features around different Bragg peaks through the dependence on the peak indices \(\left( {HKL} \right)\), as shown in "Results and Discussions" section. It is worth noting that for more convenient computing, the diffuse scattering intensity in Eq. (5) as well as the elastic energy in Eq. (3) can be expressed in terms of a three-component field variable, the tetragonality, anisotropy parameter, and Poisson’s ratio, as detailed elsewhere for brevity [29].

It is worth noting that the atomic scale of the quasi-spin Ising model is determined by the nature of the quasi-spin variable that is defined at each cubic lattice site to characterize the lattice misfit strain of the unit cell located at the lattice site. The elastic interaction is treated by continuum elasticity theory, which is a good approximation when the interaction is long-ranged (elastic interaction is of infinite range in nature). A discrete elasticity theory [3] exists that is formulated in terms of a set of Born-von Kármán force constants [30], which would improve the accuracy at the expense of computing efficiency. As shown in "Results and Discussions" section, the current model is adequate to address the diffuse scattering phenomenon, which is the focus of the current study. An atomic-scale discrete elasticity theory will be implemented in future work.

It is also worth noting that the computational grid lattice is cubic, and the lattice sites correspond to the locations of the unit cells of the cubic crystal (rather than the primitive cells). This model is applicable to cubic crystal systems (simple cubic, body-centered cubic, face-centered cubic). The lattice misfit strain at each lattice site characterizes the Bain distortion of the unit cell in ground state, which is defined for the unit cell (involving 1, 2, and 4 atoms in the respective cases of simple cubic, body-centered cubic, and face-centered cubic crystals). This treatment is of atomic scale in nature, which is different from the coarse-graining in mesoscale continuum phase field model.

High-Energy Synchrotron X-Ray Single-Crystal Diffraction

To complement the simulation study, 3D diffuse scattering from thermoelastic martensitic Ni–Mn–Ga alloy single crystals is measured experimentally at temperatures above the martensitic transformation by using in-situ high-energy synchrotron X-ray diffraction. The diffuse scattering experiments are performed at the beamline 11-ID-C of the Advanced Photon Source at Argonne National Laboratory. High-flux high-energy synchrotron X-ray single-crystal diffraction is essential for such experiments: the high energy (115 keV) provides deep penetration to probe bulk single crystal specimens, the corresponding short wavelength (0.107805 Å) enables coverage of large volume in reciprocal space to measure high-index peaks, the high flux enables short exposure time (seconds for each exposure) to measure weak diffuse scattering intensity that is critical for 3D reciprocal space mapping, and 3D single-crystal diffraction allows separation of different diffuse scattering rods associated with displacement correlations in different directions—2D measurement using rocking crystal method causes intensity overlapping problem where diffuse scattering intensities in different directions (such as along 12 < 110 > axes) are all projected into one image. A beamsize of 400 × 400 μm2 cross-section is used. Scattering intensity data are measured by 2D digital X-ray detector (PerkinElmer® XRD 1621 AN, 2048 × 2048 pixels, 200 μm pixel size, 16 bit digitization). Temperature is controlled in-situ by using cryostream (Oxford Cryosystems® 700 Plus). To scan the 3D reciprocal space, the single crystal specimen is rotated incrementally by Δω = 0.1° during each frame exposure. At each rotation position ω (defined at the midpoint of each increment Δω), the 2D detector measures the scattering intensity distribution that intersects the Ewald sphere. The 3D reciprocal space map of diffuse scattering intensity is constructed from a series of such 2D data. Details of the experiments are reported elsewhere [11, 31]. Figure 1 shows the measured 3D diffuse scattering around (800) Bragg reflection peak in Ni49.90Mn28.75Ga21.35 at 490 K and 327 K, respectively, above the martensite start temperature 323 K, as well as the measured 3D diffuse scattering around different Bragg reflection peaks \(\left( {HKL} \right)\) at 330 K.

Results and Discussions

To carry out Monte Carlo simulation based on the above formulated quasi-spin Ising model of ferroelastic phase transition and thermoelastic martensitic transformation, it is convenient to introduce normalized dimensionless temperature. To this end, an energy unit is defined as \(E_{0} = C\gamma_{0}^{2} \Delta V\), where C is an elastic constant and \(\Delta V = {V \mathord{\left/ {\vphantom {V N}} \right. \kern-0pt} N}\) is the unit cell volume of the lattice. A reference temperature is defined as \(T_{0} = {{E_{0} } \mathord{\left/ {\vphantom {{E_{0} } {k_{{\text{B}}} }}} \right. \kern-0pt} {k_{{\text{B}}} }}\), and the dimensionless temperature \(T^{*}\) is defined by normalizing the temperature T with respect to the reference temperature \(T_{0}\), i.e., \(T^{*} = {T \mathord{\left/ {\vphantom {T {T_{0} }}} \right. \kern-0pt} {T_{0} }}\). In the following, the simulation results are presented in \(T^{*}\). It is worth noting that the shear modulus \(C = C_{44}\) on (100) plane, rather than the shear modulus \(C^{\prime} = {{\left( {C_{11} - C_{12} } \right)} \mathord{\left/ {\vphantom {{\left( {C_{11} - C_{12} } \right)} 2}} \right. \kern-0pt} 2}\) on (110) plane, is used in the normalization, because \(C^{\prime}\) usually exhibits significant temperature-dependent softening in pre-martensitic austenite while \(C\) remains approximately constant. In the simulations presented in the following, 32 × 32 × 32 computational grids are employed, and isotropic elasticity is assumed, i.e., \(A = 1\), except for "Effects of Elastic Anisotropy and Softening of Shear Modulus" section which addresses the effects of elastic anisotropy on the phase transition. The simulation starts from an assumed disordered state of randomly distributed quasi-spin variables at sufficiently high temperature and continues to reach an equilibrium disordered austenite state. The temperature is gradually decreased to reach new equilibrium at each lower temperature. Thermal equilibrium is determined by the criterion of no noticeable difference in the calculated diffuse scattering (which is the focus of this study). Test runs show that 2000 lattice sweeps (Monte Carlo steps) at each temperature are sufficient according to this equilibrium criterion.

Phase Transition

Figure 2 shows simulated equilibrium states representative of disordered quasi-spin configuration at \(T^{*} = 0.22\) and ordered configuration at \(T^{*} = 0.20\), respectively. In the disordered state, the three martensite orientation variants appear to randomly distribute throughout the system, where the crystal lattice averaged over the system volume is cubic, and the local unit cell averaged over the simulation time is also cubic. In the ordered state, the orientation variants self-assemble into polytwinned plates characteristic of martensitic microstructures, which reflects development of long-range order in the martensite [32]. It is worth noting that the true ground state of the single crystal is a single-domained tetragonal martensite, which, however, can be achieved only through coarsening process driven by reduction in interfacial energy. In the current model, no interfacial energy term is introduced, thus no coarsening occurs. Since the interfacial energy term is not important for a study of disordered austenite at pre-martensitic temperature (which is the focus of this work), a simplest possible model yet able to capture the main feature of the pre-martensitic phenomenon (i.e., displacement short-range ordering and diffuse scattering) is desired. Interfacial energy term will be introduced to study the coarsened martensite morphology in future work.

Fig. 2
figure 2

Simulated equilibrium multi-variant configurations representative of a disordered austenitic state at T* = 0.22 and b ordered martensitic state at T* = 0.20. Top row illustrates spatial distributions of three orientation variants by three colors (red, green, blue). Bottom row illustrates elastic distortion in the lattice due to coherency strain (for clarity, \(\gamma_{0} = \varepsilon_{c} - \varepsilon_{a} = 0.4\) is assumed here)

According to fluctuation–dissipation theorem [33], the thermal fluctuation of the long-range order parameter (i.e., mean square fluctuation) of the phase transition is proportional to \(\left( {T - T_{{\text{C}}} } \right)^{ - 1}\) in the vicinity of the phase transition temperature \(T_{{\text{C}}}\). This relation allows a determination of \(T_{{\text{C}}}\) based on simulation results at a limited number of temperatures near \(T_{{\text{C}}}\). In ferroelastic phase transition and thermoelastic martensitic transformation, which are displacive structural transformations, strain is the long-range order parameter. Thermal fluctuations in the macroscopic strain \(\overline{\varepsilon }_{ij}^{*} = {{\overline{\varepsilon }_{ij} } \mathord{\left/ {\vphantom {{\overline{\varepsilon }_{ij} } {\gamma_{0} }}} \right. \kern-0pt} {\gamma_{0} }}\), in particular, the three principal components \(\overline{\varepsilon }_{11}^{*} \left( {t^{*} } \right)\), \(\overline{\varepsilon }_{22}^{*} \left( {t^{*} } \right)\) and \(\overline{\varepsilon }_{33}^{*} \left( {t^{*} } \right)\), are simulated at temperatures above \(T_{{\text{C}}}\), and two examples are plotted as function of the simulation time \(t^{*}\) (i.e., Monte Carlo steps) in Fig. 3(a, b). As expected, the strain fluctuations increase when T decreases towards \(T_{{\text{C}}}\). The standard deviation \(\Delta \varepsilon_{ij}^{*}\) in each strain component \(\overline{\varepsilon }_{ij}^{*} \left( {t^{*} } \right)\) is calculated with respect to its time-average value \(\left\langle {\overline{\varepsilon }_{ij}^{*} \left( {t^{*} } \right)} \right\rangle\), and the values of \(\left( {\Delta \varepsilon^{*\,} } \right)^{2} = \left( {\Delta \varepsilon_{11}^{*} } \right)^{2} + \left( {\Delta \varepsilon_{22}^{*} } \right)^{2} + \left( {\Delta \varepsilon_{33}^{*} } \right)^{2}\) and its reciprocal \(\left( {\Delta \varepsilon^{*} } \right)^{ - 2}\) are plotted as function of the temperature \(T^{*}\) in Fig. 3(c). A linear fitting of the data points in the vicinity of \(T_{{\text{C}}}\) in the functional form of \(\left( {\Delta \varepsilon } \right)^{2} \propto \left( {T - T_{{\text{C}}} } \right)^{ - 1}\) yields \(T_{{\text{C}}}^{*} = 0.209\), as shown in Fig. 3(c), which is in agreement with the simulation results presented in Fig. 2. In the discussions of the simulation results, it is more meaningful to consider the “homologous” temperature \(T_{{\text{H}}} = {{T^{*} } \mathord{\left/ {\vphantom {{T^{*} } {T_{{\text{C}}}^{*} }}} \right. \kern-0pt} {T_{{\text{C}}}^{*} }} = {T \mathord{\left/ {\vphantom {T {T_{{\text{C}}} }}} \right. \kern-0pt} {T_{{\text{C}}} }},\) which is provided in the figures and/or captions presented in this paper to relate the simulations and experiments.

Fig. 3
figure 3

Simulated fluctuations of strain components at representative temperatures above phase transition temperature: a \(T^{*} = 0.3\) and b \(T^{*} = 0.22\). c Strain fluctuation as a function of temperature. The blue line represents linear fitting in the vicinity of the phase transition temperature, which gives \(T_{{\text{C}}}^{*} = 0.209\)

It is worth noting that the above observed thermal fluctuation of the long-range order parameter (i.e., lattice strain) as a function of temperature is related to an important phonon softening anomaly in pre-martensitic phenomena. When an austenitic system approaches TC upon cooling, its lattice stability is gradually weakened, because the ground state of such a system is martensitic phase. This leads to reduced resistance to lattice vibration and thus enhanced thermal fluctuation. It is in contrast to normal systems, where thermal fluctuation decreases with decreasing temperature. The simulation results in Fig. 3 imply that the phonon softening is an intrinsic behavior (rather than an anomaly [4, 5, 10, 11]) of martensitic system.

Diffuse Scattering and Displacement Short-Range Ordering

Development of short-range order prior to phase transition producing long-range order is a phenomenon of fundamental importance. In particular, the displacement short-range ordering (i.e., spatial correlation in the heterogeneous lattice displacements) in the pre-martensitic austenite state is an interesting phenomenon. Existence of such displacement short-range ordering produces diffuse scattering in diffraction at intensities away from sharp Bragg reflection peaks, as described in "Methodology" section. Figure 4(c) shows the simulated diffuse scattering around (H00) Bragg reflection peak at T* = 0.3 (T = 1.44TC) computed using Eq. (5), i.e., \(I_{{{\text{diff}}}} \left( {{\mathbf{K}};{\mathbf{k}}} \right) \propto\)\({\left|(\mathbf{K}+\mathbf{k})\bullet \widetilde{\mathbf{v}}(\mathbf{k})\right|}^{2}\). It is worth noting that, to link the simulated diffuse scattering to the experimentally measured diffuse scattering presented in Fig. 1, a time averaging of the simulated diffuse scattering intensity over a sufficient number of Monte Carlo steps (1000 in this work) after reaching thermodynamic equilibrium is necessary, because of both the relatively small computational supercell (32 × 32 × 32 used in this work) affordable to computer simulations and the relatively long exposure time required to experimentally measure the weak diffuse scattering intensity (discussed in "Methodology" section). Figure 4(a) and (b) demonstrate the simulated instantaneous diffuse scattering, i.e., without time averaging, around (H00) Bragg reflection peak at T* = 0.3 (T = 1.44TC) and the corresponding equilibrium microstructures at different simulation time t* (i.e., Monte Carlo step), after the system equilibrates during t* = 0 to 1000. As expected, the instantaneous diffuse scattering intensity shown in Fig. 4 is very noisy due to the poor statistics associated with the small computational supercell affordable to the Monte Carlo simulations. It is observed that, more importantly, the time averaging procedure very effectively improves the statistics and the time-averaged diffuse scattering intensity shown in Fig. 4(c) achieves significantly improved signal-to-noise ratio and exhibits the major features in agreement with the experimentally measured diffuse scattering shown in Fig. 1. In this paper, simulation results of the time-averaged diffuse scattering are presented unless otherwise explicitly specified such as in Fig. 4(a, b).

Fig. 4
figure 4

Simulated diffuse scattering around (H00) Bragg reflection peak at T* = 0.3 (T = 1.44TC): a, b instantaneous scattering intensity corresponding to equilibrium microstructures at different simulation time (i.e., Monte Carlo step) t* = 1000 and t* = 1500, respectively, and c time-averaged scattering intensity

It is worth noting that diffuse scattering exists in both defect-free pre-martensite and defect-rich disordered strain glass. In the defect-free pre-martensite, the thermal lattice displacements and resultant local strain are dynamic and the state is ergodic, as shown in Fig. 4(a, b); while in the defect-rich disordered strain glass, the local strain and nanodomains are static and the state is frozen (non-ergodic). This is an important difference and the two pictures complement each other to provide a more complete insight into martensitic transformation in real materials. In particular, the existence of local defects (compositional disorder) may lift the ergodicity of the pre-martensitic state, leading to localization of the lattice dynamics.

It is also worth noting that the diffuse scattering is dominated by K and\(\widetilde{\mathbf{v}}(\mathbf{k})\), while the appearance of k in Eq. (5) accounts for the effect of local volume strain [28]. The volume strain effect is small in the case considered here and can be safely neglected, thus reducing Eq. (5) to the simplified formula \(I_{{{\text{diff}}}} \left( {{\mathbf{K}};{\mathbf{k}}} \right) \propto\)\({\left|\mathbf{K}\bullet \widetilde{\mathbf{v}}(\mathbf{k})\right|}^{2}\), which makes it easier to reveal the physical implications of the diffuse scattering. As one example, consider the diffuse scattering around a (H00) Bragg reflection peak, \(I_{{{\text{diff}}}} \propto\)\({\left|{\widetilde{v}}_{x}\left(\mathbf{k}\right)\right|}^{2}={\widetilde{v}}_{x}\left(\mathbf{k}\right){\widetilde{v}}_{x}\left(-\mathbf{k}\right)\). It is readily shown that \(I_{{{\text{diff}}}} \propto\) \(\widetilde{P}(\mathbf{k})\)\(= \int {P\left( {{\varvec{\uprho}}} \right)} \,e^{{ - i{\mathbf{k}} \cdot {{\varvec{\uprho}}}}} \,d^{3} \rho\), where

$$P\left( {{\varvec{\uprho}}} \right) = \int {v_{x} \left( {\mathbf{r}} \right)v_{x} \left( {{\mathbf{r}} + {{\varvec{\uprho}}}} \right)} \,d^{3} r$$
(6)

is the auto-correlation function of the displacement field \(v_{x} \left( {\mathbf{r}} \right)\) in the lattice. This displacement field is generated by both the spontaneous lattice misfit strain of the martensite orientation variants and the elastic strain to maintain the lattice coherency. Figure 5 shows the time-averaged auto-correlation \(\left\langle {P\left( {{\varvec{\uprho}}} \right)} \right\rangle = \int {\left\langle {v_{x} \left( {\mathbf{r}} \right)v_{x} \left( {{\mathbf{r}} + {{\varvec{\uprho}}}} \right)} \right\rangle } \,d^{3} r\), corresponding to the time-averaged diffuse scattering intensity shown in Fig. 4(c). As discussed above, the time averaging effectively improves the statistics to achieve the results equivalent to that in adequately large systems (such as in experiments). The auto-correlation function \(P\left( {{\varvec{\uprho}}} \right)\) is analogous to the Patterson function of the electron density field in a crystal structure [34] and also the pair distribution function of the atoms in a material [35]. Unlike in long-range ordered systems where sharp peaks appear in the Patterson function and pair distribution function, which correspond to well-defined interatomic distance vectors in ordered atomic arrangements, the auto-correlation function \(P\left( {{\varvec{\uprho}}} \right)\) of the heterogeneous displacement field \(v_{x} \left( {\mathbf{r}} \right)\) in the disordered austenite state exhibits significant values only at small distance \(\rho\) and vanishes rapidly with increasing distance, as shown in Fig. 5, which characterizes short-range ordering in the displacement field. The displacement short-range ordering is dictated by the elastic interactions among the orientation variants: although the thermal fluctuations at temperature above the phase transition prevent a formation of long-range ordered configurations, the orientation variants develop a spatial correlation over short distance to reduce the elastic interaction energy, despite their apparently disordered distribution. Additionally, \(P\left( {{\varvec{\uprho}}} \right)\) exhibits a pronounced directional dependence, i.e., stronger correlation in \(\left\langle {110} \right\rangle\) directions, indicating stronger displacement short-range ordering over a longer distance in \(\left\langle {110} \right\rangle\) directions. Such anisotropic features are also exhibited in the diffuse scattering shown in Fig. 4, where stronger diffuse scattering intensity is produced along \(\left\langle {110} \right\rangle\) directions. It has been shown that such \(\left\langle {110} \right\rangle\) diffuse scattering rods are attributed to the displacement waves whose wavevectors are along \(\left\langle {110} \right\rangle\) directions and displacement vectors are along \(\left\langle {1\overline{1}0} \right\rangle\) directions [28]. Indeed, the Fourier transform \(\widetilde{\mathbf{v}}(\mathbf{k})\) is the plane wave representation of the heterogeneous displacement field \({\mathbf{v}}\left( {\mathbf{r}} \right)\). It is worth noting that such displacement waves correspond to \(\left. {\left\{ {110} \right\}\,} \right|\left\langle {1\overline{1}0} \right\rangle\) shear strains and are associated with the twinning deformations on \(\left\{ {110} \right\}\) planes that transform one orientation variant into another variant of the tetragonal martensite.

Fig. 5
figure 5

Simulated time-averaged auto-correlation of displacement field \(v_{x} \left( {\mathbf{r}} \right)\). Isosurface at 10% maximum cross-correlation value is visualized within 32 × 32 × 32 computational supercell to illustrate the nature of displacement short-range ordering

Bragg Peak Dependence of Diffuse Scattering and Extinction Rule

According to Eq. (5), the diffuse scattering around a Bragg reflection peak with a general index (HKL) reflects the spatial correlation of the displacement field \({\mathbf{v}}\left( {\mathbf{r}} \right)\) involving all three displacement components, i.e., \(Hv_{x} \left( {\mathbf{r}} \right) + Kv_{y} \left( {\mathbf{r}} \right) + Lv_{z} \left( {\mathbf{r}} \right)\), which indicates a strong dependence of the diffuse scattering on the Bragg peak index (HKL) or the reciprocal lattice vector K. Figure 6 shows the simulated diffuse scattering around different Bragg reflection peaks at temperature T* = 0.25 (T = 1.20TC), which demonstrates the (HKL) or K dependence of the diffuse scattering in agreement with the experimentally measured diffuse scattering shown in Fig. 1. There are a total of 6 diffuse scattering rods (or 12 diffuse scattering spikes, i.e., every rod consists of two spikes pointing in opposite directions), each along one of \(\left\langle {110} \right\rangle\) directions. It is worth noting that not all these diffuse scattering rods and spikes appear around every Bragg peak. As a matter of fact, individual diffuse scattering rods and spikes exhibit different relative intensities around different Bragg peaks and even completely disappear around certain Bragg peaks. Such extinction rule of diffuse scattering is described by \(I_{{{\text{diff}}}} \left( {{\mathbf{K}};{\mathbf{k}}} \right) \propto\)\({\left|\mathbf{K}\bullet \widetilde{\mathbf{v}}(\mathbf{k})\right|}^{2}\). It has been shown that the heterogeneous displacement field \({\mathbf{v}}\left( {\mathbf{r}} \right)\) producing such \(\left\langle {110} \right\rangle\) diffuse scattering can be expressed as a sum of 12 branches of displacement plane waves with wavevectors along \(\left\langle {110} \right\rangle\) directions and displacement vectors along \(\left\langle {1\overline{1}0} \right\rangle\) directions [11, 28]. In particular, the extinction rule of diffuse scattering allows us to associate each diffuse scattering rod to a particular displacement plane wave branch. For example, the diffuse scattering rod along \(\left[ {110} \right]\) direction appears around (100) Bragg peak but disappears around (110) and (001) Bragg peaks, as shown in Fig. 6. Since the displacement plane wave branch producing the \(\left[ {110} \right]\) diffuse scattering rod has wavevectors along \(\left[ {110} \right]\) direction and displacement vectors along \(\left[ {1\overline{1}0} \right]\) direction, i.e., \({{\varvec{\upnu}}} = \left( {v_{x} , - v_{x} ,0} \right)\), \({\mathbf{K}} \cdot {\mathbf{v}} = 0\) at (HKL) = (110) and (001) leads to extinction of diffuse scattering around (110) and (001) Bragg peaks. Around (H00) Bragg peaks as shown in Figurers 4 and 6, the 4 diffuse scattering rods, respectively, along \(\left[ {110} \right]\), \(\left[ {1\overline{1}0} \right]\), \(\left[ {101} \right]\), and \(\left[ {10\overline{1}} \right]\) appear with equal intensity, while the 2 diffuse scattering rods, respectively, along \(\left[ {011} \right]\) and \(\left[ {01\overline{1}} \right]\) completely disappear, since the displacement vectors in the latter two branches are, respectively, along \(\left[ {01\overline{1}} \right]\) and \(\left[ {011} \right]\) directions, where \({\mathbf{K}} \cdot {\mathbf{v}} = 0\) for (H00) Bragg peaks.

Fig. 6
figure 6

Simulated diffuse scattering around different Bragg reflection peaks (HKL) at temperature T* = 0.25 (T = 1.20TC)

It is worth noting that the above observed diffuse scattering phenomenon is regarded as an important diffraction anomaly in pre-martensitic phenomena. Such characteristic diffuse scattering is not expected in cubic phase. Even though the ergodic pre-martensitic state possesses an average (both spatially and temporally) cubic structure, at every time moment the instantaneous structure deviates from the cubic structure due to the displacement short-range ordering (i.e., spatially correlated heterogeneous lattice displacements). As shown above, such heterogeneous displacement field is a sum of 12 branches of displacement plane waves with wavevectors along \(\left\langle {110} \right\rangle\) directions and displacement vectors along \(\left\langle {1\overline{1}0} \right\rangle\) directions. Such displacement plane waves appear with equal probabilities (resulting in ergodicity) and produce \(\left\langle {110} \right\rangle\) diffuse scattering because they modulate the cubic crystal lattice. Although each diffuse scattering rod exhibits K-dependent extinction rule as shown in Fig. 6, the overall diffraction consisting of the diffuse scattering around all fundamental Bragg reflection peaks (HKL) exhibit cubic symmetry, manifesting the average cubic structure of the pre-martensitic crystal. For example, the diffuse scattering around (110), (101) and (011) are related by cubic operations, and the diffuse scattering around (111) exhibits hexagonal symmetry inherent to cubic crystal.

Temperature Dependence of Diffuse Scattering

Figure 7 shows the simulated diffuse scattering around (H00) Bragg reflection peak at different temperatures above phase transition temperature. The diffuse scattering becomes stronger upon cooling towards the phase transition temperature, in agreement with the experimentally measured diffuse scattering shown in Fig. 1. As discussed above, the temperature dependence of the diffuse scattering manifests the temperature dependence of the displacement short-range ordering: with decreasing temperature, the correlation strength and length of the displacement short-range ordering increase, since the role of elastic interactions among the orientation variants characterized by \(\Delta E_{{{\text{el}}}}\) becomes more prominent with respect to the thermal energy characterized by \(k_{{\text{B}}} T\) at lower temperature. The diffuse scattering shown in Fig. 7 provides a direct means to quantitatively characterize the spatial correlation and short-range ordering of the heterogeneous displacement field at different temperatures in the pre-martensitic austenite state.

Fig. 7
figure 7

Simulated diffuse scattering around (H00) Bragg reflection peak at different temperatures above phase transition temperature \(T_{{\text{C}}}^{*} = 0.209\): a T* = 1, b T* = 0.5, c T* = 0.3, and d T* = 0.22

Effects of Elastic Anisotropy and Softening of Shear Modulus \(\user2{C^{\prime}}\)

In this section, elastic anisotropy is considered. It is worth noting that while the shear modulus \(C = C_{44}\) on (100) plane is relatively insensitive to temperature and remains approximately constant, the shear modulus \(C^{\prime} = {{\left( {C_{11} - C_{12} } \right)} \mathord{\left/ {\vphantom {{\left( {C_{11} - C_{12} } \right)} 2}} \right. \kern-0pt} 2}\) on (110) plane usually exhibits significant softening upon cooling in pre-martensitic austenite state. Thus, the anisotropy factor defined as Zener ratio between the two shear moduli, i.e., \(A = {C \mathord{\left/ {\vphantom {C {C^{\prime}}}} \right. \kern-0pt} {C^{\prime}}}\) [36], reflects the degree of elastic softening in the shear modulus \(C^{\prime}\). Based on the results from a series of Monte Carlo simulations, Fig. 8(a) shows the phase transition temperature \(T_{{\text{C}}}^{*}\) as a function of the reciprocal anisotropy factor \(A^{ - 1} = {{C^{\prime}} \mathord{\left/ {\vphantom {{C^{\prime}} C}} \right. \kern-0pt} C}\), which is proportional to the shear modulus \(C^{\prime}\). It is observed that elastic softening in the shear modulus \(C^{\prime}\) leads to decrease in the phase transition temperature. This is caused by decrease in the strength of elastic interactions among the orientation variants, resulting in more prominent role of thermal fluctuations and decrease in the phase transition temperature. Therefore, the elastic softening in the shear modulus \(C^{\prime}\) stabilizes the disordered state (austenite).

Fig. 8
figure 8

Effects of elastic anisotropy. a Simulated dependence of phase transition temperature on reciprocal anisotropy factor \(A^{ - 1}\) (proportional to shear modulus \(C^{\prime}\)), where dashed line is a guide to the eye. Simulated diffuse scattering around (H00) Bragg reflection peak at temperature T* = 0.3 with different elastic anisotropy factors: b A = 10 corresponding to \(C^{\prime}\) softening; c A = 1 corresponding to elastic isotropy; and d A = 0.6 corresponding to \(C^{\prime}\) hardening

Figure 8(b–d) show the simulated diffuse scattering around (H00) Bragg reflection peak at temperature T* = 0.3 with different elastic anisotropy factors (i.e., different degrees of elastic softening in the shear modulus \(C^{\prime}\)). It is observed that, at the same temperature T* = 0.3, softened shear modulus \(C^{\prime}\) leads to weaker diffuse scattering. This is caused by the weaker elastic interactions among the orientation variants due to decreased \(C^{\prime}\) and thus weaker spatial correlation in the heterogeneous displacement field. In addition, the same temperature T* = 0.3 corresponds to different “homologous” temperatures \(T_{{\text{H}}} = {{T^{*} } \mathord{\left/ {\vphantom {{T^{*} } {T_{{\text{C}}}^{*} }}} \right. \kern-0pt} {T_{{\text{C}}}^{*} }}\) in the systems with different phase transition temperatures \(T_{{\text{C}}}^{*}\) due to different elastic anisotropy factors (different degrees of elastic softening in the shear modulus \(C^{\prime}\)): according to \(T_{{\text{C}}}^{*}\) shown in Fig. 8(a), \(T_{{\text{H}}} \approx 7\) in the softening case shown in Fig. 8(b), and \(T_{{\text{H}}} \approx 1.1\) in the hardening case shown in Fig. 8(d), while \(T_{{\text{H}}} = 1.44\) in the isotropy case shown in Fig. 8(c). In particular, the hardening case of A = 0.6 corresponds to the lowest “homologous” temperature \(T_{{\text{H}}} \approx 1.1\) among the three cases and consequently exhibits the strongest diffuse scattering among the three cases. It is worth noting that the elastic softening in the shear modulus \(C^{\prime}\) leads to sharper diffuse scattering rods and spikes as shown in Fig. 8 (b), in better agreement with the diffuse scattering feature shown in Fig. 1 that is experimentally measured in Ni49.90Mn28.75Ga21.35 single crystal that undergoes significant \(C^{\prime}\) softening in the pre-martensitic austenite state. Softening in \(C^{\prime}\) promotes \(\left. {\left\{ {110} \right\}\,} \right|\left\langle {1\overline{1}0} \right\rangle\) displacement plane waves that correspond to shear strains on \(\left\{ {110} \right\}\) planes, resulting in sharper diffuse scattering rods and spikes confined in \(\left\langle {110} \right\rangle\) directions.

As discussed above, elastic softening tends to stabilize the high-temperature disordered austenite phase and decrease the martensitic transformation temperature. It is worth noting that such an effect can be explained in terms of increased entropy and thus decreased free energy in the disordered austenite phase. In the quasi-spin Ising model and Monte Carlo simulations presented above, the increased entropy effect is realized by more frequent switching of the orientation variants and weaker and shorter correlation of the heterogeneous displacement field resulting from the weaker elastic interactions due to the softer shear modulus \(C^{\prime}\). This effect is analogous to the stabilization of the austenite phase by a more general form of elastic softening, namely, phonon softening, where the shear modulus softening corresponds to the long-wavelength limit of transverse acoustic phonon softening. It is well known that softened phonons stabilize the phases by increasing the lattice vibrational entropy [30]. The analogy between the two mechanisms, i.e., shear modulus softening and phonon softening, goes even further in the case of cubic → tetragonal martensitic transformations. In austenitic systems that undergo martensitic transformation into tetragonal phases, the transverse acoustic TA2 phonon branches usually exhibit softening upon cooling towards the phase transition temperature. The TA2 phonons are characterized by wavevectors along \(\left\langle {110} \right\rangle\) directions and polarization (displacement) vectors along \(\left\langle {1\overline{1}0} \right\rangle\) directions [11], which are the same as the displacement plane waves discussed in "Diffuse Scattering and Displacement Short-Range Ordering" and "Bragg Peak Dependence of Diffuse Scattering and Extinction Rule" sections. Indeed, the TA2 phonons, when going to long-wavelength limit and undergoing complete softening, are believed to be responsible for the formation of tetragonal orientation variants during the martensitic transformation. Therefore, the frequent switching of the orientation variants in the Monte Carlo simulations and the associated short-range correlated evolution of the heterogeneous displacement field, which can be represented in the form of \(\left. {\left\langle {110} \right\rangle \,} \right|\left\langle {1\overline{1}0} \right\rangle\) displacement plane waves as discussed in "Diffuse Scattering and Displacement Short-Range Ordering" section, mimic the TA2 phonons in the pre-martensitic austenite state, where both stabilize the austenite phase, decrease the phase transition temperature, and condense to form tetragonal orientation variants during the martensitic transformation.

Summary

In this work, we focus on the intrinsic behaviors of a defect-free single crystal that undergoes cubic-to-tetragonal martensitic transformation to gain insight into the diffuse scattering phenomenon in the pre-martensitic austenite state. To this end, a quasi-spin Ising model of ferroelastic phase transition is developed and employed to perform atomic-scale Monte Carlo simulation of thermoelastic martensitic transformation. The quasi-spin variable associated with the lattice sites characterizes the local unit cells of the orientation variants of the ground-state martensite phase, which interact with each other through long-range elastic interactions. The atomic-scale heterogeneous lattice displacements deviating from the average cubic lattice sites are treated effectively in a manner of displacement plane waves that mimic acoustic phonons relevant to the martensitic transformation. The diffuse scattering is correlated to the short-range ordering of the atomic-scale heterogeneous lattice displacements that develop in the pre-martensitic austenite crystal lattice prior to the development of long-range order during martensitic transformation. In particular, the \(\left\langle {110} \right\rangle\) diffuse scattering rods and spikes are attributed to the displacement plane waves with wavevectors along \(\left\langle {110} \right\rangle\) directions and displacement vectors along \(\left\langle {1\overline{1}0} \right\rangle\) directions, corresponding to \(\left. {\left\{ {110} \right\}\,} \right|\left\langle {1\overline{1}0} \right\rangle\) shear strains associated with the twinning deformations on \(\left\{ {110} \right\}\) planes that transform one orientation variant into another variant of the tetragonal martensite. The effects of temperature, elastic anisotropy, and shear modulus softening on the diffuse scattering and displacement short-range ordering are investigated. The \(\left. {\left\langle {110} \right\rangle \,} \right|\left\langle {1\overline{1}0} \right\rangle\) displacement plane waves play a dominant role analogous to the transverse acoustic TA2 phonons, both of which upon elastic softening stabilize the austenite phase, decrease the phase transition temperature, and condense to form tetragonal orientation variants during the martensitic transformation. The simulated diffuse scattering is compared and agrees with the complementary synchrotron X-ray single-crystal diffuse scattering experiment.