# Purcell effect with extended sources: The role of the cross density of
states

R. Carminati<EMAIL_ADDRESS>Institut Langevin, ESPCI Paris, PSL
University, CNRS, 1 rue Jussieu, F-75005 Paris, France  M. Gurioli
<EMAIL_ADDRESS>Physics and Astronomy Department, Università degli
Studi di Firenze, via G. Sansone 1, 50019 Sesto Fiorentino, Italy

###### Abstract

We analyze the change in the spontaneous decay rate, or Purcell effect, of an
extended quantum emitter in a structured photonic environment. Based on a
simple theory, we show that the cross-density of states is the central
quantity driving interferences in the emission process. Using numerical
simulations in realistic photonic cavity geometries, we demonstrate that a
structured cross-density of states can induce subradiance or superradiance,
and change substantially the emission spectrum. Interestingly, the spectral
lineshape of the Purcell effect of an extended source cannot be predicted from
the sole knowledge of the spectral dependence of the local density of states.

##  I Introduction

In the last two decades, the ability to tailor the spontaneous decay rate of a
quantum emitter, by molding the dielectric environment at the nanoscale, has
been the driving force of a large variety of studies in cavity quantum
electrodynamics cavityQED , semiconductor exciton physics Shields2007 , and
nanophotonics Stockman-review ; Novotny-book ; Carminati-review . The
engineering of the electromagnetic local density of states (LDOS), in order to
tune the spontaneous decay rate by the Purcell effect Purcell1946 , has become
a major branch of nanophotonics and plasmonics. The standard formulation of
the Purcell effect is based on the electric dipole approximation, which
requires the electromagnetic field to be uniform over the spatial extent of
the emitter. The dipole approximation is expected to breakdown when extended
sources are used (for example extended excitons in semiconductors Andreani1999
), or when the external field is confined at the nanometer scale, or beyond
Benz2016 . A theory governing the spontaneous emission by extended excitons
has been proposed in Ref. Lodahl2012 . The theory predicts large oscillator
strengths, exceeding $10^{3}$ for quantum-dots (QDs) with diameters on the
order of $100$ nm. These ideas led to the concept of single-photon
superradiance, which has been observed in natural QDs arising in monolayer
fluctuation of GaAs quantum wells, with radiative lifetime below $100$ ps. On
a different scheme, Rydberg excitons have been observed in CuO2 bulk samples,
with a Bohr radius extending over micrometers Kazimierczuk2014 ;
Steinhauer2020 . Extended sources are also at play in microwave cavity quantum
electrodynamics with Rydberg atoms Haroche2013 , and in electron spectroscopy
where the interaction with light depends on non-local excitations Abajo2010 .
In a slightly different context, understanding and controlling the emission of
a set of single-photon sources connected through their electromagnetic
environments has become a key issue for the treatment of quantum information
in photonics Kimble2008 ; Pryde2019 , or to favor collective emission
Pustovit2009 ; Oppel2014 . This has stimulated the development of new
theoretical approaches to describe photonic states in ensembles of quantum
emitters Chang2017 ; Feist2021 . In the simple case of two (point) quantum
emitters, the role of the cross density of states (CDOS), that describes
electromagnetic spatial coherence between different points in an arbitrary
environment Caze2013 , has been put forward Canaguier2016 ; Canaguier2019 . It
was shown that the degree of quantum coherence in the emitted light is
controlled by a combination of the LDOS and CDOS, the latter describing the
coupling between the sources in the emission process.

In the case of extended quantum sources that cannot be reduced to a single
point dipole or an ensemble of independent (mutually incoherent) point
sources, it has been shown that spontaneous decay results from a non-local
interaction with the environment Lodahl2012 . Nevertheless, the exact role of
the CDOS on the Purcell effect has not been explored in its full generality.
In particular, beyond the generation of large effective oscillator strengths,
one can expect an engineering of the CDOS to affect the emission
substantially, e.g. by producing subradiance or superradiance, or by changing
the emission spectrum. The purpose of this paper is to analyze precisely the
Purcell effect for an extended quantum source. Based on simple theoretical
arguments, we show that the CDOS is the essential ingredient to understand the
change in the spontaneous emission decay rate, beyond the dipole
approximation. The theoretical results are highlighted by numerical
simulations, in a series of simple examples involving either a single high-$Q$
cavity or coupled cavities with non Lorentzian spectra. Interestingly, the
spectral emission lineshape of an extended source covering a nanostructured
environment cannot be predicted from the spectral dependence of the LDOS
alone.

##  II Spontaneous decay rate in structured environments

In this section we introduce the expressions of the spontaneous decay rate of
point and extended dipole emitters. We assume a weak-coupling regime, in which
the population of the excited state exhibits exponential decay, and the
transition dipole is not modified by the interaction with the environment. We
also restrict the discussion to electric-dipole transitions, and the theory
involves the electric part of the density of states of the electromagnetic
field. For magnetic-dipole transitions, similar expressions could be easily
deduced by using the magnetic contribution to the density of states Carminati-
review ; Joulain2003 ; Aigouy2014 .

###  II.1 Dipole emitter

In the weak-coupling regime, the spontaneous decay rate $\Gamma$ of a two-
level dipole emitter can be calculated from perturbation theory. For an
emitter with transition dipole ${\bf p}$, emission frequency $\omega$, and
located at the point ${\bf r}_{s}$, it can be shown that Agarwal1975 ;
Sipe1984

$\Gamma({\bf r}_{s},\omega)=\frac{2\mu_{0}\,\omega^{2}}{\hbar}\,|{\bf
p}|^{2}\,\mathrm{Im}\left[\hat{\bf u}\cdot{\bf G}({\bf r}_{s},{\bf
r}_{s},\omega)\hat{\bf u}\right]\,.$ (1)

In this expression ${\bf G}$ is the electric Green’s function describing the
electromagnetic response of the environment, and $\hat{\bf u}={\bf p}/|{\bf
p}|$ is the unit vector defining the orientation of the transition dipole. We
note that the internal dynamics of the emitter, described by the transition
dipole ${\bf p}$ that depends on the excited and ground state wavefunctions,
and the response of the environment factorize. In free space, the Green’s
function ${\bf G}_{0}$ is such that $\mathrm{Im}\left[\hat{\bf u}\cdot{\bf
G}_{0}({\bf r}_{s},{\bf r}_{s},\omega)\hat{\bf u}\right]=\omega/(6\pi c)$,
with $c$ the speed of light in vacuum, and the spontaneous decay rate is
$\Gamma_{0}(\omega)=\omega^{3}|{\bf p}|^{2}/(3\pi\hbar\,\epsilon_{0}\,c^{3})$.
Normalizing, we find that

$\frac{\Gamma({\bf r}_{s},\omega)}{\Gamma_{0}(\omega)}=\frac{6\pi
c}{\omega}\,\operatorname{Im}\left[\hat{\bf u}\cdot{\bf G}({\bf r}_{s},{\bf
r}_{s},\omega)\hat{\bf u}\right]\ .$ (2)

It will prove useful to recall a well-established correspondence between
classical and quantum calculations. Considering a classical electric dipole
${\bf p}$ oscillating at frequency $\omega$, the time-averaged power released
in the environment is $P({\bf r}_{s},\omega)=(\mu_{0}\omega^{3}/2)|{\bf
p}|^{2}\mathrm{Im}\left[\hat{\bf u}\cdot{\bf G}({\bf r}_{s},{\bf
r}_{s},\omega)\hat{\bf u}\right]$. Normalizing by the power $P_{0}$ radiated
in free space, we obtain

$\frac{P({\bf r}_{s},\omega)}{P_{0}(\omega)}=\frac{6\pi
c}{\omega}\,\operatorname{Im}\left[\hat{\bf u}\cdot{\bf G}({\bf r}_{s},{\bf
r}_{s},\omega)\hat{\bf u}\right]\ .$ (3)

The right-hand side is the same as in Eq. (2), showing that the change in the
spontaneous decay rate of a quantum emitter due to the electromagnetic
interaction with its environment, and the change in the power emitted by a
classical dipole, are identical Novotny-book ; Carminati-review ; Bouchet2019
.

In the particular case of a single mode cavity, the right-hand side in Eqs.
(2) and (3) defines the Purcell factor. For an emitter on resonance, it is
usually written as $3/(4\pi^{2})\lambda_{m}^{3}Q/\cal V$, with $Q$ the mode
quality factor, $\cal V$ the mode volume and $\lambda_{m}$ the resonance
wavelength Purcell1946 ; Carminati-review . Equations (2) and (3) generalize
the Purcell effect to a more complex environment, and the right-hand side can
be seen as a generalized Purcell factor. This factor actually describes the
change in the projected LDOS

$\rho({\bf
r}_{s},\omega)=\frac{2\omega}{\pi\,c^{2}}\,\mathrm{Im}\left[\hat{\bf
u}\cdot{\bf G}({\bf r}_{s},{\bf r}_{s},\omega)\hat{\bf u}\right]\,,$ (4)

which measures the coupling strength of a unit dipole $\hat{\bf u}$ to the
electromagnetic environment. The LDOS actually counts the relative
contribution of the electromagnetic modes at position ${\bf r}$ and frequency
$\omega$ Novotny-book ; Carminati-review ; Sauvan2013 . In terms of the LDOS
we find that

$\frac{\Gamma({\bf r}_{s},\omega)}{\Gamma_{0}(\omega)}=\frac{P({\bf
r}_{s},\omega)}{P_{0}(\omega)}=\frac{\rho({\bf
r}_{s},\omega)}{\rho_{0}(\omega)}\,,$ (5)

where $\rho_{0}(\omega)=\omega^{2}/(3\pi^{2}c^{3})$ is the projected LDOS in
free space.

###  II.2 Purcell effect for an extended source

A full quantum theory of spontaneous emission by an extended exciton has been
developed in Ref.Lodahl2012 , emphasizing the role of the internal dynamics of
the source. In the present work, we rather focus on the role of the
environment, and take advantage of the equivalence of the classical and
quantum descriptions to build the expression of the decay rate starting from
the expression of the emitted power. For a monochromatic classical source with
current density ${\bf j}({\bf r})$ confined within a volume $V$, and
oscillating at frequency $\omega$, the emitted power can be written as

$P=-\frac{1}{2}\mathrm{Re}\int_{V}{\bf j}^{*}({\bf r})\cdot{\bf E}({\bf
r})\,d^{3}r\,,$ (6)

where the superscript $*$ stands for complex conjugate. In terms of the
Green’s function ${\bf G}$ this becomes

$P=\frac{\mu_{0}\omega}{2}\int\mathrm{Im}\left[{\bf j}^{*}({\bf r})\cdot{\bf
G}({\bf r},{\bf r}^{\prime},\omega)\,{\bf j}({\bf
r}^{\prime})\right]\,d^{3}r\,d^{3}r^{\prime}\,,$ (7)

which, by making use of the reciprocity relation ${\bf G}({\bf r},{\bf
r}^{\prime},\omega)={\bf G}^{T}({\bf r}^{\prime},{\bf r},\omega)$, can be
rewritten as

