Deep learning for synthetic microstructure generation in a materials-by-design framework for heterogeneous energetic materials.

The sensitivity of heterogeneous energetic (HE) materials (propellants, explosives, and pyrotechnics) is critically dependent on their microstructure. Initiation of chemical reactions occurs at hot spots due to energy localization at sites of porosities and other defects. Emerging multi-scale predictive models of HE response to loads account for the physics at the meso-scale, i.e. at the scale of statistically representative clusters of particles and other features in the microstructure. Meso-scale physics is infused in machine-learned closure models informed by resolved meso-scale simulations. Since microstructures are stochastic, ensembles of meso-scale simulations are required to quantify hot spot ignition and growth and to develop models for microstructure-dependent energy deposition rates. We propose utilizing generative adversarial networks (GAN) to spawn ensembles of synthetic heterogeneous energetic material microstructures. The method generates qualitatively and quantitatively realistic microstructures by learning from images of HE microstructures. We show that the proposed GAN method also permits the generation of new morphologies, where the porosity distribution can be controlled and spatially manipulated. Such control paves the way for the design of novel microstructures to engineer HE materials for targeted performance in a materials-by-design framework.

Introduction

Propellants, explosives, and pyrotechnics (collectively termed, “energetic materials”) are key performance components in a wide range of applications, including solid rocket motors and munitions. Typically, these materials are mixtures of energetic organic crystals and materials such as plasticizers[1], metals[2], and other inclusions[3,4]. Such heterogeneous energetic (HE) materials have complex, stochastic microstructures. The sensitivity of HEs to loading, i.e. their tendency to detonate, is intimately linked to their stochastic micro-morphology[5,6]. In fact, several micro-scale mechanisms[7], such as collapse of voids/pores[8], inter-crystal friction[9], plastic deformation[10], and shear-banding[11] play important roles in determining the sensitivity of HEs. All of these mechanisms contribute toward energy localization at micro-scale heterogeneities leading to a distribution of “hot spots”. Reactions are initiated at hot spots[7] and propagate outward, culminating in full-blown detonation under suitable conditions. Naturally, designers of devices containing such energetic materials are interested in predicting and controlling the initiation sensitivity of the HE materials, from the standpoint of safety as well as to precisely control and tailor the performance of myriad systems that rely on HEs.

The frontier in computational energetic materials research is to develop predictive multi-scale models to guide the design process of novel materials with tailored performance via microstructural control[12,13]. Predictive frameworks of energetic material response to loads are a matter of concerted current development by several groups worldwide[14-20]; such capabilities are needed to establish structure–property–performance (S–P–P) linkages as necessary precursors to materials-by-design of heterogeneous materials[12,21,22]. The work presented in this paper is directed towards establishing a computational approach to relate stochastic microstructures of energetic materials to their observed performance. To achieve this overarching goal, in silico experiments on an ensemble of stochastic microstructures have been performed in the previous work[23,24], to extract quantitative data from resolved meso-scale simulations in the form of surrogate models[25,26]. These models are used to close the macro-scale governing equations in a multi-scale predictive framework. The work in this paper presents a deep learning approach to generate ensembles of synthetic microstructures that can be used for simulations; the methodology also allows for designing new microstructures, paving the way for materials-by-design of energetic materials[13,27]. Simulations on an ensemble of microstructures also facilitate uncertainty quantification (UQ) due to microstructural variability[24] i.e. aleatory uncertainties associated with the inherently stochastic microstructure. The propagation of uncertainty across scales influences the overall prediction uncertainty of HE sensitivity at the macro-scale. Therefore, ensemble simulations performed on a sufficiently large set of synthetic microstructures will be a key enabling tool for the reliability-based design[28-30] of next-generation HEs.

Need for ensemble simulations on synthetic microstructures

In recent years, multi-scale, multi-physics models that establish structure–property linkages have begun to be developed for the response of HEs to shock and impact loading[31]. A key task in the multi-scale modeling workflow is the performance of ensembles of high resolution, high-fidelity simulations of the meso-scale reactive mechanics[23]. The goal of such simulations is to capture the essential physical ingredients and to quantify energy localization at micro-scale morphological features such as pores and crystal–crystal interfaces. A predominant mechanism for energy localization in the microstructure is the creation of hot spots which are formed due to the collapse of pores[32]. Pore collapse is a well-studied problem, both theoretically[33-36] and experimentally[37]. Reactive simulations of the dynamics of pore collapse in microstructures extracted from images taken using scanning electron microscopy (SEM), X-ray computed tomography (XCT), etc have also been performed[23]. However, since the microstructure is stochastic, a large ensemble of meso-scale calculations, using a large enough number of statistically representative microstructures is needed to extract statistically meaningful information from such simulations. However, computational modelers typically do not have access to large sets of imaged data for a variety of types of energetic materials, or even for various micro-morphologies of a single type of energetic material. Image acquisition is expensive, and imaged data may be distribution-sensitive and specific to limited formulations and material types.

