# Energy storage in Madeira, Portugal: co-optimizing for arbitrage, self-
sufficiency, peak shaving and energy backup

Md Umar Hashmi  INRIA and DI ENS, ENS, CNRS, PSL University, Paris, France
Lucas Pereira  ITI, LARSyS, Ténico Lisboa, and prsma.com  Ana Bušić  INRIA and
DI ENS, ENS, CNRS, PSL University, Paris, France

###### Abstract

Energy storage applications are explored from a prosumer (consumers with
generation) perspective for the island of Madeira in Portugal. These
applications could also be relevant to other power networks. We formulate a
convex co-optimization problem for performing arbitrage under zero feed-in
tariff, increasing self-sufficiency by increasing self-consumption of locally
generated renewable energy, provide peak shaving and act as a backup power
source during anticipated and scheduled power outages. Using real data from
Madeira we perform short and long time-scale simulations in order to select
end-user contract which maximizes their gains considering storage degradation
based on operational cycles. We observe energy storage ramping capability
decides peak shaving potential, fast ramping batteries can significantly
reduce peak demand charge. The numerical experiment indicates that storage
providing backup does not significantly reduce gains performing arbitrage and
peak demand shaving. Furthermore, we also use AutoRegressive Moving Average
(ARMA) forecasting along with Model Predictive Control (MPC) for real-time
implementation of the proposed optimization problem in the presence of
uncertainty.

###### Contents

  1. 1 Introduction
  2. 2 Power system norms in Madeira
    1. 2.1 Overview of the Madeira Electric Grid
    2. 2.2 Peak Power Contracts, Tariffs and Billing Cycles
    3. 2.3 Self-Consumption and Renewables in Madeira
  3. 3 System Description
  4. 4 Co-optimizing energy storage
    1. 4.1 ToU pricing + zero feed-in-tariff + Peak-Shaving
    2. 4.2 Storage for BackUp with Arbitrage + Peak Shaving
    3. 4.3 Open Source Codes
  5. 5 Real-time Control under Uncertainty
    1. 5.1 Modeling Uncertainty: ARMA Forecasting
    2. 5.2 Model Predictive Control
  6. 6 Numerical Results
    1. 6.1 Deterministic Solution for $P_{opt}$
      1. 6.1.1 Shorter Time Scale: A day
      2. 6.1.2 Longer Time-Scale: A month
    2. 6.2 Co-optimizing with Power Backup
    3. 6.3 Real-Time Implementation (Forecast plus MPC)
  7. 7 Conclusion

##  1 Introduction

Medium sized isolated power networks often restrict the share of renewables
and enforce stringent rules necessary to maintain safety and stability of the
electrical power system. Authors in [1] indicate that a high share of wind
energy penetration could lead to large variations in active power generated
due to sudden changes in wind speed. This in effect creates a huge mismatch
between supply and demand, causing large variations in voltage and frequency
leading to hazardous operating conditions. The enforcement regulations set by
the Ministry of Industry and Energy of the Canary Government goes as far as
restricting any further increase in wind farm installations directly connected
to the power network. However, additional renewables are considered favourable
if the generated energy is self-consumed locally [2]. In [3] the use of wind-
powered hydro storage system for increasing the penetration of renewables is
proposed for the island grid of Gran Canaria (part of Canarian Archipelago).
In this solution, pumped storage acts as energy storage which is feasible due
to the geography of the island and might not be applicable to similar isolated
power networks.

The work in [4] presents a case study for Coimbra where residential solar and
energy storage is locally consumed with goal of zero energy buildings. Authors
observe that the electricity bill for the household reduces by more than 87%.
Furthermore, self-consuming intermittent renewable generation locally is
desired by the utilities, as such generation makes load balancing, frequency
and voltage regulation more challenging. Thus self-consuming renewable
generation assists the power grid to accommodate a larger share of renewables
[5]. Often utilities set the feed-in-tariff lower than the retail rate, making
self-consumption more desirable for electricity consumers. Authors in [6]
observe that by increasing self-consumption of renewables, their financial
feasibility could be achieved a decade before in some European countries.
Although self-consumption holds multiple benefits for utilities, such a
constraint creates a disparity for Distributed Generation (DG) owners, as
excess generation cannot be fed back to the grid.

In order to mitigate this disparity, these DG owners should co-optimize for
additional revenue streams. In Germany DG owners until 2012 were incentivized
for self-consumption rather than feeding power back to the grid [7]. Under the
case where DG owners are not allowed to supply power back to the grid, energy
storage can facilitate load-shifting in real-time, minimizing consumption
cost, increasing self-sufficiency [5]. Self-sufficiency is the ratio of total
energy demand met by local generation and/or storage with respect to
cumulative energy needs. Storage acts as a buffer of electrical energy and
assists in the temporal shift of energy usage. Authors in [8] describe the
various applications for energy storage in future power networks. Furthermore,
with greater integration of intermittent renewables performing arbitrage and
ancillary services will be more profitable [9]. Co-optimizing energy storage
for multiple applications has been proposed in many recent works. The authors
in [10] co-optimize for arbitrage and frequency regulation, [11] co-optimizes
peak saving and frequency regulation. In this paper, we consider the case of
Madeira Island, where utilities promote the inclusion of DG for self-
consumption only. We propose the integration of a DG source along with energy
storage battery. This battery facilitates self-consumption of locally
generated energy, assist users in selecting a lower peak power contract,
performing arbitrage and providing energy backup for instances of probable and
scheduled power outages.

The key contributions of this paper are:
$\quad\bullet$ Co-optimization: We propose a convex formulation for energy
storage control for performing arbitrage, peak demand charge saving and backup
reserve during power outages considering efficiency losses, ramping and
capacity constraints for an energy storage battery.
$\quad\bullet$ Storage profitability: The operational cycles govern storage
degradation. The performance index for monetary value per cycle introduced in
[12] is used. This index indicates that battery is financially profitable in
Madeira.
$\quad\bullet$ Real time implementation: We use Auto-Regressive Moving Average
(ARMA) processes to model temporal evolution in the MPC framework for real-
time implementation considering uncertainty, motivated by our prior work [13].

The key observations made from numerical results are:
$\quad\bullet$ Storage owners would benefit more by performing arbitrage with
contracts with more price variations.
$\quad\bullet$ Ramping capability primarily decides the ability to reduce peak
demand charge savings for consumers.
$\quad\bullet$ When DG generation is smaller than or approximately equal to
inelastic load in magnitude then self-sufficiency is governed by only the DG
generator, otherwise, storage also contributes to it by increasing self-
consumption.
$\quad\bullet$ Providing energy backup does not noticeably reduce storage
ability to perform peak demand and/or arbitrage.

The paper is organized as follows. Section 2 introduces power system norms for
consumers in Madeira. Section 3 presents the system description. Section 4
formulates the co-optimization problem of performing arbitrage, peak shaving
and energy backup. Section 5 presets the real-time control of storage under
uncertainty. Section 6 presents the numerical results. Section 7 concludes the
paper.

##  2 Power system norms in Madeira

Madeira is an archipelago in the North Atlantic Ocean, located about 1000 km
southwest of mainland Portugal. It has a population of almost 270,000. 111,000
of which living in the capital city of Funchal.

###  2.1 Overview of the Madeira Electric Grid

Madeira relies on local generation for electricity. Empresa de Eletricidade da
Madeira, S.A. (EEM) is the only DSO/TSO in Madeira, and is responsible for the
activities related to production, transport, distribution and
commercialization of electric energy. Energy generation sources used in
Madeira are: thermal energy from fossil fuels like diesel and natural gas,
hydro, wind, solid waste incineration (SWI), and photovoltaic (PV). In 2017,
Madeira consumed 800GWh energy; thermal constituted about 70% of the energy
mix, with the remaining 30% coming from renewable sources (hydro: 12.2%, wind:
9.5%, SWI: 4.8%, and PV: 3.6%) [14]. The non-domestic sector (e.g., tourism
and commerce), the domestic sector contributes 45% and 30% of total
consumption respectively. The remaining 25% are contributed by public lighting
(9%), public buildings (8%), industry (7%), and agriculture ($<$ 1%) [14].

###  2.2 Peak Power Contracts, Tariffs and Billing Cycles

As of 2018, Low Voltage (LV) customers can select between 8 levels of peak
power contract (PPC), three power tariffs, and two billing cycles [15]. Thus,
there are total 48 different contracts that users can select from.

Figure 1: Billing cycle schemes currently in place.

LV customers are subject to a maximum peak power contract (kVA), ranging from
3.45 kVA to 20.70 kVA (see Tab.1). PPC value is selected by the customer based
of their estimated electricity needs and should not be exceeded since the
supply is shut-down when that happens. The disconnection of power is done
locally and consumer can restart their energy-meter, however, sudden
interruptions should be avoided as it may damage appliances.

Table 1: Peak Power Contract for LV customer as of 2018. Peak Power Contract
(EUR/Day)
---
PPC (kVA) | Single-rate (euros) | Dual/Triple-rate (euros)
3.45 | 0.1611 | 0.1643
4.60 | 0.2096 | 0.2132
5.75 | 0.2560 | 0.2590
6.90 | 0.3040 | 0.3080
10.35 | 0.4478 | 0.4532
13.80 | 0.5902 | 0.5981
17.25 | 0.7326 | 0.7436
20.70 | 0.8751 | 0.8892

Regarding the energy tariffs, there are three options for LV customers.
Single, dual and triple-rate per kWh. In the Single-rate tariff the price is
fixed at 0.1629 euros/kWh. In the dual-rate tariff there are two time-of-use
(ToU) prices (peak at 0.1894 euros, and off-peak at 0.0982 euros), whereas in
the triple-rate tariffs have three ToU price levels (peak at 0.2153 euros,
off-peak at 0.0982 euros, and half-peak at 0.1716 euros). The off-peak, peak
and half-peak periods are defined in advance based on the notion of daily and
weekly billing cycles. In the daily cycle there are no distinctions between
weekend and workdays, whereas in the weekly cycle there are different time-of-
use periods for work-days, Saturdays and Sundays. Fig. 1 summarizes the
billing cycle schemes in practice as of 2018. For dual-rate tariff, the peak
period includes also the half-peak period.

Consumer electricity bill consists of two components: a fixed component that
depends on the contracted power (kVA), and a dynamic component governed by
actual energy consumption (kWh). Note that the fixed component is governed by
apparent power: a function of active and reactive power. In this work, we
consider PPC levels to a function of only active power as apparent power is
primarily governed by active power. Reactive power compensations as proposed
in [16] could be applied for the storage control in this work.

###  2.3 Self-Consumption and Renewables in Madeira

In Madeira island, since 2014, new mini and micro-producers are not allowed to
feed-in excess production to the local grid, thus excess generation is wasted
[14]. Due to such a constraint DGs are sized to maximize self-consumption and
minimize excess production [17]. Counterintuitively, in a period that one
should expect an explosion in the number of micro-producers leveraged by the
relatively low prices of solar PV technologies, Madeira is experiencing an
stagnation on the number of new solar PV installations due to the norms set by
EEM. The main reason for this change in the local legislation is to protect
the grid from the issues associated with the intermittent and uncertain nature
of renewable production from solar in a total energy system. In case of the
electrical LV networks in the rural areas of Madeira Island which are at the
edge of the radial distribution network, when associated with low consumption
and high production periods, it is very likely to observe the phenomena of
voltage increase [14].