$P=\frac{\mu_{0}\omega}{2}\int{\bf j}^{*}({\bf r})\cdot\mathrm{Im}{\bf G}({\bf
r},{\bf r}^{\prime},\omega)\,{\bf j}({\bf
r}^{\prime})\,d^{3}r\,d^{3}r^{\prime}\,.$ (8)

Normalizing by the power emitted in free space, we find that

$\frac{P}{P_{0}}=\frac{\int{\bf j}^{*}({\bf r})\cdot\mathrm{Im}{\bf G}({\bf
r},{\bf r}^{\prime},\omega)\,{\bf j}({\bf
r}^{\prime})\,d^{3}r\,d^{3}r^{\prime}}{\int{\bf j}^{*}({\bf
r})\cdot\mathrm{Im}{\bf G}_{0}({\bf r},{\bf r}^{\prime},\omega)\,{\bf j}({\bf
r}^{\prime})\,d^{3}r\,d^{3}r^{\prime}}\,.$ (9)

In the classical picture, the current density corresponding to an extended
exciton can be taken of the form

${\bf j}({\bf r})=-i\omega\hat{\bf u}\,\eta({\bf r},{\bf r}_{s})\,,$ (10)

where $\eta({\bf r},{\bf r}_{s})$ is the volume density of transition dipole
at the position ${\bf r}$, $\hat{\bf u}$ gives the local orientation of the
transition dipole and ${\bf r}_{s}$ defines the position of the source (the
reference point ${\bf r}_{s}$ can be chosen arbitrarily, for example as the
center of the source spatial distribution). This model allows to take into
account a possible heterogeneity of the distribution of elementary dipoles, in
amplitude and orientation, inside the extended emitter. Also note that a
single extended exciton corresponds to spatially coherent source, that can be
seen as a continuous distribution of mutually coherent elementary transition
dipoles. Inserting this expression into Eq. (9), and invoking again the
classical to quantum correspondence, we obtain

$\frac{\Gamma({\bf r}_{s},\omega)}{\Gamma_{0}(\omega)}=\frac{\int\eta^{*}({\bf
r},{\bf r}_{s})\,\eta({\bf r}^{\prime},{\bf r}_{s})\,\mathrm{Im}[\hat{\bf
u}\cdot{\bf G}({\bf r},{\bf r}^{\prime},\omega)\hat{\bf
u}^{\prime}]\,d^{3}r\,d^{3}r^{\prime}}{\int\eta^{*}({\bf r},{\bf
r}_{s})\,\eta({\bf r}^{\prime},{\bf r}_{s})\,\mathrm{Im}[\hat{\bf u}\cdot{\bf
G}_{0}({\bf r},{\bf r}^{\prime},\omega)\hat{\bf
u}^{\prime}]\,d^{3}r\,d^{3}r^{\prime}}\,.$ (11)

This expression of the normalized decay rate of an extended source accounts
for the non-local interaction with the environment, as already put forward in
Ref. Lodahl2012 . The similarity with Eq. (5) can be made more explicit by
introducing the projected CDOS Caze2013 , defined as

$\rho({\bf r},{\bf
r}^{\prime},\omega)=\frac{2\omega}{\pi\,c^{2}}\,\mathrm{Im}\left[\hat{\bf
u}\cdot{\bf G}({\bf r},{\bf r}^{\prime},\omega)\hat{\bf u}^{\prime}\right]\,.$
(12)

The projected CDOS measures the connection between two unit dipoles $\hat{\bf
u}$ and $\hat{\bf u}^{\prime}$ located at the points ${\bf r}$ and ${\bf
r}^{\prime}$. More precisely, the CDOS counts the contribution of the
electromagnetic modes connecting ${\bf r}$ and ${\bf r}^{\prime}$ at frequency
$\omega$ Caze2013 ; Carminati-review ; Sauvan2014 . It can also be understood
as a measure of the intrinsic spatial coherence between the two points ${\bf
r}$ and ${\bf r}^{\prime}$. We also note that for ${\bf r}={\bf r}^{\prime}$
the projected CDOS coincides with the projected LDOS evaluated at point ${\bf
r}$. In terms of the CDOS, we find that

$\frac{\Gamma({\bf r}_{s},\omega)}{\Gamma_{0}(\omega)}=\frac{\int\eta^{*}({\bf
r},{\bf r}_{s})\,\eta({\bf r}^{\prime},{\bf r}_{s})\,\rho({\bf r},{\bf
r}^{\prime},\omega)\,d^{3}r\,d^{3}r^{\prime}}{\int\eta^{*}({\bf r},{\bf
r}_{s})\,\eta({\bf r}^{\prime},{\bf r}_{s})\,\rho_{0}({\bf r},{\bf
r}^{\prime},\omega)\,d^{3}r\,d^{3}r^{\prime}}\,,$ (13)

where $\rho_{0}({\bf r},{\bf r}^{\prime},\omega)$ is the CDOS in free space.
This expression of the normalized decay rate (or Purcell enhancement)
generalizes Eq. (5) to extended emitters beyond the dipole approximation.
Indeed, it can be seen immediately that by using $\eta({\bf r},{\bf
r}_{s})=p\,\delta({\bf r}-{\bf r}_{s})$ in Eq. (13), which corresponds to a
point dipole source with amplitude $p$ located at the point ${\bf r}_{s}$, we
recover Eq. (5).

##  III Extended sources in photonic cavities

It is instructive to start with the simplest model of a spatially coherent
extended source, made of two mutually coherent point dipoles with amplitude
${\bf p}_{1}=p_{1}\hat{\bf u}_{1}$ and ${\bf p}_{2}=p_{2}\hat{\bf u}_{2}$, and
located at positions ${\bf r}_{1}$ and ${\bf r}_{2}$, that is

${\bf j}({\bf r})=-i\omega\hat{\bf u}_{1}p_{1}\delta({\bf r}-{\bf
r}_{1})-i\omega\hat{\bf u}_{2}p_{2}\delta({\bf r}-{\bf r}_{2})\,.$ (14)

According to Eq. (13), the spontaneous decay rate is

$\Gamma\sim|p_{1}|^{2}\rho_{11}+|p_{2}|^{2}\rho_{22}+2\mathrm{Re}(p_{1}p_{2}^{*})\rho_{12}\,,$
(15)

up to a constant prefactor (independent of the environment) that we do not
specify. Here $\rho_{ij}=2\omega/(\pi c^{2})\mathrm{Im}[\hat{\bf
u}_{i}\cdot{\bf G}({\bf r}_{i},{\bf r}_{j},\omega)\hat{\bf u}_{j}]$ is the
LDOS for $i=j$ and the CDOS for $i\neq j$. If the two dipoles have the same
amplitude $|p|$ and a phase difference $\phi$, we simply have

$\Gamma\sim|p|^{2}[\rho_{11}+\rho_{22}+2\rho_{12}\,\cos\phi]\,.$ (16)

The mutual coherence between the two sources is essential for the appearance
of the interference term involving the CDOS. This situation is different from
that analyzed in Refs. Canaguier2016 ; Canaguier2019 ; Carminati2015 , where
spontaneous emission by two independent sources (e.g., two fluorescent
molecules) in a structured environment was considered. In this case, it was
shown that intensity fluctuations are influenced by a cross term involving the
CDOS. In the situation considered here, interferences directly modify the
spontaneous decay rate. The state of interference not only depends on the
relative phase of the dipole sources, but also on the environment through the
CDOS that can take positive or negative values. The increase or decrease of
the spontaneous decay rate of an extended source due to interferences
corresponds to superradiance and subradiance Fauche2017 ; Lalanne2018 . In the
following subsections, we use numerical simulations to illustrate the
substantial role of the CDOS on the spontaneous decay rate of extended sources
in photonic cavities.

###  III.1 High-$Q$ cavity

Let us first consider a source immersed in a discrete set of high-$Q$
electromagnetic modes, whose resonances do not overlap spectrally. The Green’s
function in this case takes the following form Carminati-review

${\bf G}({\bf r},{\bf r}^{\prime},\omega)=\sum_{m}c^{2}\frac{{\bf e}_{m}({\bf
r})\otimes{\bf e}_{m}^{*}({\bf
r}^{\prime})}{\omega_{m}^{2}-\omega^{2}-i\omega\gamma_{m}}\,,$ (17)

where ${\bf e}_{m}$ are eigenmodes of the vector Helmholtz equation in absence
of losses, $\omega_{m}$ are eigenfrequencies, and $\gamma_{m}$ is a
phenomenological mode damping rate. This expression can also be seen as low-
loss limit ($\gamma_{m}\ll\omega_{m}$) of an expansion in quasi-normal modes
(QNMs) Sauvan2013 . For the CDOS we find that Carminati-review ; Sauvan2014

$\rho_{12}=\sum_{m}\frac{\gamma_{m}}{2\pi}\frac{\mathrm{Re}[(\hat{\bf
u}_{1}\cdot{\bf e}_{m}({\bf r}_{1}))\,(\hat{\bf u}_{2}\cdot{\bf
e}_{m}^{*}({\bf r}_{2}))]}{(\omega-\omega_{m})^{2}+\gamma_{m}^{2}/4}\,.$ (18)

This expression shows explicitly that the CDOS accounts for the phase
difference of each field mode ${\bf e}_{m}$ between the points ${\bf r}_{1}$
and ${\bf r}_{2}$. In the particular case of a single-mode high-$Q$ cavity,
the CDOS factorizes in the form $\rho_{12}^{2}=\rho_{11}\rho_{22}$
Canaguier2019 . For a source made of two coherent dipoles located at positions
${\bf r}_{1}$ and ${\bf r}_{2}$ in a single mode cavity, we find from Eq. (16)
that the spontaneous decay rate is

$\Gamma\sim|p|^{2}[\rho_{11}+\rho_{22}\pm
2\sqrt{\rho_{11}\rho_{22}}\,\cos\phi]\,.$ (19)

The field mode in this case is a standing wave with real amplitude, and the
plus or minus sign in the interference term corresponds to $\hat{\bf
u}_{1}\cdot{\bf e}_{m}({\bf r}_{1})$ and $\hat{\bf u}_{2}\cdot{\bf e}_{m}({\bf
r}_{2})$ having identical or opposite signs (which corresponds to a positive
or negative CDOS $\rho_{12}$).

In order to illustrate the relevance of Eq. (19) in a realistic situation, we
consider a single mode contribution in the well-known L3 photonic crystal
cavity Noda2007 . The cavity is formed by a line of three missing holes in a
triangular-lattice pattern of air holes on a slab with refractive index
$n=3.48$ and thickness $L=200$ nm. Calculated maps of the amplitudes
$E_{y}({\bf r})=\mathrm{Re}[\hat{\bf y}\cdot{\bf e}_{m}({\bf r})]$ and
$E_{x}({\bf r})=\mathrm{Re}[\hat{\bf x}\cdot{\bf e}_{m}({\bf r})]$ of the
fundamental mode $m$, obtained from a FDTD simulation, are displayed in Fig.
1(a) and 1(b).The positions of the photonic crystal holes, designing the L3
cavity, are indicated as circles. The fundamental mode $m$ is characterized by
a LDOS with a Lorentzian lineshape (not shown), with a resonance wavelength
$\lambda_{m}=1270$ nm and a quality factor $Q_{m}=\omega_{m}/\gamma_{m}\simeq
2000$. We focus on the emission by two coherent electric dipoles placed either
at the points $\alpha$ indicated in Fig. 1(a), or at the points $\beta$
indicated in Fig. 1(b). We show in Fig. 1(c) [resp. Fig. 1(d)] the Purcell
enhancement $\Gamma/\Gamma_{0}$ for two dipoles with identical amplitude and
phase difference $\phi$, and oriented along the $x$ direction [resp. $y$
direction]. The black curve corresponds to two dipoles in phase ($\phi=0$),
and the red curve to two dipoles with opposite signs ($\phi=\pi$). The decay
rate is obtained from a FDTD simulation, in which the emitted power is
calculated from the net flux of the Poynting vector through a box containing
the two dipoles, and having spectral linewidth much larger than the linewidth
of the cavity mode.

