Releasing the hidden shift current in the TTF-CA organic molecular solid via symmetry lowering

Bulk photovoltaic effect, characterized by an excitation-driven unbiased spontaneous photocurrent, has attracted substantial attention mainly due to its potential for harvesting solar energy. Here, we investigate the photovoltaic characteristics of organic molecular solids and focus on the association between the photocurrent and the crystal symmetry in the exemplary case of tetrathiafulvalene-p-chloranil. We perform comprehensive first-principles calculations, including direct evaluations of the excited-state current via real-time propagations of the time-dependent density functional theory. We find that the charge shifting in the low-temperature phase is mainly driven by the intrachain ferroelectricity, which gives rise to a photocurrent not only in the visible-light range but also near the band-edge infrared region. The shift current that is locked in the symmetry of the high-temperature phase can be released by introducing a potential asymmetry. We suggest that organic molecular solids can be exploited via appropriate engineering to lower the symmetry, aiming at room-temperature photovoltaics.


INTRODUCTION
The bulk photovoltaic effect (BPVE), which refers to the spontaneous generation of photocurrent in a non-centrosymmetric material, has attracted substantial interest for harvesting solar energy without electric bias. [1][2][3][4][5] As the main ingredient of the BPVE is rooted in a purely quantum-mechanical phenomenon, namely, the shift current mechanism, which means the spatial charge shift upon excitation from valence to conduction bands, it is not necessarily bounded by the main limitation of conventional photovoltaics, namely, the Shockley-Queisser limit (33.7%). 6 A theoretical formulation of the shift current has been established through the second-order optical response by using perturbation theory. 7 The size of the charge shift can be determined exactly from the phase evolution of the adiabatic wavefunctions, namely, the difference in the Berry connections of the occupied and unoccupied bands. 8 Thus, a precise analysis of the electronic structure is necessary for microscopic investigation of the BPVE and for design of materials for high-performance solar energy devices.
To date, the shift current and related phenomena have mainly been studied for perovskite oxide materials, such as BiFeO 3 and PbTiO 3 . 9,10 Many questions and substantial interest have now also been directed toward the usefulness of organic molecular compounds because of the small and tunable size of the bandgap, which has been desirable for harvesting solar energy in a wider range of light frequencies along with low cost and easy fabrication. 4,11 Recently, it was demonstrated that tetrathiafulvalene-p-chloranil (TTF-CA) can be utilized as an organic shift current material. 4 In the low-temperature phase (<81 K), the alternating stack of TTF and CA molecules is known to provide a large electric polarization along the stacking direction, 12 and the observed amount of the shift current provoked substantial interest toward photovoltaic systems that were derived from pure organic materials. [13][14][15] However, the limited range of the operating temperature in the previous study has casted severe doubt on whether photovoltaic systems that are based on organic materials can be exploited at elevated temperatures.
In this paper, we investigate the electronic structure and photovoltaic characteristics of the TTF-CA molecular solid. Based on the first-principles results obtained from this study, we propose methods for manipulating the electronic structures and realizing the shift current available at temperatures higher than those of ferroelectric transitions. We perform standard static density functional theory (DFT) calculations 16,17 to investigate the ground-state electronic structures and the optical response functions. Excited electronic structures are described through the principles of time-dependent density functional theory (TDDFT) 18 : to calculate the real-time profile of the excited-state current driven by the external field, we perform real-time propagations (rtp) of Kohn-Sham states (abbreviated hereafter as rtp-TDDFT), [19][20][21] and the current densities are deduced from the expectation value of the electron's velocity operator. We find that in the low-temperature (or ionic) phase, the observed shift current is mainly attributable to the intrinsic ferroelectricity in combination with the intrachain charge shift between TTF and CA molecules. A large shift current can be realized in the lowtemperature ferroelectric phase not only in the visible-light range but also in the infrared range of band-edge excitation. In the hightemperature (or neutral) phase, we find that a vertical excitation produces the interchain charge shift, of which summation over the whole Brillouin zone results in perfect cancelation, thereby leading to a vanishing shift current. We propose that the photovoltaics of the high-temperature phase can be facilitated by lowering the symmetry, namely, by introducing an asymmetric potential in the perpendicular direction.

