N Amadou

T De Resseguier

A Dragon

E Brambrink

Effects of orientation, lattice defects and temperature on plasticity and phase transition in ramp-compressed single crystal iron

Keywords: iron, plasticity, phase transition, twins, dislocations, nucleation and growth, MD simulations, ramp compression, shock wave, microvoid defects

Effects of orientation, lattice defects and temperature on plasticity and phase transition in ramp-compressed single crystal iron

Nourou Amadou, Thibaut de Resseguier, A. Dragon, E.

Introduction

The dynamic behavior of materials is an area of study at the confluence of many scientific disciplines [START_REF] Meyers | Dynamic behavior of materials[END_REF][START_REF] Remington | From microjoules to megajoules and kilobars to gigabars: Probing matter at extreme states of deformation[END_REF], motivated by both fundamental and industrial applications. Of special importance is understanding of (i) the deformation (specifically plasticity) and phase transitions mechanisms of solid matter under extreme conditions of pressure, temperature and strain rate, (ii) how deformation path and mode affect the structural phase transformations and (iii) how initial microstructure and temperature influence the overall deformation process [START_REF] Duvall | Phase transitions under shock-wave loading[END_REF][START_REF] Remington | Materials science under extreme conditions of pressure and strain rate[END_REF][START_REF] Remington | Material dynamics under extreme conditions of pressure and strain rate[END_REF][START_REF] Lorenzana | Shocked materials at the intersection of experiment and simulation[END_REF][START_REF] Rudd | Metal deformation and phase transitions at extremely high strain rates[END_REF][START_REF] Remington | From microjoules to megajoules and kilobars to gigabars: Probing matter at extreme states of deformation[END_REF].

