Innovative methods for random field establishment and statistical parameter inversion exemplified with 6082-T6 aluminum alloy

This paper aims to eliminate the disharmony between simulation and experiment, and takes the mechanical properties of 6082-T6 Al alloy as an example. In order to obtain the equivalent distribution of material properties after considering the randomicity of materials, a new inversion method combining with stochastic finite element method (SFEM) is proposed. Besides, the discrete random field in SFEM is established by an innovative method to overcome some discretization difficulties in conventional methods. In summary, the generic methods proposed in this study can give a new solution for the correlation of meso-structure and macro-performance in computational materials science.

The numerical simulation of macroscopic material properties plays an important role in material design and application 1 . Undoubtedly, the accuracy of numerical simulation is the focus of attention in computational materials science 2 , and it depends largely on the correctness of input material parameters 3 . However, many materials have the characteristics of heterogeneity, nonlinearity and discontinuity 4 , which makes it difficult to directly measure material properties. Moreover, the stress-strain curve obtained directly from a specimen experiment is not identical to the true constitutive curve, because of the uneven stress distribution caused by the boundary conditions and the large deformation (such as the necking behavior) of the specimen 5 . Therefore, it is necessary to reasonably estimate actual material parameters by inversion methods, which is precisely the parameter estimation problem in inversion theory 6 . For example, some studies combine the finite element method with simple iterative algorithm to calculate the material's actual stress-strain curve 7 , but on the assumption of homogeneous material.
Thus far, previous studies have failed to consider the randomicity and heterogeneity of material in the inversion procedure. Besides, the inversion result under the homogeneity assumption is applicable only to the particular inversion object and cannot be applied to similar experiments even made of the material of the same batch, due to the material's natural randomicity. Taking aluminum (Al) alloy as an example, its heterogeneity and randomicity caused by the stochastic meso-structure indeed should not be neglected according to the experimental results of dozens of specimens 8 . Additionally, the meso-structure of materials is commonly simulated by molecular dynamics simulation 9,10 , and it can partially describe the macro-mechanical behavior but usually on the assumption of periodic boundary, which deviates significantly from the real heterogeneous materials with complex spatial distribution 11 . Thus for the simulation of large stochastic structures, such as buildings and airplanes, producing equivalent material properties at macro scale is an effective way, just as the stochastic finite element method (SFEM) does 12 .
An essential part of SFEM is to establish discrete random fields of heterogeneous materials. The conventional means is to approximate random fields with continuous functions, and then descretize them into finite element mesh by the methods of local average, optimal linear estimation, etc 13,14 . However, this is incapable of simulating the failure behavior in high gradient areas, such as strain localization and stress concentration caused by local damage (slip, cracks, voids, etc), as a result of continuity assumption. Moreover, the determination of material's statistical parameters usually relies on rough estimation and lacks reasonable inference, which should return to the inversion method mentioned above. Experiments and constitutive model of Al alloy 6082-T6. The tensile tests of 6082-T6 alloy are designed according to ISO standard 'Metallic materials -Tensile testing' 18 , and Fig. 1(a) gives the dimensions of Al alloy specimens. The tests are carried out on ZWICK electronic universal material testing machine ( Fig. 1(b)), and the deformation rate is controlled by the feedback data of extensometer to load slowly. A total of 117 tests for tensile mechanical properties have been carried out, and Fig. 1(c) shows the ultimate failure of some specimens with different fracture locations.
Then, the stress-strain curves of Al alloy specimens in 117 tests are given in Fig. 2. Besides, the conventional method to calculate the stress σ and strain ε is expressed by the equations where F is the applied force and A 0 denotes the area of the middle cross-section. ′ l 0 and l 0 represent the effective length of the test specimens (the distance between the two fixtures) after and before stretching, respectively. It shows that the stress-strain curves of Al alloy are smooth and gradual without a yield platform like low-carbon   steel. By the way, the constitutive relation of 6082-T6 obtained by this method is not accurate enough, which will be discussed in the next section. In addition, all the strains used in this paper are of engineering strain.
In order to describe the constitutive relation of Al alloy, many constitutive (stress-strain) models have been put forward. Thereinto, Ramberg-Osgood model is the most widely used one because of its clear physical meaning and nearness to experimental stress-strain curves, and it is shown as follows: where E 0 is the elastic modulus of Al alloy, and f 0.2 is the stress value corresponding to the plastic strain of 0.2%. The exponent n represents the strain hardening degree, which controls the trend of the constitutive curve. And initially, n was determined by two parameters: = .
. in MPa, ant it has been widely used in many papers. Consequently, the simplified formula is used in this paper to reduce variables, and the constitutive curve of Al alloy can be described by only two parameters E 0 and f 0.2 (also known as yield strength f y ).
Statistical inference of the experimental results. The statistical histogram of elastic modulus E 0 and yield strength f y are illustrated in Fig. 2. The tensile specimen can be regarded as a set of parallel member bars (independent of each other). The elastic modulus E 0 of the whole specimen is the average value of all members' E 0 , so it should be of Gaussian distribution according to the central limit theorem. As for yield strength f y , a yield member will distribute its force to remaining members, and accelerate the failure of the whole specimen, thus a minimum extreme distribution is suitable for f y . Moreover, the Anderson-Darling statistical test at the significance level of 0.05 demonstrates that the material parameters E 0 and f y follow the Gaussian distribution and Gumbell distribution (extreme value distribution Type-I), respectively.
There are two statistical parameters, location parameter μ and scale parameter δ of the both types of distribution. And their probability distribution functions are shown as follows:

Methods
However, The stress-strain curves shown in Fig. 2 are the output of experimental specimens rather than the input of material, and they reflect the stress and strain level of the whole specimen instead of every point. Therefore, the stress-stain curve obtained by the conventional method should not be simply deemed to the constitutive relation. For example, a finite element model of the Al alloy specimen (element size 4 mm and total 880 elements) is established by software ANSYS as shown in Fig. 3(a). After extracting the end force and displacement of the specimen from the simulation results, the apparent material properties of the simulated specimen can be output using the same Eq. 1. When assuming the model is of homogeneous material and the parameters E 0 and f y take the mode www.nature.com/scientificreports www.nature.com/scientificreports/ values of the above statistics, the output elastic modulus ′ = .
E 83 06 GPa 0 is 17.8% greater than the input E 0 , while the output yield strength ′ f y slightly declines by 0.4 MPa. Besides, the influence of finite element discretization on output properties is also studied by changing the element size from 2 mm to 1 mm. It is found that both the ′ E 0 and ′ f y decreased by less than 0.1%, which are much smaller than the discrepancy between input and output. And this discrepancy is mainly caused by the shape and boundary conditions of specimen which makes uneven stress distribution in the specimen. Therefore, the conventional method is not accurate enough and this is also an important reason why the simulated engineering structure, which directly input the apparent stress-stain relationship, is usually more rigid than the actual structure. By the way, all the simulations in this article are of small strain analysis performed in ANSYS, that is, each element uses engineering stress/strain (rather than true stress/ strain) which is consistent with the output.
As stated in the introduction, some studies have combined the finite element method with simple iterative algorithm to calculate the actual material properties but on the assumption of homogeneous material. However, according to the failure photograph ( Fig. 1(c)) and the statistical material properties of the experimental specimens, the heterogeneous of Al alloy is indeed not negligible. Besides, the distributions in Fig. 2 also cannot actually describe the actual stochastic characteristics of material properties due to the reason proved above. Therefore, an inversion method for obtaining effective material properties needs to be investigated.
Method for establishing random field. When the material of a finite element model changes from homogeneous to heterogeneous, complex random meso-structures are formed inside the model. And an appropriate random field is necessary to describe the spatial variability of materials in SFEM. One control parameter of random field is the correlation distance 22 , and two points beyond this distance can be considered almost irrelevant. In fact, the value of this distance is often smaller than the size of elements, so all the elements can be considered mutually independent. Accordingly, a direct method is to distribute materials to each element in a purely random way, but chessboard lattice phenomenon is prone to occur in this case. For instance, a five elements model is given in the left side of Fig. 3(b), where a weak element (represents a small cavity) is surrounded by normal elements. And from the vertical displacement contours it can be seen that the cavity element of the original model is unreasonably rigidized compared with the right-side model with a refined mesh.
Therefore, it is necessary to form appropriate element clusters with at least two elements of the same material. Firstly, considering the assumption of isotropy Al alloy, Fig. 3(c) gives six types of element clusters of rotational symmetry. Next, the clusters are randomly distributed as densely as possible in the study area in a mutually exclusive manner, and then each remaining blank element is randomly merged into its adjacent cluster to make the whole study area fully filled. The final step is to assign material properties to element clusters. E 0 or f y is sampled randomly from its probability distributions. Consequently, a discrete random field is directly established with finite element mesh by this innovative method which can be called stochastic element cluster method (SEC).