From a computational mechanics standpoint, it is extremely useful to have access to an ensemble of microstructural geometries so that simulations can be performed and the statistics of micro-morphologies can be correlated with measures of sensitivity. In the absence of a large database of microstructural images, computational scientists rely on generated synthetic microstructures[19,38-43] (hereinafter abbreviated as “SynS”) as proxies to the real microstructures (hereinafter abbreviated as “RealS”); a large array of stochastic SynS that closely mimic the real sample must be created and used in in silico experiments. Generation of realistic SynS, with adequate statistical representation and ability to control global and local features in a versatile and flexible computational framework will prove to be a key component in the materials-by-design process for precise and controlled performance of energetic materials.

Previous approaches

Generating an ensemble of stochastic SynS that stand in for the RealS is a challenging task. In general, SynS can be generated using several different approaches[42,44]. Of these, approaches based on shape descriptors have been used in the past; objects are inserted into a computational domain using shape packing algorithms constrained by global shape descriptors, such as volume fractions and particle size distributions. While packing algorithms[44] can be used to generate microstructures with specified morphometric characteristics (e.g., porosity, particle size distributions), these methods are limited to regular/analytical shapes such as spheres, ellipsoids, and polygons[45]. It is also difficult to pack regular shapes for theoretical maximum densities (TMDs) significantly higher than the close-packing limit; typical pressed energetic materials of the type simulated in the current work can have TMD values greater than 90%[46]. Another popular approach, particularly to achieve high TMDs, is to start with space-filling polygons via tessellation, followed by superposition of voids or other phases in the mixture[47]. In both packing and tessellation approaches, however, the microstructure fails to mimic the morphology of real samples for a broad class of heterogeneous materials; i.e. flexibility as well as realism are lacking. In the particular context of energetic materials, both shape-packing and tessellation-based approaches[19,47] have been used to perform meso-scale simulations of the shock response of stochastic microstructures. However, the SynS generated by these previous approaches are still too “ideal”, in that the range of possible void/defect/interface shapes are not well represented in these two approaches. This shortcoming is significant for two reasons: Machine learning approaches are promising alternatives to overcome the limitations of shape-descriptor based approaches. In recent years, convolutional neural networks[50] (CNN) have been used to learn patterns and textures from RealS and generate SynS for a wide range of materials. Among a variety of architectures, Li et al.[51] proposed a method for SynS generation based on a general-purpose texture synthesis method in computer vision that uses transfer learning[52]. For a given RealS input, their method generates a SynS having the same “style” as the input RealS by minimizing the style difference, where the style of a microstructure is defined by the Gram matrix of feature maps produced by the CNN. However, the transfer learning method requires a RealS as a reference to generate SynS. Hence, from the materials-by-design standpoint, material morphology can be explored only in the “neighborhood” of the existing RealS samples. This poses a critical limitation in generating a large ensemble of microstructures that can span the space of candidate material morphologies. In another deep learning approach, Cang et al.[53] and Guo et al.[54] employed an encoder–decoder architecture to generate SynS. The encoder–decoder architecture develops a codified representation of micro-morphology by learning to compress the image pixels (encoder) and reconstruct it back to the original one (decoder). The code values learned by encoder–decoder networks parameterize micro-morphology, allowing the generation and manipulation of SynS by “turning knobs,” where the code values act as the knob control parameters. However, these networks tend to generate blurry images.

RealS possess features that contain non-ideal distributions of shape features, including outlier features such as large cracks and elongated, tortuous void structures. It has been shown experimentally that such structures possessing large surface-to-volume ratios and other shape characteristics play a significant role in energy localization and sensitivity[48].

Computational studies[34,35,49] have shown that local features in microstructures, such as inter-void distances, void shapes and orientations are in fact key aspects of sensitivity that distinguish different classes of the same material. Therefore, SynS must mimic not only global features but also local structural characteristics to be useful in building predictive models of material performance and materials-by-design frameworks.