Figure 1:  Emission by two dipoles in a L3 photonic crystal cavity on a slab.
The slab has a refractive index $n=3.48$ and a thickness $L=200$ nm. The
photonic crystal is made of holes with a diameter $d=180$ nm and a lattice
constant $a=320$ nm. (a) and (b): Field amplitudes $E_{x}({\bf
r})=\mathrm{Re}[\hat{\bf x}\cdot{\bf e}_{m}({\bf r})]$ and $E_{y}({\bf
r})=\mathrm{Re}[\hat{\bf y}\cdot{\bf e}_{m}({\bf r})]$ of the cavity mode
${\bf e}_{m}({\bf r})$. The field is calculated in a plane coinciding with the
middle of the slab. The points $\alpha$ and $\beta$ indicate the positions of
the two emitting dipoles. The locations of the photonic crystal holes are also
highlighted. (c) and (d): Spectrum of the normalized spontaneous decay rate
$\Gamma/\Gamma_{0}$ (or Purcell enhancement) of two coherent dipoles with
identical amplitudes, phase difference $\phi$, located at positions $\alpha$
and oriented along the $x$ direction (panel c), or located at positions
$\beta$ and oriented along the $y$ direction (panel d). Black line: $\phi=0$.
Red line: $\phi=\pi$. The curves corresponding to subradiant emission have
been multiplied by a factor of $10^{3}$ (panel c) and $10^{4}$ (panel d) for
the sake of visibility.

For two dipoles placed at positions $\alpha$ [Fig. 1(c)], that are connected
by a negative CDOS since the field amplitudes have opposite signs, we observe
a subradiant emission for $\phi=0$ and perfect superradiance for $\phi=\pi$
(we have verified that the emission rate in this case is twice that obtained
for a single dipole). For two dipoles placed at positions $\beta$ [Fig. 1(d)],
connected by a positive CDOS, we observe superradiance for $\phi=0$ and
subradiance for $\phi=\pi$, as expected. Note that the distances between the
two points $\alpha$ and $\beta$ are $\Delta r=320$ nm and $\Delta r=480$ nm,
respectively, that is comparable or larger than $\lambda/n$, with $\lambda$
the emission wavelength in vacuum. In bulk materials, super or subradiant
effects would be observable only for $\Delta r\ll\lambda/n$. The simulations
clearly confirm that the photonic mode imposes the constructive or destructive
interference, through the sign of the CDOS, in agreement with Eq. (19).

Next, in order to mimic the emission by a one-dimensional extended exciton (as
found for example in a carbon nanotube or a quantum wire), we consider a line
source with length $d$ (hereafter denoted by emitting nanowire). In the FDTD
simulation, the source model corresponding to the emission by a single exciton
is a linear cluster of $N$ elementary and mutually coherent dipoles with
amplitude ${\bf p}/\sqrt{N}$ (as already pointed out the source is spatially
coherent for a single exciton). In order to focus on the influence of the
environment, we consider the situation where all elementary dipoles are in
phase (but the approach allows one to use a more refined source model). The
emitting nanowire is embedded in the L3 cavity, and we can expect very
different radiative rates depending on its location and orientation. As a case
study, we consider a nanowire oriented along the $x$ direction, polarized
along the $y$ direction, and located at the center of the L3 cavity. The
emission wavelength $\lambda=1270$ nm is on resonance with the cavity mode,
and the nanowire can be seen as immersed in the field distribution shown in
Fig. 1(b). We show in Fig. 2(a) the calculated radiative rate versus the
nanowire length $d$ in the L3 cavity (green circles), and in a homogeneous
medium with the same refractive index $n=3.48$ as the cavity slab (blue
triangles). The radiative rate in the homogeneous medium increases with the
exciton size, in agreement with the result already found in Ref. Lodahl2012 ,
up to a maximum reached for a size $d=\lambda/2n$. Interestingly, a very
similar trend is observed for the nanowire in the L3 cavity.

Figure 2:  (a) Comparison of the spontaneous decay rate of an emitting
nanowire polarized along the $y$ direction and extended along the $x$
direction in the L3 cavity (green circles), and in a homogeneous medium with
refractive index $n=3.48$ (blue triangles), versus the nanowire length $d$. In
the case of the L3 cavity, the center of the nanowire coincides with the
center of the cavity. (b) Purcell enhancement $\Gamma/\Gamma_{0}$ versus the
emitting nanowire length (red dots). The field amplitude $E_{y}({\bf
r})=\mathrm{Re}[\hat{\bf y}\cdot{\bf e}_{m}({\bf r})]$ of the cavity mode at
the extremities of the emitting nanowire is also indicated (blue dashed line).
Up to a length $d\simeq 300$ nm, the center and the extremities of the
nanowire are connected by a positive CDOS (the field amplitude is positive at
the center and at the extremities). For $d\gtrsim 300$ nm, the center and the
extremities of the emitting nanowire are connected by a negative CDOS, and
destructive interferences reduce the spontaneous decay rate. The minimum is
observed for $d\simeq 560$ nm.

To address this length dependence in more details, we show in Fig. 2(b) the
Purcell enhancement $\Gamma/\Gamma_{0}$ in the L3 cavity (red dots) versus the
emitting nanowire length. In the same figure we also display the field
amplitude of the cavity mode at the extremities of the nanowire (blue dashed
line). We see that the Purcell enhancement is almost constant for $d\leq 300$
nm, i.e. as long as the nanowire is inside the central hot spot of the field
amplitude seen in Fig. 1(b). In this case, the CDOS remains positive over the
emitting nanowire extent (the field amplitude of the mode does not change sign
from the center to the extremities of the nanowire). Constructive
interferences between the emissions from different parts of the nanowire lead
to a substantial Purcell enhancement. However, for $d\gtrsim 300$ nm, the
center and the extremities of the emitting nanowire are connected by a
negative CDOS (the field amplitude of the mode is positive at the center of
the nanowire and negative at the extremities, since they overlap with the
negative spots indicated as $\beta$ in Fig. 1(b)). Destructive interferences
reduce the spontaneous decay rate, leading to a minimum of the Purcell
enhancement for $d\simeq 560$ nm. This subradiance of the emitting nanowire
arises from the dominance of negative CDOS contributions, induced by the
spatial structure of the underlying cavity mode.

###  III.2 Coupled cavities with non-Lorentzian spectrum

In addition to generating subradiance or superradiance, the CDOS contribution
in the interferometric emission process may also influence the spectral
lineshape. An interesting situation is that of two coupled photonic cavities
generating a LDOS with a non-Lorentzian spectrum, which can be achieved by
selecting two L3 cavities with unbalanced losses Gurioli2020 . This situation
goes beyond the high-$Q$ model used above, and the expansion (18) for the LDOS
and CDOS has to be replaced by a more general expansion in QNMs (for an
introduction to QNMs in photonics, see for example Lalanne2018 ). Considering
that each cavity can be described by a single QNM, the Green’s function takes
the form Sauvan2013 ; Sauvan2014

${\bf G}({\bf r},{\bf
r}^{\prime},\omega)\simeq\frac{c^{2}}{2\omega}\left[\frac{\tilde{\bf
E}_{a}({\bf r})\otimes\tilde{\bf E}_{a}({\bf
r}^{\prime})}{\tilde{\omega}_{a}-\omega}+\frac{\tilde{\bf E}_{b}({\bf
r})\otimes\tilde{\bf E}_{b}({\bf
r}^{\prime})}{\tilde{\omega}_{b}-\omega}\right]\,.$ (20)

Here $\tilde{\bf E}_{a}$ and $\tilde{\bf E}_{b}$ are QNMs corresponding to
cavity $a$ and $b$, respectively, with $\tilde{\omega}_{a}$ and
$\tilde{\omega}_{b}$ the associated complex resonance frequencies. A major
difference with Eq. (17) is that the expression above is not limited to the
low-loss (or high-$Q$) regime.

Coming back to the simple source model described by Eq. (14), with two
coherent dipoles at positions ${\bf r}_{1}$ and ${\bf r}_{2}$, we can derive
the expressions of the LDOS and CDOS entering Eq. (16). Defining
$\omega_{m}=\mathrm{Re}\tilde{\omega}_{m}$ and
$\gamma_{m}=-2\mathrm{Im}\tilde{\omega}_{m}$, with $m\in\\{a,b\\}$, the result
can be cast in a form involving Fano profiles of the form

$F(m,q,\omega)=\frac{\gamma_{m}/2}{q^{2}+1}\left[\frac{(q^{2}-1)\gamma_{m}/2+2q(\omega-\omega_{m})}{(\omega-\omega_{m})^{2}+\gamma_{m}^{2}/4}\right]\,,$
(21)

where $q$ is a parameter defining the resonance lineshape Gurioli2021 .
Introducing the phases of the QNMs at the source points
$\varphi_{i}^{m}=\arg[\hat{\bf u}_{i}\cdot\tilde{\bf E}_{m}({\bf r}_{i})]$,
with $i\in\\{1,2\\}$, and the parameters $q_{i}^{m}=\tan(\varphi_{i}^{m})$,
$q_{12}^{m}=[\tan(\varphi_{1}^{m})+\tan(\varphi_{2}^{m})]/2$, a simple
calculation leads to

$\displaystyle\rho_{11}\sim F(a,q_{1}^{a},\omega)+F(b,q_{1}^{b},\omega)$ (22)
$\displaystyle\rho_{22}\sim F(a,q_{2}^{a},\omega)+F(b,q_{2}^{b},\omega)$ (23)
$\displaystyle\rho_{12}\sim F(a,q_{12}^{a},\omega)+F(b,q_{12}^{b},\omega)\,.$
(24)

We find that the spectrum of the LDOS and CDOS is a sum of Fano lineshapes,
each arising from one QNM ($a$ and $b$), with $q$ parameters depending on the
location of the source points ($1$ and $2$). These expressions reveal the
complexity of the prediction of the spontaneous emission dynamics of an
extended source in a structured environment, that depends on the LDOS and the
CDOS, and on the spatial distribution of the phase of the QNMs. In particular,
it is interesting to note that the spectrum of the decay rate cannot be
predicted solely from the knowledge of the spectrum of the LDOS at each point.
The relative phases $\varphi_{i}^{m}$ between different source points, imposed
by the electromagnetic QNMs, play a crucial role on the lineshape of the
emission spectrum.

