Direct measurement of ferroelectric polarization in a tunable semimetal

Ferroelectricity, the electrostatic counterpart to ferromagnetism, has long been thought to be incompatible with metallicity due to screening of electric dipoles and external electric fields by itinerant charges. Recent measurements, however, demonstrated signatures of ferroelectric switching in the electrical conductance of bilayers and trilayers of WTe$_2$, a semimetallic transition metal dichalcogenide with broken inversion symmetry. An especially promising aspect of this system is that the density of electrons and holes can be continuously tuned by an external gate voltage. This degree of freedom enables investigation of the interplay between ferroelectricity and free carriers, a previously unexplored regime. Here, we employ capacitive sensing in dual-gated mesoscopic devices of bilayer WTe$_2$ to directly measure the spontaneous polarization in the metallic state and quantify the effect of free carriers on the polarization in the conduction and valence bands, separately. We compare our results to a low-energy model for the electronic bands and identify the layer-polarized states that contribute to transport and polarization simultaneously. Bilayer WTe$_2$ is thus shown to be a canonical example of a ferroelectric metal and an ideal platform for exploring polar ordering, ferroelectric transitions, and applications in the presence of free carriers.

Polar materials exhibit charge separation in the absence of an applied electric field, an effect of broken inversion symmetry and a unique polar axis in the crystal. 2,3 In certain polar systems, the charge polarization can be switched by an external electric field, an effect known as ferroelectricity. In principle, the presence or absence of ferroelectric effects depends only on the crystal class and not the details of the electronic structure. Despite this, nearly all known conventional ferroelectrics are electrically insulating. Since the first theoretical proposals for ferroelectric metals in 1965, 2 only a handful of experimental claims of ferroelectric-like phases in metallic systems have been reported, [4][5][6][7][8] and no clear case for a canonical ferroelectric metal has emerged. Many such claims fail to demonstrate two key signatures of ferroelectric behavior, direct evidence of the polarization and ferroelectric switching, due to bulk screening effects.
Here, we focus on the polar, semimetallic van der Waals crystal, T d -WTe 2 , in the limit of two atomic layers, thin enough to admit an external electric field (Fig. 1de). Few-layer crystals of WTe 2 have drawn recent interest for exhibiting a wide variety of low-temperature phases. [9][10][11][12][13] Recent transport measurements showed that bilayer (2L) and trilayer (3L) WTe 2 exhibit intrinsic, switchable electrical polarization in the conducting state, 1 and separately, surfaces of bulk WTe 2 crystals display hysteresis in piezoresponse force microscopy. 14 Subsequent first-principles calculations indicated that the net polarization points only in the out-of-plane direction, and that the underlying mechanism results from a subtle interlayer sliding between the layers in two sta-ble configurations. 15,16 These findings are exciting given the semimetallic and tunable nature of bilayer WTe 2 , which enables reaching both electron and hole bands by electrostatic gating in the ferroelectric state. While the hysteretic behavior observed in previous experiments is promising, a direct measurement of the metallic polarization is still missing. Due to methodological limitations it was previously not possible to measure the polarization while varying a pure electric field. More importantly, these limitations also prevented observing the effect of free carriers on the polarization and its dependence on carrier density, a fundamental open question for ferroelectric metals.
In this work, we directly measure the charge polarization and electronic compressibility as a function of density for electrons and holes with independent control of the electric field. We study the simplest polar WTe 2 system, a bilayer, via capacitive sensing in a dual-gated, planar capacitance device (Fig. 1b). Capacitance measures the electronic compressibility (and thus metallicity) of a 2D system. In a bilayer 2D system the top-gate and bottom-gate capacitances provide a direct measurement of the layer-specific charge distribution, 17 and thus the out-of-plane polarization. Furthermore, the parallelplate geometry enables this charge sensing with simultaneous and independent control of the vertical electric field and the carrier density in the bilayer by electrostatic gating.
Our devices each consist of a bilayer WTe 2 crystal encapsulated by two hexagonal boron nitride (hBN) dielectric layers, with metallic top and bottom gates, and con- tacts integrated into the top hBN layer 18 (Fig 1a). We measure the capacitance between the top gate and the bilayer, C t , while applying DC voltages to the top and bottom gates to tune the (nominal) total carrier density, n 0 ∝ C 0 (Fig. 1b). While the geometric contributions to the capacitance (C 0 t and C 0 b ) are constant, there are additional contributions to the measured capacitance C t from the electronic compressibility. In a bilayer system with partial electric field penetration, the layer-specific densities n i can differ between the two lay-ers (i = 1, 2) for a given total density n = n 1 + n 2 , particularly in the presence of an electric field. The electrostatic potentials of each layer φ i depend on the top and bottom gate voltages, and crucially, any built-in electric field in the bilayer. 17,19 As such, φ i in each layer can also differ, even in thermodynamic equilibrium. The electronic compressibility of a 2D bilayer is generally described by a 2 × 2 matrix, ν ij = −∂n i /∂φ j , however, in a weakly-coupled bilayer it is possible to characterize the system with only the diagonal elements, ν ii , the layer-specific compressibilities. Due to the van der Waals nature of the interlayer coupling, in bilayer WTe 2 this is indeed the case (Ext. Data Fig. 7). Subsequently, the top capacitance may be written, (see Methods for details). The second term in Eq. 1 is a quantum correction to the the geometric capacitance C 0 t , inversely proportional to the layer-specific compressibility of the top layer, ν 11 . The layer-specific densities n 1 and n 2 can be obtained by integrating the capacitance, and thus the polarization, proportional to n 1 − n 2 , can be measured.
For fixed, external field E ⊥ = 0, the measured C t as a function of the electron density exhibits a minimum near charge neutrality, n 0 ≈ 0 (Fig. 1c). The minimum in C t indicates the presence of a small band gap (incompressible state) in the WTe 2 bilayer, a feature consistent with previous observations of a sharp drop in the conductance of bilayer WTe 2 in transport. 1,20 At fixed, negative electric field the capacitance minimum becomes more prominent, suggesting that the gap is electric-field-tunable (see Ext. Data Fig. 4), similar to bilayer graphene. 19,21 However, in contrast to bilayer graphene the electric field response is robustly asymmetric around E ⊥ = 0, as shown in the full n 0 and E ⊥ dependence in Fig. 2a. This effect results from the absence of an inversion center or mirror plane between the layers of bilayer WTe 2 , implying the crystal is polar along the c-axis (Fig. 1e).
To probe the switching behavior of the polar direction, we sweep the electric field back and forth at fixed density, as shown in Fig. 2a-b. Sudden changes are observed in C t at all densities: C t jumps at a positive critical electric field value E + c when sweeping toward positive E ⊥ (Fig. 2a), whereas the critical field E − c is negative when sweeping in the negative direction (Fig. 2b), forming a hysteresis loop (Fig. 2c) at each density. Taking the difference of two representative sweeps in opposite direc- (Fig. 2d), we see that C t overlaps nearly everywhere excluding the hysteretic region between the critical fields, E ± c , where C t is multi-valued. The critical fields generally fall within |E ⊥ | 0.1 V/nm and appear to be weakly dependent on charge density (Fig. 2e). Interestingly, the sign of the hysteretic difference switches for holes (∆C t < 0) compared to electrons (∆C t > 0), and the magnitude of the difference decreases at large densities of either sign. Switching behavior is clearly present in the capacitance, but to understand how this relates to the polarization we must first recognize that the ground state structure of bilayer WTe 2 possesses two stable configurations with opposite polarization 15,16 (Fig. 1d). Similarly, from Fig. 2a-c we see that there are two stable values for C t in the central region between the two critical fields, E ± c . To reveal how the C t bistability relates to the two polarization states, we employ a self-consistent calculation of the capacitance using a k · p model for the bilayer WTe 2 bands 22 together with an exact form of Eq. 1 (see Eq. 11b in Methods). The k · p Hamiltonian (Eq. 4) provides the low-energy electronic structure, shown schematically in Fig. 3a, and layer-projected wavefunctions in each polarization state. The layer-specific densities n i (Fig. 3b) and electric potentials φ i are obtained from the bands by performing a self-consistent electrostatics calculation for the parallel-plate system (see Methods). The compressibilities ( Fig. 3c) are obtained by numerical differentiation and inserted into the general form of Eq. 1 to compute the capacitance (given by Eq. 11b). Since there are two stable configurations for the bilayer, these quantities all depend on the polarization state p = ±, yielding different values n p i , φ p i , and ν p ii in each case. Using these calculated quantities, the polarization at fixed n 0 can be determined from the difference P ± = ed i (n ± 1 −n ± 2 ) for interlayer separation d i . Evaluating this directly, we can construct a schematic polarization loop that illustrates both the dielectric polarization response as well as the spontaneous polarization from the computed layer densities (shown for n 0 ≈ 0 in Fig. 3d, with switching behavior added manually for illustration). The difference between the layer imbalance in the two polarization states yields the change in the spontaneous polarization, ∆P s ≡ P + −P − , as indicated for E ⊥ = 0 in Fig. 3d.
Whereas polarization switching reflects a reversal of the density imbalance between the layers ( Fig. 3a-b), capacitance switching reflects an inequivalence of the compressibility of the nearest layer in the two polarization states (Fig. 3c). For instance, in Fig. 3e, we show calculated curves for C t as a function of E ⊥ at a few selected electron and hole densities using the calculated compressibilites. The result of the calculation is two capacitance curves, C + t and C − t , that correspond to the polarization states P + and P − , respectively. Comparing to the experiment (Fig. 3f ), for sufficiently large positive (negative) E ⊥ , we only measure the C + t (C − t ) branch. In the hysteretic region, the polarization state P + or P − depends on the direction of the electric field sweep, and thus the measured capacitance jumps between the C + t and C − t branches at the coercive fields, E ± c . While our model does not include all the details of the electronic structure (Ext. Data Fig. 6), and thus we do not expect a perfect match with experiment, the general trends of the capacitance with E ⊥ and n 0 are captured nicely, for instance, the observation that C + t > C − t for electrons and C + t < C − t for holes. Next we will see that some of the details of the electronic structure are less relevant for difference quantities such as ∆C t and ∆P s , which leads to improved agreement between the model and experiment, as shown in Fig. 4.
Finally, we address the density dependence of the polarization to ascertain the effect of adding free charge to a polarized semimetal. Previously, we defined ∆C t as the difference between the right-sweeping and left-sweeping C t curves. Having identified the C ± t curves in the central region of each hysteresis loop, we may now equate Figure 4a shows the density dependence of ∆C t for E ⊥ = 0, equivalent to a vertical line cut from Fig. 2e. The non-monotonic curve follows the size of the C t hysteresis loop, showing a maximum and minimum on either side of charge neutrality and trending toward zero at large densities. Due to the switching process acting as a mirror operation on the WTe 2 bilayer, the difference , to the change in spontaneous polarization (see Methods). This identity enables extraction of the change in spontaneous polarization In the P − state, the colors and layer polarization would be interchanged. Each pair of valleys is centered around a point along Γ-X, labeled Q. Band parameters are exaggerated to emphasize separation of valleys. b Representative calculated layer densities, c compressibilities, and d polarization ∝ n p 1 − n p 2 in each polarization state, p = ±, from n0 ≈ 0. e Computed top capacitances, C + t and C − t in the P + and P − states, respectively, for the listed densities versus electric field, with colored symbols indicating the hysteretic path observed in experiment, while gray symbols denote portions of each capacitance branch that are inaccessible in experiment due to switching behavior (switching fields are not computed in the model; solid vertical lines are shown extending from experimental switching fields in panel f to reflect the loop observed in experiment). Capacitance is calculated using computed layer densities, potentials, and compressibilities to evaluate Eq. 11b. f Measured top capacitance hysteresis loops for matching electron and hole densities in e.
The result of this integration along E ⊥ = 0 is shown in Fig. 4b. The spontaneous polarization is strongest at charge neutrality, where there are few free carriers available to screen the built-in polarization. Adding electrons or holes to the system reveals a pronounced asymmetry in the spontaneous polarization. Holes appear to screen the polarization much more effectively than electrons for equivalent charge densities. This final effect can be understood by considering the available states in the conduction and valence bands. First-principles calculations indicate that the low-energy states are intrinsically layer-polarized (see Ext. Data Fig. 6), so despite the fact the the additional charges are free and may conduct current within each layer, they may not simply transfer freely between the layers to screen out the built-in polarization in the absence of an external field. The wavefunctions of these low-energy states posses strong layer character and thus contribute, in part, to the net polarization rather than fully screening it, contrary to a simple electrostatic picture. The degree of layer polarization and electronic compressibility differs for electrons and holes in bilayer WTe 2 (as observed in the e-h asymmetry of C t in Fig. 1c). Consequently, filling valence band states appears to suppress the net spontaneous polarization whereas filling the conduction band allows polarization to persist to relatively large densities (Fig. 4b).
Through this mechanism, it is possible to obtain a ferroelectric state that coexists with native metallicity.
In bilayer WTe 2 , the strong layer character of the wavefunctions and asymmetric density of states manifests as several distinct ferroelectric regimes with both bound and free charges contributing to the net polarization. At large hole densities the measured polarization quickly trends toward zero, entering a screened polar metal regime. In this phase, there remains underlying polarized bound charge made up of remote states, however, there are enough filled states with opposing layer character to screen this charge, resulting in a suppressed net polarization. On the other hand, at large electron densities the polarization decreases more slowly, exhibiting a persistent, metallic ferroelectric state. Near charge neutrality we observe a crossover between these p-type and n-type metallic ferroelectric states. The crossover regime is further complicated by the electric field dependence of the conduction and valence bands, as shown in Fig. 2ab. From E ⊥ = 0, the charge neutral behavior appears Fig. 4. Ferroelectric polarization in the presence of free carriers. a Measured density dependence of ∆Ct ≡ C + t − C − t at E ⊥ = 0. b Change in measured spontaneous polarization calculated by integration, ∆Ps ∝ ∆Ct dn0. Left axis provides 2D polarization units while the right scale is given in terms of charge separation between the layers, ∆ps = ∆Ps/di for WTe2 interlayer separation, di = 0.7 nm. Background shading follows the magnitude of ∆Ps, illustrating distinct regions of ferroelectric (FE) behavior. c Computed ∆Ct and d ∆Ps based on model Hamiltonian for bilayer WTe2, with inset in d showing computed ∆Ps for an extended density range (with identical units). ∆Ps is computed by two different methods: "Exact" using Eq. 18, and "Approx." using Eq. 33.
to transition from a ferroelectric semimetal with finite band overlap to a ferroelectric insulator with a maximal band gap of ∼ 3 meV from 0.3 V/nm > |E ⊥ | > 0.2 V/nm (Ext. Data Fig. 4).
These interpretations are supported by our model calculations, with both the density dependence of calculated ∆C t and the integrated ∆P s matching well with experiment (Fig 4c-d). The model also allows us to investigate the expected behavior over an extended density range, inaccessible in experiment due to dielectric breakdown of hBN. As shown in the inset of Fig. 4d, the calculated ∆P s continues to exhibit sustained ferroelectric behavior up to electron densities of at least n 0 ≈ 4 × 10 13 cm −2 , whereas the p-type ferroelectric behavior is substantially suppressed below n 0 ≈ −2 × 10 13 cm −2 .
In conclusion, we see that bilayer WTe 2 may be tuned from a sustained ferroelectric n-type metal for n 0 > 10 13 cm −2 to a variable-polarization ferroelectric p-type metal at intermediate hole densities, followed by a screened polar metal at large hole densities. In the crossover regime between the n-and p-type ferroelectric metal states, a ferroelectric insulator or semimetallic phase can be obtained depending on the electric field. The out-of-plane polarization is largest in the neighbor-hood of these latter phases, reaching a charge separation of ∼ 3.7 × 10 11 electrons/cm 2 over 0.7 nm between the layers, equivalent to a volume polarization density of ∼ 0.6 mC/m 2 . The ferroelectric behavior in each metallic phase is supported by the strong layer-polarized character of the states, a result of broken inversion and mirror symmetries combined with weak interlayer coupling. Together, these observations and the ingredients in WTe 2 that lead to them provide a recipe for engineering and measuring new ferroelectric metal systems. Looking forward, we anticipate metallic ferroelectricity will manifest in additional transition metal dichalcogenide and van der Waals structures, 23 enabling a host of new experiments in this previously unexplored regime.

