Designing polymer nanocomposites with high energy density using machine learning

Addressing microstructure-property relations of polymer nanocomposites is vital for designing advanced dielectrics for electrostatic energy storage. Here, we develop an integrated phase-field model to simulate the dielectric response, charge transport, and breakdown process of polymer nanocomposites. Subsequently, based on 6615 high-throughput calculation results, a machine learning strategy is schemed to evaluate the capability of energy storage. We find that parallel perovskite nanosheets prefer to block and then drive charges to migrate along with the interfaces in x-y plane, which could significantly improve the breakdown strength of polymer nanocomposites. To verify our predictions, we fabricate a polymer nanocomposite P(VDF-HFP)/Ca2Nb3O10, whose highest discharged energy density almost doubles to 35.9 J cm−3 compared with the pristine polymer, mainly benefit from the improved breakdown strength of 853 MV m−1. This work opens a horizon to exploit the great potential of 2D perovskite nanosheets for a wide range of applications of flexible dielectrics with the requirement of high voltage endurance.


INTRODUCTION
Dielectric materials, which control and transfer energy electrostatically, play a key role in modern electric and electronic power systems ranging from portable electronic devices to medical equipment, hybrid electric vehicles, and pulse power weapons, etc [1][2][3] . When evaluating a dielectric material, one key figure of merit is the energy density U e calculated by where E is the electric field and D is the electric displacement. Hence, both high D and high breakdown strength E b are desirable to high U e , as illustrated in Supplementary Figure 1a. Moreover, the charge-discharge efficiency η of a dielectric material, which evaluates how efficiently the dielectric can store and release charges, is another critical parameter. Especially, η could drop dramatically at high electric fields or high temperatures mainly owing to the exponential increase of conduction loss, leading to a lot of energy loss and even thermal failure 4 . Therefore, both U e and η have been regarded as two pivotal parameters of measuring a dielectric material 5 . Compared with their ceramic counterparts, polymer dielectrics are flexible, scalable, and lightweight with higher breakdown strength E b and higher reliability, and are therefore an ideal choice to meet the increasing demand of compact, bionic, and integrated systems 6,7 . However, relatively low electric displacement D of polymers limits the capability of charge storage, resulting in an unsatisfactory U e . For example, despite with high η above 90%, the commercial dielectric polymers of biaxially oriented polypropylene exhibit only D ≤ 0.012 cm −2 under the applied electric field of 600 MV m −1 and corresponding highest discharged energy density U e ≤ 4 J cm −31 .
In recent years, by introducing inorganic nanofillers with high permittivity ε into the polymer matrix to construct composition effect, polymer nanocomposites have been extensively studied and become a popular and effective method to obtain both high D and E b and thus improved U e 8,9 . Up to now, the maximal discharged energy density of state-of-the-art two-phase polymer nanocomposites with ceramic nanofillers mainly varies from~5 tõ 25 J cm −3 , whose values depend on different types of polymers and nanofillers, as partly summarized in Supplementary Table 1 and Supplementary Figure 2 [10][11][12][13][14][15][16][17][18][19][20][21][22][23][24] . Generally, nanofillers with high ε, such as BaTiO 3 and SrTiO 3 , are preferred for the improvement of D of polymer nanocomposites at the expense of possibly decreasing E b or η 25,26 . On the contrary, some low-ε nanofillers, such as Al 2 O 3 and BN, usually have greater potential for enhancing insulative and E b than their counterparts 9,14,23 . And more notably, the morphology of nanofillers also has a strong influence on the dielectric performance of polymer nanocomposites. For example, among different morphologies of Al 2 O 3 , two-dimensional (2D) nanosheets seem to be more favorable for the reduction of conduction loss and thus improvements of E b and thermal stability 27 . However, the underlying mechanisms behind those performance differences remain unclear. Therefore, building sensible microstructure-property links to gain fundamental insights on dielectric behaviors is instructive to the innovative design of two-phase polymer nanocomposites.
The phase-field method, based on thermodynamics and kinetics of inhomogeneous systems at mesoscale, has been widely developed to fundamentally understand and predict the microstructure evolution and its correlation with material behaviors under different physical stimuli 28 . Inspired by its applications including domain switching, ion transport, crack growth, and phase transformation 29,30 , the phase-field method has been extended to study the microstructure-property relations of polymer nanocomposites, especially for the dielectric breakdown process 31,32 . Furthermore, partly propelled by the Materials Genome initiative, theoretical simulation, coupled with machine learning or/and high-throughput calculations, has become a burgeoning paradigm for more rapid developments of materials [33][34][35] . Therefore, complementary to experimentations, the phase-field model combined with machine learning, could be an effective but time-and labor-saving way to discover underlying microstructure-property relations and inversely design polymer nanocomposites.
In this work, we develop a comprehensive phase-field model with three modules to investigate the nanofiller effect on the dielectric response, charge behavior, and breakdown process of polymer nanocomposites. Using five representative nanocomposites with 0D nanoparticles (np), 1D vertical nanofibers (v-nf), 1D parallel nanofibers (p-nf), 2D vertical nanosheets (v-ns), and 2D parallel nanosheets (p-ns) as examples, the underlying physical mechanisms are analyzed by simulating the redistributions of local electric field and polarization field, the transport behaviors of charges and the evolution of breakdown path. Then, by performing 6615 high-throughput phase-field calculations, the data set bridging microstructures and properties of polymer nanocomposites is constructed as the training dataset of machine learning. Using a scoring function describing the energy density, we evaluate the capability of energy storage of 2205 polymer nanocomposites and then screen some potential nanofillers using the algorithm of Back Propagation Neural Network (BPNN). Finally, we prepare P(VDF-HFP)/Ca 2 Nb 3 O 10 (poly vinylidenefluoride-hexafluoropropylene copolymer/Ca 2 Nb 3 O 10 perovskite nanosheets) nanocomposite film to validate our predictions and obtain excellent energy storage performance.