Inversion algorithm of statistical parameters.
Based on the discrete random field established above, the inverse problem of statistical parameters can be described below: Supposing there exists a functional relationship between the input and output material parameter matrices: Λ Λ = ′ g( ) , the inverse problem is to solve Λ Λ = ⁎ to make Λ ⁎ g( ) equal to the experimental output matrix Λ′ ⁎ . The dimension of the matrices is m × 2, where m denotes the number of material properties and each of them has two statistical parameters. Besides, Λ′ is calculated from the following formula: represents the estimation of two statistical parameters of each simulated property in the simulation result matrix X; F R ( ) i j is the finite element simulation of the i-th output parameter under the random field of the j-th random simulation, R j ; s is the total count of random simulations; Y k,q denotes the random variable of the q-th material property on the k-th element cluster and it subjects to the probability distribution function P x ( ) q , and n is the total number of element clusters in random specimens. In addition, f E E f f in the specific problem in this paper, and similarly μ δ μ δ , ] E E f f . Unlike conventional inverse problems, the inversion variables in this problem are statistical parameters which need a number of samples to estimate them, and the function relationship between input Λ and output Λ′ is implicit. These results in the uncertain output for a certain input. In view of this, this paper presents an effective iteration algorithm for solving statistical parameters, named as penalized regression secant (PRS), and the algorithm needs to address the following key challenges: • Function approximation. Regression method is introduced to estimate the implicit function using known iteration points. Besides, the relationship between input and output of each material parameter in this problem is monotonic in theory, so linear regression is enough to guarantee convergence. The regression line of statistical parameter λ at the t-th iteration can be expressed as λ β λ α ′ = + www.nature.com/scientificreports www.nature.com/scientificreports/ • Iteration method. The secant method is a classical and efficient iteration algorithm, so we develop it to approximate the ideal solution Λ Λ′ ⁎ ⁎ ( , ). Unlike the conventional method, the penalized regression line is applied as the secant and every regression uses all the previous iteration points.
• Simulation times. After repeated trials, it shows that 20 random simulations per iteration are enough to guarantee convergence and convergence faster. • Solution order. Since the elastic modulus unilaterally affects the yield strength, the yield strength can be calculated subsequent to the solve of E 0 finished. • Convergence criteria. Since the dependent variable is uncertain, use independent variable to define convergence criterion. For example, when the input parameter vector λ t ( ) of E 0 or f y satisfies λ where  denotes Hadamard division and x 1 is the Taxicab norm, it can be considered convergent at the t-th iteration.
In addition, the pseudo-code of the PRS algorithm for solving statistical parameters is given in the Appendix.

