An organic quantum battery

An organic quantum battery J. Q. Quach*, K. E. McGhee, L. Ganzer, D. M. Rouse, B. W. Lovett, E. M. Gauger, J. Keeling, G. Cerullo, D. G. Lidzey, T. Virgili * 1. Institute for Photonics and Advanced Sensing and School of Chemistry and Physics, The University of Adelaide, South Australia 5005, Australia 2. Istituto di Fotonica e Nanotecnologia – CNR, IFN Dipartimento di Fisica, Politecnico di Milano, Piazza Leonardo da Vinci 32, 20133 Milano, Italy 3. Department of Physics and Astronomy, University of Sheffield, Hicks Building, Hounsfield Road, Sheffield S3 7RH, U.K. 4. SUPA, School of Physics and Astronomy, University of St Andrews, St Andrews, KY16 9SS, UK 5. SUPA, Institute of Photonics and Quantum Sciences, Heriot-Watt University, EH14 4AS, UK *Corresponding authors: quach.james@gmail.com, tersilla.virgili@polimi.it

experimentally demonstrate superextensive increases in both charging power and storage capacity, in agreement with our theoretical modelling. We find that decoherence plays an important role in stabilising energy storage, analogous to the role that dissipation plays in photosynthesis 5 . This experimental proof-of-concept is a major milestone towards the practical application of quantum batteries in quantum 6 and conventional devices 7 . Our work opens new opportunities for harnessing collective effects in light-matter coupling for nanoscale energy capture, storage, and transport technologies, including the enhancement of solar cell efficiencies 8-11 . Conventional batteries operate on the basis of classical electrochemical principles developed in the eighteenth century. Quantum batteries (QB) represent a new class of energy storage devices that instead operate on distinctly quantum mechanical principles.
In particular, they are driven either by quantum entanglement, that reduces the number of traversed states in the Hilbert space compared to (classical) separable states alone 1-4,12-16 , or by cooperative behaviour that increases the effective quantum coupling between battery and source [17][18][19] .
System properties can typically be categorised as intensive (i.e. they are independent of the system size, such as density) or extensive (i.e. they grow in proportion to system size, such as mass). QBs exhibit a counter-intuitive property where the charging time is inversely related to the battery capacity. This leads to the intriguing idea that the charging power of QBs is superextensive; that is, it increases faster than the size of the battery. Typically consisting of a collection of identical quantum subsystems to which an external energy source is applied 1 , QBs have been predicted to exhibit superextensive charging rate density (charging rate per subsystem) that scales as or √ in the thermodynamic limit 4 ; therefore, the total charging rate grows faster than the system size.
If experimentally verified, this would have important implications for energy storage and capture technologies. However, there are challenges in engineering the precise environment in which such behaviour can occur, and in monitoring the ultrashort charging time scales.
Here we experimentally realise a Dicke QB 19 (DQB) using an organic semiconductor as an ensemble of two-level systems (TLSs) coupled to a confined optical mode in a microcavity. This provides a practical version of the DQB which uses the optical cavity mode to induce collective coupling between light and the molecules. In this device, constructive interference of different absorption processes leads to enhanced transition rates, analogous to superradiant emission 20 and superabsorption 21 . We demonstrate how dissipation plays a crucial role in the performance of our DQB. In a closed system, the coherent effects that lead to fast charging can also lead to subsequent fast discharging. In our open noisy system, dephasing enables a ratchet effect, where excited states are capable of absorbing but not losing energy 22 , thereby retaining the stored energy until it can be used.

Quantum battery structure
The fabricated DQBs consist of a thin layer of a low-mass molecular semiconductor dispersed into a polymer matrix that is deposited by spin-coating and positioned between two dielectric mirrors, forming a microcavity as illustrated schematically in Fig. 1a (see Methods section for fabrication details). Organic semiconductors are particularly promising for many applications as the high oscillator strength and binding energy of molecular excitons means that light can be absorbed efficiently and excitons can exist at room temperature 23 . The organic semiconductor used in this study was the dye Lumogen-F Orange (LFO), whose chemical structure is shown in Fig. 1b Fig. 1). The absorption peak at 526 nm and the emission peak at 534 nm correspond to the 0-0 transition, i.e. an electronic transition to and from the lowest vibrational state.
Operating around the 0-0 transition, the LFO molecules can reasonably be considered as a TLS. We prepared samples with 0.5%, 1%, 5%, and 10% concentrations, as these are representative of the optimal operating regimes -further increases in concentration lead to quenching, and signals from lower concentrations are indiscernible from noise. The absorption and photoluminescence spectra for the 0.5%, 5% and 10% concentrations are given in Extended Data Fig. 2.  The optical microcavities fabricated support cavity modes whose energy is determined by the optical-thickness of the LFO layer and the penetration of the optical field into the cavity mirrors 24