A. Fabrication of capacitance devices
WTe 2 is an especially air-sensitive material, quickly degrading in ambient conditions. To avoid air exposure during the fabrication process, we first integrate metal contacts into the top BN layer and then transfer this template onto the WTe 2 in a N 2 -filled glovebox. The metal contacts are prepared ahead of time by etching holes through the BN crystal and subsequently evaporating pure Au contacts to fill these holes. We pick up the top BN along with the integrated Au contacts and use this layer to pick up the WTe 2 crystal using standard dry transfer techniques. 18 We then pick up a bottom layer of BN to fully encapsulate the WTe 2 and place the stack on a graphite or PdAu bottom gate, with the WTe 2 partially overlapping the gate but with the through-hole contacts positioned as near as possible to the gated region ( Fig. 1a). At this stage the WTe 2 is sealed from the environment on all sides, allowing us to pattern leads and a top gate using standard electron-beam lithography and metal evaporation. The top gate is designed to overlap the WTe 2 crystal up to the edge of the bottom gate, aligning carefully to this edge. This ensures that the entire gated region of the crystal is dual-gated, preventing single-gate features in our capacitance measurements.

B. Capacitance measurement
Capacitance measurements were performed by connecting the bilayer WTe 2 in each device to a capacitance bridge circuit (Ext. Data Fig. 1) with a standard, known capacitor, C std . We apply an 11 kHz ac excitation δV t with an rms amplitude from 64 mV to 140 mV on the top gate of the device and a nearly of-out-phase signal δV std to C std in order to null the total ac signal at the bridge balance point. Both the amplitude and relative phase of δV std are tuned to produce a null signal. As the total capacitance of the device changes, deviations from the null voltage are amplified by a high-electron-mobility transistor (HEMT) mounted within a few millimeters of the sample in the cryostat. We determine the measured top capacitance C t from the amplified deviation from the null voltage, C t = C std (∆V null /δV std ). All measurements were performed in a dilution refrigerator with the sample between 100 mK (as in Fig. 1, Fig. 2a-d, and Fig. 3) to 30 K (Fig. 2e), though little change was observed in the capacitance features in this range.