##  3 System Description

In this work we consider a prosumer with inelastic load and rooftop solar
generation as DG and an energy storage battery. The battery serves four
purposes: (i) increase self-consumption, thus reducing waste of excess
generation if any, (ii) perform energy arbitrage, (iii) minimize the peak
demand charge and (iv) maintain battery charge level for scheduled and/or
anticipated power failures.

Notation: We consider operation over a total duration $T$, divided into $N$
steps indexed by $\\{1,...,N\\}$. The duration of each step is denoted as $h$.
Hence, $T=hN$. At time instant $i$, the information available is the end user
consumption $d_{i}$, the renewable generation ${r}_{i}$ and the storage energy
output $s_{i}$. The load without storage is denoted as $z_{i}=d_{i}-r_{i}.$
The load seen by grid is denoted as $L_{i}=d_{i}-r_{i}+s_{i}.$

Battery Model The efficiency of charging and discharging is denoted by
$\eta_{ch},\eta_{dis}\in(0,1]$, respectively. The change in energy level of
the battery is denoted as $x_{i}$=$h\delta_{i}$, where $\delta_{i}$ denotes
storage ramp rate at time instant $i$, such that
$\delta_{i}\in[\delta_{\min},\delta_{\max}],\forall i$ and $\delta_{\min}\leq
0$ and $\delta_{\max}\geq 0$ are the minimum and maximum ramp rates (kW);
$\delta_{i}>0$ implies charging and $\delta_{i}<0$ implies discharging. The
energy output of storage in the $i^{\text{th}}$ instant is given by
$s_{i}=\frac{[x_{i}]^{+}}{\eta_{ch}}-\eta_{dis}[x_{i}]^{-}$, where
$[x]^{+}$=$\max(0,x)$ and $[x]^{-}$=$-\min(0,x)$. The ramping constraint
induce limits on $s_{i}$ given by

$s_{i}\in[\delta_{\min}h\eta_{dis},\delta_{\max}h/\eta_{ch}],\quad\forall i.$
(1)

The energy stored in the battery is denoted as $b_{i}$, defined as
$b_{i}=b_{i-1}+x_{i}$. Battery capacity constraint is given as

$b_{i}\in[b_{\min},b_{\max}],\quad\forall i.$ (2)

where $b_{\min}$ and $b_{\max}$ are minimum and maximum permissible battery
charge levels respectively. We use xC-yC notation to represent the
relationship between ramp rate and battery capacity. xC-yC implies battery
takes 1/x hours to charge and 1/y hours to discharge completely.

##  4 Co-optimizing energy storage

Optimizing the energy storage is essential due to its high cost. In context of
our present work we do not evaluate the fixed electricity cost as under such a
case storage can only be used either for backing up excess generation and/or
for peak demand shaving. Backing up energy will require no look-ahead and
greedy behavior leads to optimality [18]. The optimal solution in such a case
is governed by the sign of $z_{i}$ and is given as

  * •

if ${z}_{i}\geq 0$ then battery should discharge such that
$s_{i}=\max\left\\{-{z}_{i},\delta_{\min}h\eta_{dis},({b_{i-1}-b_{\min}}){\eta_{dis}}\right\\}$

  * •

if ${z}_{i}<0$ then battery should charge such that
$s_{i}=\min\left\\{-{z}_{i},\delta_{\max}h/\eta_{ch},({b_{\max}-b_{i-1}})/{\eta_{ch}}\right\\}$.

Next we formulate co-optimization problem for storage control for arbitrage,
peak shaving and energy backup.

###  4.1 ToU pricing + zero feed-in-tariff + Peak-Shaving

The optimal arbitrage problem with battery ($P_{arb}$) is defined as the
minimization of the cost of total energy consumed denoted as
$\min\sum_{i=1}^{N}[L_{i}]^{+}p_{\text{elec}}(i)$ subject to the battery
constraints. $p_{\text{elec}}(i)$ denotes the electricity price for consuming
electricity, i.e. $L_{i}>0$, for instant $i$. Here we assume the feed-in-
tariff to be zero. This optimization framework is a special case of the
problem studied in [19], [20]. Under zero feed-in-tariff, only consumed energy
is charged and end-user gets no incentive in supplying power back to the grid
(i.e. $L_{i}<0$). This is emulated using a variable
$\theta_{i}=\max(0,L_{i})$. Keeping $L_{i}\geq 0$ will maximize the self-
consumption of renewable generation locally.

In this formulation we also consider peak demand shaving by reducing the peak
demand contract the end-user should opt for. Consider the maximum demand of
the user without storage is $P_{\max}$ then the inclusion of energy storage
would bring this maximum demand lower proportional to its ramp rate,
considering look-ahead optimization ensures battery has enough capacity to
discharge during peak demand. The end user operates the energy storage for
minimizing the cost of consumption, increasing self-consumption by reducing
the waste of excess of generation and restraining peak demand. This
optimization problem is given as

$\displaystyle(P_{opt})~{}~{}~{}$
$\displaystyle\mathmakebox[\widthof{$\underset{\displaystyle
s_{i}}{\mathrm{subject~{}to}}$}][l]{\underset{\displaystyle
s_{i}}{\mathrm{minimize}}}\quad\sum_{i=1}^{N}{p}_{\text{elec}}(i)\theta_{i}h$
(3a)
$\displaystyle\mathmakebox[\widthof{$\underset{\displaystyle\phantom{s_{i}}}{\mathrm{subject~{}to}}$}][c]{{\mathrm{subject~{}to}}}$
Ramping constraint, $\displaystyle Eq.~{}\ref{constraintramp},$ (3b) Capacity
constraint, $\displaystyle Eq.~{}\ref{constraintcapacity},$ (3c)
$\displaystyle\text{Self-sufficiency,}~{}~{}\theta_{i}$ $\displaystyle\geq 0,$
(3d) $\displaystyle\text{Arbitrage,}~{}~{}\theta_{i}$
$\displaystyle\geq[z_{i}+s_{i}],$ (3e) $\displaystyle\text{Peak
shaving,}~{}~{}[z_{i}+s_{i}]/h$ $\displaystyle\leq P_{\max}^{set}$ (3f)

The peak power threshold, $P_{\max}^{set}$, is selected close to the power
level ($P_{\max}+\delta_{\min}$), subject to
$P_{\max}^{set}\geq(P_{\max}+\delta_{\min})$. $P_{\max}^{set}$ is selected by
the electricity consumer as a PPC contract with the utility in Madeira.

###  4.2 Storage for BackUp with Arbitrage + Peak Shaving

The valuation of energy storage devices performing energy backup is hard to
quantify and often ignored in assessing the value of energy storage devices.
Regions where power network is not very reliable, consumers install energy
storage devices and a power converter to charge the battery while power is
available through the grid and immediately start discharging when the grid
supply is not active. Installing such devices provides uninterrupted power
supply to the user, enhancing the reliability of power supply at the consumer
end. Abrupt disconnection of power supply drastically affects the life of
certain appliances. Energy storage is also used for energy backup in developed
countries; loads like data centers and hospitals are critical and therefore,
local backup of energy is essential. Here we consider two types of backup
modes:
$\quad\bullet$ Pre-scheduled unavailability of power: Due to scheduled
maintenance power outages could occur. It is essential to consider such
incidents in case of Madeira as being an isolated power network it has less
inertia and less redundancy making it more prone to scheduled outages. In such
cases the time of the outage denoted as $i_{\text{incident}}$, is known a
priori and users can maintain the battery level above $b_{\text{set}}$, so as
in absence of grid supply users can meet its energy needs. This is represented
as an additional constraint denoted as $b_{i}=b_{i_{\text{incident}}}\geq
b_{\text{set}}$.
$\quad\bullet$ Probability of power loss: based on past failure incidents a
probability of power failure can be calculated; for instance, the chances of
power failure due to load shedding are much more probable during morning and
evening peak than any other time of the day. $\mathbb{P}_{i}$ denotes the
probability of power failure during the instant $i$. If $\mathbb{P}_{i}$ is
high than the user should maintain a greater charge level in the battery. The
co-optimization problem combined with planned and probable outages is given as
follows

$\displaystyle(P_{\text{madeira}})~{}~{}~{}$
$\displaystyle\mathmakebox[\widthof{$\underset{\displaystyle
s_{i}}{\mathrm{subject~{}to}}$}][l]{\underset{\displaystyle
s_{i}}{\mathrm{minimize}}}\quad\sum_{i=1}^{N}{p}_{\text{elec}}(i)\theta_{i}h-\lambda\mathbb{P}_{i}b_{i}$
(4a)
$\displaystyle\mathmakebox[\widthof{$\underset{\displaystyle\phantom{s_{i}}}{\mathrm{subject~{}to}}$}][c]{{\mathrm{subject~{}to}}}$
Ramping constraint, $\displaystyle Eq.~{}\ref{constraintramp},$ (4b) Capacity
constraint, $\displaystyle Eq.~{}\ref{constraintcapacity},$ (4c)
$\displaystyle\text{Self-sufficiency,}~{}~{}\theta_{i}$ $\displaystyle\geq 0,$
(4d) $\displaystyle\text{Arbitrage,}~{}~{}\theta_{i}$
$\displaystyle\geq[z_{i}+s_{i}],$ (4e) $\displaystyle\text{Peak
shaving,}~{}~{}[z_{i}+s_{i}]/h$ $\displaystyle\leq P_{\max}^{set},$ (4f)
$\displaystyle\text{Backup for probable outage,}~{}~{}b_{i}$
$\displaystyle=b_{i_{\text{incident}}}\geq b_{\text{set}}$ (4g)

where $\lambda$ is a scaling factor. Note $P_{opt}$ and $P_{\text{madeira}}$
are convex in nature as the objective function and associated constraints are
convex.

###  4.3 Open Source Codes

The co-optimization formulations and benchmarks presented in this chapter is
made open source. The link for the code is https://github.com/umar-
hashmi/MadeiraStorage.

##  5 Real-time Control under Uncertainty

The decision variables for the optimization problem $P_{\text{opt}}$ and
$P_{\text{madeira}}$ is the energy storage output $s_{i}$. The stochastic
variable in these settings is the the net-load excluding energy storage
output, $z_{i}$. We use AutoRegressive Moving Average (ARMA) forecasting for
modeling $z_{i}$. The details of the model is described in Section 5.1 and
model predictive control algorithm is described in Section 5.2.

###  5.1 Modeling Uncertainty: ARMA Forecasting

We define the mean behavior of past values of net load without storage at time
step $i$ as

${\bar{z}}_{i}=\frac{1}{D}\sum_{p=1}^{D}z_{i-pN}\quad\forall
i\in\\{k,...,N\\},k\geq 1,$ (5)