During dynamic compression by a planar wave front, material is uniaxially compressed which involves a large deviatoric stress [START_REF] Meyers | Defect generation in shockwave deformation[END_REF]. In practice, such compression is achieved by shock loading [9,[START_REF] Duvall | Phase transitions under shock-wave loading[END_REF],

where strain rates can reach 10 9 s -1 [START_REF] Smith | High strain-rate plastic flow in Al and Fe[END_REF]. However, shock wave is a highly dissipative process, which causes an entropy jump and high temperatures so that materials may melt [START_REF] Remington | Materials science under extreme conditions of pressure and strain rate[END_REF]. In order to preserve the solid state, the so-called "quasi-isentropic" ramp compression (where the pressure pulse is shaped to rise more progressively over time) is alternatively used [START_REF] Asay | The use of shock-structure methods for evaluating high-pressure material properties[END_REF]. In such approach, strain rates greater than 10 6 s -1 are currently available [START_REF] Edwards | Laserdriven plasma loader for shockless compression and acceleration of samples in the solid state[END_REF][START_REF] Smith | High strain-rate plastic flow in Al and Fe[END_REF][START_REF] Amadou | Direct laser-driven ramp compression studies of iron: A first step toward the reproduction of planetary core conditions[END_REF][START_REF] Smith | Time-dependence of the alpha to epsilon phase transformation in iron[END_REF][START_REF] Amadou | Kinetics of the iron α -ε phase transition at high-strain rates: Experiment and model[END_REF][START_REF] Amadou | Laser-driven ramp compression to investigate and model dynamic response of iron at high strain rates[END_REF]. They allow studying the kinetics of both plasticity and phase transitions in solid state [START_REF] Rudd | Metal deformation and phase transitions at extremely high strain rates[END_REF][START_REF] Smith | Time-dependence of the alpha to epsilon phase transformation in iron[END_REF][START_REF] Amadou | Laser-driven ramp compression to investigate and model dynamic response of iron at high strain rates[END_REF]. At moderate deformation rates, it is the emission, multiplication and propagation of dislocations that are responsible for plastic deformation where the process can be either thermally activated or limited by dissipation involving phonon emission and scattering [START_REF] Remington | Materials science under extreme conditions of pressure and strain rate[END_REF][START_REF] Remington | Material dynamics under extreme conditions of pressure and strain rate[END_REF]. At high deformation rates, it is rather twinning that dominates and accommodates the high shear stress more rapidly [START_REF] Rudd | Metal deformation and phase transitions at extremely high strain rates[END_REF]. Furthermore, although isokinetics regime was reported for the iron bcc-hcp phase transformation in a narrow range of deformation rates [START_REF] Amadou | Kinetics of the iron α -ε phase transition at high-strain rates: Experiment and model[END_REF], the kinetics of structural transformations remains largely dependent on the strain rate [START_REF] Smith | Time-dependence of the alpha to epsilon phase transformation in iron[END_REF][START_REF] Amadou | Laser-driven ramp compression to investigate and model dynamic response of iron at high strain rates[END_REF]. Understanding these various physical phenomena motivates meso-scale energy-based analyses [START_REF] Vattré | Polymorphism of iron at high pressure: A 3d phase-field model for displacive transitions with finite elastoplastic deformations[END_REF] and purposeful atomic level investigations wherein Molecular Dynamics (MD) [START_REF] Holian | Plasticity induced by shock waves in nonequilibrium moleculardynamics simulations[END_REF][START_REF] Germann | Orientation dependence in molecular dynamics simulations of shocked single crystals[END_REF][START_REF] Kadau | Microscopic view of structural phase transitions induced by shock waves[END_REF][START_REF] Bringa | Shock deformation of face-centred-cubic metals on subnanosecond timescales[END_REF].

In this context, using the recently developed modified version of the Ackland iron potential [START_REF] Gunkelmann | Polycrystalline iron under compression: Plasticity and phase transitions[END_REF], we have reported in a previous paper [START_REF] Amadou | Coupling between plasticity and phase transition in shock-and ramp-compressed single-crystal iron[END_REF] that defect-free single crystal iron of 50 K initial temperature subjected to ramp compression along the [001] direction first yields by twinning at a longitudinal pressure of ∼ 12 GPa, then upon further compression, the twins recede progressively while no dislocation activity is detected. This leads to a new increase of the deviatoric stress to values as high as 30 GPa.

Such elastic stiffening of the bcc matrix is shown to inhibit the nucleation of the hcp phase, so that the onset of the phase transformation is shifted to very high pressures (on the order of 100 GPa). In the present work, we go further to explore both deformation and bcc-hcp structural phase transition and their coupling under ramp compression in other low-index crystallographic directions, in crystals with or without microvoid defects, with different initial temperature conditions. The deformation behavior is found to be highly anisotropic in both defect-free and defective single crystal. At low temperature, plastic deformation, strongly enhanced by the presence of initial defects, is mainly dominated either by twinning alone or twinning and dislocations activities in different proportions, depending on loading direction with respect to crystal orientation. For initial temperature above 200 K, twins are observed to disappear in all cases while dislocations are not affected over the explored temperature range. The corresponding deformation anisotropy effects and initial conditions strongly affect subsequent bcc-hcp phase transition.

Method and computational details

Non Equilibrium Molecular Dynamics (NEMD) simulations were used to study deformation and bcchcp phase transition in iron crystal under planar ramp compression applied along the [001], [110] and [111] directions respectively. Samples with up to 28 millions atoms, 28.7nm × 28.7 nm (100 a 0 × 100a 0 ) cross section and length up to 402 nm (1400 a 0 ), either defect-free or containing 1.45 nm-radius microvoids (5 a 0 , with a 0 the lattice constant) regularly spaced by 28.7 nm (100 a 0 ) along the wave propagation direction were simulated. Furthermore, the initial temperature was varied from 50 K to 400K.

The simulations were realized using the Lammps molecular dynamics code [START_REF] Plimpton | Fast parallel algorithms for short-range molecular dynamics[END_REF] where interactions between atoms were modeled through the embedded atom model formalism [START_REF] Daw | Semiempirical, quantum mechanical calculation of hydrogen embrittlement in metals[END_REF][START_REF] Foiles | Embedded-atom-method functions for the fcc metals Cu, Ag, Au, Ni, Pd, Pt, and their alloys[END_REF] expressed tangibly in the modified version of the Ackland potential [START_REF] Gunkelmann | Polycrystalline iron under compression: Plasticity and phase transitions[END_REF]. This potential predicts the iron bcc to hcp phase transition at 13.75 GPa under equilibrium hydrodynamic compression [START_REF] Gunkelmann | Polycrystalline iron under compression: Plasticity and phase transitions[END_REF], in very good agreement with the experimental data [START_REF] Bancroft | Polymorphism of iron at high pressure[END_REF][START_REF] Barker | Shock wave study of the α -ε transition in iron[END_REF][START_REF] Boettger | Metastability and dynamics of the shock-induced phase transition in iron[END_REF]. It has also been successfully used in previous work to address plasticity and phase transition in [001]-oriented iron single crystal [START_REF] Amadou | Coupling between plasticity and phase transition in shock-and ramp-compressed single-crystal iron[END_REF] and polycrystalline iron [START_REF] Gunkelmann | Shock waves in polycrystalline iron: Plasticity and phase transitions[END_REF][START_REF] Gunkelmann | Interplay of plasticity and phase transformation in shock wave propagation in nanocrystalline iron[END_REF] under shock compression.

In the present work, the ramp wave is initiated by driving an effective infinite-mass wall piston with an imposed velocity v(t) along the z-direction, while periodic boundary conditions are applied in the transverse directions. The piston velocity increases linearly from 0 to Vmax within 30 ps, with Vmax ranging from 800 to 1600 m/s. The corresponding strain rates are about 10 8 -10 9 s -1 comparable with that reported under picosecond-laser dynamic compression [START_REF] Crowhurst | The α -ε phase transition in iron at strain rates up to 10 9 s -1[END_REF]. Local thermodynamic and mechanical variables such as longitudinal pressure P z , deviatoric stress, temperature, etc. are evaluated within a spatial planar bin (of 3 lattice constant width) perpendicular to the wave propagation direction in the same manner as in our previous work [START_REF] Amadou | Coupling between plasticity and phase transition in shock-and ramp-compressed single-crystal iron[END_REF]. Finally, local structural analysis is performed by adaptative Common Neighbor Analysis (CNA), Centrosymmetry and DXA, as implemented in the OVITO software [START_REF] Stukowski | Visualization and analysis of atomistic simulation data with OVITO-the open visualization tool[END_REF].

Results and Discussion

Load direction

Depending on the direction of ramp compression, defect-free iron single crystal deformation process is found to be very different leading to a directional anisotropy in the bcc-hcp structural phase transformation features. Fig. 1 shows longitudinal pressure P z and shear stress profiles at two simulation times: 30 ps when the piston reached its maximum velocity and 40 ps. Corresponding atomic configurations are shown in Fig. 2. For all propagation directions the sample first deforms elastically. The velocity of the elastic wave is found to depend on the loading direction with respect to crystal direction in consistency with previous MD simulations [START_REF] Kadau | Atomistic simulations of shockinduced transformations and their orientation dependence in bcc Fe single crystals[END_REF][START_REF] Wang | An atomic study on the shock-induced plasticity and phase transition for iron-based single crystals[END_REF]. In contrast with the [001] direction, in both [110] and [111] directions the elastic precursor is found to become steeper upon propagation thus leading to a sharp front at the foot of the compression wave. While no plasticity is clearly evident in the [111] direction before the structural phase transition, both [001] and [110] directions are found to exhibit plastic deformation.

In [001] direction, iron is found to yield via twinning at P z around 12 GPa (see Fig. 3) with subsequent relaxation of the shear stress (see Fig. 2 and Fig. 1). Then, further compression is observed and fcc atoms respectively. Under compression along the [001] direction, twins accommodate plastic deformation above the elastic limit, with no phase transition. On the contrary, compression in the [111] direction remains elastic (no twins nor dislocations) up to the bcc-hcp phase transition. Finally, in the [110] direction, yield occurs via simultaneous twins formation and hcp nucleation, then the compressive wave splits into a twinning front followed by a metastable phase mixture where dislocations develop (see Fig. 4). Next, the bcc-hcp transition evolves to completion. The elastic wave is overdriven for inhibits the nucleation of the hcp phase, so that the onset of the phase transformation in the defectfree crystal is shifted to very high pressures. Thus, for the transition to be observed a pressure as high as 100 GPa is required (see [START_REF] Amadou | Coupling between plasticity and phase transition in shock-and ramp-compressed single-crystal iron[END_REF] for more details) which has been qualitatively explained by invoking phase nucleation and growth theory [START_REF] Christian | The Theory of Transformations in Metals and Alloys[END_REF][START_REF] Smallman | Modern Physical Metallurgy[END_REF] and the fact that a transformation from the elastic state would be expected to result in a significantly slower transformation kinetics [START_REF] Knudson | Transformation kinetics for the shock wave induced phase transition in cadmium sulfide crystals[END_REF]. In the present case, where maximum pressure remains below 50 GPa, no transition is observed (Fig. 2).

On the other hand, for [110] direction, the defect-free single crystal is found to remain perfectly elastic, with the shear stress linearly increasing, until the pressure reaches about 27 GPa (Fig. 1).

Then, yielding occurs via simultaneous nucleation of twins and hcp phase (only few fcc atoms are observed at this stage) inducing a rapid relaxation in the shear stress. Upon propagation, twinning and hcp homogeneous nucleation fronts split [START_REF] Duvall | Phase transitions under shock-wave loading[END_REF], (see Fig. 2). Like in the [001] direction, further compression induces a new increase of the shear stress indicative of some hardening-like effect. During this stage where the shear stress rises, the number of hcp atoms (and fcc) in the sample is observed to remain almost constant leading to a metastable mixed regime where coexist mainly bcc and hcp phases. The rapid growth stage of the bcc-hcp phase transition is thus delayed due to the stiffening effect. Meanwhile, when the pressure reaches a value around 40 GPa, a second yielding mechanism is triggered, namely homogeneous nucleation of dislocations followed by their rapid growth (see Fig. 4) producing a new shear stress relaxation (see Fig. 1). Here, the maximum density of dislocations is about 1.6 × 10 12 cm -2 . The dislocations growth is accompanied by increasing number of fcc atoms while the number of hcp atoms decreases slightly until the shear stress attains 10 GPa. Only then, the hcp rapid growth regime is found to be triggered thus leading the structural phase transformation toward completion. The higher the piston velocity, the shorter structural phase transition time to completion, in agreement with experimental report in Cadmium sulfide (CdS) where plasticity was observed to play an important role in the structural phase transition nucleation and growth [START_REF] Knudson | Transformation kinetics for the shock wave induced phase transition in cadmium sulfide crystals[END_REF].

Thus, upon ramp loading along the [110] orientation, a double yielding is observed before the structural phase transition reaches its rapid growth regime. This behavior is consistent with the two-stage plasticity reported in other bcc single crystals dynamically compressed along their [110] direction, such as Mo in MD simulations combined with discrete to continuum method [START_REF] Kositski | Shear relaxation behind the shock front in 110 molybdenum -from the atomic scale to continuous dislocation fields[END_REF] or MgO in shock compression experiments [START_REF] Liu | Anomalous elastic-plastic transition of MgO under shock compression[END_REF].

The above complex deformation behavior contrasts with the response to compression along the [111] direction. Here, no sign of plastic deformation (twins or dislocations) is observed before the structural phase transformation, which is found to be the only way to relax shear stress in the defect-free crystal during the compression. This is consistent is with previous MD simulations of the iron single crystal response to both shock compression [START_REF] Kadau | Atomistic simulations of shockinduced transformations and their orientation dependence in bcc Fe single crystals[END_REF][START_REF] Wang | An atomic study on the shock-induced plasticity and phase transition for iron-based single crystals[END_REF] and non-equilibrium uniaxial deformation [START_REF] Wang | Molecular dynamics simulations of hcp/fcc nucleation and growth in bcc iron driven by uniaxial compression[END_REF] along [111] direction. When the longitudinal pressure P z reaches about 30 GPa some atoms are observed to start moving transversally to the compression direction, leading to a plateau in the pressure profile (see 1).

As a result, the shear stress increase also slows down. Then, the transition occurs when the pressure P z reaches 35 GPa with sharp drop in the shear stress down to around 4 GPa. The duration of this regime which could be identified with transition incubation period was observed to depend on the piston maximum velocity, i.e the maximum pressure [START_REF] Knudson | Transformation kinetics for the shock wave induced phase transition in cadmium sulfide crystals[END_REF]. The higher the piston maximum velocity the shorter the plateau duration. So, neither a hardening-like effect nor a mixed regime are observed in contrast to the other investigated directions.

Significant fractions of fcc atoms are found by the CNA analysis (see Fig. 2) in both [110] (hcp/fcc ratio ≈ 2) and [111] (hcp/fcc ratio ≈ 2.7). Although such hcp-fcc mixture was already predicted in MD simulations based on the same iron potential [START_REF] Gunkelmann | Shock waves in polycrystalline iron: Plasticity and phase transitions[END_REF], as well as by other potentials [START_REF] Kadau | Atomistic simulations of shockinduced transformations and their orientation dependence in bcc Fe single crystals[END_REF][START_REF] Wang | An atomic study on the shock-induced plasticity and phase transition for iron-based single crystals[END_REF], it has not been experimentally evidenced yet, e.g. [START_REF] Hawreliak | In situ x-ray diffraction measurements of the c/a ratio in the high-pressure phase of shock-compressed polycrystalline iron[END_REF]. Still, detecting such mixed hcp-fcc phases would require accurate probing of the local crystal structure at the atomic level over ultra-short time scales, which remains very challenging, especially since it was inferred from MD data, that the fcc phase found by CNA analysis does not result in any prominent feature in the X-ray diffraction pattern [START_REF] Gunkelmann | Shock waves in polycrystalline iron: Plasticity and phase transitions[END_REF].

Lattice defects

The presence of initial defects in a material is known to be an important factor that determines its mechanical response [START_REF] Zepeda-Ruiz | Probing the limits of metal plasticity with molecular dynamics simulations[END_REF]. Real materials contain such defects as processing and handling may intrinsically introduce inclusions, voids or imperfections even in single crystal [START_REF] Hawreliak | Analysis of the x-ray diffraction signal for the α -ε transition in shock-compressed iron: Simulation and experiment[END_REF]. Following [START_REF] Cui | Phase transformation of iron under shock compression: Effects of voids and shear stress[END_REF], in order to investigate the influence of such defective initial conditions on both plasticity and phase transition microvoids of 5a 0 -radius (with a 0 the lattice constant) were introduced every 100a 0 along Unlike in previous MD simulations of defective iron single crystal with Voter-Chen interatomic potential [START_REF] Cui | Phase transformation of iron under shock compression: Effects of voids and shear stress[END_REF][START_REF] Shao | Microscopic dynamics of structural transition in iron with a nanovoid under shock loading[END_REF], the emission of dislocations from microvoids upon further compression is observed (around a local P z of 17 GPa). This behavior is consistent with that reported for other bcc materials such as tantalum [START_REF] Ruestes | Plastic deformation of a porous bcc metal containing nanometer sized voids[END_REF]. However the growth of these dislocations is found to be very slow. Even longer lasting simulation over 200 ps showed neither accelerated growth nor homogeneous nucleation of dislocations in the [001]-oriented crystal. Nevertheless, in the presence of microvoids and subsequent emitted dislocations, no shifting of the phase transition to a high pressure was observed, which is consistent with the theory of heterogeneous nucleation [START_REF] Christian | The Theory of Transformations in Metals and Alloys[END_REF][START_REF] Smallman | Modern Physical Metallurgy[END_REF][START_REF] Porter | Phase Transformations in Metals and Alloys[END_REF]. In the context of this theory, the presence of defects allows the nuclei of new phase [START_REF] Avrami | Kinetics of phase change. I general theory[END_REF] to develop in the parent phase near defects where the work of nucleation is much reduced [START_REF] Christian | The Theory of Transformations in Metals and Alloys[END_REF]. So, in contrast with the defect-free crystal, despite some elastic stiffening of the bcc matrix, the transition is made locally possible by an important reduction in its free energy barrier. Thus, the structural phase transition is now observed to start near voids when the local P z reaches ∼ 25 GPa, which is very low compared to the 100 GPa in the case of defect-free crystal [START_REF] Amadou | Coupling between plasticity and phase transition in shock-and ramp-compressed single-crystal iron[END_REF]. Then, it grows for about 20 ps leading to a kind of plateau in the pressure profile (at the top of the plastic wave) at around 30 GPa where shear stress starts to relax (see Fig. 5). Upon completion of the structural phase transformation, the shear stress is completely relaxed and even over-relaxed [START_REF] Kadau | Atomistic simulations of shockinduced transformations and their orientation dependence in bcc Fe single crystals[END_REF][START_REF] Kadau | Microscopic view of structural phase transitions induced by shock waves[END_REF][START_REF] Wang | Molecular dynamics simulations of hcp/fcc nucleation and growth in bcc iron driven by uniaxial compression[END_REF]. This illustrates the importance of defects in the development of martensitic phase transformation.

Upon compression along the [110] direction with defective initial conditions, dislocations loops

are first observed to nucleate from microvoids when the P z reaches locally about 20 GPa. At this stage, neither twinning nor hcp nucleation occurred in the sample. Then the dislocation loops expand through the sample by stretching (see Fig. 7), occasioning a first shear stress relaxation (see Fig. 5 at 30 ps). Here, the pressure at the top of the elastic wave associated with the yielding pressure, after a slight decrease during wave propagation, is found to be established around 20 GPa, i.e at the local pressure value of dislocations emission. This value is lower than in the defect-free case (≈ 27 GPa) in agreement with the fact that without initial defect, crystal strength is over-predicted [START_REF] Zepeda-Ruiz | Probing the limits of metal plasticity with molecular dynamics simulations[END_REF]. Upon further compression, hcp phase is found to nucleate either on the dislocation loop or from the voids. While twinning generally accompanyies the hcp phase nucleation, twinning alone is also found to develop directly on the dislocations loops as it was reported in refs [START_REF] Marian | Dynamic transitions from smooth to rough to twinning in dislocation motion[END_REF][START_REF] Zepeda-Ruiz | Probing the limits of metal plasticity with molecular dynamics simulations[END_REF]. Next, the deformation process becomes similar to that of the defect-free crystal with the formation of a mixed regime and hardening like behavior until homogeneous nucleation of dislocations at pressure P z of about 38 GPa. Then, the shear stress relaxes again (see Fig. 5, profiles at 30 ps). Later when the sample is homogeneously filled by dislocations, the shear stress decreases upon further compression with no double yielding (see Fig. 5, profiles at 40 ps). The pressure at the front of the mixed regime that we consider as the onset pressure for the phase transition is around 22 GPa, which is almost the same as the yielding pressure.

Finally, in the [111] oriented crystal with defective initial configuration, dislocations emission starts above 20 GPa, which is comparable to that in the [110] direction. Then dislocations grow within less than 10 ps. An important pressure relaxation precedes the structural phase transition instead of pressure stagnation. Although twinning is found to accompany the dislocation development, no massive homogeneous nucleation of hcp phase is observed (see Fig. 7). In this case, the bcc-hcp phase transition is found to start at pressure P z of about 33 GPa, then grow whithin about 10 ps kinetics time scale up to completion around 40 GPa. So, while with defective initial conditions the deformation remains mainly dominated by twinning when compression is applied along the [001] direction, dislocations mediated deformation is the favored way of shear stress relaxation for the crystal loaded along both [110] and [111] directions. Fig. 8 shows the nucleation and growth of dislocations from microvoids for compression along the [110] direction. Mainly 1/2 < 111 > Burgers vector type dislocations are nucleated in the form of loops.

Each loop consists of two screw dislocations (with opposite Burgers vectors) moving apart, separated by a leading curved non-screw line (mixed plus edge types) that propagates away from the void, in consistence with experimental observation in pure iron under in situ straining experiments [START_REF]On the stress discrepancy at low-temperatures in pure iron[END_REF][START_REF] Caillard | Tem in situ straining experiments in fe at low temperature[END_REF].

When, the leading non-screw dislocation part meets another dislocations loop (or a microvoid), it slows down.

The interactions between different loops lead to dislocations multiplication that establishes a dense network in the sample (see Fig. 7 (b) and(c)). Here, the maximum dislocation density is found to be at around 1.4 × 10 12 cm -2 almost identic to the defect-free case. Furthermore, the, screw dislocations are found to move with average velocities around 397 m/s while the leading edge component velocity is found to be approximately 2600 m/s. Thus, the edge component moves faster than the screw one, which agrees with experimental observations [START_REF]On the stress discrepancy at low-temperatures in pure iron[END_REF][START_REF] Caillard | Tem in situ straining experiments in fe at low temperature[END_REF]. A similar behavior is observed for compression along [111]-direction. Here, dislocation density is about 4.3 × 10 12 cm -2 while the screw and edge component velocities are found to be 585 m/s and 2400 m/s , respectively. In contrast, when compression is applied along the [001] direction, dislocations are mostly of edge type (again with essentially 1/2 < 111 > Burgers vector), and their velocity is too slow to be determined (see Fig. 7(a)).

As the ramp compression is accompanied by high shear stress for both [110] and [111] directions (see Fig. 5), dislocations must flow faster in order to allow the material to deform and relieve the high shear stress on the rapid time scale of the compression [START_REF] Rudd | Metal deformation and phase transitions at extremely high strain rates[END_REF]. Thus the dislocation velocities found here are higher than that found under quasi-static loading [START_REF]On the stress discrepancy at low-temperatures in pure iron[END_REF]. However, these velocities are qualitatively consistent with previous MD simulations in iron [START_REF] Marian | Dynamic transitions from smooth to rough to twinning in dislocation motion[END_REF] and other bcc material such as tantalum in both MD [START_REF]Growth and collapse of nanovoids in tantalum monocrystals[END_REF][START_REF]Probing the character of ultra-fast dislocations[END_REF] and stochastic [START_REF]Stochastic simulation of dislocation glide in tantalum and ta-based alloys[END_REF] simulations.

Initial temperature effect

In order to investigate the effects of initial temperature on the deformation and structural phase transition, simulations (with or without defects) were done where the lattice equilibration temperature was varied from 50 to 400 K before ramp compression. Upon elevation of the initial temperature from 50 K to 300 K, the dislocations nucleation and growth are observed to remain similar, i.e. emission of mainly 1/2 < 111 > dislocation loops composed with two screw dislocation lines connected by a curved non-screw leading part for both [110] and [111] directions. However, the dislocation densities and velocities are found to increase slightly. Indeed, for compression along [110]-direction ([111]direction) with an initial temperature of 300K, the dislocation density is found to be 2.7 × 10 12 cm -2

(5.3 × 10 12 cm -2 ). On the other hand, the screw and edge velocities for [110]-direction ([111]-direction) are 506 m/s (820 m/s) and 3000 m/s (2500 m/s), respectively. Such influence of temperature on the velocity is qualitatively consistent with trends reported in the literature [START_REF]Stochastic simulation of dislocation glide in tantalum and ta-based alloys[END_REF][START_REF] Marian | Dynamic transitions from smooth to rough to twinning in dislocation motion[END_REF][START_REF]The glide of screw dislocations in bcc fe: Atomistic static and dynamic simulations[END_REF]. When compression is applied along the [001] direction, dislocations are mostly of edge type, and their velocity is still too slow to be determined. A more complete and detailed analysis of the dislocation behavior will be presented in a forthcoming paper.

In contrary, the occurrence of deformation twinning is found to be strongly temperature-dependent.

In fact, no twinning is observed during the compression in either defect-free or defective crystal when initial temperature equals or exceeds 200 K (see Fig. 9). This indicates that the deformation twinning occurring in our previous simulations is mostly a low-temperature effect of iron single crystal in consistency with the well known fact that low temperature favors twinning and with the existence of twinning temperature limit [START_REF] Meyers | Defect generation in shock wave deformation[END_REF][START_REF] Meyers | Dynamic behavior of materials[END_REF]. Thus, on varying the initial temperature both plastic deformation and structural phase transformation characteristics are found to remain almost identical to 50 K simulations for defect-free [111]-oriented already twin-free at low temperature and defective both [110] and

[111]-oriented crystals where the low temperature deformation was mainly dominated by dislocations. This behavior is coherent with previous MD results reported in iron using Voter-Chen potential where increasing the initial temperature from 50 to 300 K has no significant influence [START_REF] Kadau | Atomistic simulations of shockinduced transformations and their orientation dependence in bcc Fe single crystals[END_REF]. On the other hand, along the [110], by increasing the temperature above 200 K, the first stage of the deformation is found to become similar to that of the [111] direction in defect-free crystal at 50 K, i.e elastic deformation Atoms color codes are the same as in the Fig. 2. See text for details.

followed by pressure stagnation before the structural phase transformation. Then, the deformation proceeds by a mixed regime in which dislocations develop as in low temperature regime. So, no double yielding is observed.

Besides, along [001] direction (with or without microvoids) where twinning was the main plastic response at low temperature, the crystal response evolves from elastic-plastic to perfect elastic when the initial temperature is 400K. Here, the shear stress increases linearly with compression until the structural phase transformation proceeds in an elastic state (Fig. 10). The amplitude of atomic thermal fluctuations is the only factor that favors the structural phase transition; the transition onset pressure is observed to decrease from about 100 GPa at 50 K to about 40 GPa at 400 K. This is consistent with the martensitic character of this transition in agreement with both negative slope of the bcc-hcp phase boundary in the pressure-temperature diagram and previous experimental works on ramp compression where the transition velocity was observed to decrease on increasing the target initial temperature [START_REF] Bastea | High pressure phase transformation in iron under fast compression[END_REF].

These results show that the temperature is a key factor for understanding deformation and structural phase transition in single crystal iron. Furthermore, they also give insight into the lack of plasticity when shock wave propagates along [001] direction in iron single crystal under laser compression where initial temperature is usually about 300 K (ambient conditions), or higher (in cases of preheating due to laser-plasma interaction) [44, 15? ].

Summary and Conclusions

In summary, the effects of orientation, defects and initial temperature on the response of iron single crystals to ramp compression have been investigated using the modified version of the Ackland potential [START_REF] Gunkelmann | Polycrystalline iron under compression: Plasticity and phase transitions[END_REF].

In contrast with results obtained with the Voter Chen potential [START_REF] Kadau | Atomistic simulations of shockinduced transformations and their orientation dependence in bcc Fe single crystals[END_REF][START_REF] Wang | Molecular dynamics simulations of hcp/fcc nucleation and growth in bcc iron driven by uniaxial compression[END_REF], the elastic-plastic deformation behavior is found to strongly depend on crystal orientation, with major consequences on both the onset pressure and the kinetics of the bcc-hcp phase transition (Fig. 11). As reported previously [START_REF] Amadou | Coupling between plasticity and phase transition in shock-and ramp-compressed single-crystal iron[END_REF], ramp compression along the [001] direction from a 50 K initial temperature produces twin-mediated yielding (at Pz around 12 GPa) followed by elastic hardening which inhibits the structural phase transition (up to very high pressures on the order of 100 GPa). However, compression along other crystal orientations leads to very different responses. While perfect crystals compressed in the [111] direction remain basically elastic up to the onset of the phase transition (at Pz around 35 GPa), those compressed in their

[110] direction exhibit a complex behavior, with simultaneous twinning and massive nucleation of hcp phase (at about 27 GPa) creating a mixed-phase, metastable regime where dislocations develop (above about 40 GPa). The overall plasticity deformation behavior, i.e twinning along [001] direction and twinning followed by dislocations along [110] direction is consistent with features previously observed in defect-free bcc material such as tantalum [START_REF] Ravelo | Shock-induced plasticity in tantalum single crystals: Interatomic potentials and large-scale molecular-dynamics simulations[END_REF][START_REF] Zepeda-Ruiz | Probing the limits of metal plasticity with molecular dynamics simulations[END_REF]. They differ from those obtained with the modified analytic embedded-atom-model potential for which plasticity, fully dislocation-mediated, is only found under compression in the [110] direction above 17 GPa [START_REF] Wang | An atomic study on the shock-induced plasticity and phase transition for iron-based single crystals[END_REF].

The strong directional dependence of the bcc-hcp phase transition reported above is consistent with the theory of phase nucleation and growth [START_REF] Christian | The Theory of Transformations in Metals and Alloys[END_REF][START_REF] Smallman | Modern Physical Metallurgy[END_REF][START_REF] Porter | Phase Transformations in Metals and Alloys[END_REF]. Indeed, the thermodynamic driving force for the phase transformation to occur results from the negative bulk (chemical) free energy difference between the parent and product phases. However, the formation of nucleus of the product phase induces positive interface energy, which opposes to the transformation. Besides, in the solid state, as the product and parent phases have different crystalline structures and specific volumes, the transformation is accompanied by dimensional changes. Such changes in volume and shape cannot occur freely because of the rigidity of the surrounding matrix, and strain energy is induced, which tends to inhibit the growth of the nucleus [START_REF] Smallman | Modern Physical Metallurgy[END_REF]. Thus, a different plastic behavior along different crystallographic orientations would be expected to affect the onset of the phase transition, in consistence with our simulations and experimental observations in other shock-compressed materials such as potassium chloride [START_REF] Hayes | Polymorphic phase transformation rates in shock-loaded potassium chloride[END_REF] and cadmium sulfide [START_REF] Knudson | Transformation kinetics for the shock wave induced phase transition in cadmium sulfide crystals[END_REF].

However, this observation contrasts with plate impact experiments on iron single crystal [START_REF] Liu | Isotropic phase transition of single-crystal iron (Fe) under shock compression[END_REF] where, although elastic-plastic deformation was found to be anisotropic, the same transition onset pressure was measured for different directions of shock application. Here, several possible origins have been invoked to explain the discrepancy between the experiment and shock MD simulations among which different time scales (four magnitude order of difference) and different shear stress state during the compression [START_REF] Liu | Isotropic phase transition of single-crystal iron (Fe) under shock compression[END_REF]. Indeed, it is suggested that during the experiment the phase transition would happen under nearly hydrostatic conditions, after full completion of the plastic processes and subsequent yielding, very different from the transient conditions and high shear stresses involved in our simulations (see Fig. 1 and Fig. 5) Besides, this observation is also consistent with other experiments involving laser driven and magnetic driven ramp compression of polycrystalline iron [START_REF] Smith | Time-dependence of the alpha to epsilon phase transformation in iron[END_REF], as well as shock compression of other materials such as CdS, where the transformation kinetics was found to depend significantly on the crystallographic orientation [START_REF] Knudson | Transformation kinetics for the shock wave induced phase transition in cadmium sulfide crystals[END_REF]. Regarding high rate effects considered herein, under ramp compression involving strain rates above 10 6 s -1 , a transition is expected from thermal activation to phonon drag plasticity regime, where the growth of the new-phase nuclei within the parent matrix is limited by phonon drag mediated dislocation velocity [START_REF] Smith | Time-dependence of the alpha to epsilon phase transformation in iron[END_REF].

The introduction of initial defects in the form of microvoids favors plastic deformation, which is then observed in all three loading directions, above a yielding pressure of 12 GPa ([001] direction) or 20 GPa ([110] and [111] directions). Such anisotropic elastic limit, in violation of the Schmid law, is characteristic of bcc materials at low temperature [START_REF] Dezerald | Plastic anisotropy and dislocation trajectory in bcc metals[END_REF][START_REF] Duesbery | Plastic anisotropy in b.c.c. transition metals[END_REF], and originates from the behavior of screw dislocations. Indeed, in the presence of defects, plasticity is still dominated by twinning in the [001] case, due to a very slow dislocation kinetics in this loading direction, while it is mainly mediated by (screw) dislocations activity in the other two directions. Like in defect-free crystals, directional dependence is also found for the onset pressure of the structural phase transition and its kinetics (Fig. 11).

The values of this onset pressure (25, 22 and 33 GPa for [001], [110] and [111] directions, respectively) are lower than their counterparts in the defect-free crystal, which is consistent with the theory of heterogeneous nucleation. They are logically higher than the ∼13 GPa equilibrium value due to the very short compression times considered here, while they are in good quantitative agreement with the 25 GPa found under picosecond laser compression [START_REF] Crowhurst | The α -ε phase transition in iron at strain rates up to 10 9 s -1[END_REF] and with the 11-38 GPa range (depending on the strain rate) found under nanosecond laser ramp compression [START_REF] Smith | Time-dependence of the alpha to epsilon phase transformation in iron[END_REF][START_REF] Amadou | Kinetics of the iron α -ε phase transition at high-strain rates: Experiment and model[END_REF][START_REF] Amadou | Laser-driven ramp compression to investigate and model dynamic response of iron at high strain rates[END_REF], as well as with the 23 GPa predicted by recent MD simulations of shock waves in polycrystalline iron [START_REF] Gunkelmann | Shock waves in polycrystalline iron: Plasticity and phase transitions[END_REF].

Figure 1 :

 1 Figure 1: Wave profiles for [001] (top), [110] (middle) and [111] (bottom) crystallographic orientationsof defect-free iron samples ramp compressed (with a piston maximum velocity of 800 m/s) at 30 picosecond (thin lines) when piston reached its maximum velocity and 10 ps later (thick lines) . Green curves depict longitudinal pressure P z while red curves present shear stress.

Figure 2 :

 2 Figure 2: Atomic structure snapshots at 40 ps simulation time-step for [001] (top row), [110] (center row) and [111] (bottom row) directions in defect-free iron single crystal iron, ramp compressed to piston maximum velocity of 800 m.s -1 (corresponding to maximum P z of about 50 GPa) within 30 ps. The non hcp and fcc atoms are colored according to P z (i.e. the z-component of the atomic stress σ times the atomic volume V a ), except for atoms in twins, which are colored in brown. Violet and cyan are the hcp

Figure 3 :

 3 Figure 3: Example of atom configuration in the (110) plane under ramp compression along [001]direction to a longitudinal pressure Pz of 19 GPa. Red lines are associated with twins. See ref [23] for these twins dynamics.

Figure 4 :

 4 Figure 4: (a) Closer view of the mixed regime in a [110]-oriented defect-free crystal (Fig. 2, middle) showing particles with unknown coordination structure (by CNA analysis), identified by DXA analysis as dislocations loops. (b) Burgers vector type: dislocations are mainly 1/2 < 111 > type (green) with only few < 100 > type (violet, barely visible). (c) Dislocations local character: the red and blue lines are screw and edge dislocation types, respectively, while other colors represent mixed dislocation. Small arrows in pale blue are Burgers vectors (magnified by 10). The gray color corresponds to defect surfaces.

Figure 5 :

 5 Figure 5: Wave profiles in [001] (top), [110] (middle) and [111] (bottom) -oriented iron single crystals containing initial defects in the form of microvoids (see text for details), at 30 ps when piston reached its maximum velocity of 800 m/s and 10 ps later. Green curves depict longitudinal pressure P z while red curves represent shear stress

Figure 6 :

 6 Figure 6: Atomic structure snapshots at 40 ps for [001] (top row), [110] (center row) and [111] (bottom row) directions in defective iron single crystal ramp compressed to piston maximum velocity of 800 m.s -1 (corresponding to maximum P z of about 50 GPa) within 30 ps. Atoms color codes are the same as in Fig. 2. Deformation is characterized by dislocations emission from micro-voids (see Fig.7) before the bcc to hcp phase transition. Plasticity under compression along the [001] direction remains dominated by twins.

Figure 7 :

 7 Figure 7: Non bcc atoms structure referring Fig. 6 with the top, center and bottom images corresponding to [001], [110] and [111] directions respectively. Color codes are the same as in Fig. 2. Dislocations loops (corresponding here to atoms with non identified structure by cna analysis) are clearly evidenced around micro-voids.

Figure 8 :

 8 Figure 8: Nucleation and growth of dislocation loops between 37 ps and 40 ps in a crystal containing initial defects, ramp-compressed along its [110] direction. Burgers vectors are almost all of 1/2 < 111 > type. Loops consist of two screw dislocations with opposite Burgers vectors (in red) separated by edge dislocations (in blue), which propagate through the sample.

Figure 9 :

 9 Figure 9: Atomic structure snapshots at 20 ps simulation time-step at different equilibration-temperatures along [001] direction for defect-free iron single crystal iron, ramp compressed during 30 ps to piston maximum velocity of 1600 m.s -1 (corresponding to maximum P z of about 110 GPa). The snapshots from top to bottom correspond to 50, 100, 200, 300 and 400K initial temperature respectively.

Figure 10 :

 10 Figure 10: Profiles of shear stress for [001]-oriented single iron sample initially equilibrated at temperature of 400K. See also Fig. 9.

Figure 11 :

 11 Figure 11: Summary of the processes (with their onset pressure P z ) governing the response to ramp compression along the z-direction. The arrows present the time sequence in which those processes successively occur. The initial temperature of the crystal is 50 K. If this temperature is increased above about 200 K, twinning disappears in all cases.

Acknowledgments

Computations were performed on the supercomputer facilities of the Mésocentre de calcul de Poitou Charentes (France).

Data availability

The raw/processed data required to reproduce these findings cannot be shared at this time as the data also forms part of an ongoing study.


Florie Mesple
email: florie.mesple@cea.fr

Ahmed Missaoui

Tommaso Cea

Loic Huder

Francesco Guinea

Guy Trambly De Laissardière

Claude Chapelier

Vincent T Renard
email: vin-cent.renard@cea.fr

Heterostrain rules the flat-bands in magic-angle twisted graphene layers

The moiré of twisted graphene bilayers can generate flat bands in which charge carriers do not posses enough kinetic energy to escape Coulomb interactions with each other 1,2 leading to the formation of novel strongly correlated electronic states [START_REF] Cao | Magic-angle graphene superlattices: a new platform for unconventional superconductivity[END_REF][START_REF] Cao | Correlated insulator behaviour at half-filling in magic-angle graphene superlattices[END_REF] . This exceptionally rich physics relies on the precise arrangement between the layers. [START_REF] Carr | Pressure dependence of the magic twist angle in graphene superlattices[END_REF][START_REF] Yankowitz | Tuning superconductivity in twisted bilayer graphene[END_REF] We survey published Scanning Tunnelling Microscope (STM) measurements [START_REF] Huder | Electronic spectrum of twisted graphene layers under heterostrain[END_REF][START_REF] Kerelsky | Maximized electron interactions at the magic angle in twisted bilayer graphene[END_REF][START_REF] Choi | Electronic correlations in twisted bilayer graphene near the magic angle[END_REF][START_REF] Jiang | Charge order and broken rotational symmetry in magic-angle twisted bilayer graphene[END_REF][START_REF] Xie | Spectroscopic signatures of many-body correlations in magic-angle twisted bilayer graphene[END_REF][START_REF] Zhang | Probing the wave functions of correlated states in magic angle graphene[END_REF][START_REF] Wong | Cascade of electronic transitions in magic-angle twisted bilayer graphene[END_REF][START_REF] Choi | Tracing Correlated Chern Insulators in Magic Angle Twisted Bilayer Graphene[END_REF][START_REF] Nuckolls | Strongly Correlated Chern Insulators in Magic-Angle Twisted Bilayer Graphene[END_REF][START_REF] Zhang | Electron interactions in strain-induced zero-energy flat band in twisted bilayer graphene near the magic angle[END_REF] to prove that near the magic angle, native heterostrain, the relative deformations between the layers, dominates twist in determining the flat bands. This is demonstrated at large doping where electronic correlations have a weak effect and where we also show that tip-induced strain can have a strong influence. In the opposite situation of low doping, we find that electronic correlation further normalize the flat bands in a way that strongly depends on experimental details.

The strongly correlated electron physics recently observed in twisted graphene layers [START_REF] Cao | Magic-angle graphene superlattices: a new platform for unconventional superconductivity[END_REF][START_REF] Cao | Correlated insulator behaviour at half-filling in magic-angle graphene superlattices[END_REF] develops in flat bands 1,2 which are very sensitive to the relative arrangement between the layers. For instance, the superconducting phase has been reported to occur in a very narrow range of rotation angle between the layers around the magic angle. As another illustration, hydrostatic pressure changes the interlayer distance which also strongly influences superconductivity [START_REF] Carr | Pressure dependence of the magic twist angle in graphene superlattices[END_REF][START_REF] Yankowitz | Tuning superconductivity in twisted bilayer graphene[END_REF] . Heterostrain, the inplane deformation of one layer with respect to the other is an ubiquitous source of modification of the relative arrangement between the layers [START_REF] Huder | Electronic spectrum of twisted graphene layers under heterostrain[END_REF][START_REF] Kerelsky | Maximized electron interactions at the magic angle in twisted bilayer graphene[END_REF][START_REF] Choi | Electronic correlations in twisted bilayer graphene near the magic angle[END_REF][START_REF] Jiang | Charge order and broken rotational symmetry in magic-angle twisted bilayer graphene[END_REF][START_REF] Xie | Spectroscopic signatures of many-body correlations in magic-angle twisted bilayer graphene[END_REF][START_REF] Zhang | Probing the wave functions of correlated states in magic angle graphene[END_REF][START_REF] Wong | Cascade of electronic transitions in magic-angle twisted bilayer graphene[END_REF][START_REF] Choi | Tracing Correlated Chern Insulators in Magic Angle Twisted Bilayer Graphene[END_REF][START_REF] Nuckolls | Strongly Correlated Chern Insulators in Magic-Angle Twisted Bilayer Graphene[END_REF][START_REF] Zhang | Electron interactions in strain-induced zero-energy flat band in twisted bilayer graphene near the magic angle[END_REF] . Experiments 7 and theory [START_REF] Cosma | Moiré pattern as a magnifying glass for strain and dislocations in van der Waals heterostructures[END_REF][START_REF] Bi | Designing flat bands by strain[END_REF] have shown that heterostrain affects the flat bands which could have an impact on the strongly correlated electron physics. Bi, Yuan and Fu have even predicted that near the magic angle, the bandwith becomes insensitive to the twist angle and that the effect of heterostrain is dominant, [START_REF] Bi | Designing flat bands by strain[END_REF] calling for systematic experimental study of its effect. Such program is however difficult to implement owing to the lack of controllability of heterostrain, and its inhomogeneity previously referred to as "twist angle disorder" and which was also found to impact strongly correlated states at the macroscopic scale [START_REF] Choi | Tracing Correlated Chern Insulators in Magic Angle Twisted Bilayer Graphene[END_REF][START_REF] Lu | Superconductors, orbital magnets and correlated states in magic-angle bilayer graphene[END_REF] . In order to overcome these issues, we survey already published experimental STM data [START_REF] Huder | Electronic spectrum of twisted graphene layers under heterostrain[END_REF][START_REF] Kerelsky | Maximized electron interactions at the magic angle in twisted bilayer graphene[END_REF][START_REF] Choi | Electronic correlations in twisted bilayer graphene near the magic angle[END_REF][START_REF] Jiang | Charge order and broken rotational symmetry in magic-angle twisted bilayer graphene[END_REF][START_REF] Xie | Spectroscopic signatures of many-body correlations in magic-angle twisted bilayer graphene[END_REF][START_REF] Zhang | Probing the wave functions of correlated states in magic angle graphene[END_REF][START_REF] Wong | Cascade of electronic transitions in magic-angle twisted bilayer graphene[END_REF][START_REF] Choi | Tracing Correlated Chern Insulators in Magic Angle Twisted Bilayer Graphene[END_REF][START_REF] Nuckolls | Strongly Correlated Chern Insulators in Magic-Angle Twisted Bilayer Graphene[END_REF] in view of quantifying the effect of homogeneous heterostrain on the physics of magic-angle twisted graphene layers.

Figure 1 presents typical STM images collected from Refs. 7-9. An immediate observation is that the STM images all look very similar. This is not surprising because these samples have a twist angle very close to one another. This similarity is only apparent as evidenced by the variety of shapes of the local density of states (LDOS) measured from the dI/dV(V) spectroscopy (Fig. 1d). Despite the flat-bands should merge at the magic angle, the spectra of Ref. 8 and 9 show two van Hove singularities indicating that the flat bands are still separated. Their spacing ∆E exp is doping dependent which has been attributed to electron-electron interactions. [START_REF] Kerelsky | Maximized electron interactions at the magic angle in twisted bilayer graphene[END_REF][START_REF] Choi | Electronic correlations in twisted bilayer graphene near the magic angle[END_REF][START_REF] Jiang | Charge order and broken rotational symmetry in magic-angle twisted bilayer graphene[END_REF][START_REF] Xie | Spectroscopic signatures of many-body correlations in magic-angle twisted bilayer graphene[END_REF][START_REF] Zhang | Probing the wave functions of correlated states in magic angle graphene[END_REF][START_REF] Wong | Cascade of electronic transitions in magic-angle twisted bilayer graphene[END_REF] It can reach several tens of meV, even at large doping where correlations are not expected to renormalize strongly the bands. Still, while the two samples have a twist angle differing by only 0.1 • , doping does not explain why their ∆E exp differ by a factor 2-3. Sample to sample variation is obvious from other published data. [START_REF] Jiang | Charge order and broken rotational symmetry in magic-angle twisted bilayer graphene[END_REF][START_REF] Xie | Spectroscopic signatures of many-body correlations in magic-angle twisted bilayer graphene[END_REF][START_REF] Zhang | Probing the wave functions of correlated states in magic angle graphene[END_REF][START_REF] Wong | Cascade of electronic transitions in magic-angle twisted bilayer graphene[END_REF][START_REF] Choi | Tracing Correlated Chern Insulators in Magic Angle Twisted Bilayer Graphene[END_REF][START_REF] Nuckolls | Strongly Correlated Chern Insulators in Magic-Angle Twisted Bilayer Graphene[END_REF] Strikingly the data of Ref. 7 show a third peak at zero energy, a result which was recently reproduced [START_REF] Zhang | Electron interactions in strain-induced zero-energy flat band in twisted bilayer graphene near the magic angle[END_REF] . In order to determine whether this variety in sample behaviour can be understood within the framework of heterostrain, we use the method described in Ref. 20 to determine the precise arrangement of the layers and calculate the corresponding local density of states using a tight-binding method. 1,[START_REF] Brihuega | Unraveling the Intrinsic and Robust Nature of van Hove Singularities in Twisted Bilayer Graphene by Scanning Tunneling Microscopy and Theoretical Analysis[END_REF][START_REF] Trambly De Laissardière | Numerical studies of confined states in rotated bilayers of graphene[END_REF] Figure . 1a, b and c include our estimates of heterostrain. In all cases, biaxial heterostrain is smaller than uniaxial heterostrain which varies by a factor 3 from the smallest to the largest. Figure . 1d shows that the calculated local density of states agrees very well with experimental data. In particular the number of peaks is controlled by heterostrain which also influences their spacing pointing to a strong contribution of heterostrain to the observed phenomenology. The tight-binding calculations presented in Fig. 1f show that the van Hove singularities separate linearly with uniaxial heterostrain for all angle θ s of application of heterostrain. This angle controls the splitting of the van Hove singularities (Fig. 1g) leading to three typical behaviours reported in Fig. 1f (See supplementary information and videos for more theoretical results). These results agree with those of the continuum model [START_REF] Bi | Designing flat bands by strain[END_REF] . The more quantitative comparison of Fig. 2a shows that ∆E exp converges to the value ∆E T B predicted by our tight-binding calculation at large doping. This has to be expected because the effect of electronic correlations measured by the ratio of Coulomb to kinetic energy reduces with doping and the system evolves towards the non-interacting situation modeled by our tight binding calculations. It establishes that heterostrain controls the physics of twisted graphene layers near the magic angle at large doping. This can be viewed explicitly in Fig. 2b which shows that ∆E 0 exp , the experimental spacings of van Hove singularities at large doping, depends linearly on heterostrain as predicted by tight binding and the continuum model. On the contrary, the same data plotted as function of the twist angle do not show particular correlation with this parameter (Supplementary Fig. S5).

θ int = 1.10°ε bi = 0.13% ε uni = -0.55% θ s = 34°a θ int = 1.25°ε bi = -0.06% ε uni = 0.35% θ s = 10.6°b θ int = 1.01°ε bi = -0.0% ε uni = 0.16% θ s = -26°c f V B -V CNP (meV)
Figure 2c provides a deeper level of comparison between theory and experiment. It presents the relative difference between the ∆E T B and the experimental ∆E 0 exp as function of the tunnelling resistance R t = V b /i t (V b and i t are the tunnelling bias and current). The excellent agreement obtained for R t > 2 GΩ degrades below this value pointing to a possible influence of tip-induced strain which is known to be controlled by R t . [START_REF] Mamin | Contamination-mediated deformation of graphite by the scanning tunneling microscope[END_REF] The deformations seen in Fig 3a, b andc showing the evolution of the image of Fig. 1b for decreasing R t corroborate this interpretation (The reader may also refer to Fig. S3 Such presentation of all the data on a single graph is justified by i) the weak dependence of ∆Eexp on the twist angle near the magic angle (See Ref. [START_REF] Bi | Designing flat bands by strain[END_REF] and supplementary Fig. S5) and ii) the weak dependence of ∆E on θs if it is estimated as the spacing between the outermost singularities (Fig. 1g). This is also why the Hartree-Fock model was calculated for a twist angle of 1.1 • and θs = 30 • . On the contrary, ∆ET B were calculated with the full relative arrangement of each experimental data and show a better agreement with experiments. The twist angle has still a small influence which explains the small deviations from a purely linear strain dependence. The measurements showing a cascade of transitions are indicated by a double border. c Relative deviation of ∆E 0 exp to ∆ET B as a function of the tunnelling resistance Rt = V b /it. The figure includes new data for the sample of Fig. 1b.