To illustrate this result, we have carried out numerical simulations in the
system depicted in Fig. 3(a), made of two coupled photonic-crystal cavities
with unbalanced losses. The holes with different colors have different radii,
the green holes being designed to tailor the losses, and the pale blue hole
being designed to reduce the frequency detuning between the two cavities.
Unbalanced losses lead to highly non-Lorentzian LDOS and CDOS spectra (not
shown), resulting from the superposition of two strongly asymmetric Fano
contributions according to Eqs. (22-24). As a result, we can expect the
spectral profile of the spontaneous decay rate of a point dipole source to be
very different from that of an extended exciton. This is clearly seen in Fig.
3(b), where we show the calculated spectrum of the spontaneous decay rate
$\Gamma$ of a point dipole (black dashed line) and an emitting nanowire (red
dashed line). Both sources are polarized in the $x$ direction, and centered on
the position indicated by the red dot in Fig. 3(a). The nanowire has a length
$d=300$ nm and is oriented along the $x$ direction. We observe that the CDOS
contribution not only produces a subradiant effect reducing the decay rate to
very small values, but also induces a substantial reshaping of the spectral
profile. This example clearly demonstrates that the emission lineshape in the
Purcell effect of an extended source cannot be predicted from the sole
knowledge of the LDOS spectrum at all points.

Figure 3:  (a): Sketch of the system made of two coupled photonic-crystal
cavities on a slab, with the same thickness and refractive index as in Fig. 1.
The white holes have a diameter of 180 nm, the green holes of 140 nm and the
blue hole of 150 nm. (b): Spectrum of the decay rate for a point dipole source
(black dashed line) and an emitting nanowire (red dashed line). The red cross
in the left cavity in panel (a) indicates the location of the center of the
source, that coincides with one of the $\alpha$ points indicated in Fig. 1(a).
Both sources are polarized along the $x$ direction. The emitting nanowire has
a length $d=300$ nm and is oriented along the $x$ direction.

##  IV Conclusion

In summary, we have studied theoretically the change in the spontaneous decay
rate (or Purcell effect) of an extended emitter in a structured environment,
emphasizing the role of the CDOS in driving interferences in the emission
process. Beyond the generation of large effective oscillator strengths, that
had been put forward and exploited in previous studies, we have demonstrated
that a structured CDOS (in space and/or spectrum) may produce subradiance or
superradiance, and change substantially the emission spectrum. These results,
based on a simple but robust theory, have been highlighted by numerical
simulations in realistic photonic-crystal cavity geometries. The study clearly
shows that a knowledge of the LDOS spectrum at all points does not permit to
predict the change in the decay rate, neither in value nor in lineshape. It
also emphasizes the importance of the engineering of both the LDOS and the
CDOS for the control of spontaneous emission by extended sources.

## Acknowledgments

We thank P. Lalanne and K. Vynck for helpful discussions at the initial stage
of this study. M.G. acknowledges a visiting fellowship from Institut Langevin
and ESPCI Paris - PSL where part of this work was done.

## Funding

This work was supported by LABEX WIFI (Laboratory of Excellence within the
French Program Investments for the Future) under references ANR-10- LABX-24
and ANR-10-IDEX-0001-02 PSL*.

## References

  * (1) Cavity Quantum Electrodynamics, edited by P. Berman (Academic Press, 1994); S. Haroche and J.M. Raimond, Exploring the Quantum: Atoms, Cavities, and Photons (Oxford University Press, 2006).
  * (2) A.J. Shields, “Semiconductor quantum light sources,” Nature Phot. 1, 215-223 (2007).
  * (3) M.I. Stockman, “Nanoplasmonics: past, present, and glimpse into future,” Opt. Express 19, 22029-22106 (2011).
  * (4) L. Novotny and B. Hecht, Principles of Nano-Optics (Cambridge University Press, 2006).
  * (5) R. Carminati, A. Cazé, D. Cao, F. Peragut, V. Krachmalnicoff, R. Pierrat, and Y. De Wilde, “Electromagnetic density of states in complex plasmonic systems,” Surf. Sci. Rep. 70, 1-41 (2015).
  * (6) E. M. Purcell, “Spontaneous emission probabilities at radio frequencies,” Phys. Rev. 69, 681 (1946).
  * (7) L.C. Andreani, G. Panzarini and J.M. Gérard, “Strong-coupling regime for quantum boxes in pillar microcavities: Theory,” Phys. Rev. B 60, 13276 (1999).
  * (8) F. Benz, M.K. Schmidt, A. Dreismann, R. Chikkaraddy, Y. Zhang, A. Demetriadou, C. Carnegie, H. Ohadi, B. de Nijs, R. Esteban, J. Aizpurua and J.J. Baumberg, “Single-molecule optomechanics in picocavities,” Science 354, 726-729 (2016).
  * (9) S. Stobbe, P. T. Kristensen, J. E. Mortensen, J. M. Hvam, J. Mørk and P. Lodahl, “Spontaneous emission from large quantum dots in nanostructures: Exciton-photon interaction beyond the dipole approximation,” Phys. Rev. B 86, 085304 (2012).
  * (10) T. Kazimierczuk, D. Fröhlich, S. Scheel, H. Stolz and M. Bayer, “Giant Rydberg Excitons in Cuprous Oxide,” Nature 514, 343-347 (2014).
  * (11) S. Steinhauer, M.A.M. Versteegh, S. Gyger, A. W. Elshaari, B. Kunert, A. Mysyrowicz and V. Zwiller, “Rydberg excitons in Cu2O microcrystals grown on a silicon platform,” Commun. Materials 1, 11 (2020).
  * (12) S. Haroche, “Controlling photons in a box and exploring the quantum to classical boundary,” Rev. Mod. Phys. 85, 1083 (2013).
  * (13) F. J. García de Abajo, “Optical excitations in electron microscopy,” Rev. Mod. Phys. 82, 209 (2010).
  * (14) H. J. Kimble, “The quantum internet,” Nature (London) 453, 1023-1030 (2008).
  * (15) S. Slussarenko and G.J. Pryde, “Photonic quantum information processing: A concise review,” Appl. Phys. Reviews 6, 041303 (2019).
  * (16) V. N. Pustovit and T. V. Shahbazyan, “Cooperative emission of light by an ensemble of dipoles near a metal nanoparticle: The plasmonic Dicke effect,” Phys. Rev. Lett. 102, 077401 (2009).
  * (17) S. Oppel, R. Wiegner, G. S. Agarwal, and J. von Zanthier, “Directional Superradiant Emission from Statistically Independent Incoherent Nonclassical and Classical Sources,” Phys. Rev. Lett. 113, 263606 (2014).
  * (18) M.T. Manzoni, D.E. Chang and J.S. Douglas, “Simulating quantum light propagation through atomic ensembles using matrix product states,” Nat. Commun. 8, 1743 (2017).
  * (19) J. Feist, A. I. Fernández-Domínguez and F.J. García-Vidal, “Macroscopic QED for quantum nanophotonics: emitter-centered modes as a minimal basis for multiemitter problems,” Nanophotonics 10, 477-489 (2021).
  * (20) A. Cazé, R. Pierrat, and R. Carminati, “Spatial Coherence in Complex Photonic and Plasmonic Systems,” Phys. Rev. Lett. 110, 063903 (2013).
  * (21) A. Canaguier-Durand and R. Carminati, “Quantum coherence of light emitted by two single-photon sources in a structured environment,” Phys. Rev. A 93, 033836 (2016).
  * (22) A. Canaguier-Durand, R. Pierrat and R. Carminati, “Cross density of states and mode connectivity: Probing wave localization in complex media,” Phys. Rev. A 99, 013835 (2019).
  * (23) K. Joulain, R. Carminati, J.P. Mulet and J.-J. Greffet, “Definition and measurement of the local density of electromagnetic states close to an interface,” Phys. Rev. B 68, 245405 (2003).
  * (24) L. Aigouy, A. Cazé, P. Gredin, M. Mortier and R. Carminati, “Mapping and quantifying electric and magnetic dipole luminescence at the nanoscale,” Phys. Rev. Lett. 113, 076101 (2014).
  * (25) G. S. Agarwal, “Quantum electrodynamics in the presence of dielectrics and conductors. I. Electromagnetic-field response functions and black-body fluctuations in finite geometries,” Phys. Rev. A 11, 230 (1975).
  * (26) J. M. Wylie and J. E. Sipe, “Quantum electrodynamics near an interface,” Phys. Rev. A 30, 1185 (1984).
  * (27) D. Bouchet and R. Carminati, “Quantum dipole emitters in structured environments: a scattering approach: tutorial,” J. Opt. Soc. Am. A 36, 186-195 (2019).
  * (28) C. Sauvan, J.P. Hugonin, I. S. Maksymov and P. Lalanne, “Theory of the Spontaneous Optical Emission of Nanosize Photonic and Plasmon Resonators,” Phys. Rev. Lett. 110, 237401 (2013).
  * (29) C. Sauvan, J.P. Hugonin, R. Carminati and P. Lalanne, “Modal representation of spatial coherence in dissipative and resonant photonic systems,” Phys. Rev. A 89, 043825 (2014).
  * (30) R. Carminati, G. Cwilich, L.S. Froufe Pérez and J.J. Sáenz, “Speckle fluctuations resolve the interdistance between incoherent point sources in complex media,” Phys. Rev. A 91, 023807 (2015).
  * (31) P. Fauché, S.G. Kosionis and P. Lalanne, “Collective scattering in hybrid nanostructures with many atomic oscillators coupled to an electromagnetic resonance,” Phys. Rev. B 95, 195418 (2017).
  * (32) P. Lalanne, W. Yan, K. Vynck, C. Sauvan and J.P. Hugonin, “Light Interaction with Photonic and Plasmonic Resonances,” Laser Phot. Rev. 12, 1700113 (2018).
  * (33) S. Noda, M. Fujita and T. Asano, “Spontaneous-emission control by photonic crystals and nanocavities,” Nat. Phot. 1, 449-458 (2007).
  * (34) D. Pellegrino, D. Balestri, N. Granchi, M. Ciardi, F. Intonti, F. Pagliano, A.Y. Silov, F. W. Otten, T. Wu, K. Vynck, P. Lalanne, A. Fiore and M. Gurioli, “Non-Lorentzian Local Density of States in Coupled Photonic Crystal Cavities Probed by Near- and Far-Field Emission,” Phys. Rev. Lett. 124, 123902 (2020).
  * (35) M. Gurioli, “Vacuum Field Photonic Trap for Excitons,” Adv. Quantum Technol. 4, 2100046 (2021).


# A portfolio management of a small RES utility with a Structural Vector
Autoregressive model of German electricity markets

Katarzyna Maciejowska<EMAIL_ADDRESS>Department of
Operations Research and Business Intelligence, Wrocław University of Science
and Technology, 50-370 Wrocław, Poland

(This version: )

###### Abstract