where $D$ is the number of days in the past whose values are considered in
calculating $\bar{z}$. The forecasted net load given as
$\hat{z}_{i}=\bar{z}_{i}+\hat{X}_{i}\quad\forall i\in\\{k,...,N\\},k\geq 1,$
where $\hat{X}_{i}$ represents the forecasted difference from the mean
behavior. We define $\hat{X}_{i}\quad\forall i\in\\{k,...,N\\}$ as

$\hat{X}_{\text{k}}=\alpha_{1}\hat{X}_{\text{k}-1}+\alpha_{2}\hat{X}_{\text{k}-2}+\alpha_{3}\hat{X}_{\text{k}-3}+\beta_{1}\delta_{\text{k}}^{1}+\beta_{2}\delta_{\text{k}}^{2}+\beta_{3}\delta_{\text{k}}^{3},$
(6)

where $\delta_{k}^{m}=(z_{k-mN}-{\bar{z}}_{k-mN})$ and
$\alpha_{i},\beta_{i}\forall i\in\\{1,2,3\\}$ are constant. We use the errors
in net load without storage for the past three time steps and the error in the
same time step for past three days. At time step $i=k-1$ we calculate
$\hat{X}_{k}$ as shown in Eq 6. Using $\hat{X}$ we calculate $\hat{z}$.

###  5.2 Model Predictive Control

The vector $\hat{z}$ for instants $i$ to $N$ is fed to MPC for calculating
optimal energy storage actions for time step $i$. Similar steps are done for
$i\in\\{k+1,...,N\\}$, till the end of time horizon is reached. Real-time
algorithm is presented as ForecastPlusMPC.

Inputs:
$\eta_{\text{ch}},\eta_{\text{dis}},\delta_{\max},\delta_{\min},b_{\max},b_{\min}$,
$b_{0}$, ${p}_{elec}$, $h,N,T,i=0$

1:while $i<N$ do

2: Increment $i$=$i+1$ and forecast $\hat{z}$ from time step $i$ to $N$

3: Solve for $s^{*}=P_{\text{madeira}}({p}_{elec},\hat{z},h,N,T)$

4: $b_{i}^{*}=b_{i-1}+[s^{*}]^{+}\eta_{ch}-[s^{*}]^{-}/\eta_{dis}$ and Update
$b_{0}=b_{i}^{*}$

5:end while

Algorithm 1 ForecastPlusMPC Figure 2: Receding horizon MPC with forecasting

##  6 Numerical Results

For the numerical evaluation we use a battery with initial charge level,
$b_{0}$=1kWh, $b_{\max}$=2kWh, $b_{\min}$=0.2kWh,
$\eta_{ch}$=$\eta_{dis}$=0.95. The home has 6.25 kWp solar PV. The performance
indices used for evaluating simulations are: $\quad\bullet$ Arbitrage Gains
($G_{arb}$), $\quad\bullet$ Peak shaving gains ($G_{peak}$): difference
between nominal Peak Power Contract (PPC) and the new PPC contract after
adding storage.
$\quad\bullet$ Self-sufficiency (SS): calculated using total energy consumed,
PV generation and storage output.
$\quad\bullet$ Gains per cycle: In our prior work [12] we develop a mechanism
to measure the number of cycles of operation based on depth-of-discharge of
energy storage operational cycles. We use total gains,
$G_{T}$=$G_{arb}$+$G_{peak}$, to calculate euros/cycles gained by operating
energy storage as one of the performance index. This index puts a financial
value to operational cycles of the battery.

We perform deterministic simulations for arbitrage and peak demand shaving in
Section 6.1. Section 6.2 presents numerical results for energy storage
performing backup along with arbitrage and peak shaving. Section 6.3 compares
forecast plus MPC with results for a week with respect to the deterministic
results.

###  6.1 Deterministic Solution for $P_{opt}$

We compare various pricing contracts and propose the best contract that energy
storage owners can select for maximizing their gains. Here we consider only
the daily cycles since this is the most commonly selected option in Madeira.
Two types of simulations results are presented: simulation on a shorter-time
scale, (i.e., for a day) and for longer time-scale, (i.e., for a month). The
load data is collected from a facility in Madeira. Single-rate tariff is used
as the nominal case with respect to which profit and performance improvements
are calculated.

####  6.1.1 Shorter Time Scale: A day

Simulations using 2 and 3-level ToU price for 4 different batteries are
conducted for load data of 18th May, 2018.

Figure 3: Net Load with/without solar and with/without battery for $h=0.25$
hours; here we consider 4 different batteries

From Table 2 we can conclude that 3-level ToU provides higher gains for end-
user, fast ramping battery can increase $G_{peak}$, however, $G_{arb}$
deteriorates due to greater contribution of the battery performing peak
reduction. Fig. 3 shows the variation of net load for 3-level ToU.

Table 2: Comparison for 1 day Case | $G_{arb}$ | PPC | $G_{peak}$ | SS | $G_{T}$ | euros/cyc
---|---|---|---|---|---|---
euros | kVA | euros | % | euros |
No Battery
No PV | - | 10.35 | - | - | - | -
PV | - | 6.9 | 0.144 | 32.9 | 0.144 | -
2 Level ToU with Battery
0.25C-0.25C | 0.240 | 6.9 | 0.144 | 36.2 | 0.384 | 0.568
0.5C-0.5C | 0.235 | 6.9 | 0.144 | 36.6 | 0.379 | 0.561
1C-1C | 0.235 | 5.75 | 0.192 | 36.6 | 0.427 | 0.631
2C-2C | 0.212 | 3.45 | 0.287 | 36.3 | 0.499 | 0.363
3 Level ToU with Battery
0.25C-0.25C | 0.333 | 6.9 | 0.144 | 36.2 | 0.477 | 0.704
0.5C-0.5C | 0.366 | 6.9 | 0.144 | 36.2 | 0.510 | 0.376
1C-1C | 0.351 | 5.75 | 0.192 | 36.1 | 0.543 | 0.346
2C-2C | 0.301 | 3.45 | 0.287 | 36.1 | 0.587 | 0.379
Table 3: Comparison for Longer Time Scale Case | $G_{arb}$ | PPC | $G_{peak}$ | SS | $G_{T}$ | euros/cyc
---|---|---|---|---|---|---
euros | kVA | euros | % | euros |
No Battery
No PV | - | 17.25 | - | - | - | -
PV | - | 17.25 | 0 | 34.1 | 0 | -
2 Level ToU with Battery
0.25C-0.25C | 18.38 | 13.8 | 4.27 | 34.9 | 22.65 | 0.855
0.5C-0.5C | 18.38 | 13.8 | 4.27 | 34.9 | 22.65 | 0.855
1C-1C | 18.38 | 13.8 | 4.27 | 34.9 | 22.65 | 0.855
2C-2C | 18.37 | 10.35 | 8.54 | 34.9 | 26.92 | 0.669
3 Level ToU with Battery
0.25C-0.25C | 25.01 | 13.8 | 4.27 | 34.9 | 29.28 | 1.035
0.5C-0.5C | 25.89 | 13.8 | 4.27 | 34.7 | 30.16 | 0.640
1C-1C | 26.08 | 13.8 | 4.27 | 34.7 | 30.36 | 0.589
2C-2C | 26.06 | 10.35 | 8.54 | 34.6 | 34.60 | 0.568

Integration of storage leads to approximately 3% saving of electricity bills.
Prior to installation of storage, 5.8% of solar generation was wasted and with
addition of storage the waste is reduced to zero. The simulations also shows
that additional gains are possible by reducing the PPC shown in Fig 3 and
Table 2.

####  6.1.2 Longer Time-Scale: A month

Longer time scale simulations are conducted for the month of June, 2018. Table
3 shows that energy storage does not contribute significantly towards self-
sufficiency. Fig. 4 shows the variation of SS for each day due to PV. Long
time scale simulations also indicate that 3 level ToU is more beneficial.
Further the gains per cycle indicate that energy storage is highly profitable.
A typical LiIon battery costing around euros 500/kWh could perform around 4000
cycles at 100% DoD. Thus such a battery (of 2kWh size) to be profitable should
make more than euros 0.25/cycle. Table 2 and Table 3 show that euros/cycle is
significantly more than euros 0.25.

Figure 4: Self-Sufficiency due to Solar for June 2018

###  6.2 Co-optimizing with Power Backup

The probability of power failure used for scheduling energy storage backup is
shown in Fig. 5(a). The probability of power failure on a typical day is
primarily because of load-shedding, which happens more during peak consumption
hours.

Figure 5: (a) Probability of power failure, (b) Battery charge level

We assume that there is a scheduled power outage incident at
$i_{\text{incident}}$=6 am. During this time storage should maintain an state
of charge of 80% or higher, therefore, $b_{set}=0.8b_{\max}$.

Table 4: Comparison of gains for power backup Index | Case (i) $P_{opt}$ | Case (ii) | Case (iii) $P_{madeira}$
---|---|---|---
$G_{arb}$ | 0.3006 | 0.3006 | 0.2976
$G_{peak}$ | 0.2867 | 0.2867 | 0.2867

Fig. 5(b) shows battery charge level for three cases: (i) $P_{opt}$: arbitrage
with peak shaving, (ii) Arbitrage with peak shaving and backup for probable
outage governed by failure probability in Fig. 5(a) and (iii) $P_{madeira}$:
arbitrage, peak shaving, probable outage with scheduled outage. As evident
from Fig. 5(b) the charge level for Case (iii) maintains a higher charge level
during $i_{\text{incident}}$. For case (ii) and case (iii) maintains a high
charge level during probable outage during morning and evening peak. As shown
in Table 4, the effect on gains due to performing backup is insignificant,
less than 1% in this case. For these cases self-sufficiency remains fairly
similar.

###  6.3 Real-Time Implementation (Forecast plus MPC)

The coefficients of the ARMA forecast model is tuned using regression. For
this case, we use a battery of 1C-1C type. The comparison of total load for
deterministic and MPC simulations is shown in Fig. 6. The state of charge
(SoC) of the battery is shown in Fig. 7. It can be observed that for
stochastic simulations the battery capacity is maintained at high SoC level in
order to minimize probable outage component of the objective function.

The arbitrage gains for the deterministic case for the week is euros 5.50 and
for ARMA with MPC the gains are euros 5.01. Loss of opportunity (LoO) is
defined as = $1-$(actual arbitrage gains)/(deterministic arbitrage gains). The
LoO for this numerical experiment is 8.91%. A low value of LoO indicates the
robustness of our proposed real time framework. Peak demand shaving for this
week is euros 0.336 for both deterministic and ARMA with MPC case. Although in
this numerical experiment the peak demand is in compliance with the same
contract as for the deterministic case, however, it would be advised to select
a higher level of peak demand contract pertaining to forecast errors.

Figure 6: Net load comparison for June 1 and June 2, 2018 Figure 7: Battery
Capacity for 1st week of June, 2018

##  7 Conclusion