Experimental setup
The charging and energy storage dynamics of the DQB were measured using ultrafast transient-absorption (TA) spectroscopy 25 , allowing femtosecond charging times to be measured. In this technique, we excite the DQB with a pump pulse, and then measure the evolution of stored energy (corresponding to the number of excited molecules) with a second probe pulse, delayed by time (Fig. 1d). The probe pulse is transmitted through the top distributed Bragg reflector (DBR) of the cavity, and the reflection from the bottom DBR is measured. The differential reflectivity induced by the pump-pulse is given by where ( ) is the probe reflectivity with (without) the pump excitation. Note that control films are measured under differential transmittivity Δ ∕ .
In our experimental setup (shown schematically in Fig. 1e), TA measurements were performed in a degenerate, almost collinear configuration. Pump and probe pulses were generated by a broadband non-collinear optical parametric amplifier (NOPA) 26  is the fraction of them that actually reach the active layer of DQB. We estimate from the reflectivity data that only 6% to 8% of the initial pump excitation enters the cavity.

Results
We first show that ultrafast TA spectroscopy can monitor the population of excited molecules, even in a cavity, by comparing the control film and the DQB spectra as shown in Fig. 2a. The control film Δ ⁄ spectra are shown for various probe delays By normalising the peak value of Δ / to ℏ ↑ max / , we can directly relate Δ ( ) to the temporal behaviour of the stored energy density of the DQB: where is the molecular electronic transition frequency (taken to be equal to the cavity frequency). ↑ max is calculated through theoretical modelling. We also note that two of the DQB spectra show a negative Δ ⁄ band, which results from the change in the refractive index induced by the pump pulse 28 . Figure 2b shows the experimental values for the time-dependent stored energy density in the DQBs. In all DQBs studied, the energy density undergoes a rapid rise followed by slow decay. The timescale of the rapid rise varies with concentration. We adjust the laser power to fix photon density r across comparable DQB samples, and compare behaviour with different LFO concentrations. We found that to achieve a sufficiently high signal-to-noise ratio, it was not possible to compare all DQBs at the same r value; instead, a constant r value was maintained for matched DQB samples.
Specifically, measurements were made on DQBs with LFO concentrations of 10%, 5% and 1% with approximately constant ≃ 0.14 (respectively labelled A1, A2, and A3), and 1% and 0.5% with ≃ 2.4 (labelled B1 and B2).  Overlaying the experimental data are the corresponding theoretical predictions (see Theoretical Model section), convolved with a Gaussian response function of ~120 fs full-width-at-half-maximum, to account for both the instrument response time (~20 fs) and the cavity photon lifetime (which was extracted from the linewidth of the cavity mode). Although the signal-to-noise ratio is not sufficiently high to discern the predicted oscillatory behaviour in the A1 and A2 data, there is otherwise good agreement between the experimental data and the corresponding theoretical model.
To obtain the energetic dynamics of the DQBs, we take away the response function from the theoretical fit, as shown in Extended Data Fig. 3. Extended Data Table   1

Theoretical Model
The experimental dynamics of the DQBs can be reproduced by modelling, with the Lindblad master equation (LME), the N TLSs in an optical cavity with light-matter coupling strength g, a driving laser with a Gaussian pulse envelope and peak amplitude 0 , and three decay channels corresponding to the cavity decay ( ), TLS dephasing ( ), and TLS relaxation ( − ). To solve this many-body LME, we make use of the cumulant  Table 1.
Specifically, A1 and A2 operate in the coupling-dominated regime, where scales slightly less than −1 2 ⁄ , max scales slightly more than 0 , and max scales slightly more than 1 2 ⁄ . For the region between A2 and A3, the average scaling of falls between 0 and −1 2 ⁄ , max between 2 and 0 , and max between 2 and 1 2 ⁄ . As A2 is further in the coupling-dominated regime than A3 is in the decay-dominated regime, the average scaling values between A2 and A3 are skewed towards the coupling-dominated scalings.
the decay-dominated and coupling-dominated scalings, as reflected in Table 1.

