Probabilistic Upscaling of Material Failure Using Random Field Models – a Preliminary Investigation

Complexity of failure is reflected from sensitivity of strength to small defects and wide scatter of macroscopic behaviors. In engineering practices, spatial information of materials at fine scales can only be partially measurable. Random field (RF) models are important to address the uncertainty in spatial distribution. To transform a RF of micro-cracks into failure probability at full structural-scale crossing a number of length scales, the operator representing physics laws need be implemented in a multiscale framework, and to be realized in a stochastic setting. Multiscale stochastic modeling of materials is emerging as a new methodology at this research frontier, which provides a new multiscale thinking by upscaling fine-scale RFs. In this study, a preliminary framework of probabilistic upscaling is presented for bottom-up hierarchical modeling of failure propagation across micro-meso-macro scales. In the micro-to-meso process, the strength of stochastic representative volume element (SRVE) is probabilistically assessed by using a lattice model. A mixed Weibull-Gaussian distribution is proposed to characterize the statistical strength of SRVE, which can be used as input for the subsequent meso-to-macro upscaling process using smeared crack finite element analysis.


Introduction
Complexity of failure is reflected from sensitivity of strength to small defects and wide scatter of macroscopic behaviors.In engineering practices, spatial information of materials at fine scales, such as atomic defects, dislocations, locations and sizes of micro-particles, micro-voids, micro-cracks, grain boundaries, can only be partially measurable.Random field (RF) models are important to address the uncertainty in spatial distribution.Characteristics of a random field at a fine scale contribute to the pattern of the random field at the next coarse scale, and the contribution is governed by physical laws connecting the two scales.For instance, a failure process initiated from micro-crack propagation is critically determined by the laws of fracture mechanics and the probability distributions for the size and density of micro-cracks.To transform the random field of micro-cracks into failure probability at full structural-scale crossing a number of length scales, the operator representing physics laws need be implemented in a multiscale framework, and to be realized in a stochastic setting to deal with statistical input data.Multiscale stochastic modeling (MSM) of materials is emerging as a new methodology at this research frontier [1], which provides a new multiscale thinking by upscaling fine-scale random fields.Along this line, theoretical and computational MSM methods has been developed on stochastic homogenization [17,19,20], variational principles [18], and multiscale stochastic finite element methods (MsSFEM/MSFEM) [1,21].Other studies on multiscale stochastic models include Asokana and Zabaras [13] on random media diffusion problems, Koutsourelakis [14] on characterization of random microstructure, Chakraborty and Rahman [15] on fracture of functionally graded materials, Pugno et al [16] on strength of carbon nanotube composites, etc.
A classical example of the multiscale stochastic modeling on materials is Weibull theory that offers a closed-form solution for brittle materials to quantify evolution of failure probability passing large scale separation.Failure probability of most non-brittle materials however is not amenable to asymptotic theoretical/analytical approaches and can only be explored by means of numerical experiments.Nonlinear finite element reliability analysis, using the midpoint or average random field discretization method and first-or second-order reliability method, has been developed since 1990s to evaluate failure behaviors of quasibrittle materials [2,3].Attempts have been made to apply the spectral stochastic finite element method to plasticity and failure [4].To solve multiscale failure problems such as concrete structures, geomaterials, and earth faults with a domain size several orders of magnitude larger than the correlation length, novel multiscale schemes are desired.
In this study, a preliminary framework of probabilistic upscaling is presented for bottom-up hierarchical modeling of failure propagation across micro-meso-macro scales.In the micro-to-meso process, the strength of stochastic representative volume element (SRVE) is probabilistically assessed by using a lattice model to simulate random micro-cracking.A mixed Weibull-Gaussian distribution is proposed to characterize the statistical strength of SRVE, which can be used as input for the subsequent meso-to-macro upscaling process using smeared crack finite element analysis.

Micro-cracking in a random field
Deterministic multiscale models, hierarchical or concurrent, use representative volume element (RVE) as a vehicle to transport local effects into global solutions.In stochastic problems, the so-called stochastic RVE plays a similar role but rather than transporting deterministic or mean values, as implied in deterministic models, it transports uncertainty information, e.g.probability distribution.
In the micro-to-meso upscaling process, the RFs of SRVE strength and fracture energy can be used to account for the profusion of microscopic flaws and cracking.Local cohesive laws, including fracture surface energy and strength, at a certain fine-scale resolution are assumed to be intrinsic material properties independent of damage evolution.Let the coordinates y at micro-scale be defined as x/ε, where ε is a resolution parameter representing the ratio between macro-scale coordinates x and micro-scale coordinates y (Fig. 1 where ) , is the strain energy release rate that varies spatially dependent on local material properties.The micro-cracking in a two-dimensional RF will choose favorite initiation and propagation orientation such that the decrease of total system free energy is maximized where the path da follows the orientation  , and a  is the crack propagation incremental length chosen in a particular simulation.In a discrete model, e.g., the lattice system in Fig. 1, the path da becomes discrete, i.e. a  , and the latter corresponds to the lattice spacing between a broken bond and its nearest neighbors.
It is noted that the idea of energy minimization have appeared in a number of publications, e.g. the variational model in [12].The attempts to generalize Griffith fracture theory remain largely at mathematical formulations, and little work has been done to provide practical tools for numerical implementation.The initiation criterion (1) and the propagation criterion (2), explicitly given in our model, are directly applicable for crack simulation in continuum media [9] and discrete media as given in this study.

