Theoretical Model of Structural Phase Transitions in Al-Cu Solid Solutions under Dynamic Loading Using Machine Learning

: The development of dynamic plasticity models with accounting of interplay between several plasticity mechanisms is an urgent problem for the theoretical description of the complex dynamic loading of materials. Here, we consider dynamic plastic relaxation by means of the combined action of dislocations and phase transitions using Al-Cu solid solutions as the model materials and uniaxial compression as the model loading. We propose a simple and robust theoretical model combining molecular dynamics (MD) data, theoretical framework and machine learning (ML) methods. MD simulations of uniaxial compression of Al, Cu and Al-Cu solid solutions reveal a relaxation of shear stresses due to a combination of dislocation plasticity and phase transformations with a complete suppression of the dislocation activity for Cu concentrations in the range of 30–80%. In particular, pure Al reveals an almost complete phase transition from the FCC (face-centered cubic) to the BCC (body-centered cubic) structure at a pressure of about 36 GPa, while pure copper does not reveal it at least till 110 GPa. A theoretical model of stress relaxation is developed, taking into account the dislocation activity and phase transformations, and is applied for the description of the MD results of an Al-Cu solid solution. Arrhenius-type equations are employed to describe the rates of phase transformation. The Bayesian method is applied to identify the model parameters with fitting to MD results as the reference data. Two forward-propagation artificial neural networks (ANNs) trained by MD data for uniaxial compression and tension are used to approximate the single-valued functions being parts of constitutive relation, such as the equation of state (EOS), elastic (shear and bulk) moduli and the nucleation strain distance function describing dislocation nucleation. The developed theoretical model with machine learning can be further used for the simulation of a shock-wave structure in metastable Al-Cu solid solutions, and the developed method can be applied to other metallic systems, including high-entropy alloys.


Introduction
The theoretical description of the mechanical behavior of materials is still an urgent problem, especially under dynamic loading conditions, despite the great research activity over several decades.The gradual investigation of the dominant mechanisms leading to the evolution of material structure during loading has led to the creation of numerous models predicting material properties.Under dynamic loading, the equivalent shear stress exceeds the quasi-static yield stress due to a finite rate of defect microstructure evolution, and the purpose of a dynamic plasticity model is to take this finite rate of relaxation into account.The simplest way to describe the dynamic plasticity is to introduce a dynamic yield stress as a material parameter [1,2], which is, unfortunately, valid only for specific test conditions.Semi-empirical constitutive models [3][4][5] were created to determine the strain rate and temperature dependence of material yield behavior from the consideration of defect microstructure evolution, typically the thermally activated dislocation motion with parameters fitting to the experiments [6].In general, the dislocation-mediated plasticity is the main mechanism for most crystalline metals, which is supported by both molecular dynamics (MD) simulations [7] and experimental studies [8].The incorporation of a crystal plasticity approach with an explicit description of the dislocation activity in different slip planes is a fruitful approach for quasi-static [9] and dynamic [10][11][12] applications.This approach is successfully used to describe the shock-wave processes in solid metals [13][14][15].
In addition to the dislocation plasticity, the plastic relaxation can be caused by several other mechanisms, which can act simultaneously with the dislocation activity and with each other.One can mention the interfacial slip and diffusion processes [16], mechanical twinning [17][18][19][20] in the case of coarse-grained solids, as well as additional grain boundary migration [21][22][23] and grain boundary sliding with grain rotation [24][25][26] in the case of fine-grained solids.In addition, a structural phase transition can have a significant effect on the dynamic response of metals.The phase transition of FCC (face-centered cubic) to BCC (body-centered cubic) leads to excessive stress relaxation according to MD simulations [27] at high compressive strain rates.Shock-induced phase transitions are revealed by means of MD simulations in single-crystal aluminum [28,29], copper [30], iron [31], cerium [32], magnesium [33] and high-entropy alloys [34,35].The FCC-to-BCC phase transition in copper was experimentally revealed in the quasi-static measurements by Dewaele et al. [36] with the diamond anvil cell, as well as in dynamic measurements with shock-wave compression [37] and ramp compression [38].The large number of experimentally detected phase transitions in materials has led to the development of mechanical models that take into account phase transitions.Levitas and Javanbakht [39,40] proposed a phase field model with in-deep thermodynamics consideration for the description of martensitic transformations.Yeddu and Lookman [41] developed a variant of the phase field model based on the time-dependent Ginzburg-Landau kinetics equation for the assessment of the β(BCC)-to-ω(hexagonal) phase transformation in Zr-Nb alloys.Consideration of multiple possible variants of austenite-to-martensite transformation in a quenching and partitioning steel QP1180 allowed Yang et al. [42] to describe correctly the strain path and orientation dependencies of the transformation process.In spite of these efforts, there is still a lack of a simple and robust model of plastic relaxation with the interplay of dislocation activity and phase transitions.Combining MD data, theoretical framework and machine learning (ML) methods is a prospective approach in the field.
The equation of state (EOS) with the additional consideration of the deviatoric part of the stress tensor becomes the core part of many constitutive models for shock-wave problems.At a very high impact level, the deviatoric part becomes negligibly low [43] and the EOS completely describes the behavior of matter.In general, the construction of the EOS consists of proposing a theoretically based functional form and fitting the relevant parameters to the available experimental data and first-principles calculations [44,45].This appears to be costly for complex substances and multi-component alloys.The application of ML methods, in particular artificial neural networks (ANNs), is a perspective method of solving this problem [46].
ML is increasingly being used in material science and in the mechanics of materials [47,48].ML methods are an efficient way to account for changes in mechanical properties of materials and are quite difficult to implement by classical methods and require an analysis of a large amount of experimental data and simulation results.The problem statement such as prediction, recognition and classification determines the choice of the most effective ML algorithm such as linear and logical regression, decision trees and clustering [49][50][51].Among various ML methods, ANNs and Bayesian global optimization of model parameters are of particular interest for the present study.ANNs can be used to approximate the strain, strain rate and temperature effect on the yield surface more precisely in comparison with the known analytical dependencies [52,53], to establish structure-property relationships [48,54,55] and to construct a surrogate model of the material [56,57].Bayesian calibration of model parameters was fruitfully used to calibrate the empirical Johnson-Cook model using the data on plate impact [58] and Taylor tests [59], as well as to fit a dislocation plasticity model to the Taylor test data [60][61][62].
In this work, we simulate the behavior of pure metals, aluminum and copper, as well as Al-Cu solid solutions with different contents of elements under uniform high-strain-rate compression.The mechanical behavior of the material is described by a coupled model of phase transition and plasticity with the Arrhenius-type equations for the rates of phase transitions and a modified Maxwell-type relaxation model [63] with a variable relaxation time that takes into account the evolution of the dislocation structure.The parameters of the stress relaxation and phase fraction evolution model are fitted by the Bayesian method to the results of the MD simulations.In addition to the coupled model of phase transition and plasticity, we have developed a deep forward-propagation ANN to describe the EOS of alloys.The ANN relates pressure as the hydrostatic part of stresses, bulk modulus and temperature with density and internal energy simultaneously for different concentrations of copper in the Al-Cu solid solution ranging from 0% Cu to 100% Cu.Another forwardpropagation ANN allows the calculation of the shear modulus and dislocation nucleation threshold for the same continuous range of Cu concentration.