We present a case study for the island of Madeira and formulate a convex co-
optimization problem for performing arbitrage, peak-shaving and providing
backup during power outages. Using numerical simulations we observe that
storage owners benefit more under greater volatility, evident from higher
gains electricity customers can make with 3-level ToU compared to 2-level. We
believe an increase in storage size will make more volatile consumer contracts
more beneficial. Energy storage adds to economic value while solar PV
increases self-sufficiency for scenarios where distributed generation is lower
or comparable to the magnitude of the inelastic load. For DG generating more
than inelastic load, storage also contributes to self-sufficiency by
increasing self-consumption. We show that using storage for power backup
during probable and scheduled outages do not undermine its ability to perform
arbitrage and peak demand shaving. Considering storage operational cycle
degradation, we calculated gains per cycle indicating the storage could be
financially viable in Madeira (simple payback period $\approx$3 years).
Numerical simulation for real-time control using ARMA based forecast with MPC
shows the efficacy of the proposed scheme for storage co-optimization. Further
work is required to select an optimal battery size based on historical data.

## Acknowledgments

Funding from grants: ANR under ANR-16-CE05-0008 and EU H2020 under GA 731249
are gratefully acknowledged.

## References

  * [1] M. Papadopoulos, P. Malatestas, and N. Hatziargyriou, “Simulation and analysis of small and medium size power systems containing wind turbines,” _IEEE Transactions on Power Systems_ , vol. 6, no. 4, pp. 1453–1458, 1991.
  * [2] R. Calero and J. A. Carta, “Action plan for wind energy development in the canary islands,” _Energy Policy_ , vol. 32, no. 10, pp. 1185–1197, 2004\.
  * [3] C. Bueno and J. A. Carta, “Wind powered pumped hydro storage systems, a means of increasing the penetration of renewable energy in the canary islands,” _Renewable and Sustainable Energy Reviews_ , vol. 10, no. 4, pp. 312–340, 2006.
  * [4] F. M. Vieira, P. S. Moura, and A. T. de Almeida, “Energy storage system for self-consumption of photovoltaic energy in residential zero energy buildings,” _Renewable Energy_ , vol. 103, pp. 308–320, 2017.
  * [5] R. Luthander, J. Widén, D. Nilsson, and J. Palm, “Photovoltaic self-consumption in buildings: A review,” _Applied Energy_ , vol. 142, pp. 80–94, 2015.
  * [6] R. Bhandari and I. Stadler, “Grid parity analysis of solar photovoltaic systems in germany using experience curves,” _Solar Energy_ , vol. 83, no. 9, pp. 1634–1644, 2009.
  * [7] J. Moshövel, K.-P. Kairies, D. Magnor, M. Leuthold, M. Bost, S. Gährs, E. Szczechowicz, M. Cramer, and D. U. Sauer, “Analysis of the maximal possible grid relief from pv-peak-power impacts by using storage systems for increased self-consumption,” _Applied Energy_ , vol. 137, pp. 567–575, 2015\.
  * [8] C. A. Hill, M. C. Such, D. Chen, J. Gonzalez, and W. M. Grady, “Battery energy storage for enabling integration of distributed solar power generation,” _IEEE Transactions on smart grid_ , vol. 3, no. 2, pp. 850–857, 2012.
  * [9] M. U. Hashmi, D. Muthirayan, and A. Bušić, “Effect of real-time electricity pricing on ancillary service requirements,” in _ACM E-Energy_ , 2018, pp. 550–555.
  * [10] B. Cheng and W. Powell, “Co-optimizing battery storage for the frequency regulation and energy arbitrage using multi-scale dynamic programming,” _IEEE Trans. on Smart Grid_ , 2016.
  * [11] Y. Shi, B. Xu, D. Wang, and B. Zhang, “Using battery storage for peak shaving and frequency regulation: Joint optimization for superlinear gains,” _IEEE Transactions on Power Systems_ , vol. 33, no. 3, pp. 2882–2894, 2018\.
  * [12] M. Hashmi, W. Labidi, A. Bušic, S.-E. Elayoubi, and T. Chahed, “Long-term revenue estimation for battery performing arbitrage and ancillary services,” in _IEEE SmartGridComm_ , 2018.
  * [13] M. Hashmi, A. Mukhopadhyay, A. Busic, and J. Elias, “Storage optimal control under net metering policies,” _to be submitted IEEE Transactions on Smart Grid_.
  * [14] ACIF-CCIM, Prsma, EEM, M-ITI, and Route Monkey, “Madeira Pilot Case Study Specification and Assessment,” H2020 SMILE, Funchal, Portugal, Deliverable 4.1, Oct. 2017.
  * [15] VCharge, AAU, Route Monkey, CERTH, Prsma, and RINA-C, “Most Appropriate DR Services for each Pilot,” H2020 SMILE, Bristol, UK, Deliverable 5.1, Feb. 2018.
  * [16] M. U. Hashmi, D. Deka, A. Busic, L. Pereira, and S. Backhaus, “Arbitrage with power factor correction using energy storage,” _arXiv preprint arXiv:1903.06132_ , 2019.
  * [17] Prsma, M-ITI, EEM, and ACIF-CCIM, “Data Collection, Modelling, Simulation and Decision,” H2020 SMILE, Funchal, Portugal, Deliverable 4.3, Jun. 2018.
  * [18] H.-I. Su and A. El Gamal, “Modeling and analysis of the role of energy storage for renewable integration: Power balancing,” _IEEE Transactions on Power Systems_ , vol. 28, no. 4, pp. 4109–4117, 2013.
  * [19] M. Hashmi, A. Mukhopadhyay, A. Busic, and J. Elias, “Optimal control of storage under time varying electricity prices,” in _IEEE SmartGridCom_ , 2017\.
  * [20] J. R. Cruise, R. J. Gibbens, and S. Zachary, “Optimal control of storage for arbitrage, with applications to energy systems,” in _IEEE CISS_ , 2014, pp. 1–6.


# Probing Nonadiabatic Dynamics with Attosecond Pulse Trains and Soft X-ray
Raman Spectroscopy

Lorenzo Restaino  Department of Physics, Stockholm University, Albanova
University Centre, SE-106 91 Stockholm, Sweden  Deependra Jadoun  Department
of Physics, Stockholm University, Albanova University Centre, SE-106 91
Stockholm, Sweden  Markus Kowalewski  e-mail<EMAIL_ADDRESS>Department of Physics, Stockholm University, Albanova University Centre,
SE-106 91 Stockholm, Sweden

###### Abstract

Linear off-resonant X-ray Raman techniques are capable of detecting the
ultrafast electronic coherences generated when a photoexcited wave packet
passes through a conical intersection. A hybrid femtosecond or attosecond
probe pulse is employed to excite the system and stimulate the emission of the
signal photon, where both fields are components of a hybrid pulse scheme. In
this paper, we investigate how attosecond pulse trains, as provided by high-
harmonic generation processes, perform as probe pulses in the framework of
this spectroscopic technique, instead of single Gaussian pulses. We explore
different combination schemes for the probe pulse, as well as the impact of
parameters of the pulse trains on the signals. Furthermore, we show how Raman
selection rules and symmetry consideration affect the spectroscopic signal,
and we discuss the importance of vibrational contributions to the overall
signal. We use two different model systems, representing molecules of
different symmetry, and quantum dynamics simulations to study the difference
in the spectra. The results suggest that such pulse trains are well suited to
capture the key features associated with the electronic coherence.

##  I Introduction

Conical Intersections teller1937crossing ; yarkony1998conical ;
worth2004beyond ; domcke2011conical  (CIs) represent fast, radiationless decay
channels in electronically excited molecules (see fig. 1(a)). Virtually
present in every molecular system, CIs play a key role in charge transfer
processes balzani2001electron , reaction mechanisms zimmerman1966molecular
and in the vast majority of photochemical, photophysical and photobiological
reactions klessinger1995excited ; chergui2015ultrafast ; robb2000computational
; duan2016quantum ; domcke2004conical ; domcke2012role , such as the cis/trans
isomerization of retinal rinaldi2014comparison ; polli2010conical . At such
intersections the Born-Oppenheimer approximation breaks down causing complex
dynamics of the coupled vibronic states, which can be observed
spectroscopically cederbaum1977strong ; cederbaum1981multimode ;
young2018roadmap .

Figure 1: (a) Scheme of a CI: 1D cut of potential energy surfaces along a
generic reaction coordinate, $q$. Following photoexcitation from the ground
state (GS), the wave packet reaches the CI from the Frank-Condon (FC) region.
(b) Loop diagram loopDiagr  of the off-resonant stimulated Raman signal. The
grey area represents the preparation of the system into the excited state by
means of a pump pulse, temporally well separated from the detection process.
After photoexcitation, the system propagates freely for a delay time $T$
before being probed by the hybrid-shaped pulse. (c) Schematics of the pump
($\varepsilon_{p}$) and hybrid probe pulse
($\varepsilon_{1},\varepsilon^{*}_{0}$) setup in TRUECARS.

As the photoexcited wave packet comes closer to the conical intersection, the
energy separation between the potential energy surfaces (PESs) rapidly
decreases. Thus, detection requires an unusual combination of temporal and
spectral resolution that is not available via conventional femtosecond optical
and infrared (IR) experiments polli2010conical ; horio2009probing ;
Oliver10061 ; mcfarland2014ultrafast . However, pulses in the extreme
ultraviolet (XUV) to the short X-ray spectral region possess the required
combination to directly detect the passage through a CI galbraith2017few ;
sun2020multi ; yang2018imaging ; wolf2017probing ; chang2020revealing ;
zinchenko2021sub ; nam2021conical ; neville2018ultrafast ; rohringer2019x .
Single attosecond pulses (SAPs) have been extensively used in pump-probe
experiments. Although the availability of such pulses has seen a significant
increase thanks to high-harmonic generation (HHG) ferray1988multiple ;
varju2009physics  and free electron lasers (FELs) pellegrini2016physics
sources, the generation of SAPs still requires a complex setup. The HHG
process in gases emits a sequence of short bursts of radiation, which are
coherently driven by the generation laser, where emission events occur during
each laser half cycle. Each of these short bursts is in the attosecond regime
and their interference leads to the observation of odd harmonics. Such
attosecond pulse trains (APTs) farkas1992proposal ; harris1993atomic , unlike
isolated pulses hentschel2001attosecond ; paul2001observation , are directly
available in a HHG setup, which is a table-top source now widely found in many
laboratories. Recent theoretical developments have showed the capability of
APTs to probe the electronic coherence created by the CI in the context of
time-resolved photoelectron spectroscopy jadoun2021time  and time-resolved
electron-momentum imaging Marciniak2019 .

In this paper, we demonstrate theoretically how APTs perform as probes in the
framework of the transient redistribution of ultrafast electronic coherences
in attosecond Raman signals (TRUECARS) technique TRUECARS  to probe electronic
coherences generated by a wave packet passing through a CI. We use two
different model systems that represent molecules of different symmetry and use
quantum dynamics simulations to study the difference in the spectrum.

##  II Spectroscopic Signals and Models

###  II.1 The TRUECARS signal

The TRUECARS technique uses an off-resonant stimulated X-ray Raman process
(see fig. 1(b) and (c)) which is sensitive to coherences rather than
populations. In the X-ray Raman scheme, core-hole states are involved as
intermediates rather than common valence excited states. As shown in fig.
1(b), the two pulses making up the hybrid pulse scheme, namely
$\varepsilon_{1}$ and $\varepsilon_{0}$, drive the Raman process, which is in
turn detected by a heterodyne detection scheme, where a local oscillator is
used. The frequency and time resolved signal reads in atomic units as follows:

$\displaystyle S(\omega_{R},T)=2\Im\Bigg{\\{}\int{\mathrm{d}}t\,\,$
$\displaystyle e^{i\omega_{R}(t-T)}\varepsilon^{*}_{0}(\omega_{R})$
$\displaystyle\times$
$\displaystyle\varepsilon_{1}(t-T)\langle\hat{\alpha}(t)\rangle\Bigg{\\}}$ (1)

where $T$ is the time delay between the probe field and the preparation pulse
(see fig. 1(c)), $\langle\hat{\alpha}(t)\rangle$ is the time dependent
expectation value of the X-ray transition polarizability and $\omega_{R}$ is
the Raman shift. For the details of the signal, see Ref. TRUECARS . The
dependence on the dynamics of the system enters the TRUECARS signal through
$\langle\hat{\alpha}\rangle$. Additional details on $\hat{\alpha}$ are given
in subsection C.

The TRUECARS spectrum is characterized by an oscillating pattern of gain and
loss features in the Stokes and anti-Stokes regime, that is only visible when
a vibrational or electronic coherence is present. Assuming that both
components of the probing field have the same carrier frequency, $\omega_{X}$,
the spectrum is going to be centered at a Raman shift of $\omega_{R}=0$.

###  II.2 Modelling of the Pulse Trains

To illustrate how the pulse trains were built, we start from their definition
in the time domain:

$E_{APT}(t)=G(t)P(t)$ (2)

where $G(t)$ is a Gaussian envelope of $\sigma_{env}$ width

$\displaystyle G(t)=e^{-t^{2}/2\sigma_{env}^{2}}$ (3)

and $P(t)$ is an infinite train of pulses. The electric field of the single
pulses inside the train is defined as

$\displaystyle
E_{SAP}(t)=e^{-\left(t-\dfrac{\tau}{2}\right)^{2}/2\sigma_{SAP}^{2}}\cos\Big{[}{\omega_{X}\left(t-\dfrac{\tau}{2}\right)}\Big{]}$
(4)

with $\omega_{X}$ being the center frequency and $\tau=2\pi/\omega_{IR}$ the
period of the IR field. The expression for the electric field then reads:

$\displaystyle E_{APT}(t)=e^{-t^{2}/2\sigma_{env}^{2}}$
$\displaystyle\times\sum_{n}e^{-\left(t-n\dfrac{\tau}{2}\right)^{2}/2\sigma_{SAP}^{2}}\cos\Big{[}{\omega_{X}\left(t-n\dfrac{\tau}{2}\right)+n\pi}\Big{]}\,.$
(5)

The pulse train employed in the calculations was built according to eq. (5) by
substitution of $t$ with $t-T$ and with the following parameters:
$\sigma_{env}=2.5$ fs and $\sigma_{SAP}=0.15$ fs. For the purpose of
simulating TRUECARS spectra, the center frequency $\omega_{X}$ can assume any
arbitrary value.

To ease nomenclature, henceforth, single pulses in the femtosecond or
attosecond regime will be broadly referred to as “Gaussian pulses”. A snapshot
of the train pulse at $\omega_{IR}=1.55$ eV is displayed in fig. 2. Snapshots
at different frequencies of the IR laser are available in the SI.

Figure 2: Snapshot of the modeled APT shape at $\omega_{IR}=1.55$ eV. On the
left: envelope of the pulse train in time domain. On the right: pulse train in
frequency domain, shifted with respect to a selected central harmonic.

###  II.3 Models and Symmetry

We use group theory to identify Raman active transitions and the vanishing
integral rule to predict whether the polarizability matrix elements such as
$\langle\Psi^{\prime\prime}|\hat{\alpha}|\Psi^{\prime}\rangle$ will be zero.
In particular, they will vanish if the product of the irreducible
representations of the two relevant states and the operator does not contain
the totally symmetric representation, that is

$\Gamma(\Psi^{\prime\prime})\otimes\Gamma(\hat{O})\otimes\Gamma(\Psi^{\prime})\not\supset\Gamma^{(s)}\,.$
(6)

The two systems studied in this work belong to the $C_{i}$ and to the $C_{s}$
point groups, respectively. According to their character tables, the
polarizability tensor elements $\alpha_{xx}$, $\alpha_{yy}$ and $\alpha_{xy}$
all transform as the totally symmetric irreducible representation, i.e.
$A_{g}$ in $C_{i}$ and $A^{\prime}$ in $C_{s}$. This was taken into account
while modelling the X-ray polarizability operator, $\hat{\alpha}$, that was
here approximated to a 2x2 matrix in the electronic subspace:

$\hat{\alpha}=\begin{pmatrix}\alpha_{11}&\alpha_{12}\\\
\alpha_{21}&\alpha_{22}\end{pmatrix}\,.$ (7)

The full electronic polarizability matrix must be taken into account in order
to be basis-independent (adiabatic vs. diabatic states) and the diagonal
matrix elements cannot be neglected when transforming between representations.
The diabatic wave function, $\Psi(\bm{\mathbf{q}},t)$, is expressed in terms
of the electronic states. For a system with two electronic states, the wave
function reads:

$\Psi(\bm{\mathbf{q}},t)=\begin{pmatrix}\phi_{1}(\bm{\mathbf{q}},t)\\\
\phi_{2}(\bm{\mathbf{q}},t)\end{pmatrix}\,.$ (8)

Adiabatic and diabatic states are related by an unitary transformation, where
the advantage of the diabatic basis is the absence of derivative couplings
kowalewski2017simulating . Expanding the expectation value,
$\langle\hat{\alpha}\rangle$, yields:

$\displaystyle\langle\hat{\alpha}\rangle=$
$\displaystyle\langle\Psi|\hat{\alpha}|\Psi\rangle$ $\displaystyle=$
$\displaystyle\langle\phi_{1}|\alpha_{11}|\phi_{1}\rangle+\langle\phi_{2}|\alpha_{22}|\phi_{2}\rangle+2\Re\langle\phi_{1}|\alpha_{12}|\phi_{2}\rangle$
(9)

According to eq. (6), for a Raman active transition two conditions need to be
fulfilled here: i) the diagonal and off-diagonal polarizability matrix
elements must all transform as the totally symmetric irreducible
representation of their point group ($A_{g}$ in $C_{i}$ and $A^{\prime}$ in
$C_{s}$, respectively); ii) the electronic states $\phi_{1}$ and $\phi_{2}$
must have the same symmetry label, or else the integral
$2\Re\langle\phi_{1}|\alpha_{12}|\phi_{2}\rangle$ will vanish. The diagonal,
($\langle\phi_{1}|\alpha_{11}|\phi_{1}\rangle+\langle\phi_{2}|\alpha_{22}|\phi_{2}\rangle$),
and off-diagonal, ($2\Re\langle\phi_{1}|\alpha_{12}|\phi_{2}\rangle$), terms
refer to vibrational and electronic contributions, respectively. Each
vibrational normal mode has an associated irreducible representation, too.
Therefore, the diagonal contribution will be suppressed if the vibrational
modes transform as the wrong irreducible representation. In practice, we can
ensure the transition will be Raman active if $\phi_{1}$, $\phi_{2}$ and
$\alpha_{ij}$ all fall into $A_{g}$ for $C_{i}$ and into $A^{\prime}$ for
$C_{s}$.

The concept of point groups is based on the approximation of a rigid molecule,
that is considering the molecule as a rigid skeleton of nuclei. However, when
large amplitude motions are considered where the symmetry of the system is not
conserved, the concept is inadequate. Photoisomerization processes are well
known examples of large amplitude nuclear motion in molecules. To generate
symmetry adapted polarizability element functions in the diabatic basis, we
now introduce the complete nuclear permutation inversion (CNPI) group CNPI . A
CNPI group consists of all permutations of identical nuclei, the inversion of
all nuclear and electronic coordinates, ($E^{*}$), as well as their products.
The inversion $E^{*}$ differs from the inversion operation $i$ in point
groups. The former is an operation which results in a sign change of all
nuclear and electronic coordinates in the space-fixed coordinate system. A
similar application to electronic states and transition dipole moments can be
found in Ref. obaid2015molecular . Following these sections, details of the
modeling of the two systems as well as the polarizability matrix are shown.

####  II.3.1 $C_{i}$ Symmetry Model

This model consists of two harmonic potential wells shifted with respect to
each other (see left panel of fig. 3). It represents two excited electronic
states with a conical intersection. An example of such a molecule could be
benzene (C6H6) and its photochemistry palmer1993mc  or acetylene (C2H2)
KOPPEL2008319 .

Figure 3: Potential energy surfaces of the investigated model systems. On the
left: 1D cut of the diabatic potential energy surface along $q_{a}$, in the
$C_{i}$ system. On the right: 1D cut of the diabatic potential energy surface
along $q_{b}$, in the $C_{s}$ system. The wave packet has been previously
excited from the ground state to the Franck-Condon (FC) region on $S_{2}$ by
means of a pump pulse.

The diagonal and off-diagonal elements of the polarizability have been shaped
by a function of $q_{1}$ and $q_{2}$, which is symmetric with respect to the
two normal modes. In the framework of CNPI this translates as

$E^{*}\alpha_{ij}(q_{1},q_{2})=E^{*}\alpha_{ij}(-q_{1},-q_{2})=E^{*}\alpha_{ij}(q_{1},q_{2})\,.$

This function behaves non-linearly with respect to the nuclear coordinates in
the proximity of the CI. In this specific case, both normal modes $q_{1}$ and
$q_{2}$ transform as the totally symmetric irreducible representation $A_{g}$,
thus they are all Raman active. Plots of the symmetry adapted functions of the
polarizability matrix elements as well as the single contributions to
$\langle\hat{\alpha}\rangle$ are available in the SI.

####  II.3.2 $C_{s}$ Symmetry Model

The $C_{s}$ model consists of a Morse-like potential well and a repulsive
potential with asymptotic behavior at large values of $q$, with a conical
intersection. A 1D cut of the system is shown in the right panel of fig. 3.
Hydroxylamine (NH2OH), which has $C_{s}$ symmetry, is a well-known example for
the study of the effects of conical intersections in photodissociation
RevModPhys.68.985 . Similarly to the $C_{i}$ symmetry, the diagonal and off-
diagonal elements of the polarizability have been here shaped by a function of
$q_{1}$ and $q_{2}$. This time, the function is symmetric with respect to
$q_{2}$ and non-symmetric to $q_{1}$. A fundamental transition is Raman active
if a normal mode forms a basis for one or more components of the
polarizability. Here, the normal mode $q_{1}$ transforms as the irreducible
representation $A^{\prime\prime}$, while $q_{2}$ as $A^{\prime}$. Thus,
$q_{2}$ is Raman active. Further details, including contour plots of the
aforementioned functions, are available in the SI.

##  III Computational Details

The time evolution was simulated by solving the time-dependent (non-
relativistic) Schrödinger equation numerically with the Fourier method tannor
, where the wave function is represented on an equally spaced grid of sampling
points in coordinate space, using the in-house software QDng. In the diabatic
picture, the two-level Hamiltonian reads as:

$\hat{H}=\begin{pmatrix}\hat{V}_{1}+\hat{T}&\hat{V}_{12}\\\
\hat{V}_{12}&\hat{V}_{2}+\hat{T}\end{pmatrix}\,,$ (10)

where $\hat{T}$ is the kinetic energy operator

$\hat{T}=-\frac{\hbar^{2}}{2\mu}\sum_{i}\nabla_{q_{i}}^{2}$ (11)

while $\hat{V}_{1}$ and $\hat{V}_{2}$ are the potential energy operators,
respectively for $S_{1}$ and $S_{2}$, and $\hat{V}_{12}$ is the diabatic
coupling operator. The Arnoldi scheme arnoldi1951principle ;
friesner1989method ; saad1980variations  was employed as the propagator for
all calculations. The reduced masses, $\mu$, the time steps, $\Delta t$, and
the number of grid points employed have been all summarized in table 1. A
$\Delta t$ of 4 is equal to 96.75 as, while a reduced mass of 18000 is
approximately 10 amu. For the initial state, assumed to be the result of a
short-fs excitation pulse, the nuclear wave packet was approximated by a
Gaussian envelope on $S_{2}$. Lastly, a perfectly matched layer nissen10jcp
was employed to absorb the wave packets at the boundary and to account for the
dissociative behaviour.

|  |  $\mu$ (au) |  $\Delta t$ (au) | Grid p.
---|---|---|---|---
$C_{i}$ |  | 18000 | 4 | 256x256
$C_{s}$ |  | 30000 | 2 | 300x300
Table 1: Reduced masses, time steps and number of grid points used in the
simulations.

The propagated wave packets were used for the evaluation of the expectation
value of the polarizability, $\langle\hat{\alpha}\rangle$. Finally, the
TRUECARS spectra were calculated with eq. (II.1) at different values of the
pump-probe delay, $T$.

##  IV Results

For each system we have simulated two spectra with the conventional TRUECARS
probe scheme, containing purely electronic or vibrational signals, as
displayed in fig. 4. This was achieved by simulating the vibrational
contribution
($\langle\phi_{1}|\alpha_{11}|\phi_{1}\rangle+\langle\phi_{2}|\alpha_{22}|\phi_{2}\rangle$)
and the electronic contribution
($2\Re\langle\phi_{1}|\alpha_{12}|\phi_{2}\rangle$) separately. This
distinction will help us break down the interesting characteristics of the
spectra and highlight the differences and similarities between the probe
schemes. It should be emphasized that the distinction of vibrational and
electronic degrees of freedom near a conical intersection is fictitious, as
these degrees of freedom are mixed.

Signals in fig. 4 have been normalized with respect to the maximum value of
the vibrational contribution in the $C_{s}$ model (panel (d)) with the
following ratios: (a) $2\cdot 10^{-4}$, (b) 0.78 and (c) 0.33. By comparing
the overall intensities between the two models, we notice that the electronic
contribution in the $C_{s}$ is higher than the $C_{i}$ system. The vibrational
signals are, in both model systems, contained within the $0\leq\omega_{R}\leq
1.5$ eV range and immediately visible in the spectrum due to vibrational
coherences ($\langle\phi_{2}|\alpha_{22}|\phi_{2}\rangle$), whereas the
electronic signals spread over a broader energy range ($0\leq\omega_{R}\leq 2$
eV) and only appear after the wave packet has reached the CI. Hence, if the
energy resolution is not sufficient, the two will overlap with each other and
the electronic component may be masked by the stronger vibrational signal. A
major difference between the two components can be seen in the temporal
oscillations, where the electronic contribution to the signal oscillates much
faster than the vibrational one.

Figure 4: Comparison of electronic and vibrational contributions to the
TRUECARS signals for Gaussian pulses. (a): electronic contribution in the
$C_{i}$ model; (b): vibrational contribution in the $C_{i}$ model; (c):
electronic contribution in the $C_{s}$ model; (d): vibrational contribution in
the $C_{s}$ model. Spectra in the same row share the same Raman shift axis.
The signals have been normalized with respect to the maximum value of (d).

In contrast to the conventional hybrid probe, here pulse trains have been
employed instead of single Gaussian pulses. The following three combination
schemes have been investigated:

  * $(I)$

An APT as $\varepsilon_{1}$ and a Gaussian pulse as $\varepsilon^{*}_{0}$;

  * $(II)$

A Gaussian pulse as $\varepsilon_{1}$ and an APT as $\varepsilon^{*}_{0}$;

  * $(III)$

Two identical APTs as both fields.

In the following, we use $\sigma_{0}$ and $\sigma_{1}$ to indicate the
Gaussian width of $\varepsilon^{*}_{0}$ and $\varepsilon_{1}$, respectively.
Due to the nature of the pulse trains, multiple signals are expected in the
spectra simulated with schemes $(I)$ and $(III)$. The extra signals arise from
the side peaks of the train and are expected to be symmetric to each other,
but they are lower in intensity with respect to the central signal at
$\omega_{R}=0$. In scheme $(II)$, a Gaussian pulse in the femtosecond
timescale is employed as $\varepsilon_{1}$ and the spectrum will mostly
consist of one central signal. This is due to the limited spectral resolution
of the femtosecond narrow band pulse. Spectra similar to fig. 4 but simulated
with probe scheme $(III)$, are given in the SI.

Calculations were carried out for different values of $\omega_{IR}$ at 1.55,
0.99 and 0.83 eV (i.e., $\lambda_{IR}=$ 800, 1250, 1500 nm). The IR laser
frequency is an important parameter in these calculations because it directly
shapes the attosecond pulse train via the IR laser period, $\tau$. A smaller
$\tau$ implies more peaks in the time domain and, equivalently, less peaks in
the frequency domain. Similarly, the $\sigma_{env}$ parameter in eq. (5) can
achieve the same effect.

###  IV.1 $C_{s}$ Symmetry Model

We begin our discussion of the results starting with the $C_{s}$ model. The
time evolution of the population of the excited states is plotted in the top
panel of fig. 5. Following photoexcitation, the wave packet reaches the CI in
$\approx$12 fs with an overall population transfer of $\approx$45%. The
electronic coherence reaches a maximum of 0.15 at 15 fs, after which starts
decaying. The simulated APT TRUECARS spectra are shown in fig. 5 (c), (d) and
(e) compared to a standard single pulse TRUECARS (b). The dashed black line in
the spectra indicates the expectation value of the energy separation, $\Delta
E_{21}(t)$, between $S_{1}$ and $S_{2}$. For more details, see the S.I. of
Ref. TRUECARS . Because the spectrum is symmetric with respect to
$\omega_{R}=0$, only signals within $\omega_{R}\in[0,2.5]$ eV are shown. Extra
signals appear above 2.5 eV with schemes $(I)$ and $(III)$, however, those
only carry redundant information, as they are lower-intensity replicas of the
central peak. Initially, the vibrational coherence is the only contribution
and it is constrained within $0\leq\omega_{R}\leq 1.5$ eV, as shown in fig. 4.
Once the wave packet is in the proximity of the CI, the electronic coherence
starts to build up and becomes visible in the $0\leq\omega_{R}\leq 2$ eV
region of the spectrum. As the energy separation between the states increases
again, the oscillating pattern of gain and loss features in the Stokes and
anti-Stokes regime can be seen in the spectrum. The oscillation period
directly mirrors the coherence period: as $\Delta E_{21}(t)$ grows, the
oscillations speed up which causes the positions of the peaks in the Raman
shift $\omega_{R}$ to spread apart. Due to the shape of the potential energy
surfaces, the energy separation between $S_{2}$ and $S_{1}$ stays
approximately constant after the CI. This can be seen from the dashed black
lines in fig. 5 as well as the positions of the peaks in the Raman shift
$\omega_{R}$.

Figure 5: Comparison of pulse schemes for the $C_{s}$ model system. (a): Time
evolution of the population and the coherence of the two excited states in the
diabatic basis. The population transfer is around 45$\%$ and occurs after
about 15 fs of the wave packet propagation. The diabatic coupling is
responsible for the slight oscillation in the population between 15 and 20 fs;
(b) TRUECARS spectrum generated by an isolated Gaussian hybrid
femtosecond/attosecond probe-pulse sequence (pulse parameters $\sigma_{1}=1.5$
and $\sigma_{0}=0.15$ fs). (c) TRUECARS spectrum obtained via combination
$(I)$ with Gaussian pulse parameter $\sigma_{0}=0.15$ fs; (d) TRUECARS
spectrum obtained via combination $(II)$ with Gaussian pulse parameter
$\sigma_{1}=1.5$ fs; (e) TRUECARS spectrum obtained via combination of two
identical APTs (scheme $(III)$). Each signal has been normalized with respect
to its maximum value. The dashed black line represents the average time-
dependent separation of the adiabatic potential energy surfaces. All spectra
are simulated for $\omega_{IR}=1.55$ eV. A snapshot of the pulse train can be
found in fig. 2.

Among the three APT probe combinations displayed in fig. 5 only (c) is able to
capture the CI signature at $\omega_{IR}=1.55$ eV. By comparing panels 5(b)
and (c), we note that scheme $(I)$ can achieve similar energy resolution to
the Gaussian/Gaussian hybrid probe, while being characterized by the presence
of additional signals in the spectrum. However, probe scheme $(I)$ still
requires a single pulse in the attosecond time scale.

As the IR laser frequency decreases from 1.55 to 0.83 eV, the pulse train
peaks get closer to each other in the spectral domain and the interesting
electronic coherence fingerprint becomes visible in the spectra simulated with
scheme $(II)$, as displayed in fig. 6. This probe scheme does not require a
single pulse in the attosecond regime. Moreover, it allows to achieve a better
spectral resolution of the signals in the Raman shift than a standard TRUECARS
(fig. 5(b)). In fact, the vibrational and electronic contributions are now
sufficiently separated from each other to allow for an unambiguous assignment.
The latter is not captured by the central harmonic of the APT but by its first
side peak. This hypothesis was supported by simulations carried out varying
the $\sigma_{env}$ parameter, at the same $\omega_{IR}$ of fig. 6. The width
of the Gaussian envelope directly shapes the width of the harmonics in the HHG
spectrum. By increasing and decreasing the $\sigma_{env}$ of the pulse train
we noticed a corresponding increase and decrease in the separation between
harmonics, and therefore between the vibrational and electronic contribution
in the TRUECARS signal. The oscillation in fig. 6 appears to be slightly
shifted in the Raman shift and does not follow the dashed black line anymore.
Nevertheless, the oscillation period is the same and can still be used to
obtain information on CI. According to the calculations, for this specific
system an IR laser frequency of 0.83 eV appears to be the most suited to
resolve the time-dependent energy gap between the two PESs, as shown in fig.
6(a).

Figure 6: Simulated TRUECARS spectra for increasing values of the generating
IR laser frequency with probe scheme $(II)$, for the $C_{s}$ system. (a):
$\omega_{IR}=0.83$ eV; (b): $\omega_{IR}=0.99$ eV; (c): $\omega_{IR}=1.55$ eV.
The Gaussian pulse employed as $\varepsilon_{1}$ has a width of
$\sigma_{1}=1.5$ fs. The dashed black line represents the average time-
dependent separation of the adiabatic potential energy surfaces.