Ref. [START_REF] Choi | Electronic correlations in twisted bilayer graphene near the magic angle[END_REF] for another example).

Figure 3 shows that while spatial variations of the LDOS at high R t only reflect electronic localisation in AA regions, the response to tip-induced strain at low R t strongly depends on the position on the moiré pattern. While AB regions are weakly affected (Fig. 3e), ∆E exp is reduced by 40% at low tunnel resistance in AA regions indicating a flattening of the bands (Fig. 3d). The flattening culminates in intermediate regions between two AA regions where a 15 meV wide single LDOS peak signals that very flat bands can be engineered there (Fig. 3f). It is not surprising that these regions are the most sen-sitive to tip-induced strain because they are characterized by an excess of elastic energy [START_REF] Alden | Strain solitons and topological defects in bilayer graphene[END_REF][START_REF] Van Wijk | Relaxation of moiré patterns for slightly misaligned identical lattices: graphene on graphite[END_REF][START_REF] You | Atomic and electronic reconstruction at the van der Waals interface in twisted bilayer graphene[END_REF][START_REF] Gornostyrev | Origin of the vortex displacement field in twisted bilayer graphene[END_REF] and are therefore more easily perturbed. We note that the tip may induce a complex strain pattern including hetero-, homostrain as well as vertical displacement all of which depend on the position on the moiré. It follows that the relative arrangement between the layers can no longer be determined from STM images at low R t and that the response to tip-induced strain is largely sample dependent. This is illustrated by Fig. S5 of Ref. 9 in which, contrary to our experiment, the spacing between van Hove singularities increases with decreasing R t . We conclude that great care is needed when discussing quantitatively experimental results and that the STM tip can be used to locally engineer the relative arrangement between the layers. The very good agreement between experiments and theory obtained so far raises the question about the influence of heterostrain on electronic correlations which express most at zero doping. [START_REF] Kerelsky | Maximized electron interactions at the magic angle in twisted bilayer graphene[END_REF][START_REF] Choi | Electronic correlations in twisted bilayer graphene near the magic angle[END_REF][START_REF] Jiang | Charge order and broken rotational symmetry in magic-angle twisted bilayer graphene[END_REF][START_REF] Xie | Spectroscopic signatures of many-body correlations in magic-angle twisted bilayer graphene[END_REF][START_REF] Zhang | Probing the wave functions of correlated states in magic angle graphene[END_REF][START_REF] Wong | Cascade of electronic transitions in magic-angle twisted bilayer graphene[END_REF] Figure 2b shows that ∆E * , the experimental spacing of van Hove singularities in this regime is increased by electronic correlations and tends to increase with strain. This tendency is well pictured by our Hartree-Fock calculations presented Fig. 2b including and electronic correlations for the dielectric constant = 5 (See Ref. 28 and supplementary information for the details). These calculations however underestimate ∆E * exp . Further increasing the effect of interactions by decreasing the dielectric constant does not lead to a better agreement especially at large strains where the low energy bands become so wide that they merge into the continuum of the higher energy bands and van Hove singularities can no longer be clearly defined (Supplementary information). Also the calculations do not capture the large experimental scatter of ∆E * exp which we were not able to correlate to any experimental parameter (twist angle, heterostrain value or angle of application, bandwidths as measured by the FWHM of van Hove singularities, tunnelling resistance, temperature). The scatter also does not correlate with the appearance of polarized states leading to the cascade of transitions seen in some samples at intermediate fillings. [START_REF] Wong | Cascade of electronic transitions in magic-angle twisted bilayer graphene[END_REF][START_REF] Choi | Tracing Correlated Chern Insulators in Magic Angle Twisted Bilayer Graphene[END_REF][START_REF] Nuckolls | Strongly Correlated Chern Insulators in Magic-Angle Twisted Bilayer Graphene[END_REF] This points to a strong sensitivity of electronic correlations to some additional experimental parameters beyond those investigated here. This could be due to the substrate or more generally to the detailed electrostatic environment as suggested by several studies reviewed in Ref. 29, and to atomic lattice relaxation effects [START_REF] You | Atomic and electronic reconstruction at the van der Waals interface in twisted bilayer graphene[END_REF] calling for a systematic experimental study of the effect of those parameters.