On the other hand, another deep learning-based method, the patch-based generative adversarial networks (GAN), can generate much sharper and crisper images as we will demonstrate later. Unlike the original GAN method[55], patch-based GAN approaches evaluates the quality of generated images at different local regions, enforcing the details to be clearer and more realistic. Yang et al.[56] and Mosser et al.[57] demonstrated that GANs are not only capable of generating realistic SynS but can also be used to continuously parameterize the micromorphology; this paves the way for smoothly varying the morphology to produce new microstructures. Fokina et al.[58] also confirmed the good performance of GANs in generating SynS of ALPORAS aluminum foam. However, in these previous works, the stochastic variations in micro-morphology were not investigated or quantified. In addition, the output microstructure was fixed at a certain size and was not scalable to arbitrary sizes. Furthermore, the above-mentioned methods lacked the capacity to control micro-morphology at different local regions to produce spatially varying morphologies in a single SynS sample. In this paper, we develop a flexible and versatile algorithm for generating realistic microstructures using GAN. The new algorithm allows control of micromorphology in different regions and can be scaled to arbitrary image dimensions seamlessly.

A method for microstructure generation using GAN

Here, we employ a patch-based, fully convolutional GAN architecture as illustrated in Fig. 1. The generator takes two input vectors and defined at each of grid locations, forming an input tensor. The input tensor is then up-convolved five times, each time scaling the dimension by a factor of 2, to produce an microstructure image. The generator is trained together with a detector network symmetric to the generator, in which the upconvolution layers are replaced by regular convolution layers of the same kernel size and stride. During training, the detector network is presented with an arbitrarily chosen image, either RealS or SynS, and is tasked to determine if the presented image is real or synthetic. The determination of whether a presented image is real or synthetic occurs at each of the grid locations of the generator. The patch-wise feedback on the image quality promotes details to be more closely captured[59]. Furthermore, since the receptive fields of the detector network overlap by 32 pixels between adjacent patches, smooth and seamless connection between the patches is naturally enforced. Finally, it is worthwhile to note that the proposed architecture is fully convolutional so that arbitrary-sized images can be produced without stitching[58] by varying the size of the input tensor.

Figure 1

A schematic overview of the proposed architecture. The generator takes local stochasticity parameters and the global morphology parameters as inputs, where w and h are width and height of the spatial dimension. The discriminator takes in a shuffled distribution of the generated and real microstructures.

Two input vectors shown in Fig. 1—the local stochasticity parameters and the global morphology parameters—play a critical role in the proposed GAN model. The role of the local stochasticity parameters is the same as the “noise tensors” in the standard GAN implementation; they serve as seeds for adding stochastic variations. The global morphology parameters , on the other hand, control the overall morphological characteristics of the generated image, such as grain (or void) sizes, orientations, and aspect ratios. We achieve such a control of the global morphology by setting to be constant across different grid locations during training, while varies randomly across grid locations. Both and are uniformly distributed in the range and hence parameterize the morphological variations of a material in the domain . In the following section, we demonstrate the fidelity and versatility of the proposed GAN approach for SynS generation in comparison with the current state-of-the-art transfer learning (TL) approach presented Li et al.[51].

Results

Qualitative and quantitative comparison of real and synthetic microstructures

When mechanical pressing techniques are used to produce HE materials, defects such as voids, cracks, and inclusions are created[60]. Figure 2 (Top row labeled ‘ground truth’) shows examples of one such microstructure[61]—subsampled from an image obtained using SEM—illustrating the distribution of features. In this type of pressed HE, both inter-crystal and intra-crystal voids are nearly uniformly presented within the microstructure, i.e. crystals and voids do not appear to vary significantly in size and concentration within each sampled image nor across the different samples. As shown in previous work[23], when the pressed HE material is subject to shock loading in the range of 10–20 GPa, voids in the microstructure collapse, leading to localization of energy and formation of high temperature hot spots. The ignition and subsequent growth of these high temperature hot spots depends on the physicochemical properties of the crystalline material and on the morphology of the void and crystal distributions in the microstructure[33,48,62]. In the current paper, to produce SynS that mimic the behavior of the real (i.e. imaged) ones, we focus on comparing the distributions of key morphological metrics of void and crystal phases and the shock response of the RealS and SynS.

Figure 2

Comparison of RealS (top row) and two methods for generating SynS (middle row—TL method; bottom row—GAN). Both the benchmark TL method (Li et al.[51]) and the proposed GAN method display realistic image quality overall. However, as seen in the image callouts to the right, small artifacts and blurry boundaries are noticeable in the benchmark method while the GAN method does not produce artifacts.