MD Setup
We performed MD simulations of uniaxial dynamic compressive and tensile deformations of single crystals of Al, Cu and Al-Cu solid solutions with different component contents with an engineering strain rate of 10 9 s −1 , which is close to that achieved in state-of-art experiments on irradiation of thin metal foils with powerful sub-picosecond laser pulses [64][65][66][67]; thus, this strain rate is experimentally attainable.The single crystal is oriented so that the [100], [010] and [001] crystallographic directions are along the x-, yand z-axes (Figure 1).The MD system contains half a million atoms and is originally cubic with 50 × 50 × 50 lattice periods in size.Periodic boundary conditions are set in each of the three directions.It is known that increasing the concentration of the alloying element affects the mechanism of dislocation motion and leads to higher stresses [68][69][70][71].Usually, the proportion of the additive element in alloys is small-less than 5% Cu in Al matrix and up to 18% of Al in Cu matrix, but the alloying element can locally reach high concentrations in material.It is shown in [72] that the segregation of solid solution copper atoms reaches 15.3 wt % at grain boundaries during cryogenic surface mechanical grinding treatment.The formation of solid solution precipitates with high content of Cu in Al alloy during cyclic hardening is also possible [73]: the Cu content can reach almost 50 wt % according to the electronic microscopy images in this paper.Components of high concentrations having equal atomic fractions form high-entropy alloys with unique mechanical properties.These facts motivated us to consider copper concentrations of 10, 20, 30, 50, 70 and 80% in order to investigate theoretically the deformation behavior of the solid solution, although a part of the considered solid solutions can be thermodynamically unstable in traditional metallurgical casting.On the other hand, if they reveal unique dynamic-protecting properties, novel production techniques, such as additive manufacturing, can be applied to produce such metastable alloys; therefore, theoretical analysis of their properties is relevant.Aluminum atoms in the system were replaced by copper atoms in random order to form the solid solution.The compression and tension of the samples were modeled at constant temperature in the range of 100 to 900 K in steps of 100 K. Thus, 72 MD simulations were performed in total.The calculations were performed using the LAMMPS software package [74] with an interatomic potential of the angle-dependent potential (ADP) type [75], which is one of the most reliable classical force fields.The Nose-Hoover thermostat [76] and barostat were used to relax the stresses for 10 ps before the deformation stage.The crystal was deformed along the x-axis with fixed transverse dimensions.The deformation lasted for 300 ps, achieving an engineering strain value of 0.3.For comparison, auxiliary modeling was also carried out in the adiabatic approximation, that is, without a thermostat.Also, several MD systems were subjected to compression with other strain rates of 10 8 s −1 and 3•10 8 s −1 to reveal the rate effect.