Summary
We have demonstrated the operation of an open-system DQB driven by an external coherent light source and shown regimes of superextensive behaviour using ultrafast optical spectroscopy. We have provided direct experimental evidence of superextensive energy storage capacity and charging in the decay-and coupling-dominated regimes. Our realisation of a prototype QB highlights the fact that purely closed unitary dynamics is insufficient for realising a practical QB. The retention of energy requires finely-tuned decoherence processes, allowing the battery to charge quickly and yet discharge much more slowly. Our observation of such ratchet-type behaviour shows that realistic noisy environments are crucial for the implementation and application of useful QBs.
A prime example of how QBs could be applied to enhance existing technologies is in the field of solar energy conversion. The efficient solar energy harvesting that happens in natural photosynthesis derives from the intricate interplay of coherent and dissipative dynamics that we also exploit in our QB design, and we may learn further lessons from nature for improving QBs. Further, QBs could be made to charge more efficiently by surrounding them with antenna molecules that capture light through superabsorption 21,32,33 . Since the working principle of DQBs is closely related to superabsorbing behaviour, it may be possible to generate even faster charging by using environmental engineering and control approaches. These could keep the QB operating in the range of higher-lying energy states that are associated with maximum absorption enhancement, i.e. near the mid-point of the Dicke ladder 22 . Although our QB was charged by a coherent laser source, it opens a pathway for charging with incoherent sunlight, offering an exciting new approach to the design of solar-cell technology.

Quantum battery fabrication
The The diluted molecules are expected to be isolated at low concentration 0.1 -1%, but at higher dye concentrations, the 0-0 emission transition red-shifts by a few nm and the second peak increases in intensity due to aggregation of the dye molecules. This is evident in Extended Data Fig. 2a and b, with additional broader features observed at longer wavelengths, which we assign to intermolecular states such as excimers.
The 0.5% and 1% cavities lie in the weak-coupling regime, i.e. no polaritonic splitting could be seen in the cavity reflectivity spectrum, as shown in Extended Data Extended Data Fig. 4 shows a transfer matrix simulation of the electric field distribution of the 1% cavity (the cavities exhibit similar distributions).

Pump-probe spectroscopy
Probe and pump pulses were generated by a broadband non-collinear optical parametric amplifier (NOPA). The NOPA was pumped by a fraction (450 μJ) of the laser beam generated by a regeneratively amplified Ti:Sapphire laser (Coherent Libra) producing 100 fs pulses at 800 nm at a repetition rate of 1 kHz. A pair of chirped mirrors were placed at the output of the NOPA to compensate for temporal dispersion, and by using 7 'bounces' we were able to generate pulses with a temporal width below 20 fs. The laser beam was then split by a beam-splitter, with the probe being delayed via a translation stage and the pump being modulated mechanically using a chopper at 500 Hz.

Number of molecules in film
To The transmission of the 10% LFO concentration in film was then measured to obtain and hence n (number density of molecules in the cavity active layer), using the measured value of , with (film thickness) measured using a Bruker DektakXT profilometer. This value was then multiplied by the area of the laser beam and d to obtain . Here we assume a uniform distribution in the active layer. for other concentrations were scaled accordingly.

Lindblad master equation
The open driven system of the experiment is modelled with the Lindblad master equation, where Δ = ω − ω is the detuning of the cavity frequency from the laser driving frequency. The LFO molecules are initially in the ground state, and the laser is onresonance (Δ = 0).