RESULTS
Phase-field simulations of the nanofiller effect in polymer nanocomposites Here, we take the polymer nanocomposite P(VDF-HFP)/BaTiO 3 as an example to study how the nanofillers affect the dielectric response, charge behavior, and breakdown process by performing phase-field simulations. First, we build five representative models of polymer nanocomposites with one nanoparticle (np), one vertical nanofiber (v-nf), one parallel nanofiber (p-nf), one vertical nanosheet (v-ns), and one parallel nanosheet (p-ns), as shown in Fig. 1a−e, respectively. The simulation size is 100 nm × 100 nm × 100 nm and the applied voltage is 1 V, so the applied electric field is 10 MV m −1 . Here, the permittivities of P(VDF-HFP) matrix and BaTiO 3 nanofiller are fixed with the values of 10 and 300, respectively. To highlight the morphology features of different nanofillers, the length-diameter ratios of nanofillers are designed as shown in Fig. 1 with the same volume fraction. By performing a phase-field model, we could obtain the redistributions of local electric field and electric displacement by solving Poisson's equation. To display the local characteristics clearly, especially at the interface region, 2D cross-sectional distributions of local electric field and electric displacement are plotted in Fig. 1f-j and 1k-o, respectively. It can be seen that the local electric field could be distorted at interfacial regions owing to the large mismatch of permittivity between P(VDF-HFP) matrix and BaTiO 3 nanofiller. Note that the distribution of local electric fields is also influenced by the morphology and orientation of nanofillers, leading to different degrees of electric field accumulation. For the nanoparticle, the vertical nanofiber, and the parallel nanosheet in Fig. 1f, g, i, the local electric field around the interfacial region severely concentrates at two ends along the direction of the applied electric field. Especially for the vertical nanofiber, the highest local electric field could reach up to above 40 MV m −1 , which is about four times as much as the applied electric field. However, for the parallel nanofiber and the parallel nanosheet in Fig. 1h and 1j, the local electric field is more homogeneous than that of corresponding vertical nanofillers. More interestingly, compared to the parallel nanofiber, the parallel nanosheet with a smaller thickness along z direction is more advantageous to alleviate the accumulation of local electric field, where the highest electric field in Fig. 1j is lower than 0.4 MV m −1 . Then, the local electric displacement D induced by the electric field could be approximately calculated by the linear equation, D = ε 0 ε r E as shown in Fig.  1k−o. It is shown that the vertical nanofillers in Fig. 1l and 1n exhibit higher electric displacement than their parallel counterparts in Fig. 1m and o. For example, the electric displacement inside the vertical nanosheet is~4 × 10 −3 C m −2 in Fig. 1n, but the electric displacement in the parallel nanosheet is lower thañ 1.1 × 10 −3 C m −2 owing to the ultralow local electric field, as displayed in Fig. 1o. As a result, the effective permittivity of polymer nanocomposite with one vertical nanofiber is the highest among all five nanocomposites, as the black dotted line indicated in Fig. 1p. To study the effect of dielectric mismatch, Fig. 1p gives the predicted effective permittivity of these five polymer nanocomposites as a function of the permittivity of nanofillers from 10 0 to 10 4 . It can be seen that the vertical nanofiber is more conducive to improving the effective permittivity of polymer nanocomposites than others, which is consistent with existing experimental results 16,17 . With the increase of the permittivity of nanofillers from 10 1 to 10 4 , the enhancement of effective permittivity gradually gets saturated. For the nanoparticle and parallel nanofillers, the saturated points of effective permittivity appear early than those of vertical nanofillers, which strongly limits the contribution of high-ε nanofillers to the improvement of effective permittivity.
To clarify the underlying physical mechanism behind those differences, Fig. 1q illustrates the local polarization behavior in a simplified laminated model with one nanofiller interlayer and two bilateral polymer layers. When applying an external electric field, both polymer matrix and nanofiller could be polarized. However, if without enough compensating charges, the polarized charges would build a depolarization field, whose direction is opposite to the applied electric field. Thus, the local electric field is determined by the contributions from both the applied electric field and the reverse depolarization field. Therefore, for the parallel nanofillers shown in Fig. 1m, o, the depolarization field is almost totally applied to the whole nanofiller and thus immensely suppresses the polarization of parallel nanofillers. Furthermore, under the same polarization condition, the thinner the thickness is, the greater the permittivity is, the stronger the depolarization effect is. Thus, the parallel nanosheet exhibits lower permittivity and earlier saturation point than the parallel nanofiber as displayed in Fig. 1p. Although for vertical nanofillers in Fig. 1l and 1n, the depolarization effect could be greatly mitigated due to the large aspect ratio, resulting in higher effective permittivity as shown in Fig. 1p. Therefore, besides the intrinsic dielectric characteristics, the morphology and orientation of nanofillers are also important factors influencing the effective permittivity, which explains that why different shapes or/and orientations of nanofillers could effectively modulate dielectric performances in experiments 8,26 .
Then, we simulate the nanofiller effect on the charge transport process in polymer nanocomposites under a high electric field. As illustrated in Fig. 2a, here we consider four different states of charges which may affect the charging and discharging process: ① polarized bound charges, ②intrinsic bulk charges, ③injected charges from electrodes, and ④trapped charges at interfaces 36,37 . For example, some charges may be bound at a limited range owing to the polarization dipoles. With the increasing of external stimulus, some charges in dielectrics may be activated and then migrate or draft, forming a continuous electric current on the macro level. When at extremely high electric fields and/or temperatures, some charges accumulating at the interfaces between the electrode and the dielectric may overcome the Schottky energy barrier and then be injected into the dielectrics. Moreover, a large number of interface traps in polymer nanocomposites may influence the transport process by trapping charges. Besides, we also consider some intrinsic nature of nanofillers such as anisotropy, as more details provided in the section of Methods. For example, due to the 2D quantum confinement and surface effects, the dielectric and electrical properties of nanosheets are anisotropic 38,39 . In order to better reveal the nanofiller effect on the behaviors of the above four charges, we only display the electron density in this section as shown in Fig. 2b-p. It is shown that the maximum value of charge density could reach up to~20 C m −3 at the applied electric field of 500 MV m −1 , which depends on the applied electric field, Schottky barrier, temperature, and so on. The charge transport process depends on not only the morphology and orientation but also the intrinsic features of nanofillers. On the one hand, the vertical nanofillers are more beneficial to the charge transport than the parallel nanofillers owing to the different effective lengths in the direction of the applied electric field resulting from different aspect ratios. For example, as shown in Fig. 2e-g, charges preferentially move via the nanofiber with high electron mobility along z direction and little interface traps in x-y plane. On the other hand, for the same orientation, the charge transport process in the polymer nanocomposite with the parallel nanosheet ( Fig. 2n−p) is significantly slower than that with the parallel nanofiber ( Fig. 2h−j), mainly caused by the intrinsic ultralow out-of-plane electron mobility of nanosheets and abundant interface traps in x-y plane 40,41 . Apparently, the parallel nanosheet acts as a barrier to block charges from transferring along z direction and then bend the transport path along the interfaces/surfaces. Therefore, parallel nanosheets could be taken as the charge-bearing barrier to block the charge transport to reduce the conduction loss, while vertical nanofibers could play the opposite effect to enhance the electrical conductivity.
Next, taking the above five kinds of polymer nanocomposites as examples, we carry out the phase-field model to study the nanofiller effect on the breakdown process, as displayed in Fig. 3a-e. The breakdown process is described by two representative states with broken phase in black, as the 2D sectional views plotted in Fig. 3f−o. It can be seen that all breakdown paths start at the ends of nanofillers around the interfaces between the polymer matrix and nanofillers. For nanoparticles in Fig. 3f, vertical nanofibers in Fig. 3g, and vertical nanosheets in Fig. 3i, the breakdown phases begin to grow at the applied electric fields of 271 MV m −1 , 250 MV m −1 , and 258 MV m −1 , respectively. Although for parallel nanofibers in Fig. 3h and parallel nanosheets in Fig. 3j, the breakdown phase does not appear until the applied electric field reaches 345 MV m −1 and 500 MV m −1 , respectively. In addition to the difference in breakdown initiation, the morphology of the breakdown path is also quite different in different nanocomposites. Compared with nanoparticles and vertical nanofillers in Fig. 3k, l, and n, the breakdown paths in polymer nanocomposites with parallel nanofillers are much thicker in Fig.  3m, o. Especially for parallel nanosheets, the broken phase looks like a stout tree trunk whose thickness is almost close to the width of nanosheets. Thus, among these nanocomposites, the nanocomposite with parallel nanosheets could withstand the highest electric field of above 580 MV m −1 , which is approximately twice as many as the breakdown strength of nanocomposites with vertical nanosheets (283 MV m −1 ) or nanofibers (269 MV m −1 ). All of the above differences in the breakdown process could be elucidated from the view of the local electric field, as shown in Fig.  3p-t. For nanoparticles and vertical nanofillers, more charges and higher local electric fields concentrate at the interfaces and then interconnect with each other, easily leading to a fast-forward slender breakdown channel along the nanofillers. However, for parallel nanosheets in Fig. 3t, benefit from the homogeneous electric field and high out-of-plane voltage resistance of nanosheets, the breakdown path propagates slowly along the direction of z axis and then becomes thicker and thicker. Here, we call this phenomenon the blocking effect to describe the role of nanofillers on delaying and hindering the propagation of the breakdown path. Although for parallel nanofibers in Fig. 3c, though the shape in y-z plane is rectangular, the shape in x-z plane is still circular. Thus, the concentration of the local electric field in Fig. 3r is also serious and the breakdown path may cut corners in x-z plane, which greatly weakens the blocking effect of parallel nanofibers. Thus, the breakdown path in Fig. 3m is thinner than that in Fig. 3o. Therefore, the polymer nanocomposite with parallel nanosheets exhibits the strongest blocking effect, leading to the highest breakdown strength among those five types of polymer nanocomposites. The blocking effect could be further understood from the electron avalanche theory, as three schematic diagrams shown in Fig. 3u-w. The dielectric breakdown typically induces by the avalanche multiplication of free charge carriers, when carriers acquire sufficient kinetic energy to give a high probability of collision and ionization 42 . For nanoparticles and vertical nanofillers, impact ionization may be triggered preferentially at the interfaces with the high local electric field and then more electrons would be accelerated and collided to form a chain reaction, leading to electron avalanche breakdown immediately. However, for parallel nanosheets in Fig. 3w, the initial impact ionization could be delayed and the chain reaction is restricted in x-y plane due to insufficient kinetic energy and high insulativity perpendicular to the direction of electrostatic force. Thus, the electron avalanche would be blocked and navigated along the x-y plane of nanosheets, resulting in a tortuous breakdown path. As we discussed above, the blocking effect of parallel nanofibers is less effective than that of parallel nanosheets. Because the circular shape of nanofibers in x-z plane (Fig. 3c) would also cause severe concentration of local electric field, leading to that the dominating electron avalanche mode in x-z plane is similar to that of nanoparticles in Fig. 3u. Therefore, introducing 2D nanosheets to construct the blocking effect into nanocomposites, could effectively manipulate the behaviors of charges and breakdown path to remarkably improve the insulativity and voltage endurance of nanocomposites. In summary, the blocking effect of 2D nanosheets in polymer nanocomposites originates from two parts: one is the 2D morphology and orientation resulting in homogeneous electric field and abundant interface areas in x-y plane, the other is the anisotropy nature with high out-of-plane insulation. It is worth mentioning that the blocking effect is also influenced by some other factors such as the distribution and size of 2D nanosheets. that Next, we will use the machine learning method to find suitable nanofillers to make better use of the blocking effect to improve the dielectric properties of nanocomposite dielectrics.