C. k · p model for bilayer WTe2
To compute the layer compressibilities for bilayer WTe 2 , we adopt a k · p model for the Hamiltonian describing massive Dirac fermions in two layers with spinorbit and interlayer coupling. Beginning with a tilted massive Dirac Hamiltonian, 22 we introduce a layer index, i = 1, 2, and electrostatic potential on each layer, φ i , to include the effect of a vertical electric field as well as the polarization, H i,s = φ i +tk x,i +v(k y σ x +ηk x,i σ y )+(m/2−αk 2 )σ z , (3) along with a shifted k x -coordinate,k x,i ≡ k x + q i that accounts for the position of each layer-polarized valley (from the Q-point shown in Fig. 3c, for Q = Qk x ). The Pauli matrices, σ x,y,z , act in the orbital pseudospin space, s =↑, ↓ is the spin degree of freedom, η = ±1 is a chiral index, t tilts the Dirac cones along k x , and m gives rise to the gap at charge neutrality (αk 2 ensures convergence ofĤ i as k → ∞). For E ⊥ = 0, the electrostatic potentials take fixed values φ i (E ⊥ = 0) = φ 0 i (see Table 2) to account for the spontaneous polarization. The spin and layer degrees of freedom are coupled in the |1, ↑ , |1, ↓ , |2, ↑ , |2, ↓ basis to obtain the effective Hamiltonian, by interlayer coupling, γ, and spin-orbit coupling, where λ x , λ y govern the spin-orbit coupling strength in the x and y directions. Such a model has previously been applied to describe the Berry curvature dipole observed in bilayer WTe 2 in vertical electric fields 13 and captures the spin and shifted-valley character of WTe 2 bands determined by first-principles methods. 22,24 Here, we explicitly include a polarization state label in the Hamiltonian,Ĥ p , p = ±, to denote the two possible configurations for the WTe 2 bilayer. The Hamiltonian for the opposite polarization state,Ĥ − , is obtained by interchanging layers, 1 ↔ 2, and spin. Note that both the electrostatic potential, φ i , and shifted valley coordinate, k x,i , depend on the layer index and thus are interchanged upon a simultaneous mirror operation (M c ) and inversion of the electric field. As a result, the Hamiltonian possesses an identity, relevant for extracting the polarization from the measured capacitance, to be discussed in a subsequent section.