The RealS (top row), SynS generated by the TL method of Li et al.[51] (middle row), and SynS generated by the proposed GAN method (bottom row) are presented in Fig. 2 for qualitative comparison. Overall, the GAN generated SynS are visually more similar to the RealS, whereas artifacts and blurry crystal boundaries are observed in TL generated SynS as highlighted in closed-up images. In addition, complex geometry of natural voids and crystals in RealS are well emulated in the GAN generated SynS making them almost indistinguishable from the RealS. Quantitatively, we compare the statistical distributions of void shape descriptors such as void diameter , void aspect ratio AR, and void orientation , which are shown in Fig. 3. From each of the categories (i.e., RealS, transfer learned (TL) SynS, and GAN-SynS), 25 random samples were drawn. Morphometric analysis was performed using methods described in Roy et al.[31]. The distribution of the shape descriptors of the RealS and SynS were generally in good agreement, both between the GAN-SynS and RealS as well as between the TL-SynS and RealS. The whiskers in the figures indicate the standard deviation among the 25 microstructures while the curves correspond to the means. In general, the generated SynS have voids of sizes in the same range as the RealS; the peak of the distribution is shifted slightly in the generated microstructures. As demonstrated below, the observed differences in the size distribution have negligible effects on the computed quantity of interest, viz. the hot spot ignition and growth rate. The void aspect ratio distribution of the real and imaged microstructures is likewise in good agreement, albeit with a shift in the peak of the distribution. On the other hand, the void orientation distribution plots show good agreement; a small peak is observed at but the overall distribution of void orientations is fairly uniform, i.e. there is no strong orientational preference of the voids. As seen from the figure, transfer learning also shows overall good agreement with the RealS void distribution. For porous energetic materials, the sensitivity is strongly dependent on the characteristics of the void field in the microstructure. The void fields are interstices between packed crystals. Both the GAN and TL approaches produce crystals in the microstructure that are in good agreement with those in the real microstructures. To demonstrate this, statistical distribution of crystal shape descriptors have been obtained and are shown in the online supplementary materials (Supplementary Fig. 4).

Figure 3

Distributions of the morphometric parameters indicate the SynS are statistically similar to RealS. The probability distribution functions (PDF) in each case (RealS, TL-SynS and GAN-SynS) were computed for 25 sample images of size 25 m 25 m. The curves in different colors indicate the mean PDFs across the 25 images, while the whiskers represent the standard deviations.

In addition to the above shape descriptors, two-point correlation functions[63] were obtained to quantify the void phase morphology. The computed two-dimensional two-point correlation functions for the RealS and SynS are displayed in Fig. 4 where indicates the volume fraction of void phase at distance r. The result shows that all three microstructures have similar volume fraction of the void phase. However, unlike the transfer learning approach, the GAN generated microstructure shows a weak correlation that persists across the plot. Nevertheless, the magnitude of the correlation is small, and does not exert any discernable influence on the hot spot dynamics simulated below. Furthermore, the horizontal cross-sections of extracted from the center of the two-dimensional plots show that both the SynS are within the standard deviation of the RealS. In addition, both plots stabilize to at m, showing strong agreement of the correlation length with the RealS. Therefore, the two-point correlation function plots reveal that the generated microstructures closely resemble the statistics of the void phase in the RealS. Figure 4 indicates that the TL approach produces a 2-point correlation in better agreement with the real image that does GAN. However, we remark that the performance of TL comes at the expense of solving an optimization problem for each SynS generation, as TL generates a SynS by minimizing the style difference of the generated SynS from a reference RealS input. Furthermore, it should also be noted that a TL-SynS tends to hew close to the immediate neighborhood of the reference RealS in the morphology space due to the style-difference minimization scheme. Whereas, the proposed GAN method can explore the entire morphology space by parameterizing the space by varying the input parameter , as will be demonstrated below. On the other hand, it can be argued that TL-SynS display slightly better agreement with the RealS in morphometry due to the style difference minimization by which the TL method generates images: the TL approach generates SynS only within a narrow margin of style difference and, therefore, TL-SynS are biased towards the RealS samples; whereas the proposed GAN method generates a broader range of novel microstructures. Hence, for the challenge of generating a large ensemble of SynS from a small set of RealS images, the proposed GAN method offers greater practical benefits than the TL approach. Aside from the lack of artifacts and blurring in the GAN-SynS compared to the TL-SynS samples as shown in Fig. 2, the capacity to generate a broader range of novel microstructures using only a small number of RealS samples in the case of GAN is a meaningful advantage given the typical paucity of available microstructure images for a general HE. The GAN approach therefore provides greater flexibility and versatility in the generated microstructures and is better aligned with our goal of creating a pathway to materials-by-design.

Figure 4

(a) Two-dimensional two-point correlation functions for RealS and the two different SynS. (b) Comparison of the cross-sectional slices of the two-dimensional two-point correlation functions in (a).