Machine learning-assisted performance evaluation and targeted experiments
In this section, we first perform high-throughput phase-field simulations to build a data set as the training set of machine learning. As shown in Fig. 4, the input data set I S i ; ε j ; σ k ; V l À Á of high-throughput calculations includes four primary nanofiller features: morphology S i , permittivity ε j , electrical conductivity σ k , and volume fraction V l . By performing the phase-field simulations with three modules, we could output the corresponding data set of effective properties of nanocomposites P ε eff ijkl ; E b ijkl ; σ eff ijkl with the effective permittivity ε eff ijkl , breakdown strength E b ijkl and effective electrical conductivity σ eff ijkl . Here, we consider five morphologies, 21 permittivies from 3 to 7000, 21 electrical conductivities from 1 × 10 −16 S m −1 to 1 × 10 −6 S m −1 , and a fixed volume fraction of 1%, resulting in 6615 combinations in total. As the maps are shown in Supplementary Figure  3-5, the nanocomposites with vertical nanofillers exhibit higher effective permittivity and higher effective electrical conductivity while the nanocomposites with parallel nanofillers possess higher breakdown strength, which agrees well with the experimental results of some polymer nanocomposites with 2D insulting nanofillers, such as BN, Al 2 O 3 , and TiO 2 9,20,27 . Then, by considering the contribution of conduction loss to the discharged energy density at the high electric field with linear approximation 43 as displayed in Supplementary Figure 1b, we propose a scoring function to evaluate the capability of energy storage of different nanocomposites as follows: where Δt is the period of the applied electric field. Using the Min-Max normalization method, we could obtain the scoring data set with a value between 0 and 1, as shown in Fig. 5a and Supplementary Figure 6. The closer the scoring gets to 1, the stronger the capability of energy storage is. As a whole, nanocomposites with parallel nanosheets score higher than their counterparts, which benefits from high breakdown strength and low conduction loss. To predict performance faster and screen more potential nanofillers, we build a machine learning model with BPNN algorithm by taking the scoring dataset as the training set. Figure  5b plots the comparisons of machine learning results f ML s and phase-field scoring results f PF s of all data sets including 90% training set and 10% testing set, marked as hollow symbols. It can be seen that all hollow symbols are dispersed around the solid black line which represents f PF s ¼f ML s , indicating the good reliability of machine learning predictions. The Person correlation coefficient (PCC) and mean square error (MSE) with 0.97 and 9.01 × 10 −4 , respectively, also suggest that the accuracy of this machine learning model is very high. Then, we are dedicated to discovering more promising 2D nanofillers to guide the experimental design of polymer nanocomposites. After screening many existing 2D materials by the machine learning strategy, we have found that some 2D nanosheets consisting of perovskite building blocks of TiO 6 , NbO 6 , or TaO 6 octahedra 40 , such as Sr 2 Ta 3 O 10 , Ca 2 Nb 3 O 10 , LaNb 2 O 7 , Sr 2 Nb 3 O 10 , and Ca 2 Ta 3 O 10 , get very high scoring as solid symbols shown in Fig. 5b. All material parameters of those representative nanosheets are listed in Supplementary Table 2. Compared with other nanofillers, perovskite nanosheets could stay high out-of-plane insulating performance with low leakage current density and superior high permittivity (>100) even at the thickness of a few nanometers 38,39 , which could greatly gain the blocking effect. Therefore, 2D perovskite nanosheets with enhanced blocking effect from the superhigh capability of outof-plane insulation, high permittivity, and large blocking interface areas along the direction of the applied electric field, show great potential on high energy storage by improving the breakdown strength and decreasing conduction loss.
Among these 2D perovskite nanosheets predicted above, we choose a polymer nanocomposite P(VDF-HFP)/Ca 2 Nb 3 O 10 as the targeted materials to verify our ideas. More details about the experimental process and characterization are discussed below in Methods. As shown in Fig. 5c and Supplementary Figure 8a, after introducing Ca 2 Nb 3 O 10 nanosheets with the volume fraction of 0.1% into the P(VDF-HFP) matrix, the breakdown strength could be improved from 649 MV m −1 of pristine P(VDF-HFP) to 853 MV m −1 with an increment of 31%. As a result, the maximal discharged energy density of P(VDF-HFP)/Ca 2 Nb 3 O 10 nanocomposite reaches as high as 35.9 J cm −3 , which is 1.95 times higher than pristine P (VDF-HFP) with 18.4 J cm −3 . The charge-discharge efficiency of polymer nanocomposites with Ca 2 Nb 3 O 10 nanosheets is also promoted a little bit by decreasing the conduction loss as shown in Supplementary Figure 8b. However, due to the high intrinsic ferroelectric loss of P(VDF-HFP), the efficiency of P(VDF-HFP)/ Ca 2 Nb 3 O 10 nanocomposite is still low. Therefore, an almost double increase of energy density and slightly improved efficiency indicate that the blocking effect constructed by Ca 2 Nb 3 O 10 nanosheets does work. Therefore, we believe this strategy could be extended as a general method to substantially improve the capability of energy density, with the fast development of advanced technology of 2D materials preparation and exfoliation 39,44 . Furthermore, the phase-field model and machine learning strategy in this work could be also extended to study the energy density capability in different dielectrics, such as blending polymer, segmented copolymer, semi-crystalline polymer, and multiphase ceramics.