Returning to heterostrain, its strong impact on the flat bands of magic angle twisted graphene layers also calls for a systematic investigation of its influence on the strongly correlated phases and that of other moiré materials. In this context, it would be extremely desirable to be able to tune it. Alternatively, one could also rely of the variability in the fabrication process to generate a representative set of samples such as the one we have studied here.

Methods

Data collection

When possible, the raw data used in the survey were retrieved from the authors of the publications. Otherwise, the data were digitized from the corresponding publications. We could use the method of Ref. [START_REF] Artaud | Universal classification of twisted, strained and sheared graphene moiré superlattices[END_REF] to extract heterostrain for Refs. [START_REF] Huder | Electronic spectrum of twisted graphene layers under heterostrain[END_REF][START_REF] Kerelsky | Maximized electron interactions at the magic angle in twisted bilayer graphene[END_REF][START_REF] Choi | Electronic correlations in twisted bilayer graphene near the magic angle[END_REF] which had STM images with atomic resolution. The strength of this method is to isolate heterostrain from homostrain and possible artefacts due to the calibration of the piezoelectric tube. For Refs. (10, 11, 13a, 14, 15), we used the method described in Ref. [START_REF] Kerelsky | Maximized electron interactions at the magic angle in twisted bilayer graphene[END_REF] to determine heterostrain. Although it does not allow to determine biaxial heterostrain we found it to provide good estimates of uniaxial heterostrain for the data of Refs. [START_REF] Huder | Electronic spectrum of twisted graphene layers under heterostrain[END_REF][START_REF] Kerelsky | Maximized electron interactions at the magic angle in twisted bilayer graphene[END_REF][START_REF] Choi | Electronic correlations in twisted bilayer graphene near the magic angle[END_REF] . For Refs. (12, 13b, 13c), we did not have the data to do our own analysis and used the values of heterostrain provided by the authors. Error bars in Fig. 2 of main text present the uncertainty as estimated the full width at half maximum of a van Hove singularity.