Figure 5

Meso-scale shock simulations performed on RealS and SynS. The SynS is generated to resemble the stochastic morphology of the RealS. Note they exert similar reactive behaviors.

Simulations of reactive dynamics in real and synthetic microstructures

The void microstructure of a pressed energetic material strongly affects the meso-scale reactive dynamics of the material under shock loading. Both experiments[6,48,64] and simulations[33,34,65] have shown that crystal and void size distributions, void volume fractions, and void morphology affect the meso-scale sensitivity of the material. Here, reactive computations are used to compare the dynamics of the void collapse process and subsequent chemical decomposition of the solid HMX material for real and GAN generated microstructures. The methods for computation of void collapse process and reaction initiation from imaged microstructures are described in several previous publications[23,36,65] and only briefly outlined in the methods section of this paper.

The simulations are performed by applying shock loading from the left end of the domain at a shock pressure GPa. The size of the RealS and GAN-generated SynS are both 25 . The sample is placed in the middle of the computational domain to avoid edge effects from the domain boundaries. Extra padding regions are provided surrounding the imaged sample to account for the translation of the material during the shock passage through the material. In Fig. 5, snapshots of the evolving temperature and pressure fields are shown at several instants of time. The temperature field measures the intensity of hot spots that resulted from the process of void collapse and the pressure plot shows the blast waves that emanate from the void collapse events. As seen in the figure the temperature and pressure fields generated due to shock passage are nearly identical for the RealS and GAN-SynS. Thus, the dynamics of the void collapse and hot spot evolution is well represented by the GAN-generated microstructure. A more quantitative assessment of the physically realistic response of the SynS is obtained from Fig. 6, which shows the time evolution of a quantity of interest in meso-scale simulations, which is the total area occupied by hot spots in the domain, calculated using the approach described in the methods section. As seen in the figure the total hot spot area calculated for the RealS and SynS cases is in good agreement over the course of the evolution of the hot spots. In summary, from Fig. 5 the temperature and pressure fields as well as the quantity of interest (QoI) obtained from simulations using the synthetic microstructure captures the dynamics of real, imaged microstructures with good fidelity. Therefore the GAN-SynS can serve as suitable proxies in ensemble simulations to extract the physics of this type of pressed HMX material.

Figure 6

Hot spot area evolution of RealS and SynS with respect to time corresponding to the shock simulation in Figure 5. Notice the alignment of the trend between the RealS and the SynS.

Towards materials-by-design: controlling the micro-morphology

Controlling morphological features by changing and

As discussed above, the GAN model encodes global micro-morphology characteristics into a parameter vector . The numerical values of the elements of the vector are the “knobs,” which can be tuned to control the morphology of SynS. For example, by varying the values of microstructures with different crystal morphologies are generated in Fig. 7. The rows in the figure pertain to one set of values in the vector , i.e. to one setting of the knobs that control the microstructure. The columns from left to right, in order, show the following: the generated microstructure, the values of the elements of vector , the size distributions of voids and crystals in the microstructure plotted as PDFs, the distribution of void and crystal aspect ratios, and the distribution of void and crystal orientations. From visual observation as well as the morphometry plots shown alongside each image, it is clear that the parameter allows for subtle control of the morphology of microstructures. For example, controlled crystal size changes are shown by comparing the first and the second rows in Fig. 7. In the first row, the settings in the values of lead to a narrow range of void and crystal sizes; by changing the values of the size distribution of crystals is made wider and the void size distribution becomes strongly bimodal as seen in the second row. Meanwhile, in these two rows the void orientations and aspect ratios are largely unchanged in distribution. Controlled changes in crystal orientations are displayed in the third and fourth rows of Fig. 7. Due to the change of crystal orientation, the distribution of the void orientations is at the highest peak around in the third row, while the image in the fourth row shows void orientations predominantly at angles greater than with the highest value at . In the third and fourth rows, the distribution of void and crystal sizes and aspect ratios are seen to be similar; changes in the values of have primarily affected the orientations of the voids and crystals.

Figure 7

GAN input parameter controls the morphology of SynS. The SynS in column (a) were generated from parameters plotted in column (b). The local stochasticity parameters were randomly drawn from the uniform distribution in . Columns (c)–(e) show the morphometry of the generated images, from which the shift of distribution is noticeable as varies.

Stochastic variations in the microstructure are achieved by varying the local stochasticity parameter . Each row of Fig. 8 illustrates subtle variations in micro-morphology for a specific value of the morphology parameter , obtained by changing the values of . It is noticeable that the overall morphology looks similar across the columns, while the layout of crystals/voids and the subtle details are different. Note that significant variations in the void and crystal characteristics are achieved by varying the (comparing across rows), while subtle variation, particularly stochasticity, is achieved by varying (comparing across columns).