The changes in electricity markets expose RES producers and electricity
traders to various risks, among which the price and the volume risk play a
very important role. In this research, a portfolio building strategies are
presented, which allow to dynamically choose a proportion of electricity
traded in different electricity markets (day-ahead and intraday) and hence to
optimize the behavior of an utility. Two types of approaches are considered:
simple, assuming that the proportions are fixed, and data driven, which allows
for thier fluctuation. In order to explore the market information, Structural
Vector Autoregressive (SVAR) model is applied, which allows to estimate the
relationship between variables of interest and to simulate their future
distribution. The presented methods are evaluated with data coming from German
electricity market. The results indicate that data driven trading strategies
allow to increase the utility revenue and at the same time reduce the trading
risk, measured by the predictability of the next day income and the revenue
Value at Risk. It turns out that the approach based on Sharp Ratio provides
the most robust results.

###### keywords:

Intraday electricity market, Day ahead electricity market, Structural Vector
Autoregressive model, probabilistic forecasting, Trading strategy

[table]capposition=top

##  1 Introduction

In recent decades, thanks to development of short-term electricity markets,
new trade opportunities have opened for generation utilities and demand units.
Their operation is no longer optimized centrally and power plants managers act
in the market to maximize utility profits. Market participants can now choose,
whether to trade via the organized power exchanges, such as Nord Pool or EEX
in Europe and PJM in the USA, broker platform or to sign over-the-counter
(OTC) contracts. Although trade in power exchanges is voluntary, market
parties are encouraged to self-balance their positions. This means that no
electricity should be left intentionally for the trade on the balancing market
(Koch and Maskos (2020), Pape et al. (2016)). Self-balancing is particularly
valid for RES utilities, which generation is based on intermittent energy
sources, such as wind and solar (Gianfreda et al. (2016),Kiesel and Paraschive
(2017), Koch and Hirth (2019)).

The changes of electricity markets have been motivated and accompanied by
dynamic development of RES generation. In the year 2020, RES accounts for
50.5% of generation in Germany. The increase of RES penetration would not be
possible without a variety of support schemes introduced by European
countries. First, RES generation is granted a priority during the dispatch,
ensuring that all green electricity is efficiently traded. Second, different
financial incentives, starting with Feed-in-Tariffs (FIT) have been proposed
to increase the profitability of RES investments. In FIT, RES generators are
paid a fixed price at a guaranteed level (irrespective of the wholesale price)
for the electricity produced and fed into the grid. In Germany, the 2000
Renewable Energy Act guaranteed FIT for wind and solar generation for 20
years. For many installations, this two-decade guarantee is going to expire
soon. As the result, RES producers will need to sell their generation at
market prices either via the power exchange or bilateral contracts. At the
same time, one could observe a general change in the approach toward support
schemes and shift from FIT to Feed-in-Premium (FIP) mechanism. In case of FIP,
producers receive a premium price, which is a payment (EURO/MWh) in addition
to the wholesale price. This premium can be either fixed (Denmark, Lithuania)
or floating (Germany, Greece and Netherlands among others). It should be
noticed that the ongoing changes leads to a closer relation between revenues
of RES generators and wholesale electricity prices. They encourage financially
sustainable investments, which respond to market incentives.

As the result, RES producers become more exposed to various market risks. As
stated by Spodniak et al. (2021), the major ones are price and volume risks.
The price risk reflects the fact that electricity prices are stochastic and
depends on the unknown future levels of demand and generation structure (Weron
(2014), Uniejewski et al. (2019), Hong et al. (2020)). Additionally,
electricity producers face now a cascade of trade opportunities, which
includes different markets (bilateral, day-ahead, intrday) and contract types.
From the perspective or RES utility, the volume risk can be analyzed at two
levels: a state wise and an individual. The uncertainty about the (individual)
utility production stems from intermittent weather conditions, which change
continuously up to the delivery time. At the same time, the aggregated
generation and consumption volumes are affected by a wider range of factors,
which include both weather conditions, social events, trading strategies and
conventional power plant outages. These two types of risks, price and volume,
are closely connected with each other. In particular, the aggregated volume
risk impacts the variability of electricity prices. This property has recently
attracted attention and has been discussed in the literature (see Ketterer
(2014), Rintamäki et al. (2017), Maciejowska (2020) among others). The
research indicates that wind and solar forecast errors impacts both the
variance and the whole distribution of electricity prices and are one of the
major factors influencing the spread between the day-ahead and intraday prices
(Kiesel and Paraschive (2017), Spodniak et al. (2021)).

The exposure to price and volume risks leads to a rise of income uncertainty
and hence increase the need for appropriate risk management. The reduction of
the revenue risk can be obtained in various ways. Kath et al. (2021) show that
generators can sign a contract with a trading company, which will allow them
to sell all the produced electricity at, for example, day-ahead price and
therefore limit its trade risk. On the other hand, Maciejowska et al. (2019)
demonstrate that generators can reduce their price uncertainty by an active
trade on two markets: day-ahead and intraday. The results indicate that model
based choice of the market, which offers a higher price, can increase the
revenues of the utility and reduce its risk. Similarly, Kath and Ziel (2018)
show that the choice between different types of contracts in an intraday
market (continuous vs. auctions) can be profitable and lead to considerable
financial gains.

Although the portfolio management seems to be of great importance for
practitioners, it has not been studied much in the electricity market
literature. Most of the articles address only one source of risk, price or
volume. An exception is a paper of Faria and Fleten (2011), which proposes a
model of bidding strategy for a hydro power plant, which takes into account
the stochastic nature of both market prices and generation. Therefore, the
main goal of this research is to fill the existing literature gap. In this
paper, Structural Vector Autoregressive (SVAR) model is proposed, which allows
to analyze jointly different types of risk and hence considers the input
uncertainty (G.Corlu et al., 2020). Four sources of uncertainty are
considered: weather conditions, demand shocks and unpredictable behavior of
market participants in day-ahead and intraday markets (called speculative
shocks). It is shown that the SVAR model can be used for forecasting and
simulating of the next day revenue distribution. The outcomes are finally
employed for selection of the optimal portfolio weights. Similar to
Maciejowska et al. (2019), the resulting trade strategies are evaluated with
two types of financial measures: average revenue and associated Value at Risk
(VaR).

Vector Autoregressive (VAR) models, although widely used in time series
analysis, has not been explored much in modeling electricity markets. Silva et
al. (2012) apply Structural VAR (SVAR) to analyze the relationship between the
economic growth and electricity prices. Bernstein and Madlener (2015) used
yearly data to build a Vector Error Correction model (VECM) to asses the price
elasticity of electricity. In both articles, the macroeconomic approach is
adopted and low frequency data are analyzed. The higher frequency information,
with hourly or daily resolution, has been explored by Maciejowska (2014),
Paschen (2016), Spodniak et al. (2021). Spodniak et al. (2021) use VAR models
to asses the relevance of different short-term markets, such as day-ahead,
intraday and balancing market. They show that due to an increase of the wind
power share in the generation mix, the markets closer to the delivery are
becoming more important. In Maciejowska (2014) and Paschen (2016) the impact
of different market shocks on day-ahead electricity prices is described with
SVAR model. Paschen (2016) uses the estimates of SVAR model to obtain impulse
response functions and to analyze dynamic interrelations between spot prices
and RES power. Finally, Maciejowska (2014) shows that speculative shocks,
defined as an unpredictable behavior of electricity traders, have the largest
share in the electricity price variance. Up to my knowledge, in the literature
there are not publications exploring the potential usage of SVAR model in the
decision process of an electricity generator.

The article is structured as follows. Section 2 describes briefly the data
used in the analysis. In Section 3, a SVAR model of electricity market is
presented, which is next applied to predict a revenue distribution and to
support the decision process of a RES utility, Section 4. Section 5 presents
the results of the experiment and a statistical comparison of performance of
proposed trading strategies. Finally, Section 6 concludes.

##  2 Data

In this research, German electricity market is considered. The analysis are
based on data published by TSOs and EPEX exchange and cover the period from
01-Oct-2015 to 31-Sep-2019. Since Austria separated from the German the
bidding zone, only the data on German generation level and structure is used.
The data sets consist of day-ahead ($DA_{th}$) and intraday ($ID_{th}$)
prices, with the later one being described by an ID3 index (volume weighted
prices from the last 3 hours of trade). The electricity prices are
complemented by information on actual levels and system forecasts of
fundamental variables: the total load ($L_{th}$) and RES generation
($RES_{th}$) generation. In the remaining part of the paper, the index $h$
stands for an hour and $t$ for a day.

The properties of times series describing the German electricity market have
been extensively studied in the literature. It is well documented that both
total load and electricity prices have a strong seasonal pattern. They are on
average the highest on working days, during peak hours. The exemplary time
plots of day-ahead and intraday prices together with total load and RES
generation are presented in Fig 1. Two hours are shown, $h=4,18$, which
represent the peak and the off-peak periods of a day. They confirm the weekly
and yearly seasonal behavior of the electricity generation. Additionally, it
can be observed that day-ahead and intraday prices co-move together, with the
intraday prices being more volatile. Finally, RES generation reveals different
fluctuation pattern, with only minor differences between night and day hours
and high variability.

The statistical properties of the data are presented in Table 1. First, the
mean and the standard deviation of the variables of interest are computed
separately for each hour. Table 1 shows their average values across the day.
It could be observed that the DA prices are on average slightly lower than ID
prices, but also less volatile. Finally, the results of Augmented Dicky-Fuller
(ADF) test for the presence of unit roots is presented. The table reports the
number of hours, for which the test allows to reject the null and hence
confirms the stationarity of the series. It can be observed that majority of
analyzed data are stationary. This could be due to a relatively short period
of time – four years – used for analysis. The result is important for the
modeling approach presented below, as it supports an usage of a VAR model for
the data in levels.

Table 1: Statistical properties of the data

| RES | Load | DA | ID
---|---|---|---|---
Mean | 16.40 | 61.97 | 36.11 | 36.24
St. Dev | 8.82 | 8.18 | 15.61 | 17.06
ADF | 24 | 18 | 24 | 24

Notice: ADF indicates the number of hours for which the test rejects the null
of an unit root.

(a) Day-ahead prices

(b) Intraday prices

(c) RES

(d) Load

Figure 1: Time paths of electricity prices, RES generation and Load for hours
4:00 and 18:00

##  3 SVAR model of electricity market

In this research, VAR model is used to describe the joint behavior of
electricity prices and generation, which in turn impact the revenue of the a
small RES utility. The literature (see Weron (2014) for a review) indicates
that the electricity market has a strong daily seasonality, which impacts not
only the level of prices and generation but also its dynamics. Therefore each
hour is typically modeled separately (see Ziel and Weron (2018) for a
discussion). Here, this implies that the VAR model is built and estimated
independently for each hour, $h$.

In order to capture the major stochastic shocks influencing the variables of
interest, a four dimension model is applied:

$Y_{t,h}=A_{0,h}X_{t,h}+\sum_{p=1}^{7}A_{p,h}Y_{t-p,h}+\varepsilon_{t,h},$ (1)