STM measurements

We performed the measurements presented in 3 of the main text in order to document the effect of tipinduced strain on twisted graphene layers described in Ref. [START_REF] Huder | Electronic spectrum of twisted graphene layers under heterostrain[END_REF] . Tip-induced strain was changed by changing the bias for a given tunnel current in imaging conditions. In spectroscopic mode, the tip-sample interaction was defined by the set point prior to switching off the feedback loop and subsequent sweeping of the bias voltage between -200 mV and 200 mV. The dI/dV (V ) signal was measured using phase sensitive detection with a 2 mV oscillation at 263 Hz added to the tunnel bias. The di/dV curves were normalised to 1 at -200 mV.

  Heterostrained twisted graphene layers. Panels a, b and c present STM images of twisted graphene layers near the magic angle adapted from Refs 8 , 7 and 9 respectively. The scale bar is 10 nm in each image. Insets present the parameters describing the relative arrangement of the layers in good agreement with estimates of the original studies. 7-9 d Local density of states measured in the AA regions for each samples of panels a, b and c. e Corresponding tight-binding calculation of the LDOS including heterostrain. f Tight-binding prediction for the local density of states in AA regions as function of energy and for increasing uniaxial heterostrain. The variations are plotted for three different angles θs of application of heterostrain. g Local density of states calculated for 0.4 % of heterostrain applying along varying θs. The white dotted lines indicate the values of angles used in panel f. h Sketch of uniaxial heterostrain configuration. One layer is deformed by a uniaxial heterostrain applied along the direction defined by θs and then rotated by and angle θ with respect to the undeformed layer.