Figure 8

Effects of the local stochasticity parameter . SynS on the ith row are generated using the same morphology parameter but with different local stochasticity parameters. Notice the stochastic variations in morphology while the overall “style” is consistent.

Using to create “in between” and spatially varying microstructures

Another useful property of the GAN morphology parameter is that changes in lead to smooth and linearly proportional changes in micro-morphology, as demonstrated in Fig. 9. In the case of GAN, to span the space of morphological parameters (void/crystal size, shape, and orientation), the vector was linearly interpolated between two book-ending images to produce a smoothly varying range of microstructures. In the case of TL, there was no trivial method for linearly interpolating the microstructure, but we attempted to modify the TL algorithm by making it possible to interpolate the Gram matrices (i.e. the representation of “styles”). Figure 9 shows that an interpolation between two different microstructures can be achieved by both TL and GAN approaches. However, the TL method was not able to generate a smooth, continuously varying interpolation, while the proposed GAN method demonstrated the opposite: e.g. grain sizes GAN-SynS in Fig. 9 grow almost linearly from to . This is a clear advantage compared to the current state-of-the-art machine learning based SynS generation approaches, where such smooth and continuous morphology changes cannot be easily accomplished. The parameter linearly maps the space of micro-morphology, so that continuously varying from one value (say ) to another (say ) produces smoothly and linearly varying micro-morphologies.

Figure 9

Linear, continuous variation of micromorphology achieved via GAN. The two images on the top row display RealS samples, and the left- and the right-most images in the middle and bottom rows display SynS that mimic micro-morphology of the RealS. In between the two sides are intermediate SynS that are generated by linearly interpolating the Gram matrices (middle) and the global morphology parameters (bottom). The values of indicates the interpolation parameter . Notice that the GAN method can model “continuous growth” of the crystals that is proportional to .

The ability to generate “in-between” morphologies using the GAN approach is of key importance for the realization of a materials-by-design framework. For example, it is well-known that the sensitivity of pressed HE samples is affected by the crystal (and concomitant void) sizes[48]. However, in previous numerical simulations of the meso-scale shock response[23], it has not been possible to thoroughly quantify the behavior of a range of real, i.e. imaged microstructures of a specific material nor the change in sensitivity due to variabilities in the microstructure; this is due to the lack of sufficiently large numbers and types of imaged samples and the lack of ability to control micro-morphology in image-based shock simulations. Using the present capability, it becomes possible to conduct shock simulations with ensembles of stochastic, morphology-controlled SynS. Figure 10 shows results from shock simulations conducted on two samples of SynS with significantly different crystal/void sizes; in Fig. 10a the temperature field for a sample with small crystal/void sizes is shown on the left and for a sample with larger crystal/void sizes is shown on the right. In Fig. 10a, while the difference in the mean crystal size is rather modest (around 1 m for the left column and 2.5 m for the right column), the simulation results show significant differences in the hot spot development and therefore in the overall material sensitivity. The sample with smaller crystals/voids (left column) shows a higher density of hot spots in the control volume when compared to the sample with larger crystals/voids (right column). Figure 10b shows that there is a noticeable difference also in the pressure field between the two microstructures; similar to the hot spot temperature field in Fig. 10a, it is observed that the pressure field in the case of the small crystal size sample (left column) achieves higher overall magnitude in a larger part of the domain when compared to the pressure field for the sample with larger crystal (right column). The relatively higher sensitivity of the sample with small crystals/voids is quantitatively depicted in Fig. 11 which shows the evolution of the hot spot area with time. Clearly, the small crystal microstructure has a larger hot spot growth rate compared to the large crystal microstructure.

Figure 10

Meso-scale shock simulations performed on controlled SynS containing only small crystals and another have mostly large crystals. Since morphology of the two microstructures are different, reactive behaviors are different: SynS with small crystals tends to increase temperature and pressure quicker than SynS with large crystals.

Figure 11

Hot spot area evolution with respect to time. It is noticeable that the crystal size affects the growth rate of the hot spot area.