DISCUSSION
Here, we build an integrated phase-field model with three modules to study the relations between microstructures (0D, 1D, and 2D nanofillers) and dielectric properties (effective permittivity, breakdown strength, and electrical conductivity) of polymer nanocomposites. To understand the advantages of 2D nanosheets on improving the insulative and voltage endurance, we propose the blocking effect to describe the blocking and buffering effect of 2D nanofillers on charge transport and the breakdown process. Thus, resulting from the enhanced breakdown strength and reduced conduction loss, the discharged energy density of nanocomposites with 2D nanosheets could be substantially improved. With the machine learning strategy, we evaluate the capability of energy storage of 2205 different nanocomposites by the scoring function and screen some promising perovskite nanosheets with TiO 6 , NbO 6 , or TaO 6 octahedra. Then, we fabricate a nanocomposite film of P(VDF-HFP)/Ca 2 Nb 3 O 10 as an example to verify our predictions. Almost double increased energy density and slightly improved efficiency indicate that 2D perovskite nanosheets have great potential of tailoring the dielectric performance, which is far outperforming state-of-the-art twophase polymer nanocomposites, as shown in Supplementary  Figure 2. Therefore, in this work, we develop a simulation-aided research paradigm to address the underlying mechanism, summarize the rules and predict the dielectric performances of polymer nanocomposites and finally achieve experimental verification.