Results
As a result, Fig. 4(a) shows a set of random material distribution obtained from SEC, and it can be seen that there are no isolated elements except at the boundaries of the study area. After simulating, its contour diagram of Von Misses stress is shown in Fig. 4(b). Moreover, the failure mode is similar to the experiment results shown in Fig. 1(c), which illustrates the effectiveness of the discrete random fields established by SEC. Taking the experimental statistics result (shown in Fig. 2) as the input material distribution, the discretization errors of random element method (RE) and SEC are shown in Fig. 4(c) by changing element size from 2 mm to 1 mm. Each value in the figure is the location parameter of 50 random simulation results, μ′ f . It can be seen that the undesirable rigidification of chessboard lattice is alleviated in SEC. And the lower drop of SEC with mesh refining indicates that the assembly of element clusters can significantly reduce the discrete error.
In summary, the SEC method can obtain excellent simulation results with fewer elements, and it can manually set the shape of element cluster types to imitate specific meso-structures, such as the layered material. Considering the computational efficiency, the number of elements (average 3.5) within each cluster in this paper is sufficient for accurate simulation. Moreover, there exists a considerable gap between the input and output distribution of material properties, for example, the location parameter of ′ f y is 281 MPa (given in Fig. 4(c)) while that of input f y is 305 MPa. Therefore, the PRS inversion method is performed to obtain the equivalent input parameters of material distribution.
The inversion process of f y is shown in Fig. 5(a), and the regression lines of the first three steps are also shown in the figure Figure 5(c) shows the twenty inversion results of f y . The maximum coefficient of variation C v is less than 4%, which indicates the effectiveness of the algorithm, while that of E 0 is smaller.  www.nature.com/scientificreports www.nature.com/scientificreports/ In addition, the influence of element size to the inversion results is also studied. As listed in Table 1, with the decrease of element size, the μ ⁎ E fluctuates slightly while the other three variables have an increasing trend with varying degrees. The material distributions with smaller elements have higher dispersion, which is precisely the reverse process of homogenization from meso to macro. As a result, the element size largely affects the inversion results, which can be determined by considering computational efficiency and simulation accuracy. By the way, the size of element clusters should be greater than the correlation distance.
In summary, the SEC method can qualitatively describe the stochastic characteristic of specimens using element cluster to imitate meso-structure. Furthermore, the equivalent material statistical input parameters are obtained by the PRS inversion algorithm, which further quantitatively ensures that random simulations can accurately reflect the macro-mechanical behaviors of the specimens. In addition, the SEC method can not only simulate the stress behavior in high gradient regions, but also facilitate the processing of irregular curved surfaces and non-Gaussian fields which are hard for conventional methods 23 .

Application
Obtained from the above two methods, the equivalent stochastic constitutive relation and material spatial distribution can accurately reflect the statistical regularity of the properties of 6082-T6 material. Moreover, it can be applied to larger structures composed of the same material. For example, Fig. 6(a) shows a box beam of 6082-T6 Al alloy with the cross-section of 50 × 100 × 4 (width, height and thickness in millimeter). The beam model with a length of 1.5 m is loaded at three points and it has 100 mm margin on both sides. Then, a finite element model is built accordingly to analyze the mechanical properties of the experimental beam. Loading device of the experiment is simulated accurately by establishing a rigid plate and reasonable contact as shown in Fig. 6(b). In order to guarantee simulation accuracy and improve computational efficiency, fine grids with 2.5 × 2.5 mm 2 are divided into the middle part and 5 × 5 mm 2 in the other parts. The corresponding statistical parameters of these two sizes are listed in Table 1, respectively. In addition, the first and second order buckling modes with the amplitude of 3 mm are added to the model as initial geometric imperfections. Figure 6(b) shows the Von Mises stress contours and deformation of one random simulation. It can be seen that its local buckling is very similar to the experimental photographs of the top and front views, which shows the effectiveness of the simulation. Besides, 20 groups of random simulations are carried out and their simulated loaddisplacement curves are shown in Fig. 6(c). And it can be found that a simulated curve (red solid line) which is very close to the experimental data (black dots). Furthermore, statistical analysis of stiffness coefficient and failure load of the structure (corresponding to elastic modulus and yield strength of the constitutive curve) shows that they still maintain normal distribution and extreme value distribution respectively. As a comparison, the conventional simulation is conducted which applies average values of experimental data of specimens as the constitutive parameters. Figure 6(c) gives the simulation result and it can be seen that the elastic part is obviously stiffer than the experiment while the failure part is lower. Therefore, the correct constitutive parameters obtained from the proposed methods can solve the common problem in engineering that the simulated structure more rigid than the real structure.  www.nature.com/scientificreports www.nature.com/scientificreports/ Moreover, the accurate random simulations can replace repeated experiments, and it can greatly save experimental cost. In addition, the statistical analysis of the simulated results can provide an accurate and direct basis for the reliability design of engineering structures. Furthermore, the input parameters of other physical properties (such as thermal conductivity, electrical conductivity, etc.) can also be inverted by the methods.

Conclusion
In order to make the random simulation of specimens consistent with the experimental data, the SEC method for constructing discrete random fields and the PRS iteration algorithm for statistical parameter inversion are established. Besides, the two new methods can effectively solve the key problems in SFEM: random field discretization and parameter determination. And they together ensure that the element clusters are equivalent to the random meso-structures of specimens, both qualitatively and quantitatively. And the equivalent stochastic material distribution can be used to accurately simulate larger random structures of the same material. In addition, other physical properties of materials can also apply the methods to connect meso-structure and macro-performance.