FIG. 2 :

 2 FIG. 2: Comparison between theory and experiments. a Representative doping dependence of the spacing between van Hove singularities. The horizontal line corresponds to prediction by tight-binding calculations with heterostrain. b Experimental spacing of van Hove singularities at large doping (∆E 0 exp ), zero doping (∆E * exp ), theoretical tight-binding (∆ET B ), low-interaction Hartree-Fock (∆E =∞ HF ) and large interaction ∆E =5 HF .Such presentation of all the data on a single graph is justified by i) the weak dependence of ∆Eexp on the twist angle near the magic angle (See Ref.[START_REF] Bi | Designing flat bands by strain[END_REF] and supplementary Fig.S5) and ii) the weak dependence of ∆E on θs if it is estimated as the spacing between the outermost singularities (Fig.1g). This is also why the Hartree-Fock model was calculated for a twist angle of 1.1 • and θs = 30 • . On the contrary, ∆ET B were calculated with the full relative arrangement of each experimental data and show a better agreement with experiments. The twist angle has still a small influence which explains the small deviations from a purely linear strain dependence. The measurements showing a cascade of transitions are indicated by a double border. c Relative deviation of ∆E 0 exp to ∆ET B as a function of the tunnelling resistance Rt = V b /it. The figure includes new data for the sample of Fig.1b.