Phase-field model
In this work, the comprehensive phase-field framework mainly includes three modules, to study the dielectric response, breakdown process, and charge transport behavior of polymer nanocomposites, respectively. Besides, two basic shared modules are also built to generate different microstructures and process data.

Module of permittivity
The function of this module is to calculate the local distributions of distorted electric field and displacement and then obtain the effective permittivity in different systems. In an electrostatic equilibrium state, the local electric field in a dielectric is determined by Poisson's equation as follows: where D(r) is the electric displacement, ε 0 is the vacuum permittivity, ε(r) is the permittivity, E(r) is the electric field, and ρ(r) is the free charge density. In this module, after inputting the spatially dependent parameters at a given microstructure, the spectral iterative perturbation method is employed to solve Poisson's equation 45 . Then the effective permittivity tensor ε ij could be calculated according to where À represents the average value.

Module of breakdown
This module is used to spatiotemporally simulate the evolution of the breakdown phase under the increasing applied electric field. Based on the phase-field theory, a temporally and spatially dependent order parameter η(r,t) is introduced to describe the evolution of the breakdown phase with η(r,t) = 1 broken phase and η(r,t) = 0 unbroken phase. During the electrostatic breakdown process, the total free-energy consists of the separation energy, the gradient energy, and the electrostatic energy, i.e., Here, a double-well function is used to describe the phase separation energy as follow, where α is the energy barrier of phase separation between the broken phase and unbroken phase with a value of 10 8 J m −3 . The second term of gradient energy density could be written by where γ is the gradient energy coefficient. A value of γ = 1.6 × 10 −10 J m −1 is used to specify d 0 1 nm in the modeling. The electrostatic energy density is determined by Then, inspired by the crack propagation, the breakdown process starts when the total energy at one point is larger than its critical energy endurance. Thus, a modified Allen-Cahn equation is employed to control the dynamic evolution of the breakdown process, where L 0 is the kinetic coefficient defining the interface mobility with a value of 1.0 m 2 S −1 N −1 , H(|f elec |−|f critical |) is the Heaviside unit step function, and f critical is a position-dependent material constant representing the critical energy density. In this work, the characteristic length scale d 0 ¼ ffiffiffiffiffiffiffi ffi γ=α p and the characteristic time scale t 0 = 1/(L 0 α) can be determined by the material parameters γ, α, and L 0 . A grid size of N x Δx × N y Δx × N z Δx with grid space of Δx = d 0 is used. The simulation size is N x = N y = 1200 and N z = 1 and the time interval is Δt = 0.01t 0 .

