Elastic energy driven multivariant selection in martensites via quantum annealing

We demonstrate the use of quantum annealing for the selection of multiple martensite variants in a microstructure with long-range coherency stresses and external mechanical load. The general approach is illustrated for martensites with four different variants, based on the minimization of the linear elastic energy. The equilibrium variant distribution is then analysed under application of tensile and shear strains and for different values of the considered shear and tetragonal contributions of the different martensite variants. The interface orientations between different domains of variants can be explained using the perspective of the elastic energy anisotropy for regular stripe patterns. For random grain orientations, the response to an external elastic strain is weaker and variants changes can be interpreted based on the rotated eigenstrain tensor.


I. INTRODUCTION
The formation of microstructures is critical for many applications, as it strongly impacts the properties of materials.In many cases, simulation approaches can nowadays be used to support experimental investigations related to the understanding and prediction of microstructures.Among these approaches the phase field approach plays a central role, as it is able to predict the temporal evolution of non-equilibrium microstructures [1][2][3][4][5].To this end, continuous order parameter fields are used to discriminate between different phases, grains and in particular martensite variants, and in this way phase transformation kinetics can be simulated.Contrary to this approach, we recently presented an approach which complements the phase field picture in the sense that it separates discrete from continuous degrees of freedom [6].Specific applications are stress or strain driven transformations in shape memory alloys or martensitic transformations in general.There, we assumed for simplicity that a grain can be in one out of two martensite variant states, which differ by their eigenstrain, i.e. their stress free configuration.By application of external stresses and the presence of internal stresses, the grains can flip between the two variants, and these transitions are driven by a reduction of the overall elastic energy.
The prediction of martensite formation and variant selection has been widely discussed in the literature, see e.g.[7].In this context, the concept of a lattice correspondence between austenite and martensite gained attention, resulting e.g. in the well-known Kurdjumov-Sachs or Nishiyama-Wasserman orientation relationships.From a micromechanical perspective the strain energy minimization led to an important development, which is the phenomenological theory of martensite crystallography [8,9], which is to a large extent in agreement with experimental observations.These descriptions are different but related to global energetic considerations, which use the overall energy of the system to predict the formation of microstructures.
In the present work we focus on the long-ranged elastic effects, which are resulting from mismatches between different martensite variants.They lead to elastic misfits, which induce interactions not only between neighboring phases and grains, but also depend on and influence the entire microstructure.As discussed in [6], the proper minimization of the elastic energy requires to keep all these interactions, as an artificial cutoff leads to improper predictions of variant selections.However, this global perspective introduces an additional degree of complexity to the overall problem, as now a huge space of microstructural configurations needs to be considered.For example, for a microstructure consisting of N grains and two martensite variants, we therefore have altogether 2 N combinations, and the energy needs to be determined for each of them to find the global ground state, which corresponds to the (low temperature) thermodynamic equilibrium.It is obvious that this configuration space is too large already for moderate grain size numbers, and therefore not only a brute force energy minimization but also the use of heuristic approaches like simulated annealing are limited and require large amounts of computing times.If we consider e.g. a system consisting of just N = 20 grains, the number of configurations is 2 N ≃ O(10 6 ); for a typical FFT calculation in 2D with sufficient accuracy the single core runtime on current processors is of the order of one minute.This implies that scanning all configurations to determine the global energetic minimum requires about 17,000 core hours of computing time.The question, which arises also generally in the context of materials science, is how this optimization problem can be solved efficiently via quantum computing.
A general quantum computer with a sufficient number of qubits is not available yet, but the technology named quantum annealing (QA) [10][11][12][13][14] with several thousands of qubits and couplers is available.In contrast to general-purpose quantum computers only specific classes of problems can be handled via QA, i.e. binary quadratic optimization problems [15].The basic concept of QA is the initialization of its qubits in well defined Hamiltonians, whose ground state is unique and known [16].This Hamiltonian is then changed adiabatically under the operation at cryogenic temperatures to the desired, final Hamiltonian, what allows the conversion of the ground state to the final Hamiltonian [16,17] and therefore enables efficient global energy minimizations.The application of QA in materials science is still rare, however, few publications focus on sampling techniques via QA [18][19][20][21][22][23][24], phase transitions in the transverse field Ising model [25], critical phenomena in frustrated magnetic systems [26], energy calculations of defective graphene structures [27] and the automated design of metamaterials [28].Instead, the focus of actual research concerning quantum annealing rather lies on performance tests and benchmarking of quantum annealing against conventional approaches [29][30][31][32].
In our previous publication [6] we succeeded to map the martensite problem to a quantum annealer formulation, i.e. a binary quadratic model expressed through an Ising model or equivalently in terms of a quadratic unconstrained binary optimization [15], which allows to find the true ground state configurations even for many thousand grains in an extremely short time [6].This approach builds up on earlier works in [33][34][35][36].Although this conceptual progress demonstrates well the benefits of quantum annealing for materials science modeling, the chosen example was limited in particular by the restriction to just two variants, whereas e.g.typical shape memory alloys or martensites in steels posses many more variants [37][38][39][40].In the present paper we demonstrate how the approach can be generalized to multi-variant martensites and show how the resulting variant distributions are influenced by different tensile strains and random grain rotations.

