Dielectric breakdown driven by flexoelectric and piezoelectric charge generation as hotspot ignition mechanism in aluminized fluoropolymer films Free

The ignition of energetic materials (EMs) directly affects their safety and effectiveness in a wide range of applications, such as airbags, fire devices, and propulsion systems.¹ Numerical analyses^2–9 and experiments^10–13 have focused on thermomechanical dissipations^14–16 that drive ignition, which intimately depends on the development of localized regions of high temperatures known as hotspots. In most cases, thermomechanical dissipation associated with inelastic deformation, void collapse, fracture, and friction is the primary hotspot generation mechanism. There is a growing interest in developing alternative means for generating hotspots in EMs in a controlled manner for more precise ignition control, reaction manipulation, and novel multifunctionality. In addition to mechanical stimulus, optical,¹⁷ electromagnetic,¹⁸ and electromechanical¹⁹ stimuli have all been shown to be sound alternatives for activating ignition reactions in EMs via hotspot development.

We recently demonstrated that the ignition of poly(vinylidene fluoride-co-trifluoroethylene) binder and nano-aluminum (nAl) films [or P(VDF-TrFE)/nAl] can be influenced by their electromechanical properties²⁰ and not necessarily via thermomechanical dissipation. Specifically, experiments showed that poled films (possessing both flexoelectricity and piezoelectricity) have ignition times that are shorter than those of unpoled films (possessing flexoelectricity only) under the same load intensity, although aluminized fluoropolymers can be poled only to a relatively low attainable value of piezoelectric coefficient $(d33=−5.45pC/N)$⁠. The joint experiments and simulations also showed that flexoelectricity plays as important a role in poled films as in unpoled films, which have no piezoelectricity and yet exhibit ignition behaviors on the same order of magnitude as those of poled films. This dominance of flexoelectricity can only be explained if very significant gradients of strains occur in the materials. For the flat films under overall planar compressive loading here, the overall macroscopic, net strain gradients are relatively small. The explanation lies in the microscale heterogeneity-induced gradients. Indeed, Zaitzeff and Groven have demonstrated experimentally²¹ that the flexoelectric coefficients of PVDF/Al composites can be increased by embedded Al particles, illustrating the role microstructural heterogeneities play in the flexoelectric effects at the macroscale, owing to the heightened levels of local strain gradients within the material.

Here, we further advance the study by directly analyzing the development of hotspots during dielectric breakdown driven by flexoelectricity and piezoelectricity of the P(VDF-TrFE)/nAl composites. Since vastly different timescales and physical mechanisms are involved, a two-step sequential multi-timescale and multi-physics framework for explicit microscale simulations is developed and used to analyze both the mechanically driven electric field (E-field) development process and the subsequent dielectric breakdown process that leads to the generation of the hotspots responsible for reaction initiation. The temperature field resulting from the breakdown is analyzed to establish the hotspot conditions for the onset of self-sustained chemical reactions. Trends in the predicted ignition threshold and those measured experimentally are compared.

Experiments and computational simulations are carried out. Trends in the measured ignition sensitivity (time to observation of optical emission from the chemical reaction or “first light”) and predicted ignition sensitivity (time to sufficient hotspot generation in the material) are compared. To delineate the individual contributions of flexoelectricity and piezoelectricity, samples in the poled state and the unpoled state are subjected to loading of different intensities.

The drop-weight impact apparatus used is the BAM configuration (Bundesanstalt für Materialforschung)²² fall hammer by OZM Research to quantify impact sensitivity. This apparatus enables the application of a precisely controlled impact load based on a specified drop-height and a hammer mass. The tower has a maximum drop-height $(h0)$ of 1 m and provides frictionless guide rails for the 5 kg mass. The configuration used minimizes variability in the impactor orientation.