where $Y_{th}=(RES_{t,h},L_{t,h},DA_{t,h},ID_{t,h})^{\prime}$ is a vector of
endogenous variables. Notice that the vector $Y_{t,h}$ could be extended to
include a generation of a particular RES utility, $G_{t,h}$, if such
information is accessible. Then,
$Y_{t,h}=(G_{t,h},RES_{t,h},L_{t,h},DA_{t,h},ID_{t,h})^{\prime}$. Since in
this research no such data are available, it is assumed that an exemplary RES
producer owns a few wind and solar farms, which are spread across the Germany.
For simplicity, it is assumed that $G_{t,h}$ is proportional to the aggregated
level of RES generation, with $G_{t,h}=\rho RES_{t,h}$. As the result,
$RES_{t,h}$, represents in the model both, the total and the individual
generation of RES.

The vector $X_{t,h}$ consists of dummy variables indicating week days, TSO
forecasts of RES generation and the total load. Additionally, some nonlinear
transformation of lagged prices (similar toUniejewski et al. (2019), Ziel and
Weron (2018), Maciejowska et al. (2020)) are used: $DA_{min,t-1}$ (a minimum
price from the last day), $DA_{max,t-1}$ (a maximum price from the last day)
and $DA_{t-1,24}$ (the last known price). The parameters $A_{0,h}$ and
$A_{p,h}$ are matrices of dimensions $(4\times 12)$ and $(4\times 4)$,
respectively. They are estimated separately for each hour, $h$. The residuals,
$\varepsilon_{t,h}$, are $(4\times 1)$ random vectors with zero mean and a
variance-covariance matrix $\Sigma_{h}$. Since (1) is a reduced form of the
VAR model, then the residuals are allowed to be cross-correlated and hence the
$\Sigma_{h}$ is not diagonal.

The reduced form model, although useful for point forecasting of endogenous
variables, is not suited for risk analysis – for example impulse responses or
simulations. Therefor, a structural extension of the model (1) is used. The
SVAR model could take various forms (see Lütkepohl (2005)). Here, the B-model
is adopted, which focuses on residuals of model (1). The SVAR model assumes
that the within-sample errors, $\varepsilon_{t,h}$, are a linear
transformation of structural shocks, $u_{t,h}$, which are uncorrelated and
have a diagonal variance-covariance matrix, $\Lambda_{h}$. The relationship
between reduced form and structural shocks is described by following equation:

$\varepsilon_{t,h}=B_{h}u_{t,h}.$ (2)

The matrix $B_{h}$ is called _an instantaneous effect matrix_ , because it
describes how structural shocks affect endogenous variables in the current
time period, $t$. For example, $B_{23,h}$ describes the impact of the third
structural shock, $u_{3,th}$ on the second element of $Y_{t,h}$, which is the
total load. Notice that equation (2) implies that
$\Sigma_{h}=B_{h}\Lambda_{h}B_{h}^{\prime}$, so there is a direct relationship
between $B_{h}$ and the variance of errors $\varepsilon_{t,h}$. It is
typically assumed that either structural shocks have an identity variance-
covariance matrix, $\Lambda_{h}=I$, or the diagonal elements of $B_{h}$ are
equal to one. Here, the first approach is adopted.

Since the structural shocks are assumed to be uncorrelated, which in the
Gaussian framework implies independence, their behavior is much easier to
model and predict. They could be also simulated separately. Unfortunately, the
SVAR model cannot be directly estimated due to lack of identification. As
discussed in the literature (see Lütkepohl (2005) for a comprehensive
discussion on VAR models), the structural model requires estimation of $K^{2}$
elements of the $B$ matrix, where $K$ is the number of endogenous variables.
In the current study, $K=4$, hence the structure is defined by $16$
parameters. At the same time, the variance-covariance matrix $\Sigma$ of the
reduced form, due to its symmetry, consists of only $K(K+1)/2=10$ parameters.
This implies that there is not enough information to identify the structural
parameters. As the result, additional $K(K-1)/2$ assumptions need to be
imposed, which will restrict the parameter space. In the presented model, six
identification restrictions are need to ensure model indetifiability.

In this research, four structural shocks are considered: weather shock
$u_{1,th}$, demand shock $u_{2,th}$, day-ahead speculative shock $u_{3,th}$
and intraday speculative shock $u_{4,th}$. The energy market has its
particular features, which helps to recover the structure of SVAR model.
First, due to the dispatch priority and support schemes, RES generation does
not depend neither on demand nor on price shocks. Second, the literature
indicates the limited price elasticity of demand, because market participants
require time to adjust their production to the market situation. In
particular, the demand response to the unpredicted price innovations is
assumed to be insignificant. Finally, the spot pries are set day-ahead, before
the intraday trade rises and hence it could be assumed that they do not depend
on the intraday speculative shocks. As the result, the instantaneous effect
matrix becomes lower triangular

$B_{h}=\left[\begin{array}[]{cccc}*&0&0&0\\\ &*&0&0\\\ &*&*&0\\\ &*&*&*\\\
\end{array}\right]$ (3)

The zeros in the above matrix correspond to the ’no impact’ restrictions.
Notice that in the presented solution, the ordering of endogenous variables
and structural shocks play an important role. The model implies, for example,
that the shock $u_{1,th}$ – which is the weather shock – influences all
variables, whereas $u_{4,th}$ – which is the intraday speculative shock –
impacts only intraday prices.

##  4 Decision problem of a RES utility

In this article, SVAR model presented in Section 3 is used for designing
trading strategies of a RES utility. Additionally, it is assumed that the
generator participates in the day-ahead and the intraday market and is small
enough not to impact directly the market prices. It does not receive FIT and
hence its revenue is related to wholesale prices. The utility needs to place
an order in the day-head market at noon of the day proceeding the delivery.
Due to the stochastic nature of RES, the offered generation differs from the
actual production and therefore it needs to trade the difference between the
scheduled and the final generation in intraday market. The power plant self-
balances its position and therefore no trade is left for the balancing market.

For illustrative purposes, it is assumed that $G_{th}=\rho RES_{th}$. In order
to ensure that RES producer has a minor influence on market prices, the $\rho$
is set to equal $0.5\%$. This implies that in the year 2018, it generated
around 85 MWh per hour and accounted for 0.2% of the total electricity
production in Germany.

The utility needs to choose on the day $t-1$, which part, $g$, of the expected
generation, $\hat{G}_{th}$ is offered in the day-ahead market. The remaining
part of the production, $G_{th}-g\hat{G}_{th}$, is soled in the intraday
market. Notice that in general, the variable $g$ may change across days and
hours. Since the utility focuses on the real trade, it is assumed that it does
not speculate and hence $g\in[0,1]$. As the result, the revenue from the trade
becomes

$\pi_{t,h}(g)=g\hat{G}_{t,h}DA_{t,h}+(G_{t,h}-g\hat{G}_{t,h})ID_{t,h}.$ (4)

It should be mentioned here that values of prices and generation are not known
at the moment of taking the decision. They depend on stochastic factors, which
change throughout the day, such as the weather condition and the human
behavior. As the results, the actual level of revenue, $\pi_{t,h}(g)$ becomes
random. It is expected that its distribution is non-normal, as it includes
both the level and the product of a few random variables. Moreover, since the
generation and electricity price forecast errors are allowed to be correlated,
the expected level of revenue will generally by different than
$E(\pi_{th}(g))\neq
g\hat{G}_{th}E(DA_{th})+(E(G_{th})-g\hat{G}_{th})E(ID_{th}).$

###  4.1 Trading strategies

The utility places its order in DA market on the day preceding the deliver and
this transaction cannot be changed, as the new information on weather
conditions and prices arrive. Therefore, its revenue depends directly on the
amount of electricity soled there, $g\hat{G}_{t,h}$. As mentioned before, $g$
does not have to be constant and may adjust to the market situation. Therefore
in the remaining part of the paper, it will be indexed with a day ($t$) and an
hour ($h$): $g_{t,h}$.

####  4.1.1 Simple DA and ID strategies

First two simple approaches are considered: _day ahead strategy_ , which
assumes that $g_{t,h}=1$ and _intraday_ strategy, for which $g_{t,h}=0$. They
are two boundary approaches, which assume a fixed value of the share $g_{t,h}$
that does not depend on market conditions. In the first case, all the
predicted production is soled on DA, whereas in the second case the utility
decides to wait with the trade till the next day and leave all the generation
for the intraday market.

####  4.1.2 Data driven strategies

Next, _data driven strategies_ of choosing the level of $g_{t,h}$ are
proposed, which utilize the estimated structure of the forecast errors. Using
a bootstrap simulation, the optimal proportion of generation offered in DA
market is selected to either maximize the expected revenue ($E\pi(g_{t,h})$)
or to minimize the risk. Here, the risk is evaluated with two measures: the
Sharp Ratio ($SR(g_{t,h})$) and Value-at-Risk ($VaR(g_{t,h})$), which have
been shown relevant and useful by the financial literature.

In this research, a bootstrap procedure is proposed to generate the
distribution of $Y_{t,h}$. This allows to approximate the distribution of the
revenue for different levels of $g$ and to optimize its value. The algorithm
consists of the following steps

  1. 1.

For a selected hour $h$ and the calibration window
$\\{t_{0},t_{0}+1,...,t-1\\}$ estimate the parameters of SVAR model:
$\hat{A}_{p,h}$ for $p=0,1,...,7$ and the instantaneous effect matrix
$\hat{B}_{h}$

  2. 2.

Calculate the point forecasts of $\hat{Y}_{t,h}$ according to (1) as

$\hat{Y}_{t,h}=\hat{A}_{0,h}X_{t,h}+\sum_{p=1}^{7}\hat{A}_{p,h}Y_{t-p,h}$ (5)

and $\hat{G}_{t,h}$ as $\hat{G}_{t,h}=\rho\hat{RES}=\rho\hat{Y}_{1,th}$.

  3. 3.

Compute residuals $\hat{\varepsilon}$ of a reduced form model (1) and
corresponding structural shocks $\hat{u}=\hat{B}_{h}^{-1}\hat{\varepsilon}$

  4. 4.

Approximate the distribution of the next day $Y_{t,h}$ and profits
$\pi_{t,h}(g)$ using a bootstrap sampling of structural shocks. For each
iteration $b=1,...,B$

     * •

Pick independently a realization of each shock and obtain a $(4\times 1)$
vector $\tilde{u}^{(b)}$

     * •

Transform structural shocks into forecast errors:
$\tilde{\varepsilon}^{(b)}=\hat{B}_{h}\tilde{u}^{(b)}$

     * •

Calculate $\tilde{Y}_{t,h}^{(b)}=\hat{Y}_{t,h}+\tilde{\varepsilon}^{(b)}$,
where
$\tilde{Y}_{t,h}^{(b)}=(\tilde{RES}^{(b)}_{t,h},\tilde{L}^{(b)}_{t,h},\tilde{DA}^{(b)}_{t,h},\tilde{ID}^{(b)}_{t,h})^{\prime}$,
and $\tilde{G}_{t,h}^{(b)}=\rho\tilde{RES}^{(b)}_{t,h}$

     * •

For the chosen level of $g$ compute the corresponding profit
$\tilde{\pi}_{t,h}^{(b)}(g)$ according to (4)

  5. 5.