D. Capacitance calculation
To understand the various features of the measured capacitance, we apply this effective model to calculate the band structure and eigenstates and numerically evaluate the capacitance for a similar geometry. Starting with Eq. 4, we determine the probability densities in each layer and sum over occupied states to obtain the layer-specific densities, n 1 and n 2 for the top and bottom layers, respectively. The layer-specific compressibility elements ν ij are determined by taking partial derivatives of the densities n i with respect to the electric potentials on each layer. The resulting compressibilities in conjunction with the geometric capacitances then determine the calculated capacitances.
The eigenstates of Eq. 4 span an 8 × 8 Hilbert space of layer, spin, and orbitals. We obtain layer densities by first calculating the probability density on each layer from the eigenstates, for spin s =↑, ↓ and orbital pseudospin m. Note that here, and for the remainder of this section we have dropped explicit labels making reference to the polarization state for simplicity. Layer densities are then determined by integration, where f (p) is the Fermi-Dirac distribution. However, Eq. 8 only includes contributions from valence electrons. The contribution from ions can be determined as follows.
With the absence of an external electric field, E ⊥ = 0 and φ i = φ 0 i , and the bilayer is at charge neutrality. We then have a vanishing total charge, n 1 + n 2 + n ion 1 + n ion 2 = 0. On the other hand, the two layers have same ionic composition. As a result, From this definition, we obtain the relevant layer densities, To perform a self-consistent calculation, we treat the electrostatic potentials, φ i , as parameters and compute the layer densities, n i = n i (φ 1 , φ 2 ) for an appropriate grid of values. The compressibilities are then obtained by taking the partial derivatives ν ij = ∂n i /∂µ j = −∂n i /∂φ j over the same parameter space (where φ i = −µ i for layerspecific chemical potential µ i and with the sample connected to ground). Finally, the capacitances are defined in terms of small signal variations, C t ≡ (δn 1 + δn 2 )/δV t , C b ≡ (δn 1 + δn 2 )/δV b , and C p ≡ δn t /δV b , where C p is the penetration field capacitance, measured from the top gate to the bottom gate, for charge on the top gate n t . These expressions may be evaluated in terms of the com-pressibility elements (following Ref. 17), , (11a) , (11b) with geometric capacitances C 0 t , C 0 b , and interlayer capacitance C i , as shown in Fig. 1b, each given by C j = 0 A/d j using = 3.2 for hBN and = 1 within the bilayer. The compressibility matrix,ν, is a symmetric matrix, ν ij = ν ji , and takes the form whileĈ is a matrix of geometric capacitances, To obtain results in terms of E ⊥ and n 0 , we first apply a transformation derived from the charge balance equations, to determine V b and V t , and then calculate E ⊥ = (V t /d t − V b /d b )/2 for top and bottom hBN thicknesses, d t and d b , respectively, and n 0 = (C 0 t V t + C 0 b V b )/e. Finally, to include the effect of the two polarization states, as shown in the main text for C ± t , we compute the capacitances C p t separately for each case, p = ±, using the compressibilities, ν p ij , calculated from the eigenstates ofĤ p , derived from the appropriate form of Eq. 4.