Multiple samples of P(VDF-TrFE)/nAl²³ were analyzed and repeated for each of the drop-heights between 20 and 100 cm with equal increments of 20 cm. Each sample is a uniform square film with 10 mm in side length and 85 μm in thickness. The time of impact (defined here as $t=0$⁠) was determined based on the initial electromechanical response from the sample. The sample was placed between two metal contacts. These contacts were connected to the voltage probes of an oscilloscope, which was used to measure the voltage produced across the sample due to the impact. The ignition event was captured by measuring optical emission using a Ge switchable gain amplified detector from Thorlabs, an infrared (IR) photodiode, and a high-speed camera. The time duration between the initial impact $(t=0)$ and the onset of the ignition event (from the IR signal) is taken as the time to ignition $(tIGN)$⁠. This time is used to quantify the ignition sensitivity of the materials under different conditions. Figure 1 illustrates the experimental setup.

To explicitly analyze the electromechanical and thermal responses of the P(VDF-TrFE)/nAl subjected to impact loading, a two-dimensional microstructure model is developed based on attributes of the material used in the experiments.²⁰ The microstructure model is 8 × 3 μm² in size and comprises nAl particles embedded in a P(VDF-TrFE) polymer matrix. The particle volume fraction is $η=9%$⁠. The particles are circular core-shell structures with an outer diameter of 80 nm.²⁴ The nAl core has a radius of 36.7 nm, and the Al₂O₃ shell has a thickness of 3.3 nm, matching experimental measurements. To account for the random variations in the microstructure morphologies, which lead to statistical variations in the ignition behavior just like in experiments, a statistically equivalent microstructure sample set (SEMSS) consisting of five random but statistically similar microstructures (MS) is generated, as shown in Fig. 2.

To delineate the underlying factors influencing the development of hotspots and establish a measure for the ignition sensitivity of the P(VDF-TrFE)/nAl, the focus is on computationally predicting the temperature field and the time it takes for a sufficient portion of the microstructure to reach the critical temperature $TPIR=723K$²⁵ or when pre-ignition reactions (PIRs) would initiate near the particle–binder interface. As illustrated in Fig. 3, the simulations are carried out in two steps using an electromechanical model for resolving the microstructure-level flexoelectric and piezoelectric responses to mechanical loading (Step 1) and a thermal–electrodynamic model for analyzing the temperature rise during dielectric breakdown caused by the high E-field resulting from the flexoelectric and piezoelectric polarizations (Step 2). This is a sequential multiple timescale simulation scheme that accounts for the broad timescales and physical mechanisms in the process culminating in the onset and development of dielectric breakdown. The coupled electromechanical analysis and the thermal-electrodynamic analysis are performed successively using COMSOL Multiphysics (v5.4). The electromechanical analysis enables understanding of how the film's inherent properties lead to the generation of E-field as a function of the material deformation and applied loading and provides input for the subsequent electrothermal analysis. The thermal–electrodynamic analysis enables the temperature rises and hotspots to be quantified. The electrostatic/mechanical properties and the electrodynamic/thermal properties of the constituents in the material are listed in Tables I and II, respectively. The piezoelectric and flexoelectric properties of the P(VDF-TrFE) polymer are listed in Table III.

The P(VDF-TrFE)/nAl films are discretized using two-dimensional, isoparametric elements. Specifically, the mesh consists of quadratic, three-noded triangular, Lagrange elements approximately 10 nm in size. This ensures that there are roughly eight elements along the diameter of each nAl particle.

A coupled electromechanical model is used to solve for the E-field in the material resulting from the flexoelectric and piezoelectric responses of the microstructures to mechanical loading. As illustrated in Fig. 4, the compressive load from the drop-weight impact is modeled by prescribing a downward displacement $(uy=−t2gh0)$ along the top surface of the microstructure $(g=9.81m/s2)$⁠, while tangential displacement $(ux)$ is permitted. Here, t denotes the load duration of interest. A constant velocity is used, as the shock impedance of the impactor mass is orders of magnitude higher than that of the P(VDF-TrFE)/nAl film, and the ignition events of interest occur well before the wave reverberations cause the deformation to decelerate. This treatment is based on and consistent with the video observation of the experiment.

As outlined in Table IV, the top surface is also electrically constrained such that there is no electric displacement $(D)$ in its normal direction $(D⋅n=0)$⁠, where $n$ denotes the unit normal to a given plane. The bottom surface of the microstructure is electrically grounded $(φ=0)$ and mechanically constrained. Specifically, there is no normal displacement $(uy=0)$ perpendicular to the bottom surface, while tangential displacement $(ux)$ is permitted. The right-hand surface is traction-free and has no electric displacement gradient in its normal direction. The left-hand surface is taken as a symmetric boundary with $ux=0$ and $Dx=0$⁠.