RESULTS AND DISCUSSION
Atomic structure and electric polarization of the neutral and ionic phase TTF-CA First, we investigate the electric polarization in two phases of the TTF-CA molecular solid. 22 As depicted in Fig. 1a, b, the bulk TTF-CA molecular solid consists of 1D chains of alternating TTF (orange) and CA (blue) molecules. The neutral phase (λ = 0.0 in Fig. 1) is centrosymmetric and stabilized above T c (≅81 K), whereas the lowtemperature ionic phase (λ = 1.0 in Fig. 1) has a finite polarization in the chain direction (the x axis in Fig. 1). 23 The ferroelectric transition is mainly driven by the dimerization between TTF and CA in the chain, as depicted in Fig. 1b, which is accompanied by marginal intramolecular distortions. Detailed analysis of the geometric changes on the transition are reported in the previous literature. 24 To explicitly identify the correlation between the ferroelectricity and the structural symmetry, we linearly interpolated the structures via R λ =R n +λ(R i -R n ), where R n and R i indicate the atomic structures of the neutral and the ionic phases, respectively. Detailed coordinates for each phase are summarized in Supplementary Table 1 and Supplementary Table 2. We calculated the variation of the electric polarization by evaluating the Berry phase 25,26 in TTF-CA as a function of the phase parameter (λ), as presented in Fig. 1c, d. The obtained maximum polarization of the ionic phase (7.9 μC/cm 2 ) is comparable to the value that was measured in the previous experiment (7-9 μC/ cm 2 ). 12 Over the whole range of the structural variation (0 ≤ λ ≤ 1), the polarization is pointing exclusively in the chain direction (the x axis in Fig. 1c), and the absence of vertical polarizations (the y and z axes) is enforced by the persistent symmetry along the along the y and z axes.
We decomposed the polarization into electronic and ionic contributions, as presented in Fig. 1d, which reveal that the electronic part prevails over the compensating ionic contribution. Using the DFT electronic structures, we performed the Badercharge analysis, as presented in Supplementary Fig. 1, and found that there is no perceivable charge redistribution over the window of phase variation (0 ≤ λ ≤ 1), in contrast to previous expectations. 27 Furthermore, the electronic contribution of the polarization is not dominated by a state near the Fermi level, but is uniformly distributed over the occupied bands ( Supplementary  Fig. 2). We conjecture that this nontrivial band structure underlies the origin of the large ferroelectricity of TTF-CA, and thus, understanding the detailed correlation between the intrinsic electronic structure and its photovoltaic effect is necessary.
Electronic band structure and its optical responses of the neutral and ionic phase TTF-CA To characterize the effect of this intrinsic ferroelectricity, we calculated and analyzed the electronic structure, in the context of the photovoltaic properties, of the neutral and ionic phases in TTF-CA. As shown in Fig. 2a, b, the overall patterns of the electronic band structures appear similar; however, the bandgap of the ionic phase is substantially larger than that of the neutral phase, which is consistent with the previous theoretical result. 28 The absorption peaks near 0.46 and 0.54 eV for the neutral and ionic phases, as shown in Fig. 2c, d, respectively, correspond to the band-edge excitations of both cases. However, in the shift current spectra, as shown in Fig. 2e, f, only the ionic phase exhibits the peak structure near the 0.54 eV region, while that of the neutral phase is negligible. In the later paragraph, we explain how the photovoltaic effect is closely related to these intrinsic electronic structures.
The obtained shift current of the ionic phases is characterized by two ranges around the infrared region (0.54 eV) and the visible region (2.69 eV), as presented in Fig. 2f. A detailed analysis of the visible region is summarized in Supplementary Fig. 3. The previous experimental measurements, which were obtained by using a visible-light laser, correspond to the visible region in Fig. 2f and overlook the shift current near the infrared range. 4 However, in our first-principles computations, the excited-state current exhibits a peak in the shift current profile near 0.54 eV, which corresponds to the band-edge excitation. Although the bandgap is formed indirectly, the bottom of the conduction band is almost in parallel with the top of the valence band, as indicated by the red and blue lines in Fig. 2b, thereby providing a large electron-hole pair density of states at the band-edge excitation.
Intrachain effect of shift current in the ionic phase TTF-CA To investigate the underlying mechanism of the large shift current in the ionic phase, we separate the intrachain effect from the interchain effect. We isolated the single 1D chain of TTF-CA and performed the same calculations of the electronic structure and excited-state current. The unit-cell geometry of this 1D lattice is defined in the same way as above: each unit-cell geometry is gradually varied from λ = 0.0 (the neutral phase) to λ = 1.0 (the ionic phase). Similar to the case of the bulk molecular solid (Fig. 2a,  b), the bandgap increases with λ in the 1D chain, while the overall patterns of the electronic band structure are preserved (Fig. 3a). As shown in Fig. 3b, the highest valence bands and the lowest conduction bands mainly consist of TTF (orange color) and CA atomic orbitals (blue color), respectively. Hence, the band-edge excitation, from the valence band to the conduction band, should result in a spatial shift of the charge density from TTF to CA molecules. The intrachain shift currents for specified 1D lattices, which are defined by each value of λ, are calculated with various a ¼ x; y; z axis) for c the neutral and d ionic phases of TTF-CA. The calculated shift current spectra that correspond to the shift current density J a that is generated by a second-order nonlinear field E b E b (σ abb , a and b = x, y, z axis) for e the neutral and f the ionic phases of TTF-CA with respect to the frequency of the applied light. The peak that is indicated by the red upward arrow in (f ) originates in the band-edge transitions from the valence bands (red lines) to the conduction bands (blue lines) in (b). Fig. 3 Intrachain effect of TTF-CA. a The electronic band structure of the isolated single chains of TTF-CA whose unit-cell geometry is specified by the scaling parameter that was defined as in Fig. 1. b The partial charge densities of the conduction and valence bands in the ionic phase. c The maximum shift current of the isolated single chain of TTF-CA with respect to the phase parameter. The inset presents the shift current spectrum for each geometry defined by λ with respect to the frequency of the applied light. d Schematic illustrations of the charge shift mechanisms in the neutral and the ionic phase. The black solid arrow indicates a large charge shift that is bigger than black dotted arrows indicating a small charge shift relatively. The cases for λ=0.0, 0.3, 0.5, 0.7, and 1.0 are selected to be presented in a and the inset of c.
frequencies of the applied light, as summarized in the inset of Fig.  3c, and the maximum value is plotted in Fig. 3c with respect to the phase parameter.
When the structure is free of the ferroelectric built-in potential, the forward and backward spatial charge shifts with respect to the chain direction should be symmetric without producing a net charge flow, as illustrated in the upper panel of Fig. 3d. The spatial charge shift must be combined with the intrachain ferroelectricity to produce the net shift current, as illustrated in the lower panel of Fig. 3d. As shown in Fig. 3c, the maximum shift current of this isolated single 1D chain of the ionic phase (λ = 1.0) is approximately 85% of that of the bulk 1D chain (compare Fig. 2f and Fig. 3c). Hence, the observed shift current in the ionic bulk phase of TTF-CA is mainly attributable to the intrachain effect. To maximize the shift current along the chain, we also tested various compositions of acceptor molecules, as summarized in Supplementary Fig. 4.
Interchain effect of shift current in the neutral phase TTF-CA and its optoelectronic feasibility We discussed above that the shift current in the ionic phase originates mainly from the intrachain effect, whereas the neutral phase suffers from complete cancelation due to the isotropy in the 1D chain, as illustrated in the upper panel of Fig. 3d. For the bulk three-dimensional solid, the interchain component must be considered separately, and we assume that the interchain charge shift, if any, should be perfectly canceled in the neutral phase, as illustrated in the left panel of Fig. 4a. The momentum-resolved analysis in the neutral phase ( Supplementary Fig. 5) reveals that each Bloch state provides nonzero shift current value, which is to be canceled by each other in the overall summation (Fig. 2e). In the neutral phase, the excited carrier velocities of k and -k points are the same in magnitude but directed in opposite ways. On the ferroelectric transition, the symmetry between k and -k points is lifted, leading to the net sum of the velocity of excited carriers, which results in large photocurrents in the bulk phase even in the absence of p-n junction. Hence, the net vanishing shift current, embedded under the symmetric lattice, can be released by lowering the symmetry of the structure. The intrachain ferroelectricity depends on the atomic configuration of the molecules in the unit cell; thus, it is difficult to adjust. However, the interchain effect can be easily adapted, for example, by constructing an interface or exposing the surface. Such lowering of the symmetry might alter the cancelation mechanism, thereby possibly leading to a net shift current not only in the chain-normal but also in the chain direction, as illustrated in the right panel of Fig. 4a. To evaluate our conjecture, we varied the interchain distance (Δd) and calculated the shift current of the neutral phase via the same procedure as above, as summarized in Fig. 4b. This variation of the interchain distance renders the structure non-centrosymmetric along the perpendicular direction (y or z axis). Remarkably, the x and z directional maximum shift currents are proportional to Δd, while the y component is persistently zero over the whole considered range of Δd. The detailed energy-dependent shift current spectra upon varying Δd are summarized in Supplementary Fig. 6.
Such a symmetry lowering and the resulting imbalance in the interchain charge shift, as depicted in the right panel of Fig. 4a, can be implemented in practice by constructing an interface. For example, when the neutral phase of the TTF-CA is placed on a ferroelectric layer, such as In 2 Se 3 and BiFeO 3 , the electric field from the interface can act as a source of the asymmetric potential without requiring any external power. A stronger field can be provided by sandwiching the TTF-CA film between ferroelectric surfaces with opposite polarity. Nevertheless, it is worth noting that the TTF-CA thin films have been deposited on a substrate in recent experiments, 4 which implies the feasibility of the symmetry adjustment by a gating device or a ferroelectric substrate. To test this concept, we evaluated the excited-state current by performing the rtp-TDDFT, including a symmetry-breaking uniform field in the z axis. In this simulation, an AC electric field with the frequency of the bandgap is applied in the chain direction ðẼ AC ¼ E 0x e Àiωt ; ω ¼ E gap = hÞ to describe the light-matter interaction within the dipole approximation, 29 corresponding to the band-edge excitation. The implementation of the electric field through the vector potential is described in the "Methods" section. A uniform DC electric field is applied in the perpendicular direction ðẼ DC ¼ E DCẑ Þ to break the interchain symmetry. The amplitude of the AC field that is applied along the x axis is E 0 = 19 mV/Å and that of the DC field that is applied along the z axis is E DC = 0.36 mV/Å. We gradually increase the field over an initial ramping period of 10 fs and keep it fixed over the remainder of the calculation period. From the time-evolving Kohn-Sham states ½ψ n;k ðtÞ, the time profile of the current density is calculated by taking the velocity expectation, as expressed in Eq. (2). The timeaveraged photocurrents ½J x;avg ðtÞ ¼ ð1=tÞ R t 0 J x ðτÞdτ are presented in Fig. 4c; the oscillatory photocurrent is averaged to zero when the symmetry-breaking DC field is not applied. This real-time profile of the current obviously shows that a large chaindirectional shift current can be turned on by applying the symmetry-adjusting field in the chain-normal direction. The same Fig. 4 Interchain effect of TTF-CA. a Schematic diagrams of the interchain charge shift in the neutral phase and asymmetric-neutral phase. b The maximum shift current of the neutral phase of TTF-CA with respect to interchain distance Δd, which is depicted in the inset of the schematic configuration. c The photocurrent of chain direction, calculated by rtp-TDDFT, in the neutral phase TTF-CA with/without a uniform DC field applied in the chain-normal direction. The inset depicts the symmetry-breaking uniform field in the z direction together with an oscillating field applied in the x direction.
real-time profile of the chain-normal direction (J z,avg (t)) is plotted in Supplementary Fig. 7. In addition, we investigated the variation in the Bader charge, the polarization, and the shift current spectra on the presence of the asymmetric potential, as summarized in Supplementary Fig. 8. Surprisingly, despite no discernible change in the Bader charge and very marginal change in x axis polarization, the chain-directional shift current of the asymmetric phase is comparable to that of the ionic phase. This finding again implies that the large shift current is not attributable to the charge transfer between the TTF-CA molecules, but rather the interchain picture, shown in Fig. 4a, provides a logically justifiable explanation.
In summary, we showed that both the ionic and the neutral phase of the TTF-CA molecular solid can provide nonzero shift currents upon band-edge infrared excitation. For the lowtemperature ionic phase, which has intrinsic ferroelectric polarization, the intrachain charge shift between the TTF and CA molecules induces a substantial photocurrent. In the hightemperature neutral phase, although the interchain charge shift is intrinsically allowed for each Bloch state, the net shift suffers complete cancelation due to the lattice symmetry upon the Brillouin-zone summation. We presented the first-principles TDDFT results to support that asymmetric potentials that are introduced perpendicular to the chain direction can lift the cancelation, thereby releasing the chain-direction shift current even in the high-temperature phase. Our results can provide a deep understanding of the association of the photovoltaic characteristics with the lattice symmetry of organic solids. Thus, they can be utilized as a platform for design of materials for organic optoelectronic applications.

