INTRODUCTIONDielectric constant and polarizability are important properties of polymers, which determine the electrical behavior and properties of polymers in various applications. The accurate prediction of dielectric constant plays an important role in understanding the response of materials to electric field and designing high-performance polymers and polymer-based devices.[1] The dielectric constant can be calculated by using molecular dynamics simulations (MD) or the first-principle calculations, such as density functional theory (DFT).[2,3] MD simulations were widely used to predict the dielectric constant of poly(vinylidene fluoride) (PVDF).[4−8] At first, the dielectric constant of PVDF was calculated for crystalline and amorphous based on MSXX force field and MD simulations by Karasawa and Goddard.[4] It shows a high dielectric constant for amorphous PVDF. Yan et al. also studied the dielectric constant of amorphous PVDF at different temperature by MD simulations by applied a constant electric field.[6] The electronic contribution to the dielectric constant is hard to calculate based on MD and force field. In addition, the dielectric constant can also be calculated using DFT and perturbation theory.[9−11] The finite field method is often employed to predict polarizability of polymer based by combining with DFT.[9,10] The calculation of molecular polarizability by DFT method or finite field method is all on the basis of optimized geometry,[12] so neither of these two methods can estimate the contribution of conformation change to polarizability. Thus, the polarizability calculated by DFT is always underestimated. So far, there is no report on the calculation of polarizability considering both electronic and molecular conformational changes.PVDF is a widely used organic polymer possessing a linear chain with the repeat unit of CH2 and CF2. The dielectric constant of amorphous PVDF is about 11, which is much higher than that of many other polymer dielectric materials.[4,13,14] Its high dielectric constant mainly comes from the strong dipoles formed alternately by carbon-fluorine and carbon-hydrogen bonds in the main chain of the polymer. Due to its high dielectric constant and low dielectric loss, PVDF is an ideal material for a wide range of applications in the field of modern electronics and electrical engineering.[15] Because of its unique piezoelectric properties, PVDF is also widely used in the field of energy. It generates an electric charge in response to mechanical stress, making it a promising material for sensors, actuators, and energy harvesting devices.[1,14,16,17] PVDF is semi-crystalline and consists of both crystalline and amorphous regions, it has four different α, β, γ, δ crystal phases.[13,16,18] α, δ Crystal phases consist of trans-gauche twist (TG) conformations, β phase has all the trans zigzag structures (TTTT), and γ crystal phase is with the TTTGTTTG’ molecular sequence (as shown in Fig. 1).[16,19] The crystallinity and polymorphic phase can be controlled through various methods such as annealing, stretching.[21,22] As various of composite dielectric material with enhanced characteristics can be obtained by blending PVDF with other polymers or additives such as polyethylene, polystyrene, or polypropylene,[23−25] the dielectric property of amorphous PVDF is of great interest. The dielectric constant strongly depends on its conformation.[13] It is necessary to study the conformation of PVDF effects on dielectric property. In particular, the dielectric constant of PVDF with TTTT conformation calculated by the usual first-principles method is usually seriously underestimated.Fig 1Optimized structures of PE, PTFE, and α-, δ-, γ-PVDF at zero external field. The figure is drawn using the GaussView 6 package.[20]In this study, a new method was proposed to estimate the polarizability and dielectric constant of a polymer by combining DFT with the finite field method, which can well estimate the effect of electronic and geometric conformation changes on dielectric constant. This method was first validated by polyethylene (PE) and polytetrafluoroethylene (PTFE). Then, the effects of PVDF conformation on dielectric constant and polarizability were studied.METHODSDensity Functional Theory with Perturbation MethodBased on the perturbation theory, the energy E with respect to the homogeneously perturbational external field F is expended in a Tayler series as:[12,26]1$ E\left(\boldsymbol{F}\right)=E\left(0\right)+{\left.\frac{\partial E}{\partial \boldsymbol{F}}\right|}_{\boldsymbol{F}=0}\boldsymbol{F}+{\left.\frac{1}{2}\frac{{\partial }^{2}E}{\partial {\boldsymbol{F}}^{2}}\right|}_{\boldsymbol{F}=0}{\boldsymbol{F}}^{2}+{\left.\frac{1}{6}\frac{{\partial }^{3}E}{\partial {\boldsymbol{F}}^{3}}\right|}_{\boldsymbol{F}=0}{\boldsymbol{F}}^{3}\cdots $where the first and second derivative are the permanent dipole moment and dipole polarizability, respectively:2$ {\mu}_{0}=-{\left.\frac{\partial E}{\partial \boldsymbol{F}}\right|}_{\boldsymbol{F}=0} ; {\alpha }=-{\left.\frac{{\partial }^{2}E}{\partial {\boldsymbol{F}}^{2}}\right|}_{\boldsymbol{F}=0} $The components of the polarizability, αij where i and j indicate the Cartesian components, can be elaborated further in the context of excited electronic states under zero external field conditions (F=0) as below:[12]3$ {\alpha }_{ij}=2\sum _{n=1}\frac{\langle {{\Psi }}_{n}\left|{\widehat{\mu }}_{i}\right|{{\Psi }}_{0}\rangle \langle {{\Psi }}_{0}\left|{\widehat{\mu }}_{j}\right|{{\Psi }}_{n}\rangle }{{E}_{n}-{E}_{0}} $where $ \widehat{{\mu }} $ is the dipole operator, $ {E}_{0} $ and $ \left\{{E}_{n}\right\} $ are the energies of a ground state $ |{{\Psi }}_{0}\rangle $ and excited states $ \left\{|{{\Psi }}_{n}\rangle \right\} $, respectively. The average of the diagonal elements of the polarizability tensor is used to define the isotropic or mean polarizability αiso:4$ {\alpha }_{\mathrm{i}\mathrm{s}\mathrm{o}}=\frac{1}{3}\left({\alpha }_{xx}+{\alpha }_{yy}+{\alpha }_{zz}\right) $Density Functional Theory with Finite FieldThe finite field (FF) method employs the numerical differentiation to calculate the polarizability. Assuming that a finite electric field is applied in x direction, then the polarizability element $ {\alpha }_{xx} $ is:5$ {\alpha }_{xx}=\frac{1}{{F}_{x}^{2}}\left\{2E\left(0\right)-\left[E\left({F}_{x}\right)+E\left(-{F}_{x}\right)\right]\right\} $where $ E\left(0\right) $ and $ E\left({F}_{x}\right) $ represent the energy with a field being 0 and $ {F}_{x} $, respectively.[10] This process can be replicated for different orientations of the field to obtain all the components of polarizability. Note that the field’s intensity must be selected with caution to prevent the possibility of unphysical effects.Finite Field Combined with Geometry OptimizationInstead of relying on a static optimized geometry, the molecular conformation can change in response to an applied field, as will be discussed in the following section. This phenomenon has not been taken into account until now. To address this, we introduce a new variable, G(F), which represents the geometry that has been optimized under the influence of an applied field. Eq. (5) can be rewritten as follows:6$ {\alpha }_{xx}=\frac{1}{{F}_{x}^{2}}\left\{2E\left(0\right)-\left[E\left(\boldsymbol{G},{F}_{x}\right)+E\left(\boldsymbol{G},-{F}_{x}\right)\right]\right\} $Through this combination of the FF method and geometry optimization, both electronic and conformational changes in the molecule are taken into account.Dielectric Constants of PolymersThe dielectric constants εr of polymer materials can be calculated via the rewritten Clausius-Mossotti equation given that αiso is expressed in terms of the molecular polarizability volume as shown in Ref. [27]:7$ \frac{{\epsilon }_{{\mathrm{r}}}-1}{{\epsilon }_{{\mathrm{r}}}+2}=\frac{4{\text{π}} {N}_{{\mathrm{A}}}{\alpha }_{{\mathrm{iso}}}}{3{V}_{{\mathrm{m}}}} $where Vm is the molar volume of the materials and NA is Avogadro constant. In present work, the polarizable continuum model (PCM), a widely used implicit solvation model, is used to calculate the volume of the polymer chains.[28,29] All the calculations are conducted with the Gaussian 16 software[30] and the DFT method with the B3LYP functional at 6-31+G(d, p) basis set is employed.[31,32] Note that 6-31+G(d, p) basis set could provide a good balance between precision and computational efficiency.[2] The number of monomers for PE, PTFE, and PVDF is 6.RESULTS AND DISCUSSIONThe optimized structures of PE, PTFE, and α-, δ-, γ-PVDF at zero external field are shown in Fig. 1. The initial structures of PVDF are extracted from references.[33,34] Note that PTFE terminates with fluorine atoms instead of hydrogen to avoid the introduction of strong dipoles from carbon-fluorine and carbon-hydrogen bonds at the ends. The positive direction of the x-axis is to the right, the positive direction of the y-axis is upwards, and the positive direction of the z-axis is outwards from the surface. All the polymer chains are aligned in the direction of the x-axis for convenience.Variations of Energy with Applied Electric FieldsWe first validate the intensity of the applied field by observing the energy responses. According to Eq. (6), the variances of energy $ \Delta {E}_{i=x,y,z} $ in different directions and the isotropic or mean energy difference $ \Delta {E}_{{\mathrm{iso}}} $ are defined as below:8$ \Delta {E}_{i}=2E\left(0\right)-\left[E\left(\boldsymbol{G},{F}_{i}\right)+E\left(\boldsymbol{G},-{F}_{i}\right)\right],i=x,y,z $9$ \Delta {E}_{\mathrm{i}\mathrm{s}\mathrm{o}}=\left(\Delta {E}_{x}+\Delta {E}_{y}+\Delta {E}_{z}\right)/3 $Note that theoretically, the variances of energy are proportional to the square of the magnitude of the applied electric field $ \left|F\right| $.Fig. 2 illustrates the dependences of $ \Delta {E}_{i=x,y,z} $ and $ \Delta {E}_{\mathrm{i}\mathrm{s}\mathrm{o}} $ on the magnitude of the applied electric field $ \left|F\right| $. The legends marked without ‘opt’ represent the FF method, while those with ‘opt’ indicate a combination of the FF method and geometry optimization. Fig. 2(a) reveals the energy variances, $ \Delta E $, of PE obtained using different methods in three directions are almost falling on the top of each other. This suggests a strong resistance of the geometry conformation to the electric field due to its relatively weak polarity. However, as shown in Figs. 2(b)−2(e), the $ \Delta E $ of PTFE and PVDFs obtained from the combination method is larger than those from than FF only, indicating that the conformations of PTFE and PVDFs are more sensitive to the external field. Particularly for β-PVDF, as depicted in Fig. 2(d), $ \Delta {E}_{x,y} $ exhibits a sudden increase when |F| exceeds 0.002 au. This phenomenon could be attributed to the spatially inhomogeneous distribution of the dipole moments. As visualized in Fig. 1, the fluorine atoms are predominantly located below the main axis and the diploe moments are not mutually cancelled. Fig. 2(f) presents the isotropic energy difference $ \Delta {E}_{\mathrm{i}\mathrm{s}\mathrm{o}} $ across all systems. the $ \Delta {E}_{\mathrm{i}\mathrm{s}\mathrm{o}} $ values obtained from the FF method and the $ \Delta {E}_{\mathrm{i}\mathrm{s}\mathrm{o}} $ values of PE and PTFE obtained from FFGO do not change significantly with |F| across all field intensities. However, for PVDFs, $ \Delta {E}_{\mathrm{i}\mathrm{s}\mathrm{o}} $ surges with an increase in |F|. Notably, $ \Delta {E}_{\mathrm{i}\mathrm{s}\mathrm{o}} $ of β-PVDF increases first and most significantly due to its pronounced inhomogeneous distribution of dipole moments, just as previously discussed. Fig. 3 visually presents comparisons of the β-PVDF structure under different magnitudes of the applied field. When |F| is 0.001 au, the structure shows only minor deviations in the two terminal groups from the static structure. However, when |F| is 0.003 au, the overall structure significantly deviates from the static structure. This suggests that this intensity of the applied field may challenge the fundamental assumptions of finite field theory.Fig 2The energy differences, $ \Delta E $, in various directions (panels a−e) and the isotropic energy difference, $ \Delta {E}_{\mathrm{i}\mathrm{s}\mathrm{o}} $, (panel f) as a function of the magnitude of the applied electric field $ \left|F\right| $. Note that $ \Delta E $, $ \Delta {E}_{\mathrm{i}\mathrm{s}\mathrm{o}} $, and $ \left|F\right| $ are all expressed in the atomic units (au).Fig 3Comparisons of the β-PVDF structure under different magnitudes of |F|. In panel (a), the blue and red lines represent |F| values of 0 and 0.001 au, respectively. In panel (b), the blue and green lines represent |F| values of 0 and 0.003 au, respectively. The direction of the field aligns with the positive direction of the x-axis. The figure is drawn using the PyMOL package.[35]Variations of Geometric Conformation with Applied Electric FieldsThe changes in geometric conformation are quantitatively tracked using the root mean square displacement (RMSD), denoted as D(F), from the zero-field conformation. The RMSD is given by the equation:10$ D\left(\boldsymbol{F}\right)=\sqrt{\frac{1}{N}\sum\nolimits _{i=1}^{N}\left[{\left({x}_{i}-{{x}_{i}}'\right)}^{2}+{\left({y}_{i}-{{y}_{i}}'\right)}^{2}+{\left({z}_{i}-{{z}_{i}}'\right)}^{2}\right]} $where N is the total number of atoms in the chain, the sets of coordinates $ \left\{{x}_{1},{x}_{2}\dots {x}_{N}\right\} $ and $ \left\{{{x}_{1}}',{{x}_{2}}'\dots {{x}_{N}}'\right\} $ correspond to the coordinates in the electrical field and zero-field reference, respectively. The average RMSD is computed across all directions of the applied fields and is given by:11$ \overline{D\left(F\right)}=\frac{1}{6}\left[D\left({F}_{x}\right)+D\left({-F}_{x}\right)+D\left({F}_{y}\right)+D\left(-{F}_{y}\right)+D\left({F}_{z}\right)+D\left({-F}_{z}\right)\right] $Fig. 4 shows the dependences of D(F) and $ \overline{D\left(F\right)} $ on the magnitude of the applied electric field $ \left|F\right| $. For PE and PTFE, as shown in Figs. 4(a) and 4(b), all D(F) and $ \left|F\right| $ exhibit a linear relationship due to their homogeneous distribution and mutual cancellation of dipole moments. However, for PVDFs (shown in Figs. 4c−4e), D(F) presents varied responses to the directions and magnitudes of the applied field. As stated above, this is due to the spatially inhomogeneous distribution of dipole moments. To adhere to the fundamental hypothesis of the finite field theory, the variances of D(F) in all directions should be maintained within a reasonable range. Generally, |F| should not exceed 0.002au to prevent a sudden increase of D(F). Fig. 4(f) summarizes the average RMSD for all systems. Compared to PVDFs, the changes of $ \overline{D\left(F\right)} $ for PE and PTFE are negligible. For α and γ-PVDF, $ \overline{D\left(F\right)} $ increases linearly with |F| across all ranges. However, for β-PVDF, a linear dependence of $ \overline{D\left(F\right)} $ on $ \left|F\right| $ is observed only for small fields, with a significant conformation change occurring for larger fields. It is worth noting that the behaviors of RMSD align with those of energy.Fig 4The root mean square displacement, D(F), in various directions (panels a−e) and the average root mean square displacement, $ \overline{D\left(F\right)} $, (panel f) as a function of the magnitude of the applied electric field $ \left|F\right| $.Variations of Polarizability with Applied Electric FieldsIn Fig. 5, we illustrate the dependencies of the diagonal components of the polarizability tensor, αxx, αyy, αzz, which correspond to the ability to respond to the electric field in specific directions, and the isotropic polarizability, αiso, on the magnitude of the applied electric field |F|. Note that all polarizabilities at $ \left|F\right|=0 $ are calculated using the perturbation method and in panel d of Fig. 5 the polarizabilities are truncated at 800 au due to their high magnitude. Generally, all the polarizabilities obtained from the combination of the FF method and geometry optimization are larger than those obtained from the FF method alone. Only a small external field, e.g., 0.0001 au, is needed to cause considerable increases in polarizabilities. For the FF method, the diagonal components of the polarizability tensor and isotropic polarizability in all systems remain almost unchanged with respect to the external field in all testing strength ranges of |F|. Their values are approximately equal to those calculated using the perturbation method in present work. The αiso of PE and PTFE obtained in this work are also consistent with those reported in previous work.[2] The similarity of polarizability obtained from FF and perturbation method indicates that for the calculation of polarizability, taking only the electronic contribution into consideration may be insufficient since the resulting dielectric constants (as stated below) are lower than the values from experiments. It is worth noting that the polarizabilities are larger in longitudinal than transverse directions.Fig 5The polarizability components αxx, αyy, αzz, and isotropic polarizability αiso as functions of the magnitude of the applied electric field $ \left|F\right| $. The legends in panels (b)−(e) are the same with those in panel (a).Regarding all the polarizabilities obtained from the combination of the FF method and geometry optimization, the values of PE and PTFE remain almost unchanged with |F|, and the polarizabilities are also larger in longitudinal than transverse directions, since the changes of conformation and energy for PE and PTFE are negligible as previously mentioned. However, for PVDFs as shown in Figs. 5(c)−5(e), the polarizabilities remain unchanged with the applied field only for |F|≤0.001 au. Additionally, the polarizability components become smaller in the longitudinal direction than transverse polarizabilities. This inversion may result from the spatially inhomogeneous distribution of dipole moments and less mutual cancellation of the dipole moments in transverse directions. Therefore, the conformation is more likely to respond to transverse field (i.e., ±Fy, ±Fz). Surging trends are observed for larger |F|, which are consistent with the variances of energy and RMSD. The intensity of the applied field is too large for finite field theory, and the characteristics in this range are beyond the scope of this work.To more intuitively illustrate the variations of polarizability, we present the unit sphere representation maps of the polarizability of β-PVDF at |F|=0.001 au obtained from two methods in Fig. 6. The bold red, green, and blue axes in the center represent the diagonal components of the polarizability tensor, αxx, αyy, and αzz, respectively. The non-diagonal components are not taken into consideration. The length of the axis corresponds to the magnitude of polarizability. The color of the arrows on the sphere surface varies from blue to white to red; the shorter the arrow, the bluer it is, and the longer the arrow, the redder it is. In the FF method, as shown in Fig. 6(a), the largest polarizability aligns with the x-axis (the chain axis). However, for the combination method (as shown in Fig. 6b), the largest polarizability shifts to the y-axis. This indicates that the optimized conformation in the field has more electrons in the y-direction, resulting in significant induced dipole moments.Fig 6Unit sphere representation of polarizability tensor with only the diagonal components of β-PVDF with applied field of |F|=0.001 au. All the polarizabilities are obtained from the FF method only in panel (a), and from the combination of the FF method and geometry optimization in panel (b). The calculations for this figure were performed using the Multiwfn,[36] and the figure was rendered using VMD package.[37]From Figs. 5 and 6, the anisotropy of polarizability for different systems is clearly visible. We use the following equation to quantitatively describe the anisotropy, denoted as Δα:[38]12$ {\Delta }\alpha =\sqrt{\left[{\left({\alpha }_{xx}-{\alpha }_{yy}\right)}^{2}+{\left({\alpha }_{xx}-{\alpha }_{zz}\right)}^{2}+{\left({\alpha }_{yy}-{\alpha }_{zz}\right)}^{2}\right]/2} $Fig. 7 illustrates the dependencies of Δα on |F|. In this section, we focus only on Δα with |F| ≤ 0.001au (as shown in the inset of Fig. 7) and neglect the irregular variations of PVDFs for |F|0.001 au. In this case, Δα essentially does not change with |F|. The values of Δα are generally larger for the combination method than the FF method. However, the magnitude of increase varies considerably. In descending order, they are PE, PTFE, α-PVDF, γ-PVDF, and β-PVDF. The change in Δα for PE is negligible because of its weak dipole in chain, which consists only of repeat units of CH2. On the other hand, due to the most inhomogeneous distribution of fluorine atoms in β-PVDF (as shown in Fig. 1), its Δα exhibits the largest increase by 520%.Fig 7(a) The anisotropy of polarizability Δα as a function of the magnitude of the applied electric field $ \left|F\right| $. Δα is capped at 500 au. Δα for |F|≤0.001 au is zoomed in (b).Dielectric ConstantThe dielectric constant εr is calculated according to Eq. (7) using the isotropic polarizability obtained from both the FF and combination methods. In Fig. 8, εr is presented as a function of the magnitude of the applied electric field, |F|. The value of εr at |F|=0 is calculated using the perturbation method. A detailed view of εr for |F|≤0.001 au is provided in the inset. Similar to the polarizability discussed above, all εr obtained from the FF method are in agreement with those from perturbation method and do not change with the applied field. On the other hand, εr values obtained from the combination method are generally larger. For PE, εr shows a slight increase from 2.11 to 2.13, which are all close to the reported values.[2,39,40] For PTFE, εr increases from 1.76 to 2.05. The later value aligns well with experimental values of 2.0−2.1.[40,41] For PVDFs, εr obtained from both FF and perturbation methods are about 1.87, which are quite lower than those obtained from experiments. It is worth noting that there is a large difference in εr between amorphous (8−13) and crystalline (2.2−3.5) PVDF reported by experiments and molecular dynamics simulations.[4,42] In the combination method, εr of α-, β- and, γ-PVDF are 3.00, 4.94 and 3.02, respectively. The εr of α- and γ-PVDF in present work are well consistent with these reported values. The larger εr of β-PVDF can be attributed to the large conformation deviations between the DFT optimized structure of a single chain and the crystalline structure in a multi-chain environment. As shown in Fig. 9, the conformational discrepancy is noticeably greater in β-PVDF compared to α- and γ-PVDF. The optimized structure of a single chain exhibits somewhat helical character. However, in a crystal, there are not rotations of CF2 and CH2 groups along the main axis, and its conformation shows perfect axial symmetry. Hence, its εr is in a transitional phase between crystalline and amorphous structures. The εr of PVDFs are lower than the experimental values of the amorphous phase. For amorphous phase, the big conformation change under the applied static field can cause the larger polarizability and dielectric constant. Hence, the calculated dielectric constant for β-PVDFs can be underestimated compared with the experiment value for the amorphous phase.Fig 8(a) The dielectric constant $ {\epsilon }_{\mathrm{r}} $ as a function of the magnitude of the applied electric field $ \left|F\right| $. $ {\epsilon }_{\mathrm{r}} $ at $ \left|F\right|=0 $ are calculated using perturbation method. $ {\epsilon }_{\mathrm{r}} $ for |F|≤0.001 au is zoomed in (b).Fig 9Comparisons of PVDFs structure in crystalline phase (represented by red lines) and after DFT optimization (represented by blue lines).CONCLUSIONSIn this work, we introduced a novel approach that combines density functional theory (DFT) with the finite field method to estimate the polarizability and dielectric constant of polymers. The proposed method effectively captures the impact of electronic and geometric conformation changes on the dielectric constant. To validate the method, polyethylene (PE) and polytetrafluoroethylene (PTFE) were selected as benchmark materials. The results demonstrated that the method accurately predicted the dielectric constants of PE and PTFE, indicating its reliability and robustness. When the magnitude of the external field |F|≤0.001 au, the dielectric constants of these polymers remain constant, which are very close to the experimental values, and are generally larger than those obtained from the perturbation method and finite field method. Furthermore, we investigated the impact of conformation variations in polyvinylidene fluoride (PVDF) on its dielectric constant and polarizability. The calculated dielectric constants of α- and γ-PVDF (3.0) fall within the reasonable range of crystalline (2.2−3.5) PVDF observed in experiments. The dielectric constants of β-PVDF (4.9) corresponds to a conformation between crystalline and amorphous structures due to its helical character. These findings shed light on the relationship between PVDF’s structural properties and its electrical behavior, providing valuable insights for material design and applications. Overall, our proposed method offers a promising avenue for estimating the dielectric properties of polymers by considering the complex interplay between electronic and geometric factors. It opens up possibilities for further exploration and understanding of the dielectric behavior of various polymer systems, contributing to advancements in materials science and engineering.