Cumulant expansion
The energy density of the cavity containing identical molecules with transition energy is ( ) = This approximation is valid at large N, as corrections scale as 1⁄ . To capture the leading order effects of finite-sizes we make a second-order cumulant expansion 29-31 , i.e. we keep second-order cumulants ⟨⟨ ⟩⟩ = ⟨ ⟩ − ⟨ ⟩⟨ ⟩ and assume that the third-order cumulants vanish, which allows us to rewrite third-order moments into products of first and second-order moments 36 . In our experiments, the number of molecules in the cavity is large (>10 10 ) and we find higher order correlations are negligible. We give the equations of motion up to second order in the Supplementary Information.

Operating regimes
The decay-dominated (purple region in Fig. 3a and b) regime occurs when the collective light-matter coupling is weaker than the decay channels, √ ′ < { , , − }, where ′ = max (1, ). In this regime, the time scale of cavity dynamics is slow relative to the decay rate. Fig. 3c shows a typical time dependence of the DQB in this regime, indicating how the model parameters affect the dynamics. In this regime, the increase in the effective coupling relative to the decay strength sees an 2 superextensive scaling of the energy and power density, while rise time remains constant. Experiment A3 operates near the boundary of this regime (Fig. 3a).
In the crossover between the regimes (purple-green), the collective coupling falls between the cavity decay rate and the TLS dephasing rate, { , − } < √ ′ < . In Fig.   3a and b, − is small such that √ ′ ≫ − for all values of N, and so there is no boundary labelled for this decay rate. In this case, capacity and rise-time can simultaneously scale super-and subextensively, but at a rate slower than in the decay and coupling-dominated regimes, respectively. Experiments B1 and B2 operate in this regime (Fig. 3b).

Decay and coupling rates
The parameters needed in the theory calculations are the cavity leakage rate , the dephasing rate , the non-radiative decay rate − , the interaction strength , and the The battery capacity max , is the peak stored energy per molecule or energy density. The battery charging rate max = max ( ∕ ), is the peak charging power per molecule or charging power density.

S1. CUMULANT EQUATIONS
Here we derive the equations of motion for the expectation values using a second-order cumulant approach. As discussed in the main text, for N identical molecules of energy ω a , the energy density of the quantum battery is where σ z (t) is the expectation value of any one of the N molecules, which are all assumed identical. The system, comprising the molecules and the cavity mode, evolves under the master equation written in the main text, reproduced here ( = 1): where a (a † ) is the photon annihilation (creation) operator, σ α i for α = x, y, z are the Pauli matrices, ∆ is the energy detuning of the laser from the cavity and molecules (set to 0), g is the coupling strength of each molecule to the photon mode, η(t) is the Gaussian profile of the pump laser and κ, γ z and γ − are the cavity leakage, dephasing and non-radiative decay rates.
To determine the time evolution of σ z (t) = Tr [σ z ρ(t)] we begin by writing down the first order expectation values of the system. We adopt the notation C a (t) ≡ a(t) for photon operators, and C α=x,y,z (t) ≡ σ α (t) for spin operators, along with a similar notation for higher order expectations, e.g., C ax (t) ≡ aσ x (t) . The equations of motion for the first order expectation values are where ∂ t is short for ∂ ∂t , γ tot = 2γ z + 1 2 γ − , and for notational ease we have dropped the explicit time dependence of observables. As described in the main text, in mean field theory we would now set the second S2 order cumulants to zero. These are defined as This would result in the usual decomposition of second order expectation values into products of first order ones, C AB = C A C B , which is the assumption that molecule-molecule, molecule-photon, photon-photon and all higher order correlations are negligible. However, we instead derive equations of motion for the second order expectation values, capturing the leading order 1/N corrections to mean field theory. The second order photon correlations obey: while molecule-photon correlations follow: These now depend on third order expectation values, some of which contain multiple Pauli operators. We must note that these terms indicate Pauli operators representing different molecules and so commute -we have already taken into account the cases where the Pauli operators correspond to the same molecule by using the Pauli algebra σ α σ β = 1δ αβ + iσ γ αβγ . The molecule-molecule expectation values for the same Pauli operator acting on different molecules are Finally, the molecule-molecule expectation values for different Pauli operators acting on different molecules are In principle one can continue to write equations of motion for increasingly higher orders of expectation values, however, at large N , most essential physics is obtained at second order. We therefore truncate the cumulant expansion by setting third order cumulants to zero. These are defined as and so setting ABC = 0 closes the system of differential equations, allowing us to write ABC in terms of first and second order correlations.