FIG. 3 :

 3 FIG. 3: Tip-induced strain. a, b, c STM images measured at decreasing tunnel resistance. The tunnelling set point is it = 100 pA and the bias is V b = 175 mV in a, V b = 55 mV in b and V b = 15 mV in c. The images are 25×25 nm 2 . d, e, f Local density of states measured by STM for decreasing conductance in AA (d), AB (e) and intermediate regions (f ). The vertical dotted lines are guides for the eyes to track the position of the LDOS peaks.

Trambly de Laissardière, G., Mayou, D. & Magaud, L. Localization of Dirac Electrons in Rotated Graphene Bilayers.

Acknowledgements

The authors thank E. Andrei, S. Nadj-Perge, for sharing their data and for interesting discussions. AM is supported by Paris//Seine excellence initiative (Grant No. 2019-055-C01-A0). FG and TC acknowledge funding from the European Commision, under the Graphene Flagship, Core 3, grant number 881603, and from grants NMAT2D (Comunidad de Madrid, Spain), SprQuMat and SEV-2016-0686, (Ministerio de Ciencia e Innovación, Spain)

Author contributions:

VTR,FM and CC defined the project. FM performed the data analysis from the literature and provided the cells for tight binding calculations which were performed by AM and GTL. TC and FG performed Hartree-Fock calculations. LH, VTR and CC performed the STM experiments. VTR and FM wrote the manuscript with input from all the authors. VTR supervised the project.

Supplementary materials

Supplementary file and videos available on request to vincent.renard@cea.fr.