Figure 1. Schema of probabilistic upscaling.
A number of studies have been conducted on stochastic microcracks [5,22,23].To simulate the micro-cracking in a RF, one of the most effective numerical approaches is to use lattice models that allow disorder to be introduced naturally.Various types of discrete lattice models such as central force model, electrical fuse model, bond-bending model, and beam-type model have been used to study progressive damage of disordered quasi-brittle materials [5].In this study a random spring model is used and the network is constructed from a lattice where all the bonds between nearest-neighbor sites are Hooke's law springs with unit cross-section area.A bond i can have one of two states as follows: and the Hamiltonian of an N-bond lattice system is given as where i  , i k , and i s denote the deformation, the elastic constant, and strength of bond i, respectively.An intact bond i absorbs the driving energy i g until it reaches the rupture threshold i r .
At zero "temperature", the equilibrium of the system is obtained by minimizing the Hamiltonian in Eq. ( 4), which is a convenient postulate we use to model irreversibility of micro-cracking (see [12] and references therein).
For brittle bonds defined above, the local strength and fracture energy are linked via the local elastic constant.When local elastic constant is uniform throughout the domain, the two RFs will be completely correlated.For non-brittle bonds involving mechanisms such as plasticity, softening, etc, the RFs of strength and fracture energy are partially correlated, which will not be discussed in this study.
To simulate micro-cracking in a RF, we need to generate RF samples according to given probability distributions.An overview of the random field simulation is presented in [6].In this study a non-Gaussian RF sample is generated by point-wise translating an underlying Gaussian RF sample, which is the so-called translation method (see [6]).
Assume the marginal distribution of the microscale strength RF follows the Weibull distribution.A Weibull RF sample Y is generated from an underlying Gaussian RF sample where is the standard normal cumulative density function (cdf), and is inverse of the Weibull cdf.The correlation function of the underlying Gaussian RF is assumed to be exponential where d indicates the correlation length.The simulation of Gaussian samples follows the standard spectral representation method using Fast Fourier Transform (FFT) where the spectrum ) , ( 21  

S
is the FFT of the correlation function (8), θ are random phase angles uniformly distributed between 0 and 2, and the number of points in each dimensions in positive integer powers of 2 for FFT.On each generated Gaussian sample the point-wise translation operation (7) can be applied to obtain a corresponding RF sample following the Weibull marginal distribution.

Numerical simulation
To evaluate the meso-scale SRVE response, displacement constraints are imposed on the top and bottom of the lattice system and a secant stiffness algorithm is applied.At each step of loading, a small incremental displacement D  is applied to the SRVE, and the local stress i  in each bond is computed along with the overall stress of the lattice.The bond broken at a certain loading step is determined by adjusting the incremental displacement to where the index i labels all unbroken bonds in the lattice.With slow increase of the external load applied on the lattice, the individual bonds in the lattice will break irreversibly one after another until the system falls apart.
The effect of the SRVE size on the strain-stress curves is shown in Fig. 2. For smaller size, the area of the post-peak regime becomes larger resulting in a higher value for fracture energy.As the lattice size L increases, the softening curves show the SRVE becomes more brittle.The peak load becomes more fluctuating for a smaller size lattice, as shown in Fig. 3a and Fig. 3b.
The effect of the correlation length on the SRVE stress-strain characteristics is illustrated in Fig. 4. The SRVE strength decreases with increase of the correlation length.When the correlation length is increased, the system becomes more homogeneous that reduces the strengthening effect of heterogeneity.This phenomenon is also reflected from the snapshots at peak load in Fig. 5 where the uncorrelated RF sample exhibits more local damage dissipation than the correlated one.

Probabilistic characterization of SRVE strength
By virtue of the weakest-link or series model, strength of brittle materials can be well described by Weibull distribution.On the other end of spectrum, by using the parallel model, strength of ductile materials (plasticity) exhibit Gaussian distribution.Our lattice model shows that the failure pattern is more ductile when the size of the lattice is reduced and becomes brittle for large size, which suggests that the distribution of the SRVE strength lies between Weibull and Gaussian as noted by Bazant and Pang [7], whereas a grafted Weibull-Gaussian distribution was proposed.Hereby we further note that physical mechanisms related to the lower tail of extreme statistics are associated with high sensitivity of local defect or imperfection, which can be modeled as a weakest link system; in other words, the lower tail of strength distribution follows the Weibull distribution.Physical mechanisms for the upper tail of extreme statistics can be interpreted as flaw tolerant, i.e. a parallel system whereas the potential of the material strength is fully exploited; in other words, the upper tail follows the Gaussian distribution.Based on this understanding, we propose a mixed Weibull-Gaussian cdf with Weibull lower tail and Gaussian upper tail as and G  can be determined by using the nonlinear least-squares regression algorithm, which returns the optimized values for the parameters that minimize the sum of the squared differences between the observed responses and their fitted values.
As an illustrative example, the proposed mixed Weibull-Gaussian distribution is applied to fit the strength distribution of flawed inert silicon nitride (SNW-1000) [24].The testing data in [24] and many other publications show that strength distribution of some brittle and quasi-brittle materials deviates from the two-parameter Weibull distribution, and modified versions of Weibull distribution have been proposed.As shown in Fig 6, the sample data deviate from the linearity in the Weibull plot, and the mixed Weibull-Gaussian distribution provides a nonlinear fitting with a coefficient of correlation 9494 .0 2  R .While the example in Fig 6 is more of illustrative purpose, we emphasize the most advantage of the mixed Weibull-Gaussian distribution lies on the asymptotic matching of the lower and upper tails, and the flexibility of smooth transition between them.Further investigation is needed with respect to experimental verification, and detailed discussion of the mixed Weibull-Gaussian distribution will be provided in a separate report.
By performing Monte Carlo simulation for an ensemble of SRVE, the statistics of SRVE strength can be input into the regression model to obtain the optimized parameters of the mixed Weibull-Gaussian distribution.By using 10,000 SRVE samples, an example of the fitted mixed Gaussian-Weibull distribution is shown in Fig. 7, which demonstrates good agreement with the numerical data obtained from Monte Carlo simulation.

Meso-to-Macro Upscaling
In this part of work, the probability distribution of SRVE will be applied to individual finite elements.By using nonlinear stochastic finite elements, we can upscale the strength distribution of the mesoscale SRVE further to macroscale full structural level containing millions of microscopic bonds (Fig. 1).

Mesoscale damage model
The initial elastic moduli of the continuum are assumed to be spatially homogeneous.The isotropic damage model with an equal degradation of the bulk and shear moduli postulates the stress-strain law in the following form [8]: where  is the stress tensor,  is the strain tensor, e D is the elastic material stiffness, and  is a scalar damage variable.For quasi-static problems, the damage index ) , ( t x  is a RF varying in the fictitious time t .For an undamaged material,  is set to be zero, and the response is linear elastic.
When the continuum is progressively deformed, the initiation and propagation of micro-cracking decreases the stiffness of the meso-scale SRVE or individual finite elements, and  increases.When 1   for a particular element, the stiffness of this element vanishes, corresponding to complete damage.
By focusing on probabilistic strength RF, the following exponential softening law is used in this study (15) where e  is the largest equivalent strain reached in the previous history up to the current state, E s  0  the limit elastic strain, E Young's modulus, and f  a parameter determining the fracture energy.By employing the above softening condition and the Kuhn-Tucker complementary conditions, nonlinear finite element equations can be solved by using the Newton-Raphson algorithm where T is the tangent stiffness given by or simply the secant stiffness, and i the subscript represents the equilibrium iteration.The criterion for stopping the Newton-Raphson iteration is given by where .denotes a norm.

Numerical example
In Fig. 8, the evolution of the damage random field ) , ( t x


is shown for a sample bar subjected to longitudinal tension.As the external loading increases, the bar experiences first in elastic deformation, with the damage variable being zero for all the elements; then partly damaged and reached the limit strength point.After the limit strength point, there is a post-peak strain softening regime.The snapshots are given for damage evolution of the sample bar from Point 1 to Point 5 in the displacement-force curve.By conducting Monte Carlo simulation based on probability distribution of SRVE, probability of macro-scale strength or failure can be effectively estimated.Detailed stochastic computation will be given in a forthcoming paper on multiscale stochastic modeling of nanocomposites.
We finally note that in smeared crack modeling the problem of mesh sensitivity arises due to the invalid assumption of scale-or size-invariance of nonlinear properties for individual finite elements.In our bottom-up modeling process, the mesh sensitivity problem is circumvented by using a size matching condition, i.e. the size of the finite elements is chosen to be identical to the corresponding stochastic RVE.The nonlinear properties of finite elements in our model are size-dependent and thus consistent with physical phenomena.

Conclusion
In summary, a framework of probabilistic upscaling of failure process is described, and numerical simulations for micro-to-meso and meso-to-macro processes are conducted.Specifically, strength of Stochastic Representative Volume Element is numerically characterized by proposing a novel mixed Weibull-Gaussian distribution.
As a concluding remark, we especially note that for complex material failure problems, "computability" of deterministic approach can itself become a question [10].We expect that stochastic models can serve as an effective approach to circumvent the notorious difficulties associated with nonsmoothness, high sensitivity, etc, as the cases in chaotic systems [11].
function.The weight     s e  is in the range [0,1].When the weight is small approaching zero, the distribution is close to Gaussian; conversely when the weight increases toward the value 1, the distribution becomes more Weibull.The six parameters G b m

Figure 7 .
Figure 7. Numerical fitting of the mixed Weibull-Gaussian distribution for SRVE.

Figure 8 .
Figure 8. Stochastic finite element simulation of mesoscale damage evolution.