Module of charge
To study the charge behavior in polymer nanocomposites under high operating voltage, bipolar charge injection and transport model is built based on Schottky injection mechanism 37,46 . Here, the current density at the interface between the electrode and the dielectric is determined by where A and k B is the Richardson constant and Boltzmann constant, T is the temperature, w i is the Schottky barrier with a constant value of 1 eV in this work (i = 1 for electron and i = 2 for hole), q is the elementary electron charge, E is the electric field and ε is the permittivity. The charge behavior in polymer nanocomposites is described by transport equation, continuity equation, and Poisson's equation as following, where μ a is the apparent electron mobility, n is the charge density, D is the charge diffusion coefficient which could be calculated by the Einstein relation D/μ = k B T/q, φ is the electric potential and s is the source term.
Here, the apparent electron mobility μ a is introduced in this model to briefly consider the detrapping and trapping process at the interface regions between the matrix and nanofillers, which is expressed by where μ 0 is intrinsic electron/hole mobility with values of 1 × 10 −11 cm 2 V −1 s −1 in polymer and 0.1 cm 2 V −1 ss −1 in nanofillers, and ξ is the trap depth at the interfaces 47 . In this work, nanoparticles and nanofibers are regarded as isotropic, but the dielectric properties of nanosheets are set to strongly anisotropic with a 1000 times difference between the out-of-plane electron mobility and in-plane electron mobility.