The ground states and their optical responses
Our static DFT calculations are performed via the projected augmented plane-wave method, 30,31 as implemented in the Vienna ab initio simulation package (VASP). 32 For the exchange-correlation potential, the Heyd-Scuseria-Ernzerhof (HSE) hybrid functional is mainly used for the static calculations (Figs. 1, 2, and 3). 33 Experimentally determined structural parameters 22 are mainly used, and we confirmed that a selected set of fully relaxed structures produces essentially the same results ( Supplementary  Fig. 9). The energy cutoff for the plane-wave-basis expansion is set to 500 eV. We employ 8 × 8 × 4 and 10 × 1 × 1 k-point grids to sample the Brillouin zone for the bulk and one-dimensional (1D) cases, respectively. The shift current spectra are evaluated by self-developed code via the second-order optical response formalism. 7 A tight-binding Hamiltonian based on the maximally localized Wannier functions 34 is utilized in our calculations. For the shift current estimation, 25 × 25 × 10 and 100 × 1 × 1 k-point grids are adopted for the bulk and 1D cases, respectively. We crosschecked the accuracy of our code through comparison with the previous results. 3,10,35 Excited electronic structures To evaluate the real-time profile of the excited-state electric current, we perform the rtp-TDDFT calculations by using the plane-wave-based realtime evolution package developed by our group. [19][20][21] In this computation, the Kohn-Sham wave function, the density, and the Hamiltonian are selfconsistently evolved through the time-dependent equation as follows: i h ∂ ∂t ψ n;kr ; t ð Þ ¼ 1 2m Ài h∇ þ e cÃ ext 2 þ X λ v atomRλ t ð Þ þ V DFT ρr; t ð Þ ½ " # ψ n;kr ; t ð Þ; (1) where n is the band index andk denotes the Bloch momentum vector. The discretized time step for the time integration (Δt) is set to 2.415 s. In our calculations, the electric field is expressed by using the velocity gauge of the vector potential via the following relation: E t ð Þ ¼ À 1 c ∂Ã ext =∂t. For the initial state wavefunctions, the static ground-state band structures ½ψ n;k ðt ¼ 0Þ are used, as obtained by using the QUANTUM ESPRESSO package 36 with the Perdew, Burke, and Ernzerhof (PBE) exchange-correlation functional. 37 A 5 × 3 × 2 k-point grid is used to generate the time-propagating wave function. Using the time-evolving Bloch wavefunctions, the time profile of the current is evaluated as follows 21 : J t ð Þ ¼ À e m X n X k f n;k S hψ n;k t ð Þj b πjψ n;k t ð Þi (2) where n is the band index, f n;k is the initial occupation of the Bloch state, m is the mass of an electron, S is the area of the crystal, and the gauge-

DATA AVAILABILITY
All data generated and/or analyzed during this study are included in this article and its Supplementary Information file. The data are available from the corresponding author upon reasonable request.

CODE AVAILABILITY
The code used to evaluate shift current spectra is available from the corresponding author on reasonable request.