Results of MD Simulations: Dislocation Activity and Phase Transitions
Data from MD simulations of the uniaxial compression of crystals for a wide temperature range were analyzed using the "Polyhedral template matching" algorithm [77] integrated into the OVITO software [78].This method classifies the local structural environment of the particles.The phase fraction was calculated as the ratio of the number of atoms of a certain phase to the total number of atoms of the modeled crystal.Figure 2 shows the evolution of shear stress τ = σ xx − σ yy /2 and atomic structure for pure aluminum and copper at 300 K.Both metals are characterized by an elastic response at the beginning of deformation.The appearance of the metastable BCC phase at an engineering strain of about 0.11 leads to an elastically nonlinear behavior.The subsequent appearance of the HCP (hexagonal close-packed) phase reflects the nucleation and movement of partial dislocations (Figure 3a) along the close-packing planes and leads to a sharp drop in stresses.Thus, the HCP phase here is essentially stacking faults left as a result of the slip of partial dislocations.The shear stress level increases in copper during the following plastic flow due to strain hardening.In the case of aluminum, strain hardening also takes place, but the formation of the BCC phase is observed almost throughout the entire volume of the crystal, which leads to rapid stress relaxation; shear stress even drops down to a negative value.The pressure curve reflects the effect of polymorphic phase transition in the case of pure aluminum (Figure 4): The pressure growth slows down above 35 GPa, which is associated with a change in the structure of the loaded material.The results of the MD data give a slightly higher pressure value compared to the experimental data [36] on compression in diamond anvil cells and DFT calculations [79].In general, all considered data are consistent.Let us consider in detail the behavior of an Al-Cu solid solution by an example of the crystal with 10% Cu in aluminum at the temperature of 300 K. Figure 5 shows the atomic structure of the crystal and the evolution of phase fractions and stresses.At the beginning of the simulation, the crystal has a primarily FCC crystal structure and deforms elastically.The shear stress reaches higher values than in pure aluminum due to the strengthening by copper atoms.Then, dislocation loops are formed in the material (Figure 3b).At this stage, the fraction of the metastable BCC phase is significantly less than in pure aluminum and copper.Active multiplication of Shockley partial dislocations takes place.The shear stress relaxes down to 2 GPa due to the plastic flow at an engineering strain of 0.11.In the course of further deformation, there is a gradual increase in the fraction of the BCC phase, while the FCC and HCP phases completely disappear.The transformation of the crystal structure entails a decrease in shear stresses to negative values.After that, the shear stress increases again due to the following deformation of the stable BCC phase.For the cases of solid solutions with higher copper concentrations (Figure 6), the process of dislocation plasticity plays an insignificant role in the stress relaxation in the material.Stress relaxation occurs primarily due to the transition of the initial structure into the BCC phase and the regions of disordered atoms ("OTHER" or amorphous phase).Phase transitions completely suppress dislocation activity.For most of the modeled crystals, the process of BCC phase formation from disordered atoms is common for the considered range of engineering strain.When the Cu fraction increases up to 50%, the amplitude of the maximum shear stress decreases and the beginning of the phase transition occurs at lower strains.The distribution of phase fractions for a solid solution with 80% Cu is similar to the case with 20% Cu. Figure 3 shows an increasing trend in the dislocation density with deformation change by a decreasing trend for pure Al and all investigated Al-Cu solid solutions; this is because the phase transition decreases the volume of the FCC phase and, thus, the total length of the dislocations.The maximal dislocation density increases with the addition of copper atoms up to the Cu concentration of 20%, but decreases with a further increase in Cu concentration, see Figure 3b.Let us consider the effect of temperature on the stress-strain state of the material (Figure 7a).The material softens, and the amplitude of the shear stresses decreases with increasing temperature.The MD data analysis shows an increase in the fraction of disordered atoms and a decrease in the phase transition threshold with increasing temperature.To determine the influence of the strain rate on the phase transition process, additional modeling of solid solution compression with strain rates of 10 8 and 3•10 8 s −1 at a temperature 300 K was carried out (Figure 7b).The change in strain rate by an order of magnitude does not significantly affect the evolution of the phase composition.On the contrary, the plasticity process is sensitive to changes in the strain rate, which is evident from the strong variation in the content of the HCP phase presenting the traces of dislocation gliding, Figure 7b.Consequently, the phase transition is strain-dependent instead of strain-rate sensitive in the studied dynamic loading range.It should be noted that there is a transition from the BCC phase to the region of disordered atoms ("OTHER" phase) when the strain is greater than 0.3.The crystal is almost completely disordered at the engineering strain of 0.45.As we show further, it is possible to formulate a theoretical model of coupled phase transition and plasticity for alloys on the basis of the obtained MD data on phase transformation under dynamic compressive deformation.Such a model must take into account the formation of the BCC lattice from the FCC and "OTHER" structures, as well as the increase in the fraction of the "OTHER" structure due to the decrease in the fractions of the FCC and BCC lattices.