We can extract additional insight from the oscillation period by analyzing the
Wigner distribution keefer2020visualizing ; cavaletto2021high ;
cavaletto2021unveiling  or temporal gating spectrogram of the analytic signal
boashash1988note , which display a time-frequency map. The Wigner distribution
is defined as

$W(T,\omega)=\int_{-\infty}^{\infty}{\mathrm{d}}\tau\,\,S_{a}\left(T+\dfrac{\tau}{2}\right)S_{a}^{*}\left(T-\dfrac{\tau}{2}\right)\,e^{i\omega\tau}\;$
(12)

where $S_{a}(T)$ is the so-called analytic signal whose imaginary part is
related to the original signal, $S(T)$, by Hilbert transformation:

$S_{a}(T)=S(T)+\dfrac{i}{\pi}\bigintssss_{-\infty}^{\infty}{\mathrm{d}}s\,\,\dfrac{S(T-s)}{s}$
(13)

Here, $S(T)\equiv S(T;\omega_{R})$ is a temporal slice of the signal in fig.
6(a) at a selected Raman shift. The Wigner distribution is a quadratic
functional of the signal and so it will, in general, show interference between
the negative and positive frequency components of the signal. However, when
the analytic signal is used in the computation, no negative frequencies are
present, hence no interference will survive in the spectrogram. Figure 7(b)
and (d) show the modulus $|W(T,\omega)|$ for signal traces taken at
$\omega_{R}=1.56$ and 0.27 eV, which are interpreted as electronic and
vibrational contributions respectively. Panels 7(a) and (c) capture the
different temporal oscillations of the vibrational and electronic components
of the TRUECARS signal. The electronic Wigner distribution spans a broader
frequency window than the vibrational contribution. Such window represents the
PESs splitting in proximity of the CI. More strikingly, the energy splitting
of the main electronic feature is time-dependent and starts around 0.25 eV at
10 fs and converges to 1 eV at $\approx 25\,$fs. This is in good agreement
with the indicated splitting (black dashed line) in figs. 5 and 6.

Figure 7: Comparison between Wigner distributions of selected traces of the
TRUECARS signal. (a) and (c): Signal traces $S(T;\omega_{R})$ at
$\omega_{R}=1.56$ and 0.27 eV, respectively. (b) and (d): Normalized Wigner
spectrograms of (a) and (c). The frequency information obtained from the
spectrogram is extracted from the temporal features of the signal.

With probe scheme $(III)$ (fig. 5(e)) the characteristic features caused by
electronic coherence are concealed in the spectrum. This happens because the
small oscillation, traceable to the electronic component of the TRUECARS
signal, overlaps with the dominant vibrational contribution in the same
region.

###  IV.2 $C_{i}$ Symmetry Model

The time evolution of the population of the two excited states, $S_{2}$ and
$S_{1}$, as well as the electronic coherence magnitude, is displayed in fig.
8(a) for the $C_{i}$ model system. The wave packet takes about 15 fs to reach
the conical intersection, resulting in an overall population transfer of
$\sim$45%. The electronic coherence has a maximum of 0.0015 around 15 fs,
after which starts decaying because of the increasing energy splitting between
the surfaces.

Figure 8: Comparison of pulse schemes for the $C_{i}$ model system. (a): Time
evolution of the populations and coherence for states $S_{1}$ and $S_{2}$ in
the diabatic basis. The coherence magnitude (in green) has been magnified by a
factor of 200 for visual purposes. The population transfer ($\sim 45\%$)
occurs after about 15 fs when the wave packet reaches the CI. (b): TRUECARS
spectrum generated by an isolated Gaussian hybrid femtosecond/attosecond
probe-pulse sequence (pulse parameters $\sigma_{1}=1.5$ and $\sigma_{0}=0.15$
fs). (c): TRUECARS spectrum obtained via combination $(I)$ with Gaussian pulse
parameter $\sigma_{0}=0.15$ fs; (d): TRUECARS spectrum obtained via
combination $(II)$ with Gaussian pulse parameter $\sigma_{1}=1.5$ fs; (e):
TRUECARS spectrum obtained via combination $(III)$. The dashed black line
represents the average time-dependent separation of the adiabatic potential
energy surfaces. All are spectra are calculated at $\omega_{IR}=1.55$ eV. Each
signal has been normalized with respect to its maximum value. A snapshot of
the pulse train is displayed in fig. 2.

Following the photoexcitation on $S_{2}$, as the original wave packet
approaches the CI, the wave packet transferred on $S_{1}$ will inherit an odd
symmetry from the diabatic couplings. This generates a very weak electronic
coherence ($10^{-3}$ order of magnitude) because the integral
$2\Re\langle\phi_{1}|\phi_{2}\rangle$ gets very small. Hence, the vibronic
coherence, which TRUECARS is sensible to via the physical observable
$\langle\hat{\alpha}\rangle$, has now a dominant vibrational component that
totally conceal the interesting and characteristic features of the CI, such as
the time-dependent energy splitting. This is why we are not able to observe
them with the TRUECARS technique for this model system, not even with the
standard single Gaussian pulse probe scheme of fig. 8(b). Decreasing the IR
frequency from $\omega_{IR}=1.55$ to $0.83$ eV does not produce any
significant change. The reason why this does not occur in the $C_{s}$ model is
due to the small shift of $S_{2}$ in $q_{2}$, breaking the symmetry and making
the integral larger.

##  V Conclusions

In this paper, we tested the suitability of attosecond pulse trains as probes
for detecting electronic coherence generated at a conical intersection with
the TRUECARS technique. Model systems of different symmetry and two nuclear
degrees of freedom were used. The polarizability matrices were modelled to
obtain Raman active vibrational modes. The inclusion of the diagonal
polarizability matrix elements includes vibrational coherences, which are
inevitably a part of the signal. The full vibronic polarizability matrix must
be taken into account in order to be basis-independent and the diagonal matrix
elements cannot be neglected when transforming between adiabatic and diabatic
representations. Although this leads to a more complex signal, we could show
that it is still possible to distinguish the fast oscillating electronic
feature from the vibrational contribution by analyzing the temporal gating
spectrogram.

To gain additional insight into the TRUECARS signals, spectra originated from
purely electronic or purely vibrational contributions were simulated by
evaluating the off-diagonal and diagonal contribution separately in the time-
dependent expectation value of the polarizability operator. Due to symmetry
and Raman selection rules, the electronic coherence appears to vanish in the
$C_{i}$ model system and concealed by the much stronger vibrational
contributions.

We have explored three different schemes for the probe pulse and discussed
their features in comparison to the conventional TRUECARS scheme. We found
that, among the schemes reviewed, the combination of a Gaussian pulse as a
narrowband pulse and an attosecond pulse train as the broad band pulse proved
to be the most suitable for our purposes, offering two main advantages: first,
a more clear separation between the electronic and vibrational components of
the TRUECARS signal can be achieved by fine-tuning the IR laser frequency. The
convolution of harmonics leads to a shift of the peaks in the spectral domain,
but the energy separation can still be read off the oscillation period in the
time domain. This was corroborated by the analysis of Wigner spectrograms,
calculated for selected temporal traces of the signal. Second, single Gaussian
pulses in the attosecond timescale are no longer required. Furthermore, this
is the logical choice for the hybrid probe used in TRUECARS, since the
Gaussian pulse and the APT represent a narrowband and a broad band pulse,
respectively. The calculations showed the best results at $\omega_{IR}=0.83$
eV.

The probe scheme with an APT and a Gaussian pulse, employed as narrowband
pulse and broad band pulse respectively, also proved to be suitable to capture
the CI fingerprints. Whereas the traditional TRUECARS spectrum is composed of
one main signal centered at $\omega_{R}=0$, multiple bands are now visible.
Such bands are replicas of the main central signal with scaled intensities,
and are visible along the harmonic comb where the side peaks of the pulse
train appear. When compared, the two spectra contain similar information about
the detected vibronic coherence. Nevertheless, such combination for the probe
still requires a single attosecond pulse acting as a broadband pulse.

The use of two identical pulse trains did not reveal most features related to
the electronic coherence in the simulated spectrum. Due to the insufficient
resolution, the electronic component overlaps with the dominant vibrational
signal, resulting in a concealment of the electronic contribution.

## Supplementary Material

The supplementary material provides: i) snapshots of the modeled attosecond
pulse train shape at selected frequencies of the generating IR laser; ii)
details about modeling of the material properties as well as the potential
energy surfaces and the diabatic couplings used for the $C_{i}$ and $C_{s}$
systems; iii) TRUECARS spectra displaying the entire frequency comb of the
attosecond pulse train.

###### Acknowledgements.

This project has received funding from the European Union’s Horizon 2020
research and innovation program under the Marie Skłodowska-Curie grant
agreement No. 860553 and the Swedish Research Council (Grant No. VR
2018-05346).

## Data Availability

The data that support the findings of this study are available from the
corresponding author upon reasonable request.