S2. THEORETICAL BEHAVIOUR IN THE THERMODYNAMIC LIMIT
In Fig. S1 we present the theoretical N -dependence of the charging time τ , maximum energy density E max , and maximum power density P max over a wider range of N than shown in Fig. 3 of the main text. This shows that in addition to the decay-dominated (purple) and coupling-dominated (green) behavior described in the main text, a third region occurs at even larger N , which we discuss below.
In the main text we discussed the QB energetic dynamics around the decay-to-coupling dominated crossover regime, as this was the experimental operating region. Moving deeper into coupling-dominated regime does not necessarily improve the battery operation. This is illustrated in Fig. S1(b) which shows the simulated four points, corresponding to the circles in Fig. S1(a). Within the coupling-dominated regime, energy stored within the battery rapidly oscillates which is not a desirable feature. This occurs because the light and matter degrees of freedom hybridise to form polaritons with upper and lower branches split by Rabi frequency ±g √ N , leading to beating between these modes. These oscillations are not present in the experimentally studied crossover region. In this region, dephasing is strong enough to prevent oscillation in energy, yet weak enough to warrant superextensive charging. Therefore, this is the optimal region to produce a QB. Going deeper into the coupling-dominated regime would only be advantageous if energy was extracted from the battery on a shorter timescale than the period of oscillations, or additional mechanisms were in place to stabilise the oscillations.
At even larger N (red region) the stored energy falls with increasing N . This can be understood as arising from a condition where the polariton energy splitting exceeds the bandwidth of the pump (set by its finite pulse duration), suppressing energy absorption. Numerically, we find this occurs when N > N σ where g √ N σ = (2/5) will be explained in Section S3. To build an efficient QB in this regime, one should tune the frequency of the laser to match the polariton energies. Additionally, the time dynamics of energy absorption here change significantly, with the second half of the laser pulse causing stimulated emission, reducing the stored energy -such dynamics arises naturally from a toy model of strongly coupled modes with a splitting larger than the pulse bandwidth, and can be seen in the form of the red line in Fig. S1(b). There are two timescales in this system: the vacuum Rabi splitting (i.e. polariton detuning) g √ N and the Rabi splitting g √ rN where rN is the number of photons in the cavity. In Fig. 3(a), we show that the battery charges super-extensively once g √ N is greater than all decay channels. However, this is only true if r ≤ 1, as is true in experiments A1, A2 and A3. In Fig. 3(b), the boundaries N κ and N γ z are instead determined by g √ rN being equal to the decay rates. This is because r ≥ 1 in experiments B1 and B2. More generally, the important timescale is the larger of the polariton detuning and the Rabi splitting, and so the coupling dominant regime occurs when g Max(1, r)N is larger than all decay channels.
In Figure S2 we plot the charging time τ as a function of N and r. Here, we set κ = γ − = γ z ≡ Γ = 2 meV (note that γ z is independent of N ) so that there is only one boundary between the decay dominant and coupling dominant regimes. The green, red and dashed-black lines show the boundaries between the decay dominant and coupling dominant regimes (N = N Γ ) if g √ N , g √ rN or g Max(1, r)N are used as the relevant coupling scale respectively. Clearly, the boundary is determined by g Max(1, r)N for all values of r. We also show the boundary between the coupling dominant and non-resonant regimes (N = N σ ) as the cyan line. When r > 1, we find that N σ becomes linearly dependent on r. The prefactor (2/5) 1 4 is necessary for N σ to align with the contours of increased charging time for r > 1. All parameters are equivalent to the Q1% cavity (see main text) with the exception that the dephasing and non-radiative decay rates are equal to the cavity leakage rate (set to γ z = κ = Γ = 2 meV) and note that the dephasing rate is independent of N . Figure S3 shows how capacity, charging time and power vary as a function of laser intensity r at fixed number of molecules N . For small r, we find that the maximum energy and power densities vary linearly with r, while charging time is constant. This simply reflects the total energy in the cavity. The charging time is constant because decay channels still dominate over coherent dynamics. As r is increased beyond r = 1, the important timescale g max(1, r)N begins to scale with r, and so the boundaries separating the coupling dominant and decay dominant regions N κ and N γ z are pushed to smaller N . When these boundaries become smaller than the number of molecules in the cavity, the charging time begins to scale as 1 / √ r. Additionally, the energy density begins to saturate because there are already many more photons than there are molecules within the cavity. In Figure S3(b) we also plot the experimentally measured energy densities, and we see there is good agreement to the theoretical curve. The coloured points in Figure S3(a) indicate the charging time, maximum capacity and maximum power of the temporal dynamics of the same colour in Figure S3  As outlined in the main text, we used a reduced chi-square optimisation procedure to determine the light matter coupling g, dephasing constant γ z 0 and non-radiative decay rate γ − , as well as to estimate uncertainties on these parameters. As these quantities represent molecular properties, we would expect them to be the same in all the different experiments. For this reason, we performed a global fit rather than performing the procedure individually for each experiment. In this section we give further details on this calculation.