The Constitutive Equations in the Form of the Artificial Neural Network
Forward-propagation artificial neural networks (ANNs) are used to approximate the constitutive equations of the Al-Cu alloy with a solid solution of copper atoms with different concentrations of components.The stress-strain relation for solids in the elastic domain can be presented in the form of the generalized Hook's law [83]: in which the components of the tensor C ijpq of elastic moduli are not constants and are introduced as functions of the current strain tensor C ijpq = C ijpq (ε).We considered that the isotropic formulation and the tensor C ijpq (ε) has only two independent components, i.e., Lamme coefficients.These coefficients can be resolved through other quantities that can be calculated from numerical calculations and experiments: the shear modulus G and the bulk modulus K, which will be considered in the remaining part of the paper.The constitutive equations in this work are based on uniaxial compression data.In this case, the deformation was determined by the ratio of the deformed and initial volumes ε = V/V 0 − 1 or equivalent ε = ρ 0 /ρ − 1.Therefore, we redefined the strain dependence on ε in Equation (1) to the density ρ dependence, which will greatly simplify the use of these constitutive equations in the macroscopic simulations in both mesh-and particle-based methods.The constitutive equations in the ANN form can be used in any other numerical schemes by implementing a forward propagation on the ANN parameters (weights and biases).The first ANN, whose scheme is presented in Figure 8, is essentially a scalar EOS and describes the following vector function: where the output layer of the ANN is → y = {P, E, K}.P is the pressure; E is the specific internal energy of material, and K is the bulk modulus.The input layer of the ANN consists of the arguments of the vector function (2).T is the temperature; ρ is the density, and α CU ∈ [0, 1] is the atomic fraction of Cu atoms in the solution.The training and test data sets were constructed from MD simulation data of uniaxial compression and tension of the representative volume of the Al-Cu alloy with the solid solution of copper atoms.In the case of tension, the deformation trajectories of MD data were restricted by the beginning of the fracture.The values of density, pressure, total energy and temperature were extracted directly from the MD simulation data, and the bulk modulus was estimated from the pressure curve as follows: where i is the state point for the training and test data set, n is the averaging interval, P j is the uniaxial tension/compression pressure and ε j is the uniaxial tension/compression strain value; the step of strain in MD output data is 10 −4 .
The second ANN (Figure 9) was used to describe the second part of the constitutive equation of the Al-Cu alloy with a solid solution of copper atoms, namely, the value of the shear modulus and dislocation nucleation threshold during uniaxial deformation.The function arguments were also determined by Equation ( 2), but the function values were set as , where G is the shear modulus and Q ε is the nucleation strain distance function introduced in [65].The shear modulus at state point i was calculated as where S xx j is the component of the stress tensor deviator and i, n and j are indices, which are defined in the same way as for Equation (3).The function Q ε characterizing the nucleation threshold was constructed as a continuous one, since an ANN is capable of approximating only continuous dependencies with high precision: where ε j is the current strain points of the system under uniaxial tension/compression and ε nuc is the threshold value for the onset of dislocation nucleation under uniaxial deformation.The nucleation of dislocations under uniaxial compression at the fixed strain rate begins when the condition Q ε ≤ 0 is satisfied.The inclination of Q ε in the strain space allows one to take into account the strain-rate dependence of the nucleation threshold [84].All results obtained by MD simulations were divided into training and test data sets.The test data were used only to check the accuracy and generalizability of the ANN by the cross-validation method and did not take part in the training.The data sets contained a total of 285,715 data vectors (α CU , ρ, T, P, E and K) for the first ANN and 133,763 data vectors (α CU , ρ, T, G and Q ε ) for the second ANN.The data distribution and the data ranges are given in Table 1.All the data underwent a normalization procedure using the standard definition of MinMax normalization: where x n ij is the normalized value, index j defines the data column, e.g., P, E, K, α CU , T and ρ for the first ANN, index i defines the data row, x max j is the maximum value in the j-th column and x min j is the minimum value in the j-th column.The range of values was expanded by 20% for each column relative to its real minimum and maximum.That is, the values x max j and x min j were shifted by 10% relative to the real minimum and maximum values in the original data sets.This is because the derivative of sigmoid transfer function is small at the extreme values close to 0 and 1 and the learning process becomes slower.Another motivation is to prevent going beyond the limits.
The ANNs were constructed as classical fully connected neural networks [85] in which the values obtained through the activation function of neurons were transmitted to the next layers of the network (Figures 8 and 9 show ANNs in the form of computational graphs, where the arrows indicate the direction and sequence of the calculations).At the current layer l of the ANN, the activation function a l j depends on the input z l j to the neuron: where n is the number of layers in the ANN, j is the index by neurons in the l-th layer, k is the index by neurons in the l − 1 layer, w l jk are the weight coefficients for the j-th neuron, a l−1 k is the neuron function value from the previous layer and b l j are biases values for the j-th neuron.The ANN contains several types of layers; the first one is the input layer for which l = 0, a 0 k = (x 1 , x 2 , . . . ,x m ), where x 1 , x 2 , . . . ,x m are the arguments of the approximating function.The arguments are α CU , T and ρ for both ANNs.The second type are the hidden layers, where the main transformation of the input data takes place, for them l = 1, . . ., n − 1.The last type is the output layer, where l = n = L.The output layer gives the values of the function P, E and K for the first ANN and G and Q ε for the second ANN.
The first ANN (Figure 8) has five hidden layers of 20 neurons each and the output layer that contains 3 neurons.The artificial neurons in hidden layers are described by the PReLU activation function [86]: where a pr is the parameter of the PReLU function.The sigmoid activation function is used to describe the output layer in this ANN: Such a combination of activation functions makes it possible to describe highly nonlinear dependencies, despite the fact that the PReLU function in the positive and negative regions is separately linear [85].The second ANN (Figure 9) contains six hidden layers of 20 neurons each and the output layer that contains 2 neurons.Swish activation function is used in hidden layers [87]: where β is the parameter of Swish function.The sigmoid function in Equation ( 9) is also used as the activation function on the output layer.These ANN architectures were chosen because they show the least deviation from the MD simulation results.The next step was to train the ANNs after building the training and test data sets.The stochastic gradient descent (SGD) method was used to train the neural networks: this method minimizes the cost function that defines the accuracy of the ANN.The crossentropy cost function was chosen for both ANNs, which together with the output sigmoid activation function (9), gave the following expression: where N is the number of outputs of the ANN, y i are the real values of the function at the point and ŷi are the ANN results at the point.The cross-entropy is typically used for classification problems, but our previous studies showed its efficiency for regression problems as well [62,88].SGD is classically implemented in fully connected forwardpropagation ANNs based on a backpropagation algorithm.In this work we also used the optimization of the SGD algorithm called Adam [89].The variable learning rate is implemented in the Adam method (i.e., gradient descent steps change with training progress), which allows one to quickly "go down" to the minimum of the cost function (11).Also, the Adam method provides greater accuracy of the ANN results and, at the same time, avoids the instability of the ANN training process (exploding or vanishing gradient).Tables 1 and 2 present the main characteristics and values of the hyperparameters for the first and second ANNs.The formulas for Adam's algorithm are given in Appendix A. The ANN allowed us to predict data for continuous values of input parameters.Figure 10 shows the pressure dependencies for the intermediate values of copper concentrations at the temperature of 300 K.These data were not included in the test and training data sets.The ANN qualitatively describes both an increase in pressure during the compression of a substance and a drop in pressure during the stretching process.There is an increase in the density of the substance and a deceleration in pressure growth with an increase in the concentration of copper.Next, we compared the temperature dependencies of the bulk modulus B and shear modulus G (Figures 11 and 12) predicted by the ANNs with data from other studies in the literature.The main experimental results are given for the single-crystal aluminum [90][91][92][93] and single-crystal copper [94][95][96].A comparison of the temperature dependence of the shear modulus (Figure 12b) also contained the theoretical estimate [97] and experimental results for the polycrystal [98].The bulk modulus of copper decreases with increasing temperature by about 20 GPa at temperatures close to 900 K, both according to our results and experimental data as can be seen in Figure 11a.There is some scatter in the experimental data of the bulk modulus of pure copper from different sources, but the ANN results generally agree with the considered data and lie between the experimental points (Figure 11a).The bulk modulus of pure aluminum (Figure 11b) decreases by 5 to 30 GPa depending on the considered studies with the increasing temperature up to 900.The results of our work show that the bulk modulus of pure aluminum drops by about 10 GPa, and these results also lie between the experimental points.[96]; Chang-Himmel (1966) [95]; Gerlich-Fisher (1969) [92]; Kamm-Alers (1964) and Tallon-Wolfenden (1979) [91,93]; Sutton (1953) [90]; Wawra (1978) [98]; Raju (2002) [97].Tallon-Wolfenden (1979) [93]; Gerlich-Fisher (1969) [92]; Kamm-Alers (1964) [91].
The comparison of the obtained values of the shear modulus from the trained ANN with the experimental data for pure copper (Figure 12a) shows the underestimation of the modulus value at the initial temperature of 100 K; the difference is about 1.5 GPa.In this case, the temperature dependence of the shear modulus predicted by the ANN shows the same trend as the experimental data: There is a decrease by approximately 8 GPa at a temperature of 800 K.In the case of pure aluminum (Figure 12b), the shear modulus predicted by the ANN describes the experimental data somewhat.The deviation manifests itself at high temperatures, where experimental data show a larger drop in the values of the shear modulus.Nevertheless, the trend toward a decrease in the shear modulus is observed both in our results and in the data from other works presented in Figure 12b.The deviation of the ANN prediction may be due to the used interatomic potential.In the recent work [99], the shear modulus is calculated for various temperatures by the MD method using EAM interatomic potential.The shear modulus values calculated for <110> crystallographic orientation and at the temperature 100 K sharply drop to approximately 15 GPa and decrease with a further increase in temperature, which does not correspond to any of the data presented in Figure 12b.This stimulates the further development of interatomic interaction potentials, such as NNP (neural network potential) [100] or PINN (physics-informed neural network potential) [101].These types of interatomic potentials describe physical processes with sufficient accuracy, but at this stage of development require a lot of computational time in MD simulations.The accuracy of the constitutive relations in the form of ANNs trained on data obtained using such potentials will also increase.
The ANNs were tested by comparing the density ρ, bulk B and shear G modulus in the Al-Cu alloy with a copper solid solution, where the copper concentration of up to 6 at.% was considered (Figures 13 and 14).The results of the ANN (this work), experimental data [102] and ab initio calculations [103] were compared.A comparison of the density of the alloy (Figure 13) shows that, both in the ANN predictions and in the experiment, a linear increase in density with an increase in the concentration of copper atoms is observed.Although the experiment shows the less sharp slope of the straight line, in general, all the considered data in Figure 13 are in agreement.The bulk modulus predicted by the ANN for the temperature of 300 K (Figure 14a) is somewhat overestimated in comparison with the experimental data and ab initio calculations.In this case, the experimental data [102] are presented for the polycrystal, which may give some discrepancy.The interatomic potential used can also lead to inaccuracies in the training and test data sets built on the basis of MD simulations.The shear modulus (Figure 14b) obtained from the ANN for the temperature of 300 K shows some average behavior between the experimental results and the ab initio calculation.Nevertheless, the shear modulus for the Al-Cu alloy grows faster with an increasing copper concentration according to the ANN predictions.This may be due to the fact that the density of this material predicted by the ANN increases somewhat faster compared to the experiment (Figure 13), which can explain the faster growth of both moduli.