II. METHODS
The overall concept of the investigations is that we assume that in a microstructure, as depicted e.g. in Fig. 2, each (martensite) grain is allowed to select between different variants, such that the overall elastic energy is minimized in equilibrium.As the grains are coherently connected to each other, a variant flip as sketched in Fig. 1, leads possibly to an energetically unfavorable distortion of the surrounding grains, and therefore the discrete optimization of the variants is a nontrivial problem.
To achieve a description in the framework of continuum elasticity, we consider eigenstrain (stress free strain) expressions of the form Here, the function θ n (r) is equal to one inside grain n and vanishes everywhere else.The eigenstrain in each grain with index n can be a superposition using K "spin variables" s ν = s n,k = ±1 with the total number of spins N ν = N • K, which weight eigenstrain tensor components ϵ (0,n,k) ij . Therefore, we have altogether 2 K martensite variants in each grain by suitable selection of the parameters ϵ (0,n,k) ij .To determine the thermodynamic equilibrium state (minimum elastic energy), first the elastic problem has to be solved for a given eigenstrain.This can be expressed through the minimization of the elastic energy with respect to the displacement field u i in combination with suitable mechanical boundary conditions.In case of isotropic linear elasticity with phase and variant independent Lamé coefficient and shear modulus, the elastic energy functional reads with the system volume V and the strain tensor In the above expression, the Einstein sum convention is used.Contrary to the optimization of the displacement and strain fields as continuous variables, the variant selection becomes a discrete problem on top of the preceding step.Following the Fourier transformation approach outlined in [6], we can calculate the elastic energy E({s ν }) of the microstructure for a given "spin configuration" {s ν }.We note that the approach does not necessarily require to use Fourier transformation methods, but e.g. also finite element methods can be used, provided that the elastic interactions can be computed with sufficient accuracy, as discussed in [6].
As the expression for the elastic energy is quadratic in the framework of the linear theory of elasticity, the final expression for the elastic energy can be expressed as with a quadratic, symmetric matrix Jνη .The coefficients E 0 , Hν and Jνη can be expressed directly in terms of the Fourier transformation solution, or alternatively be calculated using any elastic solver, which delivers the total elastic energy E({s ν }) for a given spin configuration {s ν }.We note that the above form of the energy holds for arbitrary values of the spin variables and does not apply only to the special case s ν = ±1.Therefore, we obtain directly where all spins have the value 0, which we also denote as the austenite reference state.Next, we perform calculations where all spins but one vanish, and the selected one has either the value +1 or −1.From that we get Therefore, we can extract the diagonal element for the self-interactions as and consequently Finally, to calculate the interaction coefficients Jij for i ̸ = j we perform elastic calculations with two nonvanishing spins where we can directly calculate Jij using the expressions obtained above For zero stress boundary conditions, we get for the eigenstrain of the type (1) the simplification E 0 = 0 and Hi = 0 for all i.
For the use of the quantum annealer, we restrict the allowed spin values to s ν = ±1 and we also have to rewrite the expression (3) in a slightly different way, which excludes in particular self interactions of the type s i s i .The Hamiltonian reads where the offset H 0 only affects the energy value but not the spin configuration minimizing the energy.Comparison with the expression (3) above gives This general formulation has the advantage that it decouples the approach for setting up the Ising coefficients for the annealer from the specific linear elastic energy solver.Nevertheless, one has to keep in mind that the formulation (3) relies on the assumption of equal elastic constants in all phases and variants and coherent boundary conditions.However, it is valid also in situation beyond the original work [6], where we used only a single spin to discriminate between the two martensite variants in each grain.Here, the same approach also applies for multiple variants using several spins per grain.
On the quantum annealer level the determined Ising coefficients are used to define an energy landscape, where superconducting loops with clockwise or anticlockwise circulating currents define qubits with different spin states [16].These superconducting loops interact with external flux biases, which allow to control energy difference and barrier height of the constructed energy landscape [16].Starting point of the annealing process is the initialization of the system in the ground state of a known Hamiltonian H 0 ∼ − i σ x i where σ i denotes the Pauli matrices, i.e. corresponding to a strong transverse magnetic field [17,41].This Hamiltonian is turned during the annealing into the desired Ising model (10) for which an energetic minimum is sought, min {si=±1} H p [15].The Hamiltonians H 0 and H p do not commute [15], and the time taken by the initial Hamiltonian to reach the low energy state is sufficiently large to establish the validity of the adiabatic theorem of quantum mechanics [42], which postulates that a system stays in its eigenstate during adiabatic changes.In contrast to classical approaches, further quantum mechanical principles like tunneling to leave metastable regions and the usage of entangled states inside quantum annealing processors (QPU) are employed [43].Especially if energetically close low energy states exist, the QA process does not always determine the true ground state, therefore a suitable number of repetitions is performed and the lowest detected energy is chosen.Additionally, for larger systems, hybrid quantum annealing utilises classical algorithms and the interplay with quantum annealing to address areas with high computational demands using a QPU coprocessor working with generic parameters for up to 11616 spin variables on the D-Wave Advantage system [44][45][46].