A. Fitting procedure
The steps of our fitting procedure are as follows: considering the point-to-point variation. Because the uncertainty is higher near t = 0, when the pump arrives, we use different error estimates in different time windows. Specifically, we divide experiments A1 and A2 into five windows t < −300fs; −300fs < t < 300fs; 300fs < t < 700fs; 700fs < t < 1000fs; t > 1000fs. For experiments A3, B1 and B2 we found that four windows t < −300fs; −300fs < t < 300fs; 300fs < t < 1000fs; t > 1000fs was sufficient. In each window, the uncertainty estimate for each experiment is taken from the variance over a narrow range of points (typically 150fs) where there is no strong time dependence.
3. For each set of parameters, we performed an "internal" chi-square minimisation to find the optimal scaling factor S between the battery energy E(t) and the measured differential reflectivity ∆R/R, and a time shift between the theory and experiment T 0 . That is, we minimise with respect to S and T 0 . We treat the result of this minimisation as the chi-squared value which we use in the following steps to estimate the meaningful parameters g, γ z , γ − and their uncertainties.
Estimating the scaling factor S from first principles is difficult because of reflections by the cavity mirror, hence this factor is found by the best fit value. The time shift reflects uncertainty of delays in the optics, so that it is not a-prior clear when the peak of the pump pulse arrives. After this shift, we define t = 0 as the moment the pump arrives. This is important when calculating the charging time τ , which we defined as the time from the arrival of the pump until the battery reaches half maximum energy.
4. We then use the chi-square value described above, and divide by the total number of degrees of freedom k eff = k−3 (where k is the total number of data points), to arrive at the final reduced chi-squareχ 2 map. A slice of this three dimensional reduced chi-square map is shown in Figure S4 for γ − = 0.0263 meV, which is the optimal non-radiative decay rate given in the main text. The optimal parameter set used in the main text that optimisesχ 2 is shown as the red point in Figure S4. We findχ 2 min = 2.76, suggesting our estimated measurement uncertainties onχ are reasonable, but likely underestimates.
5. Finally, the 68% confidence interval for each parameter was estimated by considering the contour for whichχ 2 =χ 2 min + 1 k eff ∆ * where ∆ * = 3.51 is extracted from the reduced chi-square distribution for 3 parameters and error tolerance (68%), see [1]. In the right panel of Fig. S4 we show the contour as a white line, and the actual parameter values which lie within this 68% contour as black points.  FIG. S4. An illustration of the reduced-chi square optimisation procedure to find the optimal parameters for the theoretical model and their 68% confidence intervals, where ∆ = ∆ * /k eff . The chi-square contour plots shown in this figure are slices of the full three-dimensional map at the optimum non-radiative decay rate γ − = 0.0263 meV used in the main text. In the yellow region of the bottom right cornerχ 2 > 10, which we do not show to emphasise smaller variations inχ 2 .

B. Residuals of the best fit
To check whether systematic errors arise from our fitting procedure, Fig. S5 shows the residual errorsi.e. difference between the theoretical curves and the experimental data-for the five experiments shown in Figure 2 of the main text. As is clear, there are no discernible features in the residuals that are consistently present across the different experiments. This indicates that the theoretical curves account for the essential characteristics of the data.