Estimate the expected value and the variance of incomes, $E\pi_{t,h}(g)$ and
$\sigma_{t,h}^{2}(g)$ as the mean and the mean squared deviation of
$\tilde{\pi}_{t,h}^{(b)}(g)$ across $b$. Using this information, compute the
Sharp ratio as $SR_{t,h}(g)=E\pi_{t,h}(g)/\sigma_{t,h}(g)$ and approximate
$VaR_{t,h}^{\tau}(p)$ by a $\tau$ quantile of $\tilde{\pi}^{(b)}_{t,h}(g)$

The selection of the level of $g_{t,h}$ depends on the optimality condition
considered. The first approach picks $g_{t,h}^{*}$, which maximizes the
expected revenue $E\pi_{t,h}(g)$. This method, although profitable, may result
in an increase of the transaction risk. Therefore, alternative trading
strategies are also examined, which aim at minimizing the risk by choosing
$g_{t,h}^{*}=argmax_{g}{SR_{t,h}(g)}$ or
$g_{t,h}^{*}=argmax_{g}{VaR_{t,h}^{\tau}(g)}$. The brief description of the
approaches together with their notation are presented in Table 2.

Table 2: Trading strategies Name | Description | Optimality criteria | $g$
---|---|---|---
$S_{DA}$ | Sell the expected generation on DA market | - | 1
$S_{ID}$ | Sell all the generation on Intraday market | - | 0
$S_{E\pi}$ | Maximize expected profit | $\max E\pi(g)$ | $g^{*}_{t,h}$
$S_{SR}$ | Maximize Sharp ratio | $\max SR(g)$ | $g^{*}_{t,h}$
$S_{VaR}$ | Maximize Value of Risk for the $5\%$ quantile | $\max VaR_{0.05}(g)$ | $g^{*}_{t,h}$

###  4.2 Evaluation of trading strategies

The trading strategies presented above can be compared according to various
dimensions. Here, three features are used to evaluate their performance: the
level of revenue, its predictability and variability.

####  4.2.1 The level of income

In this article, the average hourly revenue is used to analyze the performance
of presented approaches. It is computed as follows:

$\bar{\pi}=\frac{1}{T_{eval}\times 24}\sum_{t,h}\pi_{t,h}(g_{t,h}),$ (6)

where $T_{eval}$ is the number of days used for evaluation. In order to
verify, if a chosen strategy, $i$, yields a higher average income than the
starategy $j$, a new variable is defined

$d_{t,h}=\pi_{t,h}(g^{(i)}_{t,h})-\pi_{t,h}(g^{(j)}_{t,h}).$

If both strategies are characterized by the same expected value then
$Ed_{t,h}=0$. Hence, the natural hypothesis in this setup are
$H_{0}:Ed_{t,h}=0$ and $H_{1}:Ed_{t,h}>0$. Under the null, both strategies
provide the same average revenue, whereas under the alternative, the strategy
$i$ is more profitable than $j$. In order to verify the hypothesis, a Diebold
and Mariano (1995) type testing procedure is applied to the mean daily level
of $d_{t,h}$: $d_{t}=1/24\sum_{h=1}^{24}d_{t,h}$. The variable $d_{t}$ can be
viewed as a counterpart of the loss differential, for example a difference of
squared forecast errors. The statistic takes the form

$DM=\frac{\bar{d}}{\sqrt{2\pi\hat{f}_{d}(0)/T_{eval}}},$ (7)

where $\bar{d}$ is the average value of $d_{t}$ and $\hat{f}_{d}(0)$ is an
estimator of spectral density of $d_{t}$ at frequency 0. As
$2\pi\hat{f}_{d}(0)/T_{eval}$ is an estimator of the variance of $\bar{d}$,
$DM$ converges asymptotically to a standard normal distribution, $N(0,1)$.

####  4.2.2 Risk

The risk associated with a strategy is evaluated according to two features:
possibility to predict accurately the next day revenue and the variability of
the income. In order to asses the revenue forecast quality, we use the
outcomes of a SVAR model. The revenue predictions are calculated as the
average of $\tilde{\pi}^{(b)}(g)$ across bootstrap iterations:

$\hat{\pi}_{t,h}=\frac{1}{B}\sum_{b=1}^{B}\tilde{\pi}^{(b)}_{t,h}(g_{t,h})$
(8)

First, it is checked, how much the predicted revenue differs from its actual
values, $\pi_{t,h}(g_{t,h})$. In this research, similar to other EFP papers,
the forecast quality is evaluated with Root Mean Squared Error (RMSE) and Mean
Absolute Error (MAE)

$RMSE=\sqrt{\frac{1}{T_{eval}\times
24}\sum_{t,h}(\pi_{t,h}-\hat{\pi}_{t,h})^{2}}.$ (9)

$MAE=\frac{1}{T_{eval}\times 24}\sum_{t,h}|\pi_{t,h}-\hat{\pi}_{t,h}|.$ (10)

In order to statistically verify if forecasts stemming from two different
strategies, $i$ nad $j$, are equally accurate, the DM test is applied to the
following loss differentials:

$d_{t}=\frac{1}{24}\sum_{h=1}^{24}{(u_{t,h}^{2}(i)-u_{t,h}^{2}(j))}$

for RMSE and

$d_{t}=\frac{1}{24}\sum_{h=1}^{24}{|u_{t,h}(i)-u_{t,h}(j)|}$

for MAE measure. The $u_{t,h}(i)$ and $u_{t,h}(j)$ are forecast errors of
analyzed strategies:
$u_{t,h}(i)=\pi_{t,h}(g^{(i)}_{t,h})-\hat{\pi}_{t,h}(g^{(i)}_{t,h})$ and
$u_{t,h}(j)=\pi_{t,h}(g^{(j)}_{t,h})-\hat{\pi}_{t,h}(g^{(j)}_{t,h})$. The test
statistics have a form analogous to (7).

Finally, the risk associated with the variability of income is measured by the
Value-at-Risk of revenues for a given hour. In order to aggregate the results,
the average VaR is used to compare the outcomes of different strategies.

##  5 Results

###  5.1 Experiment design

In order to asses different trading strategies, a forecasting experiment is
run, in which a moving window methodology is applied. In the experiment, the
calibration window includes 731 observations and the results are assessed with
the last two years from 01-Oct-2017 to 31-Sep-2019 ($T_{eval}=730$
observations). Moreover, in order to bring the experiment as close as possible
to the empirical problem, different information sets are used for hours before
and after the time of taking the decision. In this work it is assumed that the
generator places on order on DA market at 12:00. Therefore, when performing
predictions it knows the actual generation and intraday prices only for hours
from midnight till 10:00. For the remaining hours, the actual values
($L_{t-1,h}$, $RES_{t-1,h}$) are not known and are replaced by their TSO
predicted levels (see Maciejowska et al. (2021) for more details). In case of
$ID_{t-1,h}$, which is also unavailable, it is assumed that is has no impact
on endogenous variable and hence the last column of $A_{1,h}$ matrix is set to
equal zeor. Finally, following the energy forecasting literature (Ziel (2016),
Uniejewski et al. (2019), Maciejowska et al. (2020)), three lags are selected
in the VAR model: $p=1,2,7$.

###  5.2 Comparison of trading strategies

Let us first look at the aggregated results for each trading strategy. The
outcomes are reported in Table 3, in which the second column presents the
average hourly revenue, $\bar{\pi}$. The remaining three columns focus on risk
that is measured by $RMSE$ and $MAE$ of revenue forecasts and $VaR_{1\%}$ and
$VaR_{5\%}$ of income. It should be noticed that the benchmark results of
$S_{DA}$ strategy are expressed in nominal values. The remaining outcomes are
presented relative to the benchmark, as the percentage difference
($\%\Delta$). Hence for all the strategies, apart from $S_{DA}$, values lower
than zero indicate that the given strategy reduces the indices and larger then
zero prove the opposite.

The analysis of simple strategies, which assume a constant level of $g_{t,h}$,
confirm that the transactions made in ID market are slightly more profitable
then in DA market. However, the rise of the average revenue in $S_{ID}$ is
achieved at the cost of a higher risk. The revenue forecasts suffer from the
loss of accuracy (RMSE increases by 9% and MAE by 7.5% as compared to
$S_{DA}$). Moreover, the income itself is substantially more volatile – the
$VaR_{1\%}$ drops by almost then 40%.

At the same time, the data driven approaches provide results characterized by
a higher income and lower risk than the benchmark. As expected, the $S_{E\pi}$
strategy allows to earn the largest revenue among all the strategies. It
reduces RMSE but its impact on risk is mixed: it lowers $VaR_{1\%}$ and
increases $VaR_{5\%}$. The $S_{VaR}$ strategy, which aims at maximizing the
$VaR_{5\%}$ of revenue, leads to a rise of income by less than 1%, which is
the weakest result among the data driven approaches. It provides more accurate
predictions than the simple $DA$ strategy and reduces RMSE and MAE by more
than 6%. Similar to $S_{E\pi}$ case, its impact on risk is ambiguous.

Finally, the strategy based on the Sharp Ratio provides the most promising
results. It allows to increase the revenue and at the same time reduces the
risk measure by both RMSE and VaR. The decrease of uncertainty is substantial,
as it reduces RMSE and MAE by 8% and 6.7%, respectively. Moreover, it
increases 1% and 5% VaR by 19.5% and 10.3%.

Fig. 2 presents the asymptotic $p$-values of test statistics used to compare
the profitability and Fig. 3 depicts comparison of RMSE (left panel) and MAE
(right panel) forecast accuracy measures. It should be recalled that if the
$p$-value of $DM_{i,j}$ statistic is smaller than 10%, then it implied that
the strategy $i$ provides a higher revenue or more accurate predictions than
the strategy $j$. The results confirm previous outcomes that all strategies
provide higher income than the benchmark, $S_{DA}$. Moreover, $S_{E\pi}$
(denoted on plots as Profit), brings the largest revenue among the presented
approaches. When the RMSE and MAE of forecast errors are considered, the
outcomes show that $S_{DA}$ allows to predict the next day revenue more
accurately than $S_{ID}$, which is the worst strategy according to this
measure. Finally, the data driven methods aiming at minimizing the risk leads
to a significant reduction of both forecast accuracu measures.

To sum up, the aggregate results indicate that using bootstrap method based on
SVAR model for forecasting generation and constructing trading strategies
could bring additional profit and at the same time reduce the risk. Hence, it
it is preferable from the perspective of a small RES utility. Among presented
approaches, $S_{SR}$ is the most attractive one. It does not only increase
substantially VaR but brings on average 41.16 EURO more revenue per hour than
the benchmark, which is equivalent to an increase of profit within two years
by 721 117.59 EURO.

Table 3: The average hourly revenue and risk measures Strategy | $\bar{\pi}$ | $RMSE$ | $MAE$ | $VaR_{1\%}$ | $VaR_{5\%}$
---|---|---|---|---|---
$S_{DA}$ | 3215.6 | 1196.0 | 635.5 | -1803.1 | 758.7
| $\%\Delta$
$S_{ID}$ | 1.20 | 8.95 | 7.46 | -39.96 | -1.04
$S_{E\pi}$ | 1.84 | -3.24 | 0.66 | -4.13 | 7.07
$S_{SR}$ | 1.28 | -8.00 | -6.67 | 19.51 | 10.29
$S_{VaR}$ | 0.91 | -6.50 | -6.29 | -4.44 | 6.09
Figure 2: Results of equal profitability test: $p$-values