Shear Stress Relaxations Due to Dislocation Activity and Phase Transitions 4.1. Theoretical Model of Dislocation Plasticity and Phase Transitions
The results of the MD simulations show the relaxation of shear stresses in the process of uniaxial compression of single crystals of the Al-Cu solid solution.The presence of a stacking fault and the formation of dislocations indicate the plastic response of the material mediated by dislocation activity, while the evolution of phase fractions accompanies phase transitions during deformation.The propagation of a plane shock wave caused, for example, by detonation or a high-velocity impact, increases pressure up to tens of gigapascals in the condition of uniaxial deformation.Such severe deformation should be analyzed in terms of finite strains.Following the idea of multiplicative decomposition [40,[104][105][106], the macroscopic strain gradient F is represented as follows: where F e is the elastic part of the deformation to be defined for the stress calculation, while F p is the contribution of plastic (irreversible) deformation due to the dislocation activity and F tr is the contribution of phase transitions.In the considered case of uniaxial deformation, the macroscopic strain gradient can be calculated as follows: where F 11 is the ratio of the current length L x along the deformation axis to the initial length L 0 .The atom rearrangement during phase transition is equivalent in a certain sense to an imposed deformation.Therefore, the contribution of phase transition can be represented as a linear combination of the atomic fractions of the emerging phases (ω BCC for BCC and ω OTHER for disordered atoms-"OTHER" phase): Here, we assumed F tr 11 F tr 22 F tr 33 = 1 and did not consider the volume change at the phase transitions, because this volume change influencing the pressure was already taken into account in the ANN-based EOS.In Equation ( 14), n 1 and n 2 are the model parameters to be fitted to the MD data.
The next step is the equations of balance among the atomic fractions of the FCC, BCC and OTHER (disordered or amorphous) structures.The fraction of the FCC structure decreased due to the formation of the BCC and OTHER structures.An analysis of the MD data showed that the mutual transformation of the BCC and OTHER structures is observed during compression.Consequently, the rate of change in phase fractions in the material during deformation can be written as follows: Here, ω FCC represents the sum of atomic fractions of both the FCC and HCP phases, because the HCP phase arises as a result of dislocation activity inside the FCC phase, but not due to a phase transition in the classical sense.The parameters in the expression for the coefficients in linear combinations depend on the formation of a certain structure in the phase transition.Based on the Arrhenius equation, it can be written as follows: 1, at the transition FCC → BCC 2, at the transition FCC → OTHER 3, at the transition BCC → OTHER 4, at the transition OTHER → BCC , (16) where K 0β is the phase formation rate parameter, U β is the threshold energy, V 0 = 10 −30 m 3 is the activation volume, σ 11 and σ 22 are the components of the Cauchy stress tensor and ε 1β and ε 2β are the characteristic strains of phase transition along the loading direction and transverse direction, respectively, such that V 0 ε 1β σ 11 + 2ε 2β σ 22 is a positive work of phase transition.The parameters of Equation ( 16) are also fitted to the MD data.
The components of the plastic strain gradient tensor are defined through the total plastic strain w: where 1 holds because of the volume conservation during the dislocation slip.Dislocation-mediated plastic deformation occurs, when shear stress exceeds the threshold yield stress Y b /2, but the inertness of the dislocation system holds the shear stress level above this threshold under dynamic loading.This inertness can be taken into account through the relaxation time approach [107].Due to the evolution of the dislocation structure, the characteristic time of stress relaxation changes with the dislocation density [108], such that the following modified Maxwell model can be implemented [63,109]: where S is the stress deviator, H(•) is the Heaviside step function, ρ D is the scalar density of dislocations and b is the modulus of Burgers vector.According to the results of the MD simulations [110], the dynamic drag coefficient of a dislocation B is directly proportional to temperature; however, following [111], a nonlinear growth of the friction coefficient is observed for copper and aluminum when approaching melting temperatures.Now we considered the evolution equation for the dislocation density.The first process to be considered was the dislocation nucleation.The formation of a dislocation loop of radius R and length 2πR in each elementary volume R 3 of a metal is a random process occurring in a single crystal with frequency c t /R due to thermal fluctuations, and its rate per unit volume can be represented as follows [79,84]: where Q ε is the nucleation strain distance function [84], which shows the threshold of dislocation nucleation as discussed in Section 3, c t = G/ρ is the transverse speed of sound, G is the shear modulus, ρ is the material density and k B is the Boltzmann constant.Following [62], the nucleation energy barrier is presented as W = k n Gb 3 , where k n is a parameter.Expressing the concentration of atoms through the density ρ of the solid solution and the average mass of one atom m 1 , we wrote down the dislocation nucleation rate as follows: The nucleation threshold and shear modulus in our work were calculated using an ANN (see Section 3).
A certain part of plastically dissipated power was spent on the formation of new dislocations.This idea was used in [112] to formulate the dislocation multiplication rate, and we applied a similar approach: where ε D is an efficient energy of the dislocation multiplication per length b of the new dislocation line.The annihilation rate of dislocations was written down in the standard form: where K a is the annihilation coefficient.As a result, the expression of dislocation kinetics took the following form: .
The static yield strength was calculated using Taylor's law of hardening: where A is the hardening parameter.The stresses can be represented through the spherical and deviatoric stress components Usually, the pressure is calculated from the equation of state.In our work, we used an ANN (see Section 3).The input parameters of the ANN were the atomic concentration of copper in solid solution, temperature and current density of the substance: From Equation ( 12), we expressed the elastic part of the strain gradient and found the Cauchy-Green strain tensor: Thus, we obtained expressions for the stress deviator and shear stresses using Hook's law [113]: In the formulated model of shear stress relaxation, it was necessary to determine a set of parameters: K 0β , U β , V 0 , ε 1β and ε 2β in the model of phase structure evolution (Equation ( 16)), as well as n 1 and n 2 in the model of phase transition-induced stress relaxation (Equation ( 14)) and B, Y 0 , A, k n , ε D and K a in the model of dislocation plasticity.

Optimization and Results of Phase Transition Model
The parameter identification method based on the Bayesian approach is a powerful tool of machine learning [85,114].It involves calculating the probability for many sets of randomly chosen parameter values.The probability characterizes the degree of agreement of the calculated values with the reference data.The probability is taken to be equal to one at the beginning of the modeling.In the process of modeling, the difference between the compared values is accumulated and the probability decreases correspondingly.The parameters are selected in such a way that the resulting probability is maximal for the best parameters of the model.In this case, the modeled curve of the compared quantity will be as close as possible to the reference curve.For the phase composition evolution model, the probability characterizes the degree of correspondence of the BCC and OTHER (disordered or amorphous) phase fractions calculated from the system of Equations ( 15) and ( 16) with MD data.The probability was calculated as the product of the probabilities at each modeling step i: where ω i1 and ω i2 are the atomic fractions of the BCC and OTHER phases, ∆ω MD j are the ranges of values of fractions in MD simulations and k is the normalization factor.We considered the HCP phase as a defective structure in the FCC lattice.Distribution maps were constructed to determine the accuracy, range and step of parameter variation.Figure 15 shows the parameter and probability distributions in a solid solution with 70% Cu at 300 K.It can be seen from the graphs that the set of parameters with maximum likelihood is not unambiguously defined.That is, insignificant fluctuations in the value of parameters will not significantly affect the shape of the curves of the distribution of the phase fractions.The set of parameters corresponding to the highest probability was further used to solve Equations ( 15) and (16).The calculated values of the phase fractions during uniaxial compression for the alloy with 70% Cu are in agreement with the MD modeling data at 300 K (Figure 16).The model describes the transition from the FCC to the BCC structure and regions with disordered atoms, as well as the further disordering of the BCC phase.

Optimization and Results of Stress Relaxation Model
The degree of agreement of shear stresses and pressures with the results of the MD calculations was estimated in order to determine the contribution of phase transitions to the stress relaxation process (parameters n 1 and n 2 ) and the parameters of the dislocation plasticity model.The shear stresses were calculated from the stress relaxation model and the pressure from the equation of state in the form of the ANN.The expression for calculating the probability has the following form: where ∆τ MD and ∆P MD are the ranges of shear stress and pressure values, respectively, in MD simulations.Figure 17 shows the probability maps for different parameter ranges for an alloy with 10% Cu at 300 K.The variation in values of yield strength Y 0 and dislocation friction coefficient B almost does not affect the probability value.A similar low sensitivity of the results of the dislocation-based model of pore growth to the dislocation friction coefficient was previously reported in [115] in the case of iron.On the contrary, favorable values of the parameter n 1 are localized in a narrow range.For the hardening coefficient A, nucleation parameter k n and dislocation loop formation energy ε D , there are regions for which the probability takes zero value.In this case, the shear stress plot is far from the reference curve.Let us consider the results of the shear stress relaxation model with the fitted parameters by the Bayesian method by an example with a copper concentration of 10% presented in Figure 18.The model qualitatively well describes the stresses in the material under uniaxial compression (Figure 18a).The difference is in the amplitude of the maximum shear stresses.The model also indicates hardening in the material at strains 0.1-0.2 and a sharper drop in stresses at the initiation of dislocation activity (at a strain of about 0.09) and phase transitions (at a strain of about 0.17).The ANN describes the pressure in the material with good accuracy (Figure 18b).Another ANN calculates the shear modulus and nucleation threshold of dislocations.After the nucleation of dislocations, there is an active multiplication of dislocations at a rate of the order of 10 27 s −1 m −2 (Figure 19).Due to the dynamic loading with an ultra-high strain rate applied to the material, the dislocation density exceeds 10 16 m −2 , which is in the range of the experimental data for Cu [116] and MD modeling [117,118].Our model value of the dislocation density in Figure 19a is an order of magnitude lower than that determined by the MD in Figure 3, but this discrepancy is explained by the fact that we took into account only mobile dislocations in Equation ( 18) of the relaxation plasticity model, while most of the dislocation segments in the MD are immobilized due to intersections with other dislocations.More complex dislocation kinetics with an explicit separate accounting of mobile and immobilized dislocations similar to [63] can be applied in future work.The plastic strain rate increases as the dislocation density increases dramatically, but after an initial spike at a strain of about 0.1, it drops by an order of magnitude.

Conclusions
The MD study of uniaxial compression of Al, Cu and Al-Cu solid solutions was performed.The stress relaxation process is provided by (a) dislocation plasticity and (b) phase transition.There is a sharp drop in stresses in pure metals due to the initiation of plastic deformation induced by dislocation nucleation and motion.The formation of a metastable BCC phase is also observed, which precedes the onset of the plasticity process.In contrast to copper, aluminum reveals a phase transition from the FCC to the BCC structure almost throughout the entire volume of the simulated crystal at pressures of about 36 GPa.The phase transition is accompanied by a decrease in shear stresses to negative values.A similar behavior is observed for Al-Cu solid solutions, but there is a difference in that the phase transition becomes the dominant stress relaxation mechanism.The phase transition completely suppresses the dislocation activity for alloys with Cu concentrations in the range 30-80%, while the dislocation plasticity prevails again for pure Cu.An increase in temperature leads to a softening of the material and a decrease in the amplitude of the shear stresses.
A theoretical model of stress relaxation was developed, taking into account the main mechanisms of plasticity, the dislocation activity and phase transitions, by an example of an Al-Cu solid solution.The Bayesian method was applied to identify the model parameters.Two forward-propagation ANNs trained by MD data for uniaxial compression and tension were used to approximate the constitutive relations.The ANNs were used to describe the equation of state and determine the beginning of plasticity and the elastic modulus.
Further development of the present work is connected with using the developed model for the simulation of a shock-wave structure in an Al-Cu solid solution, as well as with the application of the model to other metallic systems, including high-entropy alloys.Most of the solid solutions considered here are thermodynamically unstable in traditional metallurgical casting, and experimental data are still lacking, but if they reveal unique dynamic-protecting properties, novel production techniques, such as additive manufacturing, can be applied to produce such metastable alloys.Therefore, a theoretical analysis of their properties is relevant.
Author Contributions: N.G.: methodology, software, validation, formal analysis, investigation, writing-original draft preparation and visualization; E.F.: methodology, software, validation, formal analysis, investigation, writing-original draft preparation and visualization; A.M.: conceptualization, methodology, resources, data curation, writing-review and editing, supervision, project administration and funding acquisition.All authors have read and agreed to the published version of the manuscript.

Figure 1 .
Figure 1.Schematic representation of crystal loading in MD modeling.The compression direction [100] is shown by arrows.

Figure 2 .
Figure 2. Shear stresses and evolution of phase fractions for pure aluminum (a) and pure copper (b) during compression.

Figure 3 .
Figure 3. Dislocation density evolution for (a) pure aluminum and copper and (b) Al-Cu solidsolutions during compression at the temperature of 300 K obtained from the MD.The dislocation length in FCC phase was divided by the total volume of the system; therefore, it tends to zero after the phase transition from FCC to BCC and "Other" structures.

Figure 5 .
Figure 5. Evolution of phase fractions and shear stresses in a solid solution with 10% Cu in aluminum.The atoms are colored according to the legend in the graph; the atoms forming FCC structure are not shown.Adapted from ref. [80].

Figure 7 .
Figure 7. (a) Distribution of shear stresses at different temperatures in solid solution with 20% Cu and (b) evolution of crystal structure of solid solution with 10% Cu at different strain rates.

Figure 8 .
Figure 8.The first part of constitutive equation of Al-Cu solid solution in the form of fully connected forward-propagation ANN that calculates the pressure, total energy and bulk modulus.

Figure 9 .
Figure 9.The second part of constitutive equation of Al-Cu solid solution in the form of fully connected forward-propagation ANN that calculates the shear modulus and nucleation strain distance function.

Figure 10 .
Figure 10.The dependence of pressure on density for intermediate values of copper concentrations predicted by the ANN.

Figure 13 .
Figure 13.The comparison of density of the Al-Cu alloy with solid solution of copper atoms obtained in this work and the experiment of Senoo-Hayashi (1988) [102].

Figure 14 .
Figure 14.Bulk modulus (a) and shear modulus (b) of Al-Cu solid solution for low concentrations of copper atoms in the aluminum matrix predicted by the ANN and compared with the experimental data at room temperature fromSenoo-Hayashi (1988) [102] and ab initio calculations from Ma (2012)[103].

Figure 15 .
Figure 15.Probability distribution in the parameter space in the case of 70% Cu for the following couples of parameters: (a) ε 23 −ε 12 and (b) ε 11 − U 1 .

Figure 16 .
Figure16.Comparison of the evolution of phase fractions with parameters fitted by Bayesian method with the results of MD calculations for an alloy with 70% Cu, where the sum of FCC and HCP phases is in green, BCC phase is in blue and disordered (OTHER) structure is in grey.

Figure 17 .
Figure 17.Probability distribution in the parameter space in the case of 10% copper for the following couples of parameters: (a) dislocation formation energy-hardening coefficient, (b) dislocation formation energy-initial yield strength, (c) dislocation friction coefficient-dislocation nucleation parameter and (d) dislocation friction coefficient-parameter of shear stress relaxation due to FCC-to-BCC phase transition.

Figure 18 .
Figure 18.Comparison of (a) shear stresses and (b) pressures from MD simulations with the results from a stress relaxation model with parameters fitted by Bayesian method for a copper concentration of 10% at 300 K.

Figure 19 .
Figure 19.(a) Plastic strain rate and dislocation density and (b) nucleation, dislocation multiplication and annihilation rates during uniaxial compression of Al-Cu solid solution with 10% Cu at 300 K.The dislocation density ρ D in FCC phase is multiplied by the phase fraction ω FCC in (a) in order to obtain an averaged density through the whole material volume.

Table 1 .
Main characteristics and hyperparameters of the ANNs and the distribution of training and test data sets.

Table 2 .
The mean absolute percentage error (MAPE) in the final configuration of both ANNs.