The microstructures are initially stress-free and charge-free. This is a two-dimensional model, and the plane-strain condition is assumed. This treatment allows the essential microstructural effects to be captured. A full three-dimensional model may be used in the future.

The governing equations and constitutive relations used to model the coupled electromechanical behavior are described here. A quasi-static model is assumed because the load duration is orders of magnitude longer than the time for the stress wave to traverse the sample. The electrical response is governed by Gauss's law and Faraday's law, while the mechanical response is governed by the conservation of momentum. The governing equations in the reference configuration are

where $D$⁠, $E$⁠, and $ρq$ represent the electric displacement vector, the electric field vectors, and the volumetric free-charge density, respectively; $S$ is the second Piola–Kirchhoff stress tensor. The deformation gradient tensor $(F)$ can be expressed as $Fij=δij+∂ui/∂xj$ with $δij$ denoting the Kronecker delta and $ui$ being the components of the displacement $u=u(x)$⁠. Finite strain theory is used to account for the large deformation of the film caused by the external impact.

The electric displacement comprises a dielectric polarization part caused by the induced dipole moment, a piezoelectric polarization part caused by stress, and a flexoelectric polarization part caused by strain gradients. This electromechanical constitutive relation^26–32 under finite deformation is

where $J≡det(F)$⁠, $Ej$⁠, and $εjk,l=∂εjk/∂xl$ represent the Jacobian, electric field, and the gradient of the Green–Lagrange strain $εij=12(ui,j+uj,i+uk,iuk,j)$⁠, respectively; $κ$⁠, $dijk$⁠, and $μijkl$ are the isotropic dielectric permittivity (absolute), the piezoelectric coefficient, and flexoelectric coefficient tensors, respectively. The electric field vector is equal to $Ek=−∂φ/∂xk$ with $φ$ being the scalar electric potential. The values of $κ$ are given in Table I in terms of the vacuum permittivity $(ϵ0=8.854×10−12F/m)$⁠.

Since the specimen is poled in the thickness (vertical, impact, or “2”) direction in the global frame in Fig. 4, the non-zero components of the piezoelectric coefficient tensor $dijk$ are d₂₂₂= d₃₃, d₂₁₁= d₃₁, d₂₃₃= d₃₂, d₁₁₂= d₁₅, and d₃₃₂= d₂₄, where $dIJ$ (I = 1, 2, 3; J = 1, 2, …, 6 = 11, 22,  … , 12) are the commonly used reduced coefficients in the material frame with “3” indicating the poling orientation. The reduced coefficients are listed in Table III.

The flexoelectric response is also isotropic such that the fourth-order flexoelectric coefficient tensor can be fully determined using two independent constants, just like the isotropic linear elastic modulus tensor [see Eq. (6)]. Specifically,

where $μTR=μijji$ (i, j = 1, 2, 3; no summation is implied over repeated indexes i and j) denotes the transverse flexoelectric coefficient whose value is listed in Table III; and $λ≡μSH/μTR$ is the ratio of the shear flexoelectric coefficient $μSH=μiijj=μijij$ (i, j = 1, 2, 3; no summation is implied over the repeated indices i and j) to the transverse flexoelectric coefficient $μTR$⁠. The longitudinal coefficient is equal to $μLT=μiiii=(1+2λ)μTR$ (i = 1, 2, 3; no summation is implied over the repeated index). Since it is difficult to measure $μSH$ and $μLT$ directly and since no values are reported for them in the literature for P(VDF-TrFE), a parametric study is carried out on a similar material (THV/Al samples⁴⁰) by varying $λ$ over a range to determine $μSH$ and $μLT$ via best fit to experimental data. This fit uses the measured macroscopic transverse flexoelectric constant values of the THV/Al composites with different Al concentrations. It is found that $λ=−0.625$ provides the best description of the experimental results. The flexoelectric part of the polarization in Eq. (2) can be expressed as