E. Simplifications and approximations
In practice, quantum capacitances of 2D electron gases tend to dominate over geometric capacitances from external electrodes and gates due to the small energy spacing between levels. For instance, the density of states (and thus the non-interacting electronic compressibility) in monolayer graphene is of order ρ MLG (0) ∼ 10 13 eV −1 cm −2 , yielding a quantum capacitance of order C q ∼ 20 fF/µm 2 . In comparison, a metallic gate separated by 10 nm of hBN from a parallel plate generates a geometric capacitance one order of magnitude smaller, C g ∼ 2 fF/µm 2 . Typical thicknesses for hBN gate dielectrics are even larger, such as the device shown in the main text with d t ≈ 15 nm and d b ≈ 19 nm, resulting in even smaller geometric capacitances. Additionally, in each layer of bilayer WTe 2 , screening from carriers in the nearby layer and weak interlayer coupling significantly reduces the off-diagonal compressibilities relative to the diagonal terms in Eq. 12. Together, these considerations imply a hierarchy, |ν 11 |, |ν 22 | C 0 t , C 0 b , C i |ν 12 |, |ν 21 | (ignoring factors of e 2 preceding the ν ij terms). As a result, the diagonal compressibilities ν ii dominate the behavior of C t and C b , enabling an approximate form for Eqs. 11, as shown in Eq. 1 in the main text. While the numerical computations shown in the paper make use of the full Eq. 11b, deviations of the approximate formulas from the exact quantities are negligible and thus the simpler form is given in main text Eq. 1 to facilitate understanding.