In addition to its ability to create “in between” microstructures by smoothly interpolating between two microstructures that define a range of possible morphologies, the proposed GAN method can exercise spatial control of the micro-morphology, providing the ability to spatially grade morphologies within a sample. This is achieved by spatially varying the morphology parameter . Although is constrained to be constant across the spatial grid during training, at the inference time, different locations can have different values of . Through this, for example, one can generate a smoothly graded SynS by allowing to vary continuously from one location to the other, as in Fig. 12a. Here, we show that the generated SynS has crystal and void sizes that vary spatially, with small voids/crystals on the left end of the sample transitioning smoothly to larger voids/crystals on the right end. This type of spatial variation in the microstructure can allow for more precise control on the rate of burning of the energetic material. Furthermore, one can “paint” micro-morphologies by placing different “shades” of , as demonstrated in Fig. 12b, where a “Hawkeye” pattern shown on the left is imprinted into the crystal size distribution, resulting in a microstructure on the right with a morphology that reflects this pattern. In both examples, it is noteworthy that there are no observable stitching boundaries or awkward image transitions. Using this facility, combined with the on-going micro-scale additive manufacturing (or micro 3D printing) techniques[13,21,27,66], the present controlled micro-morphology generation technique can lead to the realization of microstructure-engineered materials, which is being pursued by the authors in ongoing work.

Figure 12

(a) A graded microstructure generated by linearly interpolating two global morphology parameters. A global morphology parameter corresponding to small crystal sizes and a parameter corresponding to large crystals were obtained manually. The image was then generated by spatially grading and linearly: , where is the horizontal position on the image as a fraction of the image width. (b) A layout of global morphology parameters (left) and the corresponding GAN-SynS. Foreground regions were “painted” with in panel (a) while the background was painted with .

Discussion

A novel GAN model for SynS generation has been developed in this work. The method is shown to produce microstructures that qualitatively and quantitatively replicate real microstructures obtained from images. The approach also provides control over the generated morphology, including the ability to spatially vary the morphological properties seamlessly. Compared to the recent deep learning based material synthesis methods[51,53,54,56-58], the advantage of the present method is the ability to scale to arbitrary size without stitching or quilting, to produce a linear and continuous control of morphology, and an ability to generate a microstructure with spatially controlled morphology.

There are several avenues for further development of the current approach. For example, the current morphology parameter is at a single scale, and cannot yet capture a large variation of crystal sizes, which occurs frequently in natural images. As described in depth in the methods section, an input grid location of the generator controls pixel area of the generated image. We found that the GAN was not able to generate crystals bigger than 125 pixels in size, which currently limits our capacity to cover the full range of morphometry. To this end, multi-scale GAN architectures recently proposed in GAN literature (e.g., Karras et al.[67]), in which the input parameters are injected at different scales of convolution are being explored. Another interesting direction of study being currently pursued is to unravel the intuitive, semantic meaning of the global morphology parameters . This can potentially be achieved via correlation analyses between each dimension of the morphology parameters and the conventional morphometric measures, such as void/crystal diameters, aspect ratios, etc. Our preliminary study, in which we conducted a multivariate linear regression to test the correlation between GAN parameter and morphometric measures (Supplementary Figs. 5, 6 and Supplementary Tables 1, 2), revealed statistically significant correlations between morphometry parameters and GAN parameters. Hence, a further investigation may allow the GAN be reparameterized with morphometry parameters, permitting more intuitive handles to the users. It is also remains to establish the optimal dimension of . Unfortunately, theoretical works on this issue are lacking in the literature, other than naïve search methods for finding the minimal loss by varying the dimensions. Therefore, more fundamental research should be directed on this avenue.

With regard to design and analysis of energetic materials, an immediate extension of this work would be to use the GAN method to extensively study the sensitivity of heterogeneous energetic materials to applied loads; for this purpose we will characterize the effects of micromorphology on different physical QoIs, such as hot spot ignition and growth rates[26]. This can be accomplished by performing computational simulations on a large ensemble of morphology-controlled microstructures, as briefly demonstrated in Figs. 10 and 11. This will lead to a comprehensive quantitative characterization of the detonation behavior of HE materials as a function of micro-morphology, and will provide surrogate models for bridging meso-scale behaviors to macro-scale simulations[16]. Further development will lead to a materials-by-design framework where a desired performance or property of a HE material can be engineered through optimization of material morphology, analogous to topology optimization in mechanical component design. The ability of the proposed GAN method to spatially control micro-morphology will play a key role in the realization of such a framework.

Methods

Data

We used a raw image of a class V cyclotetramethylene-tetranitramine (HMX) pressed energetic material[60,61] obtained using SEM for the training of the GAN model. The original image was in size, from which we sampled image patches by cropping the image at a random position. The physical resolution of the image is 52 nm per pixel, or 19 pixels per 1 m. The sample patches, therefore, correspond to in physical dimension. For training, a total of 12,500 of such samples were generated and the pixel values were re-scaled from [0, 255] to [− 1, 1].