Machine learning
The training data set of machine learning is generated by performing highthroughput phase-field calculations of three modules with the same grid size of 1200 × 1200 × 1. All data in machine learning is normalized between 0 and 1 by Min-Max normalization as follows: x ¼ where x min and x max are the minimum and maximum of the data. Then, we propose a scoring function of Eq. (2) to quantify the capability of energy storage based on the linear approximation as illustrated in Supplementary  Figure 1b. Here, we assume that the measured electric displacement at zero field D r mainly comes from the leakage current without considering the ferroelectric loss. Then, the relationship between the effective Z.-H. Shen et al.
conductivity σ eff and D r can be derived as follows 43 : (16) where f is the frequency of the applied electric field. Here, Δt ¼ 1=f ¼ 1=10 ¼ 0:1s is used in this approximation. Then, the energy loss caused by electrical conduction could be calculated as follows: Therefore, the maximal discharged energy density could be given under the linear approximation: Thus, the scoring function is closely related to the energy density, whose values could reflect the capability of energy storage of polymer nanocomposites. Then, we normalize the scoring data by the Min-Max method and then perform the machine learning model via the BPNN algorithm with 90% of the scoring data as training set and 10% of the scoring data as testing set.
To verify the reliability of the machine learning model, PCC and MSE are used to evaluate the performance of the machine learning model. Generally, PCC between two variables x and y is defined as follows: where N is the number of pairs of scores and ∑ is the sum notation. PCC returns a value between −1 and 1 to reflect that how strong the relationship between the two variables is. A higher absolute value indicates a stronger relationship between variables. MSE is usually used to relaying the concepts of precision, bias, and accuracy during the prediction process by machine learning. It could be calculated by taking the average of the square of the difference between the original and predicted values of the data as follows: where m is the number of pairs, y i and y i are the original and predicted data, respectively. The smaller the MSE, the closer the predicted data is to the actual data.