The second Piola–Kirchhoff stress $(S)$ and the Green–Lagrange strain $(ε)$ are used to write the coupled electromechanical constitutive law, which accounts for converse piezoelectricity and flexoelectricity as

where the elastic stiffness tensor can be expressed as

with E and $ν$ being the Young's modulus and Poisson's ratio (see Table I), respectively.

The electrostatic results from the coupled electromechanical analysis (Step 1) serve as part of the initial condition for the thermal–electrodynamic analysis, which focuses on the modeling of the dielectric breakdown, the transient evolution of the E-field, the energy dissipation due to electrical conduction, and the temperature field within the material. The interest is in the extent of dissipative heating along the interfaces between the polymer and the particles. The minimum time required for a sufficient portion of the microstructure to attain the critical temperature for reaction initiation $TPIR=723K$²⁵ is determined and taken as a measure for the ignition sensitivity of the material.

Dielectric breakdown is explicitly modeled in the P(VDF-TrFE) polymer only since the aluminum core of the particles is electrically conductive and the aluminum oxide (Al₂O₃) shell has a very high breakdown strength that effectively precludes its breakdown under the conditions analyzed. In addition to the continuity of the E-field, conservation of energy is maintained between the two models at the beginning of the thermal–electrodynamic simulation at each material point. This treatment is a phenomenological approach to account for the electromechanical coupling during the breakdown simulation, which does not explicitly track the mechanical deformation, and the piezoelectric and flexoelectric processes due to the vastly different timescales of the processes.

The electrodynamic breakdown simulation is governed by Gauss's law, Faraday's law, and Ampère's law (with Maxwell's correction) in the form of

where the inductive effect is assumed to be negligible and $J$ is the total current density whose constitutive relation can be expressed as

with $σ$ denoting the electrical conductivity. The first and second terms represent the conduction current density induced by the electric field and the displacement current density, respectively.

Dielectric breakdown is modeled explicitly using a conditional function

where $σ$ represents the current conductivity, $σ0$ is the initial electrical conductivity prior to breakdown, $tbd$ is the specified time of breakdown initiation, $σbd$ is the post-breakdown conductivity $(σbd≫σ0)$⁠, and $Ebd$ is the dielectric breakdown strength of the constituent. This criterion states that the local electrical conductivity at a material point changes irreversibly to a higher value $σbd$ once the local E-field exceeds the breakdown strength.

Many external factors, such as temperature and surface roughness, may influence the post-breakdown conductivity of materials. Due to a lack of experimental data on the post-breakdown conductivity of materials in local areas, a constant value is assumed here, as is often the case in simulations.⁴¹ 

The thermal analysis as part of the breakdown simulation is governed by the conservation of energy as

where $ρ$⁠, $cp$⁠, T, and k denote mass density, specific heat, temperature, and thermal conductivity, respectively. The heat source term on the right-hand side accounts for the resistive heating (Joule) caused by the dielectric breakdown that drives the heating.

The phase transformation of P(VDF-TrFE) from a semi-crystalline solid state to a molten state affects the post-breakdown temperature in the material and is phenomenologically considered in the thermal analysis. This transition occurs at a melting temperature of $Tm=446K$⁠. The effect of this endothermic melting process on temperature is modeled using the apparent specific heat⁴² in the form of

which accounts for the energy absorption due to latent heat. In the above expression, $cp1$ and $cp2$ denote the specific heat capacities of the solid and molten states at constant pressure, respectively, and $ϕ$ and L are the volume fraction of the solid phase and the latent heat, respectively. The mass fraction $(πm)$ can be expressed as $πm=0.5−ϕ$⁠.

The electromechanical response of the P(VDF-TrFE)/nAl is first analyzed. This analysis focuses on the flexoelectric and piezoelectric effects under the impact load shown in Fig. 4. The distributions of induced electric potential and the E-field are shown in Fig. 5 for a film in the poled and unpoled states at $t=140μs$ after impact, and the drop-height is $h0=20cm$⁠. The poled case [Fig. 5(a)] shows a potential field with a significant gradient in the thickness (vertical) direction owing to the piezoelectric contribution, while the unpoled case [Fig. 5(c)] exhibits an insignificant overall potential development. However, both poled and unpoled cases show a significant E-field development in the material [Figs. 5(b) and 5(d)] near the nAl particles, with the poled case having higher field levels compared with the unpoled case. The local E-field levels are the highest near the particles in both cases, with the field intensities in the unpoled case in such areas on the same order as those in the poled case. Specifically, the interfacial E-field intensities in poled samples are only ∼11% higher on average than those in the unpoled samples, suggesting that the flexoelectric effect plays a key role in the interfacial areas.