GAN

GAN[55] is a generative model capitalizes an adversarial training of two or more competing neural networks. Typically, an image generator network and a real/synthetic discriminator are trained together, which poses a minimax game where tries to maximize the misdetection rate while tries to minimize it:where is the distribution of real images and z is the input parameter to . The convergence is reached at the equilibrium, where the discriminator is no longer able to distinguish synthetic images generated by from real images drawn from the data set.

In our work, the GAN loss in Eq. (1) is evaluated multiple times at grid patches and averaged over grid locations, resulting the following loss function:where is the prediction of at jth row and i-th column on the grid in Fig. 1.

The input tensor is composed of the local stochasticity parameter and the global morphology parameter vector defined at each grid location (j, i). During the training, we set the global morphology parameter for all i, j for some constant randomly drawn from the uniform distribution in [− 1, 1], while the local stochasticity parameter is independently drawn from the uniform random distribution in [− 1, 1]. At the inference time, can vary across the grid, as demonstrated in Fig. 12 and controlled by the user.

The architectures of the generator and the discriminator are symmetric. The discriminator contains a stack of five convolution layers of kernel size and stride 2. The generator has the mirrored architecture, where the convolutions are replaced by upconvolutions of the same kernel size and the stride 1/2. Such architecture results in a receptive field, which determines the window size that one grid parameter controls.

For the training, the adaptive moment estimation (ADAM) optimizer with the learning rate of 0.0002 and the batch size of 10 was used. We set and as the dimensions of and , respectively. During the training, was used.

Morphometry

Morphometry refers to the quantification of the microstructural features using shape descriptors or correlation functions. The shape descriptors calculated in this paper are used to quantify the size, shape and orientation of voids and crystals, the two primary features in the microstructure of pressed HMX. In the present work, the shape descriptors are obtained using level-set based morphometry techniques described in detail in Roy et al.[31] Briefly, the gray-scale microstructural images are segmented in the image (pixel) space using active contouring to obtain level-set fields[68]. With the zero level set value defining the crystal shapes, various approaches have been developed in Roy et al.[31] to calculate the sizes of voids and crystals, their aspect ratios and their orientations. The orientation of voids () are defined as the mean of the orientations of the voids relative to the direction of the shock propagation (from left to right across this paper). It is calculated by averaging over the local orientations within each void, where the local orientation is calculated at each Cartesian grid point in the interior of the voids using the volume orientation method[69]. These three quantities are then computed for each void/crystal and histograms (pdfs) are developed to obtain the distribution of void sizes, aspect ratios and orientations of all voids/crystals in the domain. The second approach to morphometry employed here is the use of correlation functions; n-point correlation functions provide rich information on the structures embedded in the domain. Typically, in the interest of computational cost. Two-point correlation functions[63] determine the probability of finding 2 random points with position vectors p and in a given phase. Geometrically, the two-point correlation function can be interpreted as the probability of having both ends of the chord vector r in the same phase. The two-dimensional two-point correlation function is calculated as the following probability function: where, I(p) and are the indicator function values at locations p and respectively, is the magnitude of the length of the chord vector in the normal direction . The indicator function I(p) is used to identify phase at a particular point in the image space. In this work since the void phase is of interest, I(p) is equal to 1 if the position vector p indicates void phase and 0 otherwise. In the present work, the two-point correlations are obtained using the open source software developed by de Geus et al.[70].

Meso-scale simulations of reactive mechanics

Meso-scale simulations are performed to simulate the collapse of voids in the microstructure due to the passage of a shock wave. The shock is parameterized by the pressure . Reactive calculations are performed using methods discussed extensively in previous publications[65]. Several validation exercises have been demonstrated in previous work, providing high confidence in the physical correctness of the shock computations[36]. By performing the meso-scale simulations, temperature, pressure and species field data are utilized to quantify the response of the pressed material to the imposed shock, as shown in previous work[23]. In the present context, the QoIs used to quantify the effect of microstructure on the sensitivity of the material, are calculated by following the evolution of the temperature field and the reaction product mass fraction in the sample. The temperature field T(x, t) in the domain measures the intensity of a hot spot resulted from the process of void collapse. Higher temperature hot spots formed due the collapse of voids in the material lead to higher chemical decomposition rates. The reaction zone defines the hot spot in the domain, which is defined as the region where the temperature of the material exceeds the value of the temperature ()reached after the passage of a planar shock wave. The hot spot area is a significant quantity of interest for determining sensitivity and is calculated as the area of the domain where the temperature . The hot spot area is recorded throughout the simulation to track the evolution of the hot spots.

Supplementary Information.