Fabrication and characterization of P(VDF-HFP)/Ca 2 Nb 3 O 10
First, the solid powders of layered perovskite KCa 2 Nb 3 O 10 (KCNO) were prepared by a solid-state reaction method using K 2 CO 3 , CaCO 3, and Nb 2 O 5 powders. Then, KCNO powders were stirred in a 6 mol L −1 HNO 3 aqueous solution at 60°C for 24 h. Next, the solid sediments were filtered, washed by deionized water, and then dried in the air, and become HCa 2 Nb 3 O 10 ·1.5H 2 O. Then, mixing with TBAOH aqueous solution with a molar ratio of 1:1, the mixture was ultrasonicated at 60°C for 40 h to obtain a colloidal suspension. Next, PEI aqueous solution was introduced dropwise to produced white precipitate with the sandwiched structure of Ca 2 Nb 3 O 10 and PEI. Finally, the precipitate was washed and centrifuged with deionized water and anhydrous alcohol, and then dispersed into DMF. For nanocomposite films, P(VDF-HFP) powders were dissolved into DMF and stirred for 10 h, and then mixed with Ca 2 Nb 3 O 10 nanosheets. Then, the suspension was cast onto the quartz glass and dried at 110°C for 10 h in air. Finally, the film was uniaxially stretched at 80°C. For energy storage characterization, the ferroelectric loops were measured at 10 Hz by Radiant Technologies Precision Premier II (Radiant Tech.) equipped with a high voltage amplifier (TREK MODEL 609B). The breakdown strength is described by a two-parameter Weibull distribution function.

DATA AVAILABILITY
All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.