F. Capacitance relation to polarization
To extract the spontaneous polarization from the measured capacitance, we seek the change in the polarization, obtained by integrating the difference of polar response functions, defined up to a constant of integration, ∆P 0 , which we take to be zero at large hole densities, where the density dependent spontaneous polarization appears to saturate. The polar response in each state, ±, is defined by Eq. 18 together with Eq. 19 is employed to compute the "Exact" spontaneous polarization in the model calculations shown in Fig. 4d. In the experiment, we do not have direct access to the partial derivatives in Eq. 19, however the symmetry of our device allows an approximate form of Eq. 18 to be related to ∆P s . Following Ref. 19, we express the polar response in Eq. 19 as, where we have dropped the label indexing the polarization state. We then employ to obtain the partials ∂n i /∂V t,b , and finally rewrite the polar response in terms of the average geometric capacitance, C = (C 0 b + C 0 t )/2, and the asymmetry, (22) For the device shown in the main text (Device A in Table 1), δ ≈ 0.11 and C /2C i ≈ 0.067. Keeping only the leading order terms in these two small parameters, we obtain Therefore, we have To evaluate this further, we first rearrange Eq. 1 and explicitly show the E ⊥ dependence (ignoring factors of e 2 ), which implies .
(26) Similarly, we have Using the fact that switching the polarization state is equivalent to a mirror operation between the layers, we find a symmetry that is evident in the identity relating the Hamiltonian in the two polarization states, Eq. 6, and the subsequent compressibilities derived from each. Physically, this implies that measuring the capacitance from one side of the device in the two polarization states is related to a measurement of the capacitance from both sides of the device by geometrical factors. Thus, we invoke Eq. 28 to relate Eq. 26 and Eq. 27, or equivalently for hBN dielectrics on both sides of the device. Finally, we insert this identity into Eq. 24, solely in terms of measured ∆C t (E ⊥ ). Integrating this expression from the largest hole density measured, n l , up to n 0 , we obtain Evaluating at E ⊥ = 0, we arrive at the approximate integral used to determine the measured polarization in Fig. 4b, n l ∆C t (0, n 0 ) dn 0 .
(33) The same expression is employed in Fig. 4d to calculate the "Approx." curve based on computed capacitance data, illustrating the small deviation from the exact polarization introduced by making the approximations outlined in this section. A small constant of integration is included in the total measured spontaneous polarization shown in Fig. 4b, P 0 = 0.42 e/cm, inferred by matching the magnitude of the measured ∆P s (n 0 ) with the calculated curve. The latter curve exhibits saturation of ∆P s at large hole densities (taken to be zero) and thus offers a lower bound on the constant of integration and ultimately the maximum value of the polarization near charge neutrality.