## References

  * (1) E. Teller.  The crossing of potential surfaces.  J. Phys. Chem., 41(1):109–116, 1937.
  * (2) D. R. Yarkony.  Conical intersections: Diabolical and often misunderstood.  Acc. Chem. Res., 31(8):511–518, 1998.
  * (3) G. A. Worth and L. S. Cederbaum.  Beyond Born-Oppenheimer: molecular dynamics through a conical intersection.  Annu. Rev. Phys. Chem., 55:127–158, 2004.
  * (4) W. Domcke, D. R. Yarkony, and H. Köppel.  Conical intersections: theory, computation and experiment, volume 17.  World Scientific, 2011.
  * (5) V. Balzani, P. Piotrowiak, M. A. J. Rodgers, J. Mattay, D. Astruc, et al.  Electron transfer in chemistry, volume 1.  Wiley-VCH Weinheim, 2001.
  * (6) H. E. Zimmerman.  Molecular orbital correlation diagrams, Mobius systems, and factors controlling ground-and excited-state reactions. ii.  J. Am. Chem. Soc., 88(7):1566–1567, 1966.
  * (7) M. Klessinger and J. Michl.  Excited states and photochemistry of organic molecules.  VCH publishers, 1995.
  * (8) M. Chergui.  Ultrafast photophysics of transition metal complexes.  Acc. Chem. Res., 48(3):801–808, 2015.
  * (9) M. A. Robb, M. Garavelli, M. Olivucci, and F. Bernardi.  A computational strategy for organic photochemistry.  Rev. Comput. Chem., 15:87–146, 2000.
  * (10) H.-G. Duan and M. Thorwart.  Quantum mechanical wave packet dynamics at a conical intersection with strong vibrational dissipation.  J. Phys. Chem. Lett., 7(3):382–386, 2016.
  * (11) W. Domcke, D. Yarkony, and H. Köppel.  Conical intersections: electronic structure, dynamics $\&$ spectroscopy, volume 15.  World Scientific, 2004.
  * (12) W. Domcke and D. R. Yarkony.  Role of conical intersections in molecular spectroscopy and photoinduced chemical dynamics.  Annu. Rev. Phys. Chem., 63:325–352, 2012.
  * (13) S. Rinaldi, F. Melaccio, S. Gozem, F. Fanelli, and M. Olivucci.  Comparison of the isomerization mechanisms of human melanopsin and invertebrate and vertebrate rhodopsins.  Proc. Natl. Acad. Sci. U.S.A., 111(5):1714–1719, 2014.
  * (14) D. Polli, P. Altoe, O. Weingart, K. M. Spillane, C. Manzoni, D. Brida, G. Tomasello, G. Orlandi, P. Kukura, R. A. Mathies, et al.  Conical intersection dynamics of the primary photoisomerization event in vision.  Nature, 467(7314):440–443, 2010.
  * (15) L. S. Cederbaum, W. Domcke, H. Köppel, and W. Von Niessen.  Strong vibronic coupling effects in ionization spectra: The “mystery band” of butatriene.  Chem. Phys., 26(2):169–177, 1977.
  * (16) L. S. Cederbaum, H. Köppel, and W. Domcke.  Multimode vibronic coupling effects in molecules.  Int. J. Quantum Chem., 20(S15):251–267, 1981.
  * (17) L. Young, K. Ueda, M. Gühr, P. H. Bucksbaum, M. Simon, S. Mukamel, N. Rohringer, K. C. Prince, C. Masciovecchio, M. Meyer, et al.  Roadmap of ultrafast X-ray atomic and molecular physics.  J. Phys. B, 51(3):032003, 2018.
  * (18) S. Rahav and S. Mukamel.  Ultrafast nonlinear optical signals viewed from the molecules perspective: Kramers-Heisenberg transition amplitudes vs. susceptibilities.  Adv. At. Mol. Opt. Phys., 59:223, 2010.
  * (19) T. Horio, T. Fuji, Y. Suzuki, and T. Suzuki.  Probing ultrafast internal conversion through conical intersection via time-energy map of photoelectron angular anisotropy.  J. Am. Chem. Soc., 131(30):10392–10393, 2009.
  * (20) T. A. A. Oliver, N. H. C. Lewis, and G. R. Fleming.  Correlating the motion of electrons and nuclei with two-dimensional electronic–vibrational spectroscopy.  Proc. Natl. Acad. Sci. U.S.A., 111(28):10061–10066, 2014.
  * (21) B. K. McFarland, J. P. Farrell, S. Miyabe, F. Tarantelli, A. Aguilar, N. Berrah, C. Bostedt, J. D. Bozek, P. H. Bucksbaum, J. C. Castagna, et al.  Ultrafast X-ray Auger probing of photoexcited molecular dynamics.  Nat. Commun., 5(1):1–7, 2014.
  * (22) M. C. Galbraith, S. Scheit, N. V. Golubev, G. Reitsma, N. Zhavoronkov, V. Despré, F. Lepine, A. I. Kuleff, M. J. J. Vrakking, O. Kornilov, et al.  Few-femtosecond passage of conical intersections in the benzene cation.  Nat. Commun., 8(1):1–7, 2017.
  * (23) K. Sun, W. Xie, L. Chen, W. Domcke, and M. F. Gelin.  Multi-faceted spectroscopic mapping of ultrafast nonadiabatic dynamics near conical intersections: A computational study.  J. Chem. Phys., 153(17):174111, 2020.
  * (24) J. Yang, X. Zhu, T. J. A. Wolf, Z. Li, J. P. F. Nunes, R. Coffee, J. P. Cryan, M. Gühr, K. Hegazy, T. F. Heinz, et al.  Imaging CF3I conical intersection and photodissociation dynamics with ultrafast electron diffraction.  Science, 361(6397):64–67, 2018.
  * (25) T. J. A. Wolf, R. H. Myhre, J. P. Cryan, S. Coriani, R. J. Squibb, A. Battistoni, N. Berrah, C. Bostedt, P. Bucksbaum, G. Coslovich, et al.  Probing ultrafast $\pi\pi$*/n$\pi$* internal conversion in organic chromophores via K-edge resonant absorption.  Nat. Commun., 8(1):1–7, 2017.
  * (26) K. F. Chang, M. Reduzzi, H. Wang, S. M. Poullain, Y. Kobayashi, L. Barreau, D. Prendergast, D. M. Neumark, and S. R. Leone.  Revealing electronic state-switching at conical intersections in alkyl iodides by ultrafast XUV transient absorption spectroscopy.  Nat. Commun., 11(1):1–7, 2020.
  * (27) K. S. Zinchenko, F. Ardana-Lamas, I. Seidu, S. P. Neville, J. van der Veen, V. U. Lanfaloni, M. S. Schuurman, and H. J. Wörner.  Sub-7-femtosecond conical-intersection dynamics probed at the carbon K-edge.  Science, 371(6528):489–494, 2021.
  * (28) Y. Nam, D. Keefer, A.r Nenov, I. Conti, F. Aleotti, F. Segatta, J. Y. Lee, M. Garavelli, and S. Mukamel.  Conical intersection passages of molecules probed by X-ray diffraction and stimulated Raman spectroscopy.  J. Phys. Chem. Lett., 12:12300–12309, 2021.
  * (29) S. P. Neville, M. Chergui, A. Stolow, and M. S. Schuurman.  Ultrafast X-ray spectroscopy of conical intersections.  Phys. Rev. Lett, 120(24):243001, 2018.
  * (30) N. Rohringer.  X-ray Raman scattering: a building block for nonlinear spectroscopy.  Philos. Trans. R. Soc. A, 377(2145):20170471, 2019.
  * (31) M. Ferray, A. L’Huillier, X. F. Li, L. A. Lompre, G. Mainfray, and C. Manus.  Multiple-harmonic conversion of 1064 nm radiation in rare gases.  J. Phys. B, 21(3):L31, 1988.
  * (32) K. Varju, P. Johnsson, J. Mauritsson, A. L’Huillier, and R. Lopez-Martens.  Physics of attosecond pulses produced via high harmonic generation.  Am. J. Phys., 77(5):389–395, 2009.
  * (33) C. Pellegrini, A. Marinelli, and S. Reiche.  The physics of X-ray free-electron lasers.  Rev. Mod. Phys., 88(1):015006, 2016.
  * (34) G. Farkas and C. Tóth.  Proposal for attosecond light pulse generation using laser induced multiple-harmonic conversion processes in rare gases.  Phys. Lett. A, 168(5-6):447–450, 1992.
  * (35) S. E. Harris, J. J. Macklin, and T. W. Hänsch.  Atomic scale temporal structure inherent to high-order harmonic generation.  Opt. Commun., 100(5-6):487–490, 1993.
  * (36) M. Hentschel, R. Kienberger, C. Spielmann, G. A. Reider, N. Milosevic, T. Brabec, P. Corkum, U. Heinzmann, M. Drescher, and F. Krausz.  Attosecond metrology.  Nature, 414(6863):509–513, 2001.
  * (37) P. Paul, E. S. Toma, P. Breger, G. Mullot, F. Augé, P. Balcou, H. G. Muller, and P. Agostini.  Observation of a train of attosecond pulses from high harmonic generation.  Science, 292(5522):1689–1692, 2001.
  * (38) D. Jadoun and M. Kowalewski.  Time-resolved photoelectron spectroscopy of conical intersections with attosecond pulse trains.  J. Phys. Chem. Lett., 12:8103–8108, 2021.
  * (39) A. Marciniak, V. Despré, V. Loriot, G. Karras, M. Hervé, L. Quintard, F. Catoire, C. Joblin, E. Constant, A. I. Kuleff, and F. Lépine.  Electron correlation driven non-adiabatic relaxation in molecules excited by an ultrashort extreme ultraviolet pulse.  Nature Commun., 10(1):337, 2019.
  * (40) M. Kowalewski, K. Bennett, K. E. Dorfman, and S. Mukamel.  Catching conical intersections in the act: Monitoring transient electronic coherences by attosecond stimulated X-ray Raman signals.  Phys. Rev. Lett., 115:193003, 2015.
  * (41) M. Kowalewski, B. P. Fingerhut, K. E. Dorfman, K. Bennett, and S. Mukamel.  Simulating coherent multidimensional spectroscopy of nonadiabatic molecular processes: From the infrared to the X-ray regime.  Chem. Rev., 117(19):12165–12226, 2017.
  * (42) H. C. Longuet-Higgins.  The symmetry groups of non-rigid molecules.  Mol. Phys., 6(5):445–460, 1963.
  * (43) R. Obaid and M. Leibscher.  A molecular symmetry analysis of the electronic states and transition dipole moments for molecules with two torsional degrees of freedom.  J. Chem. Phys., 142(6):064315, 2015.
  * (44) I. J. Palmer, I. N. Ragazos, F. Bernardi, M. Olivucci, and M. A. Robb.  An MC-SCF study of the S1 and S2 photochemical reactions of benzene.  J. Am. Chem. Soc., 115(2):673–682, 1993.
  * (45) H. Köppel, B. Schubert, and H. Lischka.  Conical intersections and strong nonadiabatic coupling effects in singlet-excited acetylene: An ab initio quantum dynamical study.  Chem. Phys., 343(2):319–328, 2008.
  * (46) D. R. Yarkony.  Diabolical conical intersections.  Rev. Mod. Phys., 68:985–1013, Oct 1996.
  * (47) D. J. Tannor.  Introduction to quantum mechanics: a time-dependent perspective.  2007\.
  * (48) W. E. Arnoldi.  The principle of minimized iterations in the solution of the matrix eigenvalue problem.  Q. Appl. Math., 9(1):17–29, 1951.
  * (49) R. A. Friesner, L. S. Tuckerman, B. C. Dornblaser, and T. V. Russo.  A method for exponential propagation of large systems of stiff nonlinear differential equations.  J. Sci. Comput., 4(4):327–354, 1989.
  * (50) Y. Saad.  Variations on Arnoldi’s method for computing eigenelements of large unsymmetric matrices.  Linear Algebra Appl., 34:269–295, 1980.
  * (51) A. Nissen, H. O. Karlsson, and G. Kreiss.  A perfectly matched layer applied to a reactive scattering problem.  J. Chem. Phys., 133(5):054306, 2010.
  * (52) D. Keefer, T. Schnappinger, and S. de Vivie-Riedle, R.and Mukamel.  Visualizing conical intersection passages via vibronic coherence maps generated by stimulated ultrafast X-ray Raman signals.  Proc. Natl. Acad. Sci. U.S.A., 117(39):24069–24075, 2020.
  * (53) S. M. Cavaletto, D. Keefer, and S. Mukamel.  High temporal and spectral resolution of stimulated X-ray Raman signals with stochastic free-electron-laser pulses.  Phys. Rev. X, 11(1):011029, 2021.
  * (54) S. M. Cavaletto, D. Keefer, J. R. Rouxel, F. Aleotti, F. Segatta, M. Garavelli, and S. Mukamel.  Unveiling the spatial distribution of molecular coherences at conical intersections by covariance x-ray diffraction signals.  Proc. Natl. Acad. Sci. U.S.A., 118(22), 2021.
  * (55) B. Boashash.  Note on the use of the wigner distribution for time-frequency signal analysis.  IEEE Trans. Signal Process., 36(9):1518–1521, 1988.