III. RESULTS
To illustrate the workflow, we consider a case with four variants per grain, i.e. two spins per grain.Similar to the previous work [6] we use specific cases of shear and tetragonal distortions, which are now combined and lead to eigenstrains where n numerates the grains.For the first application we assume that the individual tensors of the shear (involving the strain parameter ϵ 0 ) and tetragonal deformations (associated with ϵ 1 ) are the same in all grains; later we will also discuss the case of grain rotations, where the tensors differ from grain to grain.Fig. 1 visualizes the different strain variants and illustrates the color codings used for the following plots.The grain structure is generated using a Voronoi tesselation in two dimensions, and we assume plane strain conditions, using isotropic linear theory of elasticity with periodic boundary conditions.According to the above approach, first the self energies related to each spin are calculated, followed by a computation of spin-spin interactions, from which the Ising vector components H n and interaction parameters J ij are extracted.
Apart from the self generated internal stresses due to the different variants, also a homogeneous external (average) strain field ⟨ϵ ij ⟩ can be superimposed, which leads to an additional contribution to the magnetic field vector H n , which couples the external strain ⟨ϵ ij ⟩ to the k = 0 mode of the Fourier transform of the eigenstrain tensor.This means that a change of the external strain leaves the coupling constants J ij invariant and only modifies a contribution to H n , which can be treated analytically, therefore minimizing the computational cost for repeated calculations with different external boundary conditions.Fig. 2 shows a sequence of equilibrium microstructures and their response to a tensile strain in horizontal (x) direction for ϵ 1 = ϵ 0 .Here, the computer generated, random 2D microstructure consists of 200 grains.Each grain selects its optimum variant state (s n,1 , s n,2 ) to minimize the total energy.In panel (a) the fixed external strain vanishes, which implies that in the ground state the volume fractions of the blue and green (resp.red and grey) variants should be close to each other, as then the positive and negative contributions to the shear and tetragonal (eigen)strain cancel each other to minimize the elastic energy.The arrangement of the differently coloured patches is then a result of the mutual spin-spin interactions.We can observe in panel (a) a structure of horizontal stripes with alternating layers of shear strain (blue-grey and red-green).Such stripes with (100) orientation have been discussed already in [6], as the formation of these structures minimizes the overall shear deforma- tion.Additionally, tilted stripe pairs red-blue and and green-grey correspond to the same tetragonal strain and annihilating shear strain, similar to the patterns for pure tetragonal eigenstrain discussed in [6].
The patterns formed by several variants can be better understood by analyzing the elastic anisotropy of the transformation.Specifically, the normal to the interface between two variants must be close to a direction that minimizes the elastic kernel of the transformation from one variant to the neighboring one.This can be demonstrated for plate-like domains [47], but is also very commonly observed in complex microstructures whose anisotropy is of elastic origin, such as cuboidal microstructures in superalloys [48], Widmanstätten structures [49] or chessboard structures in Co-Pt alloys [50].
The elastic kernel can be obtained analytically in Fourier space [47] or can be calculated in real space, using a regular stripe pattern of both variants.Indeed, in this geometry, only Fourier modes perpendicular to the stripes remain, and the elastic energy is proportional to the elastic kernel.In the following, the orientation dependence of the pattern is analyzed by considering successively each pair of variants.We assume fixed vanishing average strain and equal volume fractions of the two involved variants.The angular dependence of relevant combinations of stripe pairs is shown in Fig. 3.In these pole figures the angle 0 • corresponds to vertical stripes, and positive angles to anticlockwise rotations of the pattern.From this figure, one can conclude that both green-grey and red-blue interfaces are preferred either at 0 • or 90 • , and this is in agreement with the observed interface orientation between these variant pairs in Fig. 2(a).
Both green-red and grey-blue pairs have preferred orientations around ±38 • .Due to constraints of the chosen system sizes and the used periodic boundary conditions, typically deviations from the analytical prediction for infinite system sizes are expected, as discussed in detail in [6].In the equilibrium microstructure in Fig. 2(a) we indeed find inclined interfaces between green-red and grey-blue pairs, which are in line with the theoretical expectation.
Finally, green-blue and red-grey interfaces show a very pronounced anisotropy with a strong cusp at either 17 • and −62 • (green-blue) or 62 • and −17 • (red-grey).Due to this strong orientation dependence these interfaces either have to appear at these preferred angles or they are strongly supressed.In fact, the corresponding interface lengths in Fig. 2(a) are much shorter than the other interfaces.
We note for the subsequent panels of Fig. 2 that the elaborated anisotropy of the interactions is not affected by a homogeneous external strain, and therefore the polar plots in Fig. 3 remain valid also in these cases.In these pictures, the applied tensile external strain in horizontal direction is incrementally increased.This favors the variants with a positive "tetragonal spin value", and consequently the fraction of green and grey grains increases.They form more and more pronounced bands in the (100) directions, as explained above, partly interrupted by disappearing inclined red/blue bands with unfavorable tetragonal strain contribution.Finally, only grains with s n,2 = +1 remain, and the volume fractions with positive and negative shear are essentially the same, as much as the irregular grain sizes permit.
It is useful to emphasize the periodicity of the patterns by looking at a 2 × 2 supercell of the simulation domain, as shown in Fig. 4. The first striking observation is that the pattern in panel (c) for ⟨ϵ xx ⟩/ϵ 0 = 0.5 differs significantly from the other cases.These latter cases show both horizontal and tilted stripes, as analyzed above.There, grey-green and red-blue interfaces appear at 0 • or 90 • , and green-blue as well as red-grey interfaces play at most a minor role, as discussed above.The red-green and grey-blue interfaces are roughly at ±45 • as a compromise between the preferred analytical predictions for infinite systems and the constraints by the periodic boundary conditions.
The special case ⟨ϵ xx ⟩/ϵ 0 = 0.5 shows a remarkable topological change, as now interfaces between green-blue, green-red, grey-blue and grey-red appear, whereas the previously important red-blue and green-grey interfaces are absent.Moreover, the newly appearing interfaces are at their preferred orientations.For example, the marked 24 • angle inside the red lozenge is the difference between the 62 • minimum energy orientation for the red-grey interface and the 38 • orientation for the red-green pattern.This leads to the geometry sketched in Fig. 4(c), from which the overall periodic pattern is constructed.According to the equilibrium angles, this demands that the joint volume fraction of the diamond-shaped green and grey domains is sin(80 • )/[sin(24 • ) + sin(80 • )] ≈ 0.7, and therefore this exceptional pattern is only observed close to a specific external strain value.We note that this chessboard structure is very similar to the one observed in Co-Pt and (CuAu) 1−x -Pt x systems [50].However, the symmetry is lower in the present case because the green and grey regions are not tetragonal.
The same microstructure as in Fig. 2 is subjected to a shear strain ⟨ϵ xy ⟩ in (110) direction, and the results of the energy minimization are shown in Fig. 5.The starting configuration without external strain (panel (a)) is therefore identical to Fig. 2(a).An increase of the shear strain favors variants with positive s n,1 , as they accommodate the given external strain to a large amount.Therefore, for large positive shear strain ⟨ϵ xy ⟩/ϵ 0 values only grey and blue patches remain.As they differ by the tetragonal eigenstrain only, they align as inclined bands to minimize the total elastic energy, as discussed in detail in [6].In particular, grey-green and red-blue interfaces are at 0 • or 90 • , as predicted in Fig. 3. Similarly, blue-green and red-grey interfaces have orientations close to the expected ones (±17 • , ±62 • ).For blue-grey and red-green interfaces the periodicity constraint enforces orientations close to 45 • .
An interesting case is the three-variant configuration in Fig. 5(d) (apart from a negligible fraction of the red variant).A 2 × 2 supercell of the equilibrium microstructure is shown in Fig. 6.The overall 45 • alignment is imposed by periodic boundary conditions.The green domain slightly changes the orientation of the blue-grey interfaces so that they are rotated closer to the optimal orientation 90 • − 38 • = 52 • .This configuration would not be favorable at −45 • (instead of the 45 • ) because the blue-grey interfaces would be rotated away from the optimal orientation.Instead, for −45 • , red domains appear in the blue platelets (see Fig. 5(e)).
For the present four variant setup, the response of the microstructure does not only depend on the ratio ⟨ϵ ij ⟩/ϵ 0 , but also on the relative relevance of shear and tetragonal deformations, as expressed through the dimensionless parameter ϵ 1 /ϵ 0 .Whereas in the above investigations we have fixed this ratio to ϵ 1 /ϵ 0 = 1, Fig. 7 shows the influence of this parameter for fixed external strain.Here, we note that a change of ϵ 0 → aϵ 0 and ϵ 1 → bϵ 1 does not require to repeat the entire computations for the interaction parameters J ij and the external magnetic field H i .Instead, it follows from the theory of linear elasticity that the coupling constants J ij are quadratic in the rescaling factors a, b, and H i scales linearly with them.In detail, for a spin pair (i, j), where both of them refer to shear transformations, the coupling constants are rescaled according to J ij → a 2 J ij , and similarly for two "tetragonal spins" as J ij → b 2 J ij .For mixed interactions between a tetragonal and a shear degree of freedom, the rescaling obeys J ij → abJ ij .The "magnetic" terms H n scale linearly with either a or b, depending on the spin type.
The starting case (center panel in Fig. 7) is again the already discussed reference situation with ϵ 1 = ϵ 0 .In all three pictures, the boundary conditions are fixed to vanishing average strain, ⟨ϵ ij ⟩ = 0.The right panel shows the equilibrium microstructure for ϵ 1 /ϵ 0 = 5, which means that the weight of the tetragonal deformations is significantly larger than that of the shear.Here, the variants select essentially diagonal bands, as discussed in [6], in order to match the periodic boundary conditions.The blue and green variants have cancelling tetragonal and shear eigenstrains, such that for equal volume fractions the zero strain boundary conditions are met best.We can expect that also the opposite colour pair (red-grey) can lead to an energetically equivalent solution, and therefore a weak, locally "noisy" pattern, which mixes these two configurations, may appear.We expect this effect to be more pronounced in larger systems.One may expect that the opposite case with ϵ 1 /ϵ 0 = 0.2, which is dominated by the shear eigenstrain, should similarly lead to bands with alternating shear in ⟨100⟩ direction, as discussed for the two-variant case in [6].However, as can be seen from the left panel of Fig. 7, this is not the case.Instead, each horizontal stripe consists of blue and grey (or red and green) patches, which share the same shear eigenstrain but have opposite tetragonal strain.The reason for this assymetry -as compared to the tetragonal dominated case discussed above -is that a small, but nonvanishing tetragonal eigenstrain leads also to a tangential mismatch strain at the interfaces between the ⟨100⟩ bands.Therefore, a purely blue horizontal band would lead to a contraction in tangential direction, whereas the compensating green horizontal band would expand in the same direction, hence leading to large mismatch coherency stresses.Therefore, it is energetically favorable to generate ⟨100⟩ structures with an additional zigzag substructure, such that all four variants are present with equal volume fraction, and only the irregularity of the microstructure adds some minor noise.
It is again instructive to look at the angular depen-  dence of the elastic energy for infinite stripe patterns (again for zero mean strain boundary conditions), as it depends on the ratio ϵ 1 /ϵ 0 , see Fig. 8. Here, one can clearly see that the interaction becomes very anisotropic for green-grey interfaces at low ratios ϵ 1 /ϵ 0 and essentially isotropic at high values, with minimum energetic cost at 0 • and 90 • .For the green-red interfaces, the trend is opposite.Green-red (and similarly grey-blue) interfaces have preferred interface orientations ±38 • .The green-blue and similarly grey-red interfaces have the strongest dependence on ϵ 1 /ϵ 0 and therefore the overall energy minimization clearly favors the optimal orientations for these interfaces.From this perspective, Fig. 7(a) is particularly interesting, and its zigzag structure becomes obvious in the 2 × 2 supercell representation in Fig. 9.The white lines emphasize the structure and illustrate that for this microstructure all interface orientations are optimal and not affected by periodicity constraints.
The last example concerns a case where the grains are oriented randomly, i.e. rotated around the [001] axis, as depicted in the left panel of Fig. 10.Like in the previous publication [6], the response to a uniaxial strain ⟨ϵ xx ⟩ is less pronounced than in the cases with uniform grain orientation above.The center panel shows the equilibrium microstructure for vanishing mean strain, whereas the right panel is for a rather high tensile strain ⟨ϵ xx ⟩/ϵ 1 = 2.5.
Therefore, for the chosen example, a rotation by π/4 converts the tetragonal deformation to a shear transformation (and vice versa), apart from the eigenstrain contribution in z direction, which only affects the energetic cost but not the strain distribution.Therefore, we expect for a tensile load in [100] direction, that grains with orientation close to 0 or π (the response is π-periodic), to have a tetragonal spin preferentially being in state +1; for a rotation near 90 • , the tetragonal contribution prefers the s = −1 state.For both of these orientations (and nearby values), the shear contributions are indifferent, as they have only a vanishing or small contribution ϵ (0) ′ xx , which is sensitive to the tensile strain.As explained above, for rotations by approximately π/4 (or 3π/4) the role of shear and tetragonal deformations is exchanged, but the original tetragonal transformation has an energy penalty due to the contribution in z direction for the used plane strain setup.
For some selected grains in Fig. 10 the change of the variant states is discussed in the following.Grain 1 has an orientation close to 0 (or π).In the zero strain state, it is marked in red, hence the shear spin equals −1 and the tetragonal component is also −1.Under stronger tensile load it turns grey, hence both spins flipped their sign, leading to a stronger alignment with the external field.Grain 2 has an orientation close to 90 • , is originally red (shear: -1, tetragonal -1) and becomes blue in the right panel.In this case, the dominant effect is to have a negative tetragonal spin values, as discussed above.Finally, for grain 3, which has a mixed orientation, we have effec-tively a transition from blue to grey, where the positive shear spin is relevant due to the orientation near π/4, whereas the flipping tetragonal spin component is less relevant.