(a) RMSE

(b) MAE

Figure 3: Results of equal forecast accuracy test for two loss functions(RMSE,
MAE): $p$-values

###  5.3 Distribution of $g$

Table 4: The averaged proportion of generation offered on DA market: % points Strategy | $\bar{g}$ | $g=0$ | $0<g<1$ | $g=1$
---|---|---|---|---
$S_{DA}$ | 100.0 | 0.000 | 0.000 | 100.0
$S_{ID}$ | 0.000 | 100.0 | 0.000 | 0.000
$S_{E\pi}$ | 49.05 | 50.95 | 0.000 | 49.05
$S_{SR}$ | 51.36 | 12.47 | 84.46 | 3.08
$S_{VaR}$ | 52.72 | 9.68 | 85.84 | 4.48

Let us now look more into details of the presented outcomes. First, the data
driven approaches allow to choose a proportion of predicted generation, which
is offered in DA market. Table 4 shows the average level of optimal $g_{t,h}$
for different strategies ($\bar{g}$) together with the proportion of its
values in three groups: $g=0$, $0<g<1$ and $g=1$. All the quantities are
expressed in % points.

It is clear that simple strategies provide two boundary outcomes with
$\bar{g}=1$ or $\bar{g}=0$. This is due to the fact that for these two
approaches all the predicted generation is offered either in DA or ID market.
The more complex results are observed for the data driven strategies, for
which the average level of $g$ oscillates around 50%. First, it can be noticed
that in case of $S_{E\pi}$ strategy $g_{t,h}$ takes only the extreme values:
$0$ or $1$. It selects $ID$ market in 50.95%, which confirms the previous
findings that indicate larger profitability of this market.

When the strategies aiming at minimizing the risk, $S_{SR}$ and $S_{VaR}$, are
considered, it can be observed that the algorithm selects $g_{t,h}=0$ and
$g_{t,h}=1$ in around 11% and 4% cases, respectively. Hence, for around 85% of
hours, the proportion of generation soled in DA market falls into the interval
$(0,1)$. The distribution of $0<g_{t,h}<1$ for $S_{SR}$ is depicted in Fig. 4.
The plot shows that in majority of cases, $g_{t,h}$ takes the value between
$0.5$ and $0.8$. Moreover, its distribution is skewed to the right indicating
that the DA market is more attractive than the ID one.

The average levels of $g_{t,h}$ across 24 hours is presented in Fig. 5. The
most risky, $S_{E\pi}$ strategy is characterized by lower values of $g$ and
hence tends to sell a larger share of the generation in ID market than
$S_{VaR}$ and $S_{SR}$. When these two are considered, it can be observed that
the highest average value of $g$ is obtained for hours 15-17, when on average
more than 70% of predicted generation is offered in DA market. On the
contrary, during the evening hours, from 19-24, these strategies suggest
selling majority of the production in the ID market. For night and early
morning hours the results are mixed, as the average value of $g$ is close to
50%.

Figure 4: Frequences of $g$ values such that $0<g_{t,h}<1$ for strategy
$S_{SR}$ Figure 5: Average value of $g_{t,h}$ across hours

##  6 Conclusions

The changes in the electricity markets expose RES generators to various risks,
among which the price and the volume risk play a very important role. RES
generators, which revenue depends on the market prices and the offered
quantity, can now actively build a portfolio from different types of
contracts. In this research, it is assumed that it trades produced electricity
either in day-ahead or intraday market. It can be noticed that due to an
intermittent generation and stochastic electricity prices, the entrepreneur
acts under strong uncertainty.

In this research, it is assumed that the trading portfolio constructed by the
utility depends on two major factors: predicted level of production and the
chosen share of generation, $g_{t,h}$, offered in the DA market. Since the
utility does not speculate, it is assume that $g_{t,h}$ belongs to an interval
$(0,1)$. As a result, two types of trading strategies are considered: simple
strategies – which assume a fixed value of $g_{t,h}$ – and data driven
strategies. In the later one, the SVAR model is used to predict future level
of generation and to select an optimal level of $g_{t,h}$ in order to either
maximize the revenue or minimize the transaction risk.

The performance of presented trading strategies is next compared using the
data from German electricity market. The results of the research indicate that
the transactions in ID market are on average more profitable than in DA
market. They are however burden with a significantly larger risk. Second, the
data driven strategies provide revenues larger than the benchmark, $S_{DA}$,
strategy and at the same time allow to reduce risk measured by the
predictability (RMSE, MAE) of a next day revenue and its variability ($VaR$).
Among the proposed approaches, the strategy maximizing the Sharp Ratio is the
most promising one, as it provides the most robust outcomes.

Finally, the selected shares of predicted generation offered in DA market,
$g_{t,h}$, are analyzed. The results show that the strategy aiming at
maximizing the profit, chooses only the two boundary values of $g_{t,h}$ and
hence offers all the forecasted production in either DA or ID market. On the
contrary, approaches minimizing the risk, are more prone to build the
portfolios from both markets at the same time.

This research provides a promising approach for constructing trading portfolio
of a small RES producer. This comprehensive method allows to explore
information on different aspects of the market: the structure of generation
and electricity prices. The results indicate that using such diversified
information set can lead to a significant increase of profits and a reduction
of a transaction risk. Finally, the proposed approach can be further extended
to allow for other types of contracts or trading strategies.

## Acknowledgments

This work was partially supported through SONATA BIS grant no.
2019/34/E/HS4/00060 and SONATA grant no. 2015/17/D/HS4/00716.

## References

  * Bernstein and Madlener (2015) Bernstein, R., Madlener, R., 2015\.  Short- and long-run electricity demand elasticities at the subsectoral level: A cointegration analysis for German manufacturing industries.  Energy Economics 48, 178–187.
  * Diebold and Mariano (1995) Diebold, F.X., Mariano, R.S., 1995\.  Comparing predictive accuracy.  Journal of Business & Economic Statistics 13, 253–263.
  * Faria and Fleten (2011) Faria, E., Fleten, S.E., 2011\.  Day-ahead market bidding for a Nordic hydropower producer: taking the Elbas market into account.  Computational Management Science 8, 75–101.
  * G.Corlu et al. (2020) G.Corlu, C., Akcay, A., Xiec, W., 2020.  Stochastic simulation under input uncertainty: A review.  Operations Research Perspectives 7, 100162.
  * Gianfreda et al. (2016) Gianfreda, A., Parisio, L., Pelagatti, M., 2016.  The Impact of RES in the Italian Day–Ahead and Balancing Markets.  The Energy Journal 37.
  * Hong et al. (2020) Hong, T., Pinson, P., Wang, Y., Weron, R., Yang, D., Zareipour, H., 2020\.  Energy forecasting: A review and outlook.  IEEE Open Access Journal of Power and Energy 7, 376–388.
  * Kath et al. (2021) Kath, C., Nitka, W., Serafin, T., Weron, T., Zaleski, P., Weron, R., 2021\.  Balancing generation from renewable energy sources: Profitability of an energy trader.  Energies 13, 205\.
  * Kath and Ziel (2018) Kath, C., Ziel, F., 2018.  The value of forecasts: quantifying the economic gains of accurate quarter-hourly electricity price forecasts.  Energy Economics 76, 411–423.
  * Ketterer (2014) Ketterer, J., 2014.  The impact of wind power generation on the electricity price in germany.  Energy Economics 44, 270–280.
  * Kiesel and Paraschive (2017) Kiesel, R., Paraschive, F., 2017\.  Econometric analysis of 15-minute intraday electricity prices.  Energy Economics 64, 77–90.
  * Koch and Hirth (2019) Koch, C., Hirth, L., 2019.  Short-term electricity trading for system balancing: An empirical analysis of the role of intraday trading in balancing Germany’s electricity system.  Renewable and Sustainable Energy Reviews 113\.
  * Koch and Maskos (2020) Koch, C., Maskos, P., 2020.  Passive balancing through intraday trading: Whether interactions between short-term trading and balancing stabilize Germany’s electricity systems.  International Journal of Energy Economics and Policy 10, 101–112.
  * Lütkepohl (2005) Lütkepohl, H., 2005.  New Introduction to Multiple Time Series Analysis.  Springer-Verlag, Berlin.
  * Maciejowska (2014) Maciejowska, K., 2014.  Fundamental and speculative shocks, what drives electricity prices?  2014 11th International Conference on the European Energy Market (EEM).
  * Maciejowska (2020) Maciejowska, K., 2020.  Assessing the impact of renewable energy sources on the electricity price level and variability - a quantile regression approach.  Energy Economics 85, 1–22.
  * Maciejowska et al. (2019) Maciejowska, K., Nitka, W., Weron, T., 2019.  Day-ahead vs. intraday – forecasting the price spread to maximize economic benefits.  Energies 12, 631\.
  * Maciejowska et al. (2021) Maciejowska, K., Nitka, W., Weron, T., 2021.  Enhancing load, wind and solar generation for day-ahead forecasting of electricity prices.  Energy Economics 99.
  * Maciejowska et al. (2020) Maciejowska, K., Uniejewski, B., Serafin, T., 2020.  PCA forecast averaging – predicting day-ahead and intraday electricity prices.  Energies 12, 3530\.
  * Pape et al. (2016) Pape, C., Hagemann, S., Weber, C., 2016.  Are fundamentals enough? Explaining price variations in the German day-ahead and intraday power markets.  Energy Economics 54, 376–387.
  * Paschen (2016) Paschen, M., 2016.  Dynamic analysis of the german day-ahead electricity spot market.  Energy Economics 59, 118–128.
  * Rintamäki et al. (2017) Rintamäki, T., Siddiqui, A.S., Salo, A., 2017.  Does renewable energy generation decrease the volatility of electricity prices? An analysis of Denmark and Germany.  Energy Economics 62, 270–282.
  * Silva et al. (2012) Silva, S., Soares, I., Pinho, C., 2012.  The impact of renewable energy sources on economic growth and CO2 emissions - a SVAR approach.  European Research Studies Journal 15, 133–144.
  * Spodniak et al. (2021) Spodniak, P., Ollikka, K., Honkapuro, S., 2021.  The impact of wind power and electricity demand on the relevance of different short-term electricity markets: The Nordic case.  Applied Energ doi:https://doi.org/10.1016/j.apenergy.2020.116063.
  * Uniejewski et al. (2019) Uniejewski, B., Marcjasz, G., Weron, R., 2019.  Understanding intraday electricity markets: Variable selection and very short-term price forecasting using LASSO.  International Journal of Forecasting 35, 1533 – 1547.
  * Weron (2014) Weron, R., 2014.  Electricity price forecasting: A review of the state-of-the-art with a look into the future.  International Journal of Forecasting 30, 1030–1081.
  * Ziel (2016) Ziel, F., 2016.  Forecasting electricity spot prices using LASSO: On capturing the autoregressive intraday structure.  IEEE Transactions on Power Systems 31, 4977–4987.
  * Ziel and Weron (2018) Ziel, F., Weron, R., 2018.  Day-ahead electricity price forecasting with high-dimensional structures: Univariate vs. multivariate modeling frameworks.  Energy Economics 70, 396–420.