G. First-principles calculations
As the essential ingredient to interpret the capacitance in experiments, the low-energy model in Eq. 4 is well supported by first-principles calculations. The most essential element for the capacitance is the polarization or the occupation difference between the top and bottom layer. In Ext. Data Fig. 6, we show such layer character in the band structure near one group of valleys from the first-principles calculations compared with those from the low-energy model. By tracing the valence band, we find that bands from the low-energy model and firstprinciples calculation show similar tilting, together with similar layer-occupancy variation. This striking similarity demonstrates that the low-energy model is indeed capable of capturing the essential physics for the capacitance calculation.
The supporting density functional theory calculations shown in Ext. Data Fig. 6b are performed with the Vienna Ab initio Simulation Package (VASP) using the PBEsol functional. The plane wave basis cutoff is 300 eV. The atomic structures are relaxed until forces on every atom are smaller than 0.01 eV/ A and the vacuum layer is larger than 15 A. The van der Waals (vdW) functional used here is the zero damping DFT-D3 method with small modifications: the vdW correction is only applied to atomic pairs from different layers, since monolayers can be calculated well without vdW correction. This modification ensures that the structure of the sub-layer of the bilayer will not be affected by an inappropriate vdW correction, leading to the same lattice constant for the monolayer and bilayer system.       Table 2. Parameters employed in k · p calculation, following Eqs. 3-5. These parameters are obtained by fitting the polarization from the low-energy model to Fig. 4a, and hence are different from those in Ref. 22.