IV. CONCLUSIONS
From a materials science perspective, quantum quantum is only starting to influence the opportunities and methodologies of modeling approaches.It can be expected, that in future general purpose quantum computing will have the potential to strongly accelerate materials science related simulations and to enable completely new possibilities.As of today, sufficiently large quantum computers are not yet available, and also suitable algorithms have not yet been developed.Nevertheless, quantum annealing, as a specific type of adiabatic quantum computing, is commercially available with large numbers of qubits already today.The type of problems which can be addressed is limited, as it is required to map the materials science related problem of interest to the minimization of an Ising Hamiltonian.
In the present work we have demonstrated, that the selection of multiple martensite variants can be mapped to an Ising problem for specific conditions and have demonstrated how the energy minimization can be performed efficiently via quantum annealing.The superiority of QA compared to classical computing has been demonstrated in our previous publication [6].The multiplicity of the domains together with the long-range, anisotropic and shape-dependant character of the elastic interactions make the determination of the ground state a particularly difficult problem, and therefore QA opens completely new possibilities.The central result of the present publication is the generalization to multi-variant systems and the analysis of the arising microstructures.Also, we demonstrate how the proper coefficients of the Ising model can be obtained from general linear elastic calculations.To illustrate the approach, we use a four variant per grain case, where each grain is represented by two Ising spins.We use a linear combination of shear and tetrag-onal deformations to demonstrate the methodology.As a result, the equilibrium microstructure depends on the applied external strains (or stresses), as well as on the strength of the distortions induced by the different variants.The observed patterns can largely be explained by computations of the anisotropy of laminar arrangements of the different variants, showing that the mismatch of neighboring domains dominates the overall elastic energy.
For many applications in science, the formation of microstructures is critical and efficient simulation approaches to support experimental findings are highly desired.Our developed quantum annealing approach presents an efficient new technological opportunity for the determination of equilibrium microstructures with long-range elastic interactions, where multiple martensite variants need to be considered.

FIG. 1 .
FIG.1.Visualization of the variants used in the present study.Spin s1 corresponds to the shear deformation, s2 to the tetragonal deformation, see Eq. (14).

FIG. 3 .
FIG. 3. Orientation dependence of the normalized elastic energy for equal volume fraction stripes consisting of pairs of variants for ⟨ϵij⟩ = 0 and ϵ1 = ϵ0.

FIG. 6 . 2 ×
FIG. 6. 2 × 2 supercell of the simulation domain for external shear ⟨ϵxy⟩/ϵ0 = 0.75 from Fig. 5(d) with ϵ1 = ϵ0.Apart from a negligible amount of the red variant three phases are present here.The sketch shows the geometrical explanation of the equilibrium pattern.

FIG. 8 .
FIG. 8. Polar representation of the orientation dependence of the elastic energy (scaled to the maximum value) for mismatching pairs of variants.The calculations are for fixed strain boundary conditions ⟨ϵij⟩ = 0.