Because the reaction initiates in these areas, the flexoelectric effect plays a dominant role in determining the hotspot formation, and consequently, the ignition and the piezoelectric response play only a secondary, modulating role for the materials and the particle size considered. As a result, only small differences in ignition behavior are observed in experiments and the simulations to be discussed later in this paper. Figure 5(d) also demonstrates that the flexoelectric effect alone (without a piezoelectric contribution) can cause sufficiently high E-field levels along the particle interfaces. Ultimately, dielectric breakdown can occur when the interfacial E-field reaches or exceeds the breakdown strength of the P(VDF-TrFE) binder under sufficient deformation.

Figure 6 shows the evolution of the E-field in MS #1 shown in Fig. 4 as the loading progresses. Increasing deformation leads to higher stress and strain gradients, causing the E-field to increase and eventually reach the breakdown strength $(Ebd)$ of the material at approximately 140 μs in areas near the particles.

Figure 7 shows the mechanical and electrical fields near the nAl particles of the material in the poled state in Fig. 6 at $t=140μs$ for $h0=20cm$⁠. The property mismatch between the particles and the matrix causes the stress and strain gradients in the areas around the particles to be much higher than in areas away from the particles. This local magnification of stress and strain gradients enhances the development of the piezoelectric and flexoelectric polarization, respectively, leading to high local electric displacement and E-field levels. This enhancement of the interfacial electromechanical responses is primarily due to the high elastic stiffness mismatch between the nAl particles and the polymer binder (Young's modulus ratio is ∼46). To put this in perspective, the average E-field near the interface of each particle is calculated. The interfacial E-field is approximately twice the average E-field level in the entire binder. This is primarily due to flexoelectricity: the ratio of interfacial flexoelectric polarization to interfacial piezoelectric polarization is ∼2.4 for the poled case. Here, the converse flexoelectric effect is practically negligible, as the event is mechanically driven.

The results of the electromechanical simulations at different times are used as input for the thermal–electrodynamic simulations, allowing the trend in the relation between loading and hotspot development to be analyzed. Figures 8–10 show the breakdown behavior for onset at 140 μs after impact loading $(tbd=140μs)$ for the case in Figs. 5–7 (poled; MS #1 in Fig. 2; $h0=20cm$⁠).

In Fig. 8, the interfacial E-fields at the onset of dielectric breakdown $(tbd=140μs)$ and at 20 ns after breakdown initiation $(tbd+20ns)$ are shown. Here, 20 ns is chosen to illustrate the attenuation of the E-field due to dissipation during breakdown, as further evolution of the fields does not lead to a significant further temperature increase in the material. At this time, the average E-field in the microstructure is ∼48% of the value at the onset of breakdown. In contrast, the average E-field along the particle–binder interfaces at this time is ∼19% of the corresponding value at the onset of breakdown. This difference indicates that the E-field intensity attenuates much more quickly in areas around the particles (which are the sites of breakdown) than in areas away from the particles. Still, the decrease in the E-field in areas away from the particles is due to the fact that energy stored in the fields away from the particles is delivered to the breakdown sites via the propagation of electromagnetic waves and dissipated to cause a temperature increase.

The E-field decrease and energy flow are due to the electric current in areas of breakdown that have high electrical conductivity. Figure 9(a) shows the distribution of the electrical conductivity in the binder. Red indicates the high post-breakdown conductivity $σbd$⁠, and blue indicates very low conductivity $σ0$ before breakdown. Note that breakdown is an irreversible process; therefore, the breakdown areas only grow with time but stabilize quickly after the initial onset. Subsequently, energy flow from the breakdown sites continues, causing temperatures in the breakdown sites to increase and eventually plateau. As shown in Fig. 9(b), the electric energy density is substantially lower in regions of breakdown, as much of the initial energy has been dissipated as heat throughout the breakdown process.

Figure 10 provides another perspective on the breakdown process and the state of the material at 20 ns after breakdown initiation at $tbd=140μs$ for the case of Figs. 5–9. In Fig. 10(a), the areas of the binder that have melted are shown in red (melting temperature $Tm=446K$⁠). In Fig. 10(b), the temperature distribution in the binder is shown; note that the temperature for reaction initiation or temperature for pre-ignition reaction $TPIR=723K$⁠. It can be seen that essentially, all interfaces between the particles and binder have temperatures exceeding the melting temperature $Tm$ and a significant fraction of the interfaces have temperatures at or above the reaction initiation temperature for $TPIR$⁠. These areas are the “hotspots” sufficient for triggering pre-ignition reactions in the material. It is the evolution of this hotspot field that will be analyzed systematically later in order to establish a macroscopic condition or time for a pre-ignition reaction.

To assess how the time for the initial onset of breakdown affects the outcome, a parametric study is carried by assuming that breakdown starts at different times in the impact loading process. For the material case in Figs. 5–10 (poled; MS #1; $h0=20cm$⁠), Figs. 11 and 12 show the evolutions of the E-field intensity and temperature, assuming that the breakdown initiation time is $tbd=130μs$ and $tbd=150μs$⁠, respectively. Just like for $tbd=140μs$ (Figs. 8–10), the E-field intensity decreases with time after breakdown initiation and the temperature increases, both approaching a steady state. Overall, the temperature is higher for higher t_bd, with the highest temperatures equal to 698 K for $tbd=130μs$⁠, 763 K for $tbd=140μs$⁠, and 821 K for $tbd=150μs$⁠. This is due to the fact that the initial E-field is higher at later $tbd$ times. For example, the average interfacial E-field at $tbd=150μs$ is ∼15% higher than that at 130 μs. This trend will be systematically quantified later as well.

Figure 13 shows the time histories of the E-field intensity $|E|$⁠, the current density $|J|$⁠, and the temperature rise $ΔT$ along the particle–binder interfaces for the material case shown in Figs. 5–10 (poled; MS #1 in Fig. 2; $tbd=140μs$⁠; $h0=20cm$⁠). Here, the interfacial results are averaged over the entire set of particle interfaces and normalized by their respective maximum value (the highest values for the E-field intensity $|E|max$ and the current density $|J|max$ occur at the beginning of the breakdown process, and the highest value for temperature rise $ΔTmax$ is the steady state limit temperature in late stages of the breakdown process). All three quantities approach a plateau approximately 20 ns after breakdown initiation $(tbd)$⁠. Beyond ∼20 ns, the temperature is expected to gradually decrease due to the effect of thermal conduction and since the exothermic chemical reaction is not modeled here. For these reasons, the peak temperature distribution at 20 ns after $tbd$ is used for a hotspot condition trend assessment.

There is an interest in establishing a threshold condition associated with the ignition of the material. In the experiment, the ignition time is taken as the threshold, which is the minimum time it takes for a visible light emission (or flash) to be observed from the ignition reaction. Since the chemical reaction is not modeled, it is not possible to calculate this time in the current study. However, because the full temperature field is available, analyses can be carried out on the overall conditions of heating relative to the temperature at which reaction initiates $(TPIR=723K)$ near the particle interface. The interfacial temperature is computed for each nAl particle by averaging the local temperature throughout a thin outer ring (with a thickness of 1 nm) around the particle–binder interface. The time $(tPIR)$ at which a pre-specified fraction $(θ)$ of the particle–binder interfaces reaches this pre-ignition reaction temperature is taken as a measure for the ignition threshold under a given load intensity. The attainment of this critical condition (e.g., $θ=5%$²⁰) activates the pre-ignition reactions only in the local hotspots, which precedes the ignition time $(tIGN)$ associated with light emission observed in experiments (i.e., $tPIR<tIGN$⁠). Nonetheless, the overall trends in these two time measures are similar. This is the basis of our analyses and comparison here.

In Fig. 14, the average steady state temperature along the interface of each particle is shown for MS #1 in Fig. 2 in both poled and unpoled states under an impact loading of $h0=20cm$⁠. The particles are colored by their respective average interfacial temperature. Particles whose interfacial temperature is at or above the $TPIR$ can lead to reaction initiation, as they are associated with hotspots that are considered “critical.” The poled case has more critical hotspots or more particles with average interfacial temperatures at or above $TPIR$⁠. Obviously, an increased time to breakdown initiation leads to a higher interfacial temperature due to the increase in the deformation level from impact loading. For the poled case, the percentage of particles whose interfacial temperature is at least $TPIR$ is 2.5% for $tbd=130μs$⁠, 11.4% for $tbd=140μs$⁠, and 26.5% for $tbd=150μs$⁠. For an unpoled microstructure, the percentage is 0.2% for $tbd=130μs$⁠, 1.3% for $tbd=140μs$⁠, and 6.1% for $tbd=150μs$⁠.

For each of the five microstructures in the SEMSS in Fig. 2, the probability distribution of the average interfacial temperature is calculated and the data are fitted to the lognormal probability density function (PDF). Figure 15 shows the results for the five microstructures (MS #1–5) for $tbd=150μs$ and $h0=20cm$⁠. The corresponding cumulative distribution function (CDF) of the average interfacial temperature is obtained by integration and is shown in Fig. 16. Significant variations in the PDF and CDF are present among the five random microstructures in Figs. 15 and 16. Taken together, this set of results allows the statistical variations in the behavior to be quantified. The red-colored portion of each plot in Figs. 15 and 16 denotes the fraction of the particles with average interfacial temperatures at or above the ignition temperature $TPIR$ for $tbd=140μs$⁠. This fraction is slightly different among the five microstructures due to statistical fluctuations. The purple dashed line denotes the temperature reached or exceeded by 5% of the particles in each microstructure. The critical time at which point the initial breakdown occurs $(tbd)$ that subsequently leads to the attainment of 5% of the particles²⁰ to reach the ignition temperature $TPIR$ is taken as the pre-ignition reaction time $(tPIR)$ for the respective microstructure and the load intensity (i.e., $h0$⁠) examined, shown as

where $Np$ is the total number of nAl particles in a given microstructure; $Hi[Tiint−TPIR]$ is a Heaviside function whose value is equal to unity if $Tiint≥TPIR$ and zero if otherwise.

This process, as illustrated in Fig. 18, involves establishing the 5% threshold temperature for a series of $tbd$ values for each microstructure. The results for the five microstructures allow for a statistical assessment of the pre-ignition reaction time $tPIR$ for the entire microstructure sample set.

Figure 17 shows the range and average value $t¯PIR$ of the predicted times to pre-ignition reaction $(tPIR)$⁠. For $h0=20cm$⁠, $t¯PIR$ are 130.5 and 144.9 μs for poled and unpoled samples, respectively. Both have the same coefficient of variation of ∼0.018. This indicates that poled P(VDF-TrFE)/nAl undergoes a pre-ignition reaction at loading durations that are ∼10% shorter than those for unpoled films. The difference is due to the effect of piezoelectricity. Obviously, the primary driving force for breakdown is due to flexoelectricity.

Finally, the predicted $tPIR$ as a function of drop-heights for poled and unpoled microstructures is shown in Fig. 18 (dashed lines). For a comparison of trend, the experimentally measured time to ignition $tIGN$ is also shown (solid lines). The $tPIR$ is shorter than $tIGN$⁠, with the difference between the two being the time for reaction in hotspots to grow and propagate in the material leading to macroscopically observable light emission and thermal flash. Note that $tIGN$⁠, as measured in the experiment, is based on the infrared detector measurement of flash, an event that occurs after the start of the pre-ignition reaction at the local hotspot level. Both the simulation and experiment show that the difference between the times for poled P(VDF-TrFE)/nAl is ∼10% shorter than those for unpoled films. Although the determination of $tPIR$ from the simulation is based on a somewhat arbitrary threshold of 5% particles reaching the ignition temperature, any threshold level would produce the same trend that is consistent with the trend in the experimental measurement. Consequently, the analysis allows the relative importance of piezoelectricity and flexoelectricity to be established. Note also that both sets of data can be reasonably described by an empirical relation of the form $t=(α+βh0)−1/2$⁠.

A sequential multi-timescale and multi-physics framework for explicit microscale computational simulations of the electromechanical and thermal–electrodynamic responses of P(VDF-TrFE)/nAl composites to impact loading is developed. The framework integrates multi-physics models that explicitly resolve the microstructure and coupled mechanical–electrical–thermal processes underlying the charge accumulation, deformation, dielectric breakdown, and temperature rise in the material. The sequential models are integrated such that they utilize the same set of microstructures and allow for the information to pass from the first model (longer μs timescale, coupled electromechanical) to the next (shorter ns timescale, thermal–electrodynamic). Integration of the models allows the effects of mechanical loading, flexoelectricity, piezoelectricity, dielectric breakdown, Joule heating, thermal conduction, and most importantly, the development of critical hotspots to be tracked.

Analysis shows that dielectric breakdown due to charge accumulation from flexoelectricity and piezoelectricity can lead to sufficient hotspot development and initiation of chemical reactions within the P(VDF-TrFE)/nAl films. The electric field development (E-field) under mechanical loading and the consequent dielectric breakdown in the material are predominantly driven by the flexoelectric response and further modulated by the piezoelectric response of the binder.

The temperature field resulting from the breakdown analysis is used to establish the hotspot conditions for the onset of pre-ignition reactions. The results demonstrate that temperatures well above the ignition temperatures can be generated. Simulations show that flexoelectricity plays a primary role, and piezoelectricity plays a secondary role in the processes leading to reaction initiation, consistent with the experimental observation. In particular, the time to ignition (experiment) and the time to pre-ignition reaction (simulation) of poled microstructures (possessing both piezoelectricity and flexoelectricity) are ∼10% shorter, respectively, than those of unpoled films (possessing only flexoelectricity).

It is worthwhile to point out that, because the exothermic chemical reaction is not explicitly considered, the computational analyses entail assessing the critical condition for the microstructure to reach an arbitrary threshold of temperature rise (here, taken as 5% of the nAl particles reaching or exceeding the ignition temperature of the material) in order to estimate the time to pre-ignition reaction. Although this has allowed trends to be delineated, only the specific time to reaction initiation can be estimated. A more accurate prediction without the use of an arbitrary threshold requires an explicit account of the exothermic ignition reaction. This should be and is the focus of a separate, ongoing study whose outcome will be reported in the future.

Finally, although the calculations here focused on one specific material, the developed framework can be used for other heterogeneous materials as well.

This research was sponsored by the Air Force Office of Scientific Research MURI (Award No. FA9550-19-1-0008) (Program Manager: Dr. Mitat Birkan). Any opinions, findings, conclusions, or recommendations expressed in the article are those of the authors and do not necessarily reflect the views of the United States Air Force.

The authors have no conflicts to disclose.

Ju Hwan (Jay) Shin: Conceptualization (lead); Data curation (lead); Formal analysis (lead); Investigation (lead); Methodology (lead); Resources (equal); Validation (lead); Visualization (lead); Writing – original draft (lead); Writing – review & editing (lead). Derek K. Messer: Conceptualization (supporting); Data curation (supporting); Formal analysis (supporting); Investigation (supporting); Methodology (supporting); Resources (supporting); Validation (supporting); Visualization (supporting); Writing – original draft (supporting); Writing – review & editing (supporting). Metin Örnek: Conceptualization (supporting); Data curation (supporting); Formal analysis (supporting); Investigation (supporting); Methodology (supporting); Resources (supporting); Validation (supporting); Visualization (supporting); Writing – original draft (supporting); Writing – review & editing (supporting). Steven F. Son: Conceptualization (equal); Data curation (supporting); Formal analysis (supporting); Funding acquisition (equal); Investigation (supporting); Methodology (supporting); Project administration (equal); Resources (equal); Supervision (lead); Validation (supporting); Writing – original draft (supporting); Writing – review & editing (equal). Min Zhou: Conceptualization (lead); Data curation (equal); Formal analysis (equal); Funding acquisition (lead); Investigation (equal); Methodology (equal); Project administration (lead); Resources (equal); Supervision (lead); Validation (equal); Visualization (supporting); Writing – original draft (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.