# Photovoltaic power forecasting using quantum machine learning

Asel Sagingalieva  Terra Quantum AG, 9000 St. Gallen, Switzerland  Stefan
Komornyik  HAKOM Time Series GmbH, 1230 Vienna, Austria  Arsenii Senokosov
Terra Quantum AG, 9000 St. Gallen, Switzerland  Ayush Joshi  Terra Quantum AG,
9000 St. Gallen, Switzerland  Alexander Sedykh  Terra Quantum AG, 9000 St.
Gallen, Switzerland  Christopher Mansell  Terra Quantum AG, 9000 St. Gallen,
Switzerland  Olga Tsurkan  Terra Quantum AG, 9000 St. Gallen, Switzerland
Karan Pinto  Terra Quantum AG, 9000 St. Gallen, Switzerland  Markus Pflitsch
Terra Quantum AG, 9000 St. Gallen, Switzerland  Alexey Melnikov  Terra Quantum
AG, 9000 St. Gallen, Switzerland

###### Abstract

Predicting solar panel power output is crucial for advancing the energy
transition but is complicated by the variable and non-linear nature of solar
energy. This is influenced by numerous meteorological factors, geographical
positioning, and photovoltaic cell properties, posing significant challenges
to forecasting accuracy and grid stability. Our study introduces a suite of
solutions centered around hybrid quantum neural networks designed to tackle
these complexities. The first proposed model, the Hybrid Quantum Long Short-
Term Memory, surpasses all tested models by over $40\%$ lower mean absolute
and mean squared errors. The second proposed model, Hybrid Quantum Sequence-
to-Sequence neural network, once trained, predicts photovoltaic power with
$16\%$ lower mean absolute error for arbitrary time intervals without the need
for prior meteorological data, highlighting its versatility. Moreover, our
hybrid models perform better even when trained on limited datasets,
underlining their potential utility in data-scarce scenarios. These findings
represent a stride towards resolving time series prediction challenges in
energy power forecasting through hybrid quantum models, showcasing the
transformative potential of quantum machine learning in catalyzing the
renewable energy transition.

##  I Introduction

Electricity generation prediction, especially for photovoltaic (PV) systems,
is a crucial tool for renewable energy adoption sabri2022novel ; li2020review
. The global economy must radically reduce emissions to stay within the
$1.5\degree$C pathway (Paris Agreement) and the transition to renewable energy
sources is necessary to achieve these objectives ellabban2014renewable .
According to the IEA, solar PV’s installed power capacity is poised to surpass
that of coal by 2027, becoming the largest in the world.

Accurate PV power forecasts are vital for multiple facets of the energy
industry such as long-term investment planning, regulatory compliance for
avoiding penalties, and renewable energy management across storage,
transmission, and distribution activities. Several studies show that an
increase in forecasting accuracy reduces electricity generation from
conventional sources. Increased accuracy also reduces operating costs of
systems through reducing the uncertainty of PV power generation anido2016value
. They support improving the stability and sustainability of the power grid
through optimizing power flow and counteracting solar power’s intermittent
nature lee2021pv . Such predictions are foundational in increasing the
economic viability and improving the adoption of solar energy as they inform
pricing and economic dispatch strategies, bolster competitiveness and over
time reduce reliance on reserve power. Additionally, they assist in managing
energy storage effectively and integrating PV systems into the power grid
bozorgavari2019two , which is essential for the enduring success of renewable
energy solutions diagne2013review .

Traditional methods for predicting PV power have primarily relied on
statistical models, machine learning algorithms, or a blend of both
ahmed2019review . These approaches encompass a diverse toolkit, ranging from
time series forecasting and artificial neural networks sabri2022novel ;
gomez2020photovoltaic ; kaloop2021novel , to support vector machines
kulkarni2018output ; shi2012forecasting , k-nearest neighbor methods
rocha2022computational , and random forest models inman2013solar . However,
the intermittent and non-linear nature of solar power generation, influenced
by a wide range of meteorological factors, poses a significant challenge to
the performance of these conventional models blaga2019current .

In light of these challenges, quantum machine learning (QML) emerges as a
promising avenue. This rapidly evolving field, which melds the principles of
quantum mechanics with classical machine learning Schuld_2014 ;
biamonte2017quantum ; Dunjko.Briegel.2018 ; qml_review_2023 , can offer
enhanced capabilities for improving the forecasting accuracy of time series
tasks emmanoulopoulos2022quantum , including PV power generation
sushmit2023forecasting . QML’s potential arises from quantum features like
superposition and entanglement, promising exponential speedups in certain
tasks schuld2018implementing . Moreover, QML algorithms produce inherently
probabilistic results, aptly suited for prediction tasks, and they may
potentially function within an exponentially larger search space, amplifying
their efficacy Lloyd2013 ; 365700 ; Lloyd1996UniversalQS ; schuld2020circuit .
Nonetheless, implementing quantum algorithms bears its own set of challenges,
such as the need for error correction and sensitivity to external interference
Aaronson2015 . Yet, in spite of these challenges, hybrid quantum-classical
models, especially hybrid quantum neural networks (HQNNs), have showcased
their potential in diverse industrial realms, including healthcare
jain2022hybrid ; pharma ; steatosis_sagingalieva , energy
kurkin2023forecasting ; sushmit2023forecasting , aerospace
rainjonneau2023quantum , logistics haboury2023supervised  and automotive
sagingalieva2023hyperparameter  industries.

Figure 1: The input to a hybrid quantum model is presented as a chronological
data table, documenting hourly meteorological parameters including ambient
temperature ($T_{a}$), module temperature ($T_{m}$), and solar irradiance
($I_{3}$, $I_{15}$), alongside the mean PV power output ($P$). The model is
designed to leverage this data to generate predictions of PV power output for
a short-term forecast aimed at near-future output, typically the next hour,
and a long-term forecast that extends to a broader temporal horizon.

In this article, we present three types of hybrid quantum models as potential
solutions for PV power forecasting. We assess the performance of our proposed
models using a publicly accessible dataset, encompassing a comprehensive array
of meteorological variables as well as hourly mean PV power measurements
spanning a $21$-month period. This dataset, along with the data preprocessing
and analytical methodologies employed, is described in detail in Section II.1.

Our first proposed HQNN architecture, articulated in Section II.2,
incorporates classical fully connected layers with a vanilla variational
repetitive quantum layer (VVRQ). Our second model, delineated in Section II.3,
constitutes a hybrid quantum adaptation of the classical recurrent neural
network, termed the Hybrid Quantum Long Short-Term Memory with quantum depth-
infused layer (HQLSTM). While the first two models can predict the power for a
certain hour ahead, the third model, presented in Section II.4, a Hybrid
Quantum Sequence-to-Sequence Neural Network with quantum depth-infused layer,
HQSeq2Seq, after training is capable of forecasting PV power for arbitrary
time intervals without requiring prior meteorological data.

Remarkably, despite having fewer parameters, our hybrid quantum models
outperform their classical counterparts in terms of more accurate predictions,
including trained on a reduced dataset. We summarize our conclusions and
outline future research directions in Section III.

##  II Results

The application of HQNNs in addressing time series prediction challenges,
specifically in forecasting PV power output offers several advantages.
Primarily, their capability to operate within an exponentially larger
computational search space enables them to efficiently capture intricate data
patterns and relationships kordzanganeh2023exponentially . This feature not
only enhances forecast accuracy biamonte2017quantum ; senokosov2023quantum
but also streamlines the learning process, requiring fewer iterations for
model optimization asel1 . Furthermore, the inherent capacity of quantum
technologies to manage the uncertainty and noise ubiquitous in data offers
more resilient and trustworthy predictions schuld2018implementing . This is
particularly pertinent to power forecasting, given the inherent noise in
meteorological data. Recent research also suggests that quantum models can be
represented as partial Fourier series, positioning them as potential universal
function approximators schuld2021effect , thereby broadening their
applicability and efficacy in predictive tasks.

In terms of architecture, an HQNN is an amalgamation of classical and quantum
components. The classical segments may consist of fully connected layers,
convolutional layers, or recurrent layers, while the quantum segments are
typically represented by variational quantum circuits (VQCs) or their
contemporary modifications kordzanganeh2023parallel ; sedykh2023quantum .

Figure 2: (a) Mean and standard deviation of the PV power value for each hour
of the day. (b) Mean and standard deviation of the PV power value for each
month of the year. The plot shows the PV power reaching maximal values in June
and July. (c) Correlation matrix of input features. (d) Joint distribution of
features.

###  II.1 Dataset

To underscore the advantages of hybrid quantum models using empirical
evidence, we selected a publicly accessible dataset malvoni2016data  from a
conventional generation plant situated in the Mediterranean region. This
dataset not only provides comprehensive data but also allows benchmarking with
results from various algorithms available in literature. A comparative
analysis of our model’s predictions and those from the study by
kaloop2021novel  is provided in Section II.5. The dataset, presented as a
numerical table showcased in Fig. 1, encompasses variables like hourly mean
ambient temperature ($T_{a}$), hourly mean module temperature ($T_{m}$),
hourly mean solar irradiance recorded on two tilted planes with tilt angles of
3 and 15 degrees ($I_{3},\,I_{15}$), and hourly mean PV power ($P$) spanning
$21$ months, accounting for more than $500$ days.

Beyond the scope of constructing models for predicting the output of PV
panels, this dataset’s utility extends to other applications. It aids in
planning distributed battery energy storage systems aghaei2019flexibility ,
devising novel energy collection systems kulkarni2018performance , and
researching the degradation patterns of photovoltaic panels malvoni2017study .
The dataset’s multifaceted applicability emphasizes its significance.

To ensure the validity and precision of the data, meticulous preprocessing and
analysis were undertaken. We discovered approximately 20 anomalies in the
original dataset. To maintain a continuous timeline, missing data points were
replaced with the arithmetic mean of the preceding and succeeding day’s
values. Additionally, data corresponding to the date “12/31/13” was excluded
as it contained all-zero values, suggesting an error in data collection. As a
result, we obtained an uninterrupted dataset ranging from 4:00 AM on “3/5/12”
to 12:00 AM on “12/30/13”.

Additional in-depth analysis of the dataset was also conducted for a more
nuanced understanding. Fig. 2(a) delineates the hourly distribution of PV
power across the entire period. As expected, peak PV power values occur during
midday, whereas night time values plummet to zero. Fig. 2(b) portrays monthly
PV power fluctuations, which are more volatile compared to daily patterns,
likely attributable to the limited number of full-year periods in the dataset.
Fig. 2(c) presents a correlation matrix for the dataset features, identifying
solar irradiances $I_{3}$ and $I_{15}$ as the features most correlated with PV
power. Finally, the joint distribution of dataset features depicted in Fig.
2(c) further confirms that solar intensity is the feature most highly
correlated with PV power.

Figure 3: The architectures of: (a) Hybrid Quantum Neural Network with VVRQ
layer, (b) Hybrid Quantum Long Short-Term memory, (c) QDI layer used in HQLSTM
model, (d) Hybrid Quantum Seq2Seq with QDI layer.

###  II.2 HQNN

This section introduces our first proposed model, referred to as the HQNN. As
illustrated in Fig. 1, the model accepts weather data spanning $24$
consecutive hours as its input. The output is a prediction of the PV power for
the upcoming $25$th hour. The HQNN presented at the Fig. 3(a) is a combination
of classical fully-connected layers, in our case with $120$, $17$ and $8$
neurons, and a VVRQ layer, which is a VQC, consisting of $q$ qubits and $d$
repetitions of variational layers, each distinguished by unique weights. The
choice of $120$ neurons is methodical: the model ingests $5$ distinct features
for each of the $24$ hours, resulting in a total of $120=5*24$. The
determination of the remaining parameters stemmed from an extensive
hyperparameter optimization process, detailed in the subsequent sections.

Initially, every qubit in the VVRQ layer is set to the state
$\left|0\right\rangle$. We subsequently encode the classical data by
converting it into rotation angles around one of the $X$, $Y$, $Z$ axes using
$R_{x}$, $R_{y}$, $R_{z}$ gates respectively. This conversion employs the
angle embedding technique PhysRev.70.460 . For each qubit, the rotation angle,
denoted by $x_{j}$, is determined by the $j$-th component of the input vector.

Following this, the variational layer is applied, which can either utilize
“basic” or “strong” entanglements. For the “basic” entanglement, each qubit
undergoes a rotation by an angle $w_{j}^{i}$ around the $X$ axis, subsequently
followed by a layer of CNOT gates Barenco_1995 . Conversely, for the “strong”
entanglement, each qubit is sequentially rotated by the angles
$w^{(Z_{1})}_{ji}$, $w^{(Y_{2})}_{ji}$, and $w^{(Z_{3})}_{ji}$ around the $Z$,
$Y$, and $Z$ axes, respectively. This sequence is then followed by a layer of
CNOT gates. In both cases, the variables $i$ and $j$ play crucial roles in
determining the operations. The variable $i$ signifies the particular wire to
which the operation is applied, and it takes values from the set
${1,2,\dots,q}$. Meanwhile, the variable $j$ represents the number of
variational layers and ranges from ${1,2,\dots,d}$.

Lastly, all qubits are measured in Pauli-$Z$ basis, yielding the classical
vector $v\in\mathbb{R}^{q}$. This output serves as input for a subsequent
classical fully-connected layer. This layer processes information from $q$
neurons into $1$ neuron that predicts the power value.

The proposed HQNN model will be compared with its classical analog – a
Multilayer Perceptron (MLP) that consists of $4$ fully connected layers with
$120$, $32$, $3$, $3$, and $1$ neurons. The number of neurons in each layer
was selected by a hyperparameter optimization procedure, detailed in the
subsequent sections.

###  II.3 HQLSTM

This section presents a description of our second hybrid model – HQLSTM, which
is a hybrid analog of the classical LSTM model hochreiter1997long , with which
predictions will be compared in the following sections. LSTM architectures
have garnered significant attention in the realm of time series forecasting,
including in predicting PV power gao2023short ; sharadga2020time ; tovar2020pv
.

HQLSTM models have proven themselves well for solving problems from various
fields. Examples of successful use of this model are the tasks of natural
language processing wang2022adverse , the detection of software
vulnerabilities akter2022automated , and predicting solar radiation
yu2023prediction .

In this proposed model we added a quantum layer to each of the LSTM gates
chen2020quantum . Let’s take a closer look at our implementation, depicted in
Fig 3(b). The input to the model is:

  1. 1.

The current step information, represented by a green circle, $x(t)$. This is a
tensor of size $5$, reflecting the five features for an hour, which include
meteorological data and the PV power itself.

  2. 2.

The information from the previous step, denoted by a purple circle, $h(t-1)$.
It consists of a tensor of size $h_{dim}$. For the initial step, this is
simply a zero vector.

These inputs are processed through classical fully-connected layers to yield
vectors with a uniform dimension of $4n_{q}$. These vectors are then
concatenated through a bitwise addition operation.

Subsequently, this concatenated vector is partitioned into four distinct
groups, for the four gates of the LSTM cell. As illustrated in Fig. 3(b-c),
each group is directed to the input of its corresponding quantum layer,
symbolized by the QDI square.

The outputs from QDI layers are transformed via classical fully-connected
layers to standardize their dimensions to $h_{dim}$. Following this,
activation functions together with appropriate for each of the $4$ gates
transformations, similar to the classical LSTM, are applied to the outputs
originating from the quantum layers. This processing culminates in the
derivation of the new cell state $C(t)$ and the hidden state $h(t)$ vectors.

The process operates in a cyclical manner. For each iteration, the vector from
the current time step, $x(t)$, and the hidden vector from the previous step,
$h(t-1)$, serve as inputs to the HQLSTM. This iterative process is executed as
many times as the input width; in our case input width equals $24$.
Subsequently, all the hidden vectors are concatenated to produce a single
composite vector. This vector is then processed through a fully-connected
layer consisting of a single neuron, which outputs a value that predicts the
PV power.

In our first proposed architecture, the HQNN, the quantum layer functioned as
a vanilla layer, where variational layers were sequentially placed after the
encoding layer. In contrast, in the HQLSTM approach we used a QDI layer pharma
as depicted in Fig. 3(b). Here, variational layers are positioned multiple
times before the encoding layer (green rectangular) and additionally (purple
rectangular), within each encoding layer (blue rectangular) to increase the
layer’s expressivity.

###  II.4 HQSeq2Seq

Here we present a hybrid version of the Sequence-to-Sequence (Seq2Seq) model,
first introduced in sutskever2014sequence . Seq2seq models are widely used in
natural language processing tasks kalchbrenner2013recurrent , where the length
of input and output sequence is not pre-determined and can be variable. We can
also apply the principle of Seq2Seq models to the power prediction task
mu2023improved . That means we can feed the neural network with time series
with arbitrary length and prompt it to give us the forecast for any hours
ahead. In this problem setting, the longer the input time series is, the
better the model prediction is. The same applies to the required output
length: the shorter it is, the easier it is for the model to generate the
forecast.

The seq2Seq model is a type of encoder-decoder model. The encoder is given the
entire input sequence, which it uses to generate a context vector. This vector
is used as an input hidden state for the decoder, so it literally provides it
with “context”, according to which the decoder will generate the forecast.
Thereby, the hidden dimensions of the encoder and the decoder must match. The
decoder creates the output sequence step by step. It starts with only one
entry: the one which is the last known. Based on this entry and the context
vector, the decoder generates the second entry and appends it to the existing
one. Now, the obtained two-entry sequence is once again fed into the decoder
to generate the third entry. Then, the cycle repeats until the length of the
generated sequence matches the length requested by the user.

We create and compare two models with Seq2Seq architecture: the classical
Seq2Seq and the hybrid model called HQSeq2Seq. Both of these models have
identical LSTMs acting as encoders and decoders. In the classical model, the
decoder’s hidden output vector is mapped to the “Power” value with a single
linear layer, while in HQSeq2Seq it is processed by a QDI layer pharma .

In the QDI layer, instead of attempting to use a qubit for each feature asel1
, we employed the data re-uploading technique schuld2021effect ;
data_reuploading . Specifically, we work with $4$ qubits and structure them
into a lattice of depth $4$ (depicted as a blue big rectangular in Fig. 3(d)).
Each of our $16$ input features leading to the quantum layer is intricately
encoded within this lattice. The first four features are mapped onto the
initial depth, followed by the subsequent features in blocks of four. Encoding
these classical features into the quantum domain, we adopt the “angle
embedding” using $R_{z}$ gate. This operation effectively translates the input
vector into a quantum state that symbolizes the preceding classical layer’s
data.

Entangling variational layers, signified by purple squares, are interposed
between every encoding layer, ensuring optimal Fourier accessibility. Each
variational layer has two components: rotations governed by trainable
parameters, and sequential CNOT gates. The rotations are implemented by
quantum gates that metamorphose the encoded input in line with the variational
parameters, while the CNOT operations handle the entanglement of the qubits,
facilitating quantum superposition.

Each lattice depth, represented by each blue square, encompasses a variational
layer (purple square). Moreover, prior to all encoding layers, we introduce a
variational layer (designated by a green square) for enhanced model
representation. Consequently, the total weight count in the quantum segment of
our network is $20$. In the measurement phase, except for the first qubit, all
qubits execute a CNOT operation targeting the first qubit, ensuring the
$Y$-measurement spans all qubits. Therefore, the quantum layer’s output serves
as the power value prediction for a specific hour.

The input size of the encoder and decoder can differ, which is a substantial
benefit. For instance, we can use all of the $5$ features to create a context
vector, but request to generate the forecast for only $1$ feature. Exploiting
this advantage, we will feed the Seq2Seq model with a window of all known
features and demand the forecast only for the “Power” one.

For simplicity’s sake, we will train both models with fixed input and output
length of $96$ hours and then try to vary the length in the testing stage.

Figure 4: Results of training and testing of HQNN & MLP and HQLSTM & LSTM
models. (a-b) Training and testing history of models. Solid and dotted lines
represent training and testing of models respectively. Filled space shows
standard deviation of models that were averaged over different testing
subsets. Histogram for reduced training dataset. (c) Histogram shows mean MAE
and RMSE metrics for our models on testing subset, averaged over different
testing subsets. (f-g) Train and test learning curves for the Seq2Seq and
HQSeq2Seq models. (h) Example of the classical and hybrid Seq2Seq models
inference on the testing data. The models get $124$ hours of data as an input
(before the dashed line) and gives forecast of “Power” into $137$ hours ahead
(after the dashed line). The solid black line represents the ground truth
value of the “Power” feature. MSE/MAE errors in this particular example are
$0.0052/0.0349$ for Seq2Seq and $0.0040/0.0292$ for HQSeq2Seq.

###  II.5 Training and results

In the study, six distinct models were employed for PV power prediction based
on weather features: HQNN, MLP, HQLSTM, LSTM, HQSeq2Seq, and Seq2Seq. To train
the models, the mean square error (MSE) was chosen as the loss function:

$\text{MSE}=\frac{1}{N}\sum\limits_{n=1}^{N}\left(x_{n}-y_{n}\right)^{2},$

where $N$ is number of predictions, $x$ denotes the predicted PV power, and
$y$ represents the actual PV power value.

To test the models, in addition to the MSE loss metric, we also used the mean
absolute error (MAE), root mean squared error (RMSE) and variance account
factor (VAF).

Here $\vec{x}=\left(x_{1},x_{2},\dots,x_{N}\right)$ and
$\vec{y}=\left(y_{1},y_{2},\dots,y_{N}\right)$ represent vectors of predicted
and target PV power values respectively, where $N$ is the number of predicted
values.

All the machine learning simulations for this study were conducted on CPUs, on
the QMware cloud qmw_qmw_2022 ; kordzanganeh2023benchmarking  device. The
classical part of our modeling was structured using the PyTorch library
PyTorch , while the quantum part was implemented using the PennyLane
framework. Notably, PennyLane provides an assortment of qubit devices. For our
requirements, we selected the lightning.qubit device which is a custom backend
for simulating quantum state-vector evolution. To compute the gradients of the
loss function relative to each parameter: for the classical components of our
hybrid models, we employed the widely-recognized backpropagation algorithm
backprop ; and for the quantum part, we used the adjoint method as highlighted
in Refs. Luo2020yaojlextensible ; efficient_calculation .

Hyperparameters | Range | Best value
---|---|---
HQNN
| number of neurons in
---
the second layer
$8-128$ | 17
number of qubits | $2-10$ | 8
number of variational layers | $1-10$ | 7
embedding |  $R_{x}$, $R_{y}$, $R_{z}$ | $R_{x}$
measurement |  $X$, $Y$, $Z$ | $Z$
variational part | basic, strongly | basic
initial learning rate | $1-1000\times 10^{-4}$ | $3\times 10^{-2}$
MLP
| number of neurons in
---
the first layer
$8-128$ | 32
| number of neurons in
---
the second layer
$8-128$ | 3
| number of neurons in
---
the third layer
$8-128$ | 3
initial learning rate | $1-1000\times 10^{-4}$ | $1\times 10^{-2}$
HQLSTM
| number of neurons in
---
the second layer
$3-25$ | 20
dropout range | $0-99$ | 0.239
number of qubits | $2-5$ | 4
number of variational layers | $1-5$ | 1
number of quantum layers | $1-5$ | 3
initial learning rate | $1-100\times 10^{-3}$ | $0.52\times 10^{-2}$
LSTM
| number of neurons in
---
the second layer
$3-25$ | 21
dropout range | $0-99$ | 0.158
initial learning rate | $1-100\times 10^{-3}$ | $0.5\times 10^{-2}$
Table 1: The table shows which hyperparameters are being optimized, the limits of change, and the best values found during hyperparameter optimization. Model | HQNN | MLP | HQLSTM | LSTM
---|---|---|---|---
train, MAE | $0.0458$ | $0.0651$ | $0.0382$ | $0.0454$
test, MAE | $0.0659$ | $0.0887$ | $0.0343$ | $0.0570$
train, MSE | $0.0059$ | $0.0144$ | $0.0056$ | $0.0054$
test, MSE | $0.0120$ | $0.0204$ | $0.0058$ | $0.0096$
test, RMSE | $0.1097$ | $0.1428$ | $0.0743$ | $0.0937$
Table 2: Summary of the results for the proposed models. In a direct
comparison, the HQNN outperforms the MLP in both training and testing losses,
as evidenced across three critical metrics: RMSE, MAE, and MSE. Notably, the
HQNN’s MSE is $41\%$ less than MLP’s one, which means that HQNN has more
reliable forecasts, despite the fact that it has $1.8$ times fewer parameters.
At the same time, the quality of HQLSTM’s on the test dataset is $40\%$ better
compared to LTSM’s one on MAE and MSE metrics, and is $21\%$ better on RMSE
metric, although the first one has less than half weights. In a broader
comparison encompassing all four models, HQLSTM emerges as the most precise
model on all metrics, namely on $52\%$ more precise than HQNN having two times
fewer trainable parameters.

####  II.5.1 HQNN & MLP and HQLSTM & LSTM

Both the HQNN and its classical analog, the MLP, were trained for $20$ epochs.
In contrast, the HQLSTM and its classical counterpart, the LSTM, were trained
for over $50$ epochs. The Adam optimiser kingma2017adam  from the PyTorch
framework was used to update the parameters of the models in order to minimize
their loss functions. The comprehensive training process, accompanied by the
results, is delineated in Fig. 4 and in Table 2.

In this study, we employed cross-validation as a fundamental technique to
assess the performance of our models across distinct testing subsets. The
application of cross-validation is pivotal to safeguard against potential data
leakage from the training dataset into the testing dataset. To achieve this, a
rigorous approach was adopted wherein a 24-hour time window, from each side of
the subsets, was systematically excluded from the dataset.

Furthermore, we opted to partition the dataset into training and testing sets
in a $4:1$ ratio. This strategy was implemented to promote a comprehensive
evaluation of our models, as we carried out model training and assessment on
five distinct data splits. Subsequently, a meticulous averaging process was
employed to consolidate the results obtained from these splits, and the mean
values thus derived served as the primary metrics for inter-model comparisons.

The utilization of cross-validation techniques in our methodology
significantly bolsters the robustness and reliability of our results, as they
diminish the reliance on specific train-test partitioning, thereby enhancing
the credibility of our findings.

In a head-to-head comparison between HQNN and MLP, the former exhibits
superior performance regarding training and testing losses across three key
metrics: RMSE, MAE, and MSE. Specifically, HQNN surpasses MLP’s power
prediction accuracy by $41\%$ estimated by MSE loss function and by $26\%$ by
MAE loss, all the while boasting $1.8$ times fewer parameters ($2266\ \&\
3987$).

On juxtaposing HQLSTM with LSTM, the former outperforms in training and
testing loss across all three aforementioned metrics. Remarkably, HQLSTM’s has
better predictive ability (on $40\%$ better) than LSTM assessed by the MSE and
MAE metrics, and it achieves this with less than half the number of parameters
($1109\ \&\ 2857$). Moreover, HQLSTMs are more resistant to overfitting, while
classical LSTM suffers from it.

In a broader comparison encompassing all four models, HQLSTM emerges as the
most precise model on all metrics, namely on $52\%$ more precise than HQNN,
having two times fewer trainable parameters.

It is worth noting that we performed hyperparameter optimization technique
using the Optuna optimizer optuna_link . The set of optimized parameters,
limits of their variation, and best sets of hyperparameters for all our $4$
models are presented in Table 1.

Moreover, we scanned external articles that refer to this dataset and found
only one article that solves $1$ hour ahead PV power prediction using neural
networks. In comparison, our HQNN model is better than the model from the
external article according to VAF metrics by $40\%$ ($91\,\&\,65$).

Further, to confirm that hybrid models train better including on a smaller
dataset, an additional experiment was conducted, wherein the volume of
training data was intentionally reduced. The results are shown in Figure 4.
The hybrid models performed better with less data, having fewer losses and
better prediction capabilities than classical models.

####  II.5.2 HQSeq2Seq & Seq2Seq

After preprocessing the dataset, which is described in the Section II.1, it
spanned $12775$ hours from 3/5/12 4:55 AM to 8/19/13 10:00 AM for training,
and $3194$ hours from 8/19/13 11:00 PM to 12/30/13 00:00 AM for testing.

Although models are capable of being trained on sequences of arbitrary
lengths, we chose to use sequences of fixed $96$ hours for simplicity. In this
case, the encoder gets $96$ hours of all available features, while the decoder
is asked to extrapolate only the “Power” feature of the data $96$ hours ahead.
Training for $15$ epochs with the Adam optimizer (learning rate $0.001$)
proved to be enough for the models to converge (Fig. 4 (f-g)). As an example
of inference, we pass the time series of the length different from $96$ into
both models and prompt them to give us a forecast for $137$ hours ahead (Fig.
4 (h)). We can conclude that both models transition from the fixed sequence
length to an arbitrary one quite well. It may even be possible to improve
these results by introducing variable-length sequences into the training
stage. We also measured the dependency of test loss on the size of the
training dataset for Seq2Seq and HQSeq2Seq, shown in Fig. 4 (e). As one can
see, the test RMSE loss of the hybrid model is less for any size of training
data, which proves that the hybrid model shows an advantage over the classical
model, including on a trimmed dataset.

##  III Discussion

In this work we introduced three hybrid quantum approaches to time series
prediction task. The first two models allow one to predict the power of solar
panels for $1$ hour ahead, using weather features for the previous $24$ hours.
The third model allows to predict a longer-term user-defined forecast,
showcasing the versatility of our models for various planning tasks.

The first approach is the HQNN, a combination of classical fully connected
layers and quantum layer, which is a VVRQ, analog of classical fully-connected
layer. We compared this hybrid model with its classical counterpart, MLP, and
demonstrated that, even though HQNN has $1.8$ times less variational
parameters, it has $41\%$ better predictive ability, estimated by MSE error.

The second approach is HQLSTM, a hybrid quantum analogue of classical LSTM.
Here, QDI layer is inserted into each gate of the LSTM cell. This approach
provides a $40\%$ improvement in prediction using the MAE and MSE metrics
comraped to its classical counterpart. Our proposed architecture is a unique
combination of classical and quantum layers, which we believe to be a
breakthrough in solving time series prediction tasks. Comparing HQNN and
HQLSTM models, the second one was better on $52\%$ than HQNN having two times
fewer weights.

The third approach is hybrid Seq2Seq model, a classical Seq2Seq model,
consisting of $2$ LSTMs with quantum layer at the end. This approach allows
one to predict the PV power not only for an hour ahead, but for any number of
hours ahead, without knowing the weather features in advance. The addition of
the proposed QDI layer improves the accuracy of the predictions, reducing the
MAE error by $16\%$ compared to a purely classical Seq2Seq model. Our proposed
architecture is a unique combination of classical and quantum layers, which we
believe to be a breakthrough in solving time series prediction tasks.

Also, for all models, we conducted an additional experiment in which our
models were trained on a reduced dataset, and confirmed that hybrid models
have better learning capabilities and have less loss trained on any amount of
dataset compared to their classical counterparts. This confirmation can serve
as an excellent motivation to use hybrid networks for applications where data
collection is a complex task.

It is worth noting that all parameters are trainable in our layers; the
architecture and hyperparameters were selected by the Optuna optimizer.
Moreover, we compare our models to a paper that solved the same problem using
the same dataset, and demonstrate that our best HQLSTM is $40\%$ more accurate
in predicting power using VAF metric.

To fully unlock the potential of HQNNs in time series prediction problem,
further research and testing of models on other datasets is necessary. Also,
the development of more efficient optimization VQC training and implementation
methods, larger-scale quantum hardware could lead to even more significant
performance improvements.

Furthermore, while this work was done on a public dataset with an emphasis on
hybrid quantum models for better forecasting performance, the quality and
source of data plays a crucial role in overall effectiveness in the real
world, especially considering weather data. Accurate weather forecast is a
crucial input into any high performing and useful PV prediction given its
dynamism and influence on PV output. An interesting area of research is cloud
prediction using satellite and weather data for geo-locations, directly
impacting solar irradiance and therefore PV output. The added complexity could
enhance the need for hybrid quantum models to increase computational
efficiency and higher quality forecasts.

To summarize, our developments provide three hybrid quantum approaches for
time series problems that demonstrate the possibility of combining classical
and quantum methods. Our proposed models show improved performance compared to
classical models with similar architecture when using fewer variation
parameters. We believe that these results pave the way for further research in
developing hybrid models that leverage the strengths of both classical and
quantum computing.

## References

  * (1) M. Sabri and M. El Hassouni.  A novel deep learning approach for short term photovoltaic power forecasting based on gru-cnn model.  E3S Web of Conferences, 336:00064, 2022.
  * (2) B. Li and J. Zhang.  A review on the integration of probabilistic solar forecasting in power systems.  Solar Energy, 210:68–86, 2020.
  * (3) O. Ellabban, H. Abu-Rub, and F. Blaabjerg.  Renewable energy resources: Current status, future prospects and their enabling technology.  Renewable and Sustainable Energy Reviews, 39:748–764, 2014.
  * (4) Anthony R. Florita a Carlo Brancucci Martinez-Anido a, Benjamin Botor a.  The value of day-ahead solar power forecasting improvement.  Science, 129:192–203, 2016.
  * (5) D. Lee and K. Kim.  PV power prediction in a peak zone using recurrent neural networks in the absence of future meteorological information.  Renewable Energy, 173:1098–1110, 2021.
  * (6) S.A. Bozorgavari, J. Aghaei, S. Pirouzi, V. Vahidinasab, H. Farahmand, and M. Korpås.  Two-stage hybrid stochastic/robust optimal coordination of distributed battery storage planning and flexible energy management in smart distribution network.  Journal of Energy Storage, 26:100970, 2019.
  * (7) M. Diagne, M. David, P. Lauret, J. Boland, and N. Schmutz.  Review of solar irradiance forecasting methods and a proposition for small-scale insular grids.  Renewable and Sustainable Energy Reviews, 27:65–76, 2013.
  * (8) A. Ahmed and M. Khalid.  A review on the selected applications of forecasting models in renewable power systems.  Renewable and Sustainable Energy Reviews, 100:9–21, 2019.
  * (9) J. López Gómez, A. Ogando Martínez, F. Troncoso Pastoriza, L. Febrero Garrido, E. Granada Álvarez, and J. A. Orosa García.  Photovoltaic Power Prediction Using Artificial Neural Networks and Numerical Weather Data.  Sustainability, 12:10295, 2020.
  * (10) M. R. Kaloop, A. Bardhan, N. Kardani, P. Samui, J. W. Hu, and A. Ramzy.  Novel application of adaptive swarm intelligence techniques coupled with adaptive network-based fuzzy inference system in predicting photovoltaic power.  Renewable and Sustainable Energy Reviews, 148:111315, 2021.
  * (11) R. S. Kulkarni, D. B. Talange, and N. V. Mate.  Output Estimation of Solar Photovoltaic (PV) System.  In 2018 International Symposium on Advanced Electrical and Communication Technologies (ISAECT), pages 1–6. IEEE, 2018.
  * (12) J. Shi, W.J. Lee, Y. Liu, Y. Yang, and P. Wang.  Forecasting Power Output of Photovoltaic Systems Based on Weather Classification and Support Vector Machines.  IEEE Transactions on Industry Applications, 48(3):1064–1069, 2012\.
  * (13) O. P. Rocha, A.M. da Silva1, and A. Á. B. Santos.  Computational Model for Photovoltaic Solar Energy Forecasting Based on the K-Nearest Neighbor Method.  Journal of Bioengineering, Technologies and Healt, 5(3):168–172, 2022.
  * (14) R. H. Inman, H. T. C. Pedro, and C. F. M. Coimbra.  Solar forecasting methods for renewable energy integration.  Progress in Energy and Combustion Science, 39(6):535–576, 2013\.
  * (15) R. Blaga, A. Sabadus, N. Stefu, C. Dughir, M. Paulescu, and V. Badescu.  A current perspective on the accuracy of incoming solar energy forecasting.  Progress in energy and combustion science, 70:119–144, 2019.
  * (16) M. Schuld, I. Sinayskiy, and F. Petruccione.  An introduction to quantum machine learning.  Contemporary Physics, 56(2):172–185, 2015.
  * (17) J. Biamonte, P. Wittek, N. Pancotti, P. Rebentrost, N. Wiebe, and S. Lloyd.  Quantum machine learning.  Nature, 549(7671):195–202, 2017.
  * (18) V. Dunjko and H. J Briegel.  Machine learning & artificial intelligence in the quantum domain: a review of recent progress.  Reports on Progress in Physics, 81(7):074001, 2018.
  * (19) Alexey Melnikov, Mohammad Kordzanganeh, Alexander Alodjants, and Ray-Kuang Lee.  Quantum machine learning: from physics to software engineering.  Advances in Physics: X, 8(1):2165452, 2023.
  * (20) D. Emmanoulopoulos and S. Dimoska.  Quantum Machine Learning in Finance: Time Series Forecasting.  arXiv preprint arXiv:2202.00599, 2019.
  * (21) M. M. Sushmit and I. M. Mahbubul.  Forecasting solar irradiance with hybrid classical–quantum models: A comprehensive evaluation of deep learning and quantum-enhanced techniques.  Energy Conversion and Management, 294:117555, 2023.
  * (22) M. Schuld, M. Fingerhuth, and F. Petruccione.  Implementing a distance-based classifier with a quantum interference circuit.  Europhysics Letters, 119(6):60002, 2018.
  * (23) S. Lloyd, M. Mohseni, and P. Rebentrost.  Quantum algorithms for supervised and unsupervised machine learning.  arXiv preprint arXiv:1307.0411, 2013.
  * (24) P.W. Shor.  Algorithms for quantum computation: discrete logarithms and factoring.  Proceedings 35th Annual Symposium on Foundations of Computer Science, pages 124–134, 1994.
  * (25) S. Lloyd.  Universal quantum simulators.  Science, 273:1073–1078, 1996.
  * (26) M. Schuld, A. Bocharov, K. M Svore, and N. Wiebe.  Circuit-centric quantum classifiers.  Physical Review A, 101(3):032308, 2020.
  * (27) S. Aaronson and L. Chen.  Complexity-Theoretic Foundations of Quantum Supremacy Experiments.  arXiv preprint arXiv:1612.05903, 2016.
  * (28) P. Jain and S. Ganguly.  Hybrid quantum generative adversarial networks for molecular simulation and drug discovery.  arXiv preprint arXiv:2212.07826, 2022.
  * (29) Asel Sagingalieva, Mohammad Kordzanganeh, Nurbolat Kenbayev, Daria Kosichkina, Tatiana Tomashuk, and Alexey Melnikov.  Hybrid quantum neural network for drug response prediction.  Cancers, 15(10):2705, 2023.
  * (30) Luca Lusnig, Asel Sagingalieva, Mikhail Surmach, Tatjana Protasevich, Ovidiu Michiu, Joseph McLoughlin, Christopher Mansell, Graziano de’ Petris, Deborah Bonazza, Fabrizio Zanconati, Alexey Melnikov, and Fabio Cavalli.  Hybrid quantum image classification and federated learning for hepatic steatosis diagnosis.  arXiv preprint arXiv:2311.02402, 2023.
  * (31) A. Kurkin, J. Hegemann, M. Kordzanganeh, and A. Melnikov.  Forecasting the steam mass flow in a powerplant using the parallel hybrid network.  arXiv preprint arXiv:2307.09483, 2023.
  * (32) Serge Rainjonneau, Igor Tokarev, Sergei Iudin, Saaketh Rayaprolu, Karan Pinto, Daria Lemtiuzhnikova, Miras Koblan, Egor Barashov, Mo Kordzanganeh, Markus Pflitsch, and Alexey Melnikov.  Quantum algorithms applied to satellite mission planning for Earth observation.  IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 16:7062–7075, 2023.
  * (33) Nathan Haboury, Mo Kordzanganeh, Sebastian Schmitt, Ayush Joshi, Igor Tokarev, Lukas Abdallah, Andrii Kurkin, Basil Kyriacou, and Alexey Melnikov.  A supervised hybrid quantum machine learning solution to the emergency escape routing problem.  arXiv preprint arXiv:2307.15682, 2023.
  * (34) Asel Sagingalieva, Andrii Kurkin, Artem Melnikov, Daniil Kuhmistrov, et al.  Hybrid quantum ResNet for car classification and its hyperparameter optimization.  Quantum Machine Intelligence, 5(2):38, 2023.
  * (35) Mo Kordzanganeh, Pavel Sekatski, Leonid Fedichkin, and Alexey Melnikov.  An exponentially-growing family of universal quantum circuits.  Machine Learning: Science and Technology, 4(3):035036, 2023.
  * (36) Arsenii Senokosov, Alexander Sedykh, Asel Sagingalieva, and Alexey Melnikov.  Quantum machine learning for image classification.  arXiv preprint arXiv:2304.09224, 2023.
  * (37) Michael Perelshtein, Asel Sagingalieva, Karan Pinto, Vishal Shete, Alexey Pakhomchik, Artem Melnikov, et al.  Practical application-specific advantage through hybrid quantum computing.  arXiv preprint arXiv:2205.04858, 2022.
  * (38) M. Schuld, R. Sweke, and J. J. Meyer.  Effect of data encoding on the expressive power of variational quantum-machine-learning models.  Physical Review A, 103(3):032430, 2021.
  * (39) Mo Kordzanganeh, Daria Kosichkina, and Alexey Melnikov.  Parallel hybrid networks: an interplay between quantum and classical neural networks.  Intelligent Computing, 2:0028, 2023.
  * (40) Alexandr Sedykh, Maninadh Podapaka, Asel Sagingalieva, Nikita Smertyak, Karan Pinto, Markus Pflitsch, and Alexey Melnikov.  Quantum physics-informed neural networks for simulating computational fluid dynamics in complex shapes.  arXiv preprint arXiv:2304.11247, 2023.
  * (41) M Malvoni, M.G. De Giorgi, and P.M. Congedo.  Data on photovoltaic power forecasting models for Mediterranean climate.  Data in Brief, 7:1639–1642, 2016.
  * (42) J. Aghaei, S.A. Bozorgavari, S. Pirouzi, H. Farahmand, and M. Korpås.  Flexibility planning of distributed battery energy storage systems in smart distribution networks.  Iranian Journal of Science and Technology, Transactions of Electrical Engineering, 44(3):1105–1121, 2019.
  * (43) Rohan S. Kulkarni, Rajani B. Shinde, and Dhananjay B. Talange.  Performance Evaluation and a New Thermal Model for a Photovoltaic-Thermal Water Collector System .  In 2018 International Conference on Smart Grid and Clean Energy Technologies (ICSGCE), pages 106–111. IEEE, 2018.
  * (44) M. Malvoni, M. G. De Giorgi, and P. M. Congedo.  Study of degradation of a grid connected photovoltaic system.  Energy Procedia, 126:644–650, 2017.
  * (45) F. Bloch.  Nuclear Induction.  Physical Review, 70(7-8):460–474, 1946.
  * (46) Adriano Barenco, Charles H. Bennett, Richard Cleve, David P. DiVincenzo, Norman Margolus, Peter Shor, Tycho Sleator, John A. Smolin, and Harald Weinfurter.  Elementary gates for quantum computation.  Physical Review A, 52(5):3457–3467, 1995.
  * (47) S. Hochreiter and J. Schmidhuber.  Long short-term memory.  Neural computation, 9(8):1735–1780, 1997.
  * (48) H. Gao, S. Qiu, J. Fang, N. Ma, J. Wang, K. Cheng, et al.  Short-Term Prediction of PV Power Based on Combined Modal Decomposition and NARX-LSTM-LightGBM.  Sustainability, 15(10):8266, 2023.
  * (49) H. Sharadga, S. Hajimirza, and R. S. Balog.  Time series forecasting of solar power generation for large-scale photovoltaic plants.  Renewable Energy, 150:797–807, 2020.
  * (50) M. Tovar, M. Robles, and F. Rashid.  PV Power Prediction, Using CNN-LSTM Hybrid Neural Network Model. Case of Study: Temixco-Morelos, México.  Energies, 13(24):6512, 2020.
  * (51) X. Wang, X. Wang, and S. Zhang.  Adverse Drug Reaction Detection from Social Media Based on Quantum Bi-LSTM with Attention.  IEEE Access, 11:16194–16202, 2022.
  * (52) M. S. Akter, H. Shahriar, and Z. A. Bhuiya.  Automated vulnerability detection in source code using quantum natural language processing.  In Inernational Conference on Ubiquitous Security, pages 83–102. Springer, 2022.
  * (53) Y. Yu, G. Hu, C. Liu, J. Xiong, and Z. Wu.  Prediction of solar irradiance one hour ahead based on quantum long short-term memory network.  IEEE Transactions on Quantum Engineering, 2023.
  * (54) S. Y.C. Chen, S. Yoo, and Y. L. L. Fang.  Quantum long short-term memory.  arXiv preprint arXiv:2009.01783, 2020.
  * (55) I. Sutskever, O. Vinyals, and Q. V Le.  Sequence to sequence learning with neural networks.  Advances in neural information processing systems, 27, 2014.
  * (56) N. Kalchbrenner and P. Blunsom.  Recurrent continuous translation models.  In Proceedings of the 2013 conference on empirical methods in natural language processing, pages 1700–1709, 2013.
  * (57) Y. Mu, M. Wang, X. Zheng, and H. Gao.  An improved LSTM-Seq2Seq-based forecasting method for electricity load.  Frontiers in Energy Research, 10:1093667, 2023.
  * (58) Adriá n Pérez-Salinas, Alba Cervera-Lierta, Elies Gil-Fuster, and José I. Latorre.  Data re-uploading for a universal quantum classifier.  Quantum, 4:226, 2020.
  * (59) QMware.  QMware — The first global quantum cloud.  https://qm-ware.com/, 2022.
  * (60) Mohammad Kordzanganeh, Markus Buchberger, Basil Kyriacou, Maxim Povolotskii, Wilhelm Fischer, Andrii Kurkin, et al.  Benchmarking simulated and physical quantum processing units using quantum and hybrid algorithms.  Advanced Quantum Technologies, 6(8):2300043, 2023.
  * (61) PyTorch — Machine Learning framework.  https://pytorch.org/, 2022.
  * (62) D. E. Rumelhart, E. Hinton, G, and R. J. Williams.  Learning representations by back-propagating errors.  Nature, 323(6088):533–536, 10 1986.
  * (63) X. Z. Luo, J. G. Liu, P. Zhang, and L. Wang.  Yao.jl: Extensible, Efficient Framework for Quantum Algorithm Design.  Quantum, 4:341, 2020.
  * (64) T. Jones and J. Gacon.  Efficient calculation of gradients in classical simulations of variational quantum algorithms.  arXiv preprint arXiv:2009.02823, 2020.
  * (65) Diederik P. Kingma and Jimmy Ba.  Adam: A method for stochastic optimization.  arXiv preprint arXiv:1412.6980, 2017.
  * (66) Optuna — A hyperparameter optimization framework.  https://optuna.org/, 2023.


# Refined Ehrhart series and bigraded rings

Praise Adeyemo  Faculty of Science, University of Ibadan, Ibadan, Oyo State,
Nigeria<EMAIL_ADDRESS>and  Balázs Szendrői  Mathematical
Institute, University of Oxford, Oxford OX2 6GG, United Kingdom
<EMAIL_ADDRESS>
###### Abstract.

We study a natural set of refinements of the Ehrhart series of a closed
polytope, first considered by Chapoton. We compute the refined series in full
generality for a simplex of dimension $d$, a cross-polytope of dimension $d$,
respectively a hypercube of dimension $d\leq 3$, using commutative algebra. We
deduce summation formulae for products of $q$-integers with different
arguments, generalizing a classical identity due to MacMahon and Carlitz. We
also present a characterisation of a certain refined Eulerian polynomial in
algebraic terms.

## Introduction

Given a $d$-dimensional lattice polytope $\Sigma$, the number of lattice
points in its $r$-fold dilation is given by a polynomial in $r$, its Ehrhart
polynomial. The generating function of the number of lattice points as $r$
varies gives rise to a related notion, the Ehrhart series of $\Sigma$. On the
other hand, the polytope $\Sigma$ can also be used to define a projective
toric variety $X_{\Sigma}$ which, under a further assumption, comes with a
natural embedding into projective space. The Ehrhart polynomial and Ehrhart
series of $\Sigma$ get identified with the Hilbert polynomial and Hilbert
series of this embedding.

Recently Chapoton [4] generalised Ehrhart theory by counting lattice points in
slices defined by hyperplanes orthogonal to a fixed direction, giving a one-
parameter refinement of the Ehrhart series of $\Sigma$. In this short paper,
we observe that Chapoton’s refinement corresponds on the geometric side to a
natural bigrading on the projective coordinate ring of $X_{\Sigma}$. We use
this observation to calculate the refined Ehrhart series, for all possible
choices of slicing direction, for the standard simplex and cross-polytope of
arbitrary dimension $d$, as well as for a hypercube of dimension $d\leq 3$.
Our main results are Propositions 2.1, 2.2, 2.3 and 2.5. We explain why a
general formula for the refined Ehrhart series of a hypercube of dimension
$d>3$ is not really of interest. We also use these ideas to formulate in
Proposition 2.7 an interpretation of a refined Eulerian polynomial, due to
MacMahon and Carlitz, in terms of commutative algebra.

The projective coordinate ring of $X_{\Sigma}$ is in fact naturally graded by
$\mathbb{Z}\oplus\mathbb{Z}^{d}$, compare e.g. [7, Sect.13.3]. It is
straightforward to write out generalizations of the formulae of this note in
the multi-graded case. We chose to focus on the bigraded case, as that
corresponds to Chapoton’s refinement, and also that is the case that
(sometimes) gives simplified expressions, such as the one corresponding to the
refined Eulerian polynomial.

### Notation

We work in the fixed lattice $\mathbb{Z}^{d}$, with standard basis vectors
${\bf e}_{1},\ldots,{\bf e}_{d}$ orthonormal with respect to the inner product
$({\bf v},{\bf w})\mapsto{\bf v}\cdot{\bf w}$. A $d$-dimensional lattice
polytope $\Sigma\subset\mathbb{R}^{d}$ is a compact convex polytope with
vertices contained in our lattice $\mathbb{Z}^{d}\subset\mathbb{R}^{d}$.

### Acknowledgements

We would like to thank Miles Reid, John Shareshian, Anna Seigal, Bernd
Sturmfels and especially Mateusz Michalek for useful correspondence and
comments. The ideas in this paper were formulated during a visit of P. A. to
the University of Oxford, supported by EPSRC GCRF grant EP/T001968/1, part of
the Abram Gannibal Project.

##  1\. Standard definitions and examples

###  1.1. Polytopes and Ehrhart theory

Let $\Sigma\subset\mathbb{R}^{d}$ be a $d$-dimensional lattice polytope. For a
non-negative integer $r$, let $r\Sigma\subset\mathbb{R}^{d}$ denote the
$r$-fold dilation of $\Sigma$. Then it is known that there exists a
polynomial, the Ehrhart polynomial $L_{\Sigma}\in\mathbb{Q}[x]$, with the
property that for all non-negative integers $r$, we have

$\\#\\{r\Sigma\cap\mathbb{Z}^{d}\\}=L_{\Sigma}(r).$

We can also consider the Ehrhart series

$E_{\Sigma}(t)=\sum_{r=0}^{\infty}L_{\Sigma}(r)t^{r}.$

Then it is known that

$E_{\Sigma}(t)=\frac{H_{\Sigma}(t)}{(1-t)^{d+1}},$

where $H_{\Sigma}(t)$ is a polynomial of degree $d$ with non-negative integer
coefficients.

######  Example 1.1.

Let $\Sigma=\Delta_{d}\subset\mathbb{R}^{d}$ be the unit $d$-simplex spanned
by the standard basis vectors of $\mathbb{R}^{d}$. The Ehrhart polynomial is

$L_{\Delta_{d}}(r)={r+d\choose d},$

whereas the Ehrhart series is

$E_{\Delta_{d}}(t)=\frac{1}{(1-t)^{d+1}}.$

######  Example 1.2.

Let $\Sigma=\Xi_{d}$ be the $d$-dimensional cross-polytope, the convex hull of
all unit coordinate vectors and their opposites. Then (see for example [2,
Thm.2.7]), the Ehrhart polynomial of $\Xi_{d}$ is

$L_{\Xi_{d}}(r)=\sum_{k=0}^{d}{d\choose k}{r-k+d\choose d}.$

The corresponding Ehrhart series is

(1.1) $E_{\Xi_{d}}(r)=\frac{(1+t)^{d}}{(1-t)^{d+1}}.$

######  Example 1.3.

Let $\Sigma=C_{d}$ be the unit hypercube of dimension $d$. The Ehrhart
polynomial is clearly

$L_{C_{d}}(r)=(1+r)^{d}.$

Thus the Ehrhart series is

(1.2)
$E_{C_{d}}(t)=\sum_{r=0}^{\infty}(1+r)^{d}t^{r}=\frac{A_{d}(t)}{(1-t)^{d+1}},$

where $A_{d}(t)$ is the $d$-th Eulerian polynomial, introduced by Euler to sum
this very series. The first few examples are

$A_{2}(t)=1+t,\ A_{3}(t)=1+4t+t^{2},\ A_{4}(t)=1+11t+11t^{2}+t^{3},$

and

$A_{5}(t)=1+26t+66t^{2}+26t^{3}+t^{4}.$

All nontrivial coefficients of the Eulerian polynomial $A_{d}(t)$ are positive
integers. In fact, as is well known,

$A_{d}(t)=\sum_{\sigma\in S_{d}}t^{\mathrm{desc}(\sigma)},$

where for a permutation $\sigma\in S_{d}$ of $\\{1,\ldots,d\\}$,
$\mathrm{desc}(\sigma)\in\mathbb{N}$ denotes the number of descents of
$\sigma$, the number of positions $i<d$ with $\sigma(i)>\sigma(i+1)$.

###  1.2. Toric geometry of polytopes

Given a convex lattice polytope $\Sigma\subset\mathbb{R}^{d}$, consider the
$\mathbb{C}$-vector space

$R_{\Sigma}:=\mathbb{C}\langle x_{{\bf v},r}\,|\,r\geq 0,{\bf v}\in
r\Sigma\cap\mathbb{Z}^{d}\rangle.$

This becomes a $\mathbb{C}$-algebra by the multiplication rule

$x_{{\bf v},r}\cdot x_{{\bf v}^{\prime},r^{\prime}}:=x_{{\bf v}+{\bf
v}^{\prime},r+r^{\prime}}.$

It is moreover graded, with $r$-th graded piece $R_{\Sigma,r}$ generated by
the variables $x_{{\bf v},r}$ for ${\bf v}\in r\Sigma\cap\mathbb{Z}^{d}$; in
particular $R_{\Sigma,0}\cong\mathbb{C}$. By definition, $E_{\Sigma}(t)$ is
the Hilbert series of the graded algebra $R_{\Sigma}$.

Recall that the lattice polytope $\Sigma$ is called normal, if for any $r>1$,
every lattice point ${\bf v}\in r\Sigma\cap\mathbb{Z}^{d}$ can be written as a
vector sum ${\bf v}=\sum_{i=1}^{r}{\bf v}_{i}$ of $r$ lattice points ${\bf
v}_{i}\in\Sigma\cap\mathbb{Z}^{d}$. If $\Sigma$ is normal, then $R_{\Sigma}$
is generated as a $\mathbb{C}$-algebra by its grade-$1$ piece $R_{\Sigma,1}$.
Let

$S_{\Sigma}:=\mathbb{C}[x_{{\bf v},1}\,|\,{\bf v}\in\Sigma\cap\mathbb{Z}^{d}]$

be the free graded algebra on the lattice points in $\Sigma$.

Consider $X_{\Sigma}=\operatorname{Proj}R_{\Sigma}$, a $d$-dimensional
projective toric variety. If $\Sigma$ is normal, we have a toric embedding
(compare [10, Chap. 14])

$\nu_{\Sigma}\colon
X_{\Sigma}=\operatorname{Proj}R_{\Sigma}\hookrightarrow\operatorname{Proj}S_{\Sigma}\cong\mathbb{P}^{\\#\\{\Sigma\cap\mathbb{Z}^{d}\\}-1}.$

In this case, the Ehrhart polynomial $L_{\Sigma}$ can be identified with the
Hilbert polynomial of $X_{\Sigma}$ in this embedding, and thus its Hilbert
series gets identified with the Ehrhart series of $\Sigma$.

######  Example 1.4.

For $\Sigma=\Delta_{d}\subset\mathbb{R}^{d}$ the standard $d$-simplex with
$d+1$ vertices and no interior lattice point, the graded algebra
$R_{\Delta_{d}}$ is the free polynomial algebra

$R_{\Delta_{d}}=\mathbb{C}[x_{{\bf v},1}\mid{\bf
v}\in\Delta_{d}\cap\mathbb{Z}^{d}]\cong\mathbb{C}[x_{0},\ldots,x_{d}]$

generated by variables corresponding to the origin, respectively the other
vertices. The attached projective embedding is $\nu_{\Delta_{d}}\colon
X_{\Delta_{d}}\cong\mathbb{P}^{d}$.

######  Example 1.5.

Consider the unit hypercube $\Sigma=C_{d}$ of dimension $d$, with $2^{d}$
vertices and no interior lattice point. As this is also normal, the graded
algebra $R_{C_{d}}$ is still generated by degree one variables. It is known
that the surjection $S_{C_{d}}\rightarrow R_{C_{d}}$ has a kernel

$I_{C_{d}}=\langle x_{{\bf v},1}x_{{\bf v}^{\prime},1}-x_{{\bf w},1}x_{{\bf
w}^{\prime},1}:{\bf v}+{\bf w}^{\prime}={\bf w}+{\bf v}^{\prime}\rangle\lhd
S_{C_{d}}$

generated by quadratic relations. We get the embedding

$\nu_{d}\colon
X_{C_{d}}\cong(\mathbb{P}^{1})^{d}\hookrightarrow\operatorname{Proj}S_{d}\cong\mathbb{P}^{2^{d}-1}$

of the projective toric variety $(\mathbb{P}^{1})^{d}$ by its complete linear
series $\mathcal{O}(1,\ldots,1)$.

It will be important below to spell out some special cases. For $d=2$, we have
the unit square as our polytope. Re-labelling degree-$1$ variables for easier
readability, the algebra $R_{C_{2}}$ is

$R_{C_{2}}\cong\mathbb{C}[x_{00},x_{01},x_{10},x_{11}]/\langle
x_{00}x_{11}-x_{01}x_{10}\rangle,$

corresponding to the Segre embedding
$\nu_{C_{2}}\colon(\mathbb{P}^{1})^{2}\hookrightarrow\mathbb{P}^{3}$, with
Ehrhard–Hilbert series

$E_{C_{2}}(t)=\frac{1+t}{(1-t)^{3}}=\frac{1-t^{2}}{(1-t)^{4}}.$

For $d=3$, and $\Sigma=C_{3}$, the unit cube, the algebra has a more
interesting structure. Again, re-labelling the variables in an obvious way, we
have

$R_{C_{3}}\cong\mathbb{C}[x_{000},x_{001},x_{010},x_{011},x_{100},x_{101},x_{110},x_{111}]/I_{C_{3}},$

where $I_{C_{3}}$ a codimension $4$ ideal, the ideal of the embedding
$\nu_{C_{3}}\colon(\mathbb{P}^{1})^{3}\hookrightarrow\mathbb{P}^{7}$. It has a
$9\times 16$ Gorenstein resolution: the $8$ generators satisfy $9$ quadric
relations related by $16$ cubic syzygies, as evidenced by the Ehrhart–Hilbert
series

(1.3)
$E_{C_{3}}(t)=\frac{1+4t+t^{2}}{(1-t)^{4}}=\frac{1-9t^{2}+16t^{3}-9t^{4}+t^{6}}{(1-t)^{8}}.$

######  Example 1.6.

Finally consider $\Sigma=\Xi_{d}$ be the $d$-dimensional cross-polytope, once
again a normal polytope, with $2d$ vertices and one interior lattice point,
the origin. Then it is clear from basic linear algebra that all relations
between the degree-$1$ variables are generated by the relations

$x_{{\bf e}_{j},1}x_{-{\bf e}_{j},1}=x_{{\mathbf{0}},1}^{2},$

one for each basis vector ${\bf e}_{j}\in\mathbb{R}^{d}$. Re-labellling the
degree-$1$ variables as $x_{0}$ for the origin, $x_{01},\ldots,x_{0d}$ for the
positive basis vectors, and $x_{11},\ldots,x_{1d}$ for the negative ones, we
thus have

$R_{\Xi_{d}}\cong\mathbb{C}[x_{0},x_{01},\ldots,x_{0d},x_{11},\ldots
x_{1d}]/\langle x_{0j}x_{1j}-x_{0}^{2}:1\leq j\leq d\rangle.$

The corresponding toric embedding is

$\nu_{\Xi_{d}}\colon
X_{\Xi_{d}}\hookrightarrow\operatorname{Proj}\mathbb{C}[x_{0},x_{01},\ldots,x_{0d},x_{11},\ldots
x_{1d}]\cong\mathbb{P}^{2d}$

presenting $X_{\Xi_{d}}$ as a complete intersection of $d$ quadric
hypersurfaces in $\mathbb{P}^{2d}$. The Ehrhart–Hilbert series is now

$E_{\Xi_{d}}(r)=\frac{(1+t)^{d}}{(1-t)^{d+1}}=\frac{(1-t^{2})^{d}}{(1-t)^{2d+1}}.$

In the special case $d=2$, we get

$\nu_{\Xi_{2}}\colon
X_{\Xi_{2}}\cong\left\\{x_{01}x_{11}=x_{02}x_{12}=x_{0}^{2}\right\\}\subset\mathbb{P}^{4},$

a canonical del Pezzo surface of degree $4$ with four $A_{1}$ singularities.

##  2\. The refinement

###  2.1. The definition

We follow Chapoton [4]. Fix a $d$-tuple of integer weighs ${\bf
a}=(a_{1},\ldots,a_{d})$. For integers $k$, consider the quantity

$L^{\bf a}_{\Sigma}(r,k)=\\#\\{v\in r\Sigma\cap\mathbb{Z}^{d}:{\bf
v}\cdot{{\bf a}}=k\\}.$

In other words, we are slicing the polytope with hyperplanes perpendicular to
the direction ${\bf a}$, and counting lattice points in the different layers.
We can then set

$L^{\bf a}_{\Sigma,r}(q)=\sum_{k\in\mathbb{Z}}L_{\Sigma}(r,k)q^{k},$

a Laurent polynomial in a new variable $q$, and consider the series

$E^{\bf
a}_{\Sigma}(t,q)=\sum_{r=0}^{\infty}\sum_{k\in\mathbb{Z}}L_{\Sigma}(r,k)\,q^{k}t^{r}.$

Setting $q=1$, these specialise to the value at $r$ of the Ehrhart polynomial
and the Ehrhart series, respectively. Chapoton [4] proves that under a non-
negativity assumption as well as genericity conditions on the weights, the
refined Ehrhard series is a rational function in $(t,q)$ with a controlled
denominator. We will not need these non-negativity and genericity assumptions
in what follows.

From the point of view of the algebra $R_{\Sigma}$, the weights furnish it
with a bigrading, with the generator $x_{\bf v}$ for ${\bf v}\in
r\Sigma\cap\mathbb{Z}^{d}$ having degree $(r,{\bf v}\cdot{\bf a})$. By
definition, the series $E^{\bf a}_{\Sigma}(t,q)$ is the Hilbert series of this
bigraded algebra.

###  2.2. The simplex

Return to $\Sigma=\Delta_{d}$, the unit $d$-simplex. Choose first ${\bf
a}=(1,\ldots,1)$. This is a non-generic choice in the sense of [4]: we are
slicing the $d$-simplex into parallel regular $(d-1)$-simplices. We clearly
have

$L^{(1,\ldots,1)}_{\Delta_{d}}(r,k)={k+d-1\choose d-1}\mbox{ for }0\leq k\leq
r.$

The polynomial

$L^{(1,\ldots,1)}_{\Delta_{d},r}(q)=\sum_{k=0}^{r}{k+d-1\choose d-1}q^{k}$

specialises at $q=1$ to the binomial coefficient ${r+d\choose d}$ by standard
rules of the Pascal triangle. It easy to compute [1], by elementary methods,
the corresponding refined Ehrhart series

$E^{(1,\ldots,1)}_{\Delta_{d}}(t,q)=\frac{1}{(1-t)(1-qt)^{d}}.$

Next choose ${\bf a}=(1,2,\ldots,d)$, a generic choice of slicing. There is no
short expression for the quantity $L^{(1,2,\ldots,d)}_{\Delta_{d}}(r,m)$. On
the other hand, the polynomial

(2.1) $L^{(1,2,\ldots,d)}_{{\Delta_{d}},r}(q)=\left[r+d\choose
d\right]_{q}\mbox{ for }0\leq k\leq r$

is known to be [8, 1] a $q$-binomial coefficient. So in this case, the
substitution $q=1$ specialises the $q$-binomial coefficient to the standard
one.

It is easy to write down the refined Eckhart series of the $d$-simplex in full
generality.

######  Proposition 2.1.

For an arbitrary $d$-tuple of integers ${\bf a}=(a_{1},\ldots,a_{d})$, we have

$E^{\bf a}_{\Delta_{d}}(t,q)=(1-t)^{-1}\prod_{j=1}^{d}(1-q^{a_{j}}t)^{-1}.$

###### Proof.

As noted before, the $\mathbb{C}$-algebra $R_{\Delta_{d}}$ is freely generated
by the variables $x_{0},x_{1},\ldots,x_{d}$ corresponding to the lattice
points ${\mathbf{0}},{\bf e}_{1},\ldots,{\bf e}_{d}$ of $\Delta_{d}$. The
bigraded degree of these generators is $(1,0)$, respectively $(1,a_{j})$ for
$j=1,\ldots,d$. The statement follows. ∎

For the specific weights ${\bf a}=(1,2,\ldots,d)$, we get

$E^{(1,2,\ldots,d)}_{\Delta_{d}}(t,q)=\prod_{j=0}^{d}(1-q^{j}t)^{-1}=\sum_{r=0}^{\infty}\left[k+d-1\choose
d-1\right]_{q}t^{r}$

where the second equality is one form of the $q$-binomial theorem, agreeing
with (2.1).

###  2.3. The cross-polytope

Let $\Sigma=\Xi_{d}$, the $d$-dimensional cross-polytope. The next result is
also straightforward.

######  Proposition 2.2.

For an arbitrary $d$-tuple of integers ${\bf a}=(a_{1},\ldots,a_{d})$, we have

(2.2) $E^{\bf
a}_{\Xi_{d}}(t,q)=\frac{(1+t)^{d}(1-t)^{d-1}}{\displaystyle\prod_{j=1}^{d}(1-q^{a_{j}}t)\prod_{j=1}^{d}(1-q^{-a_{j}}t)}.$

###### Proof.

As discussed in Example 1.2, the $\mathbb{C}$-algebra $R_{\Xi_{d}}$ is
generated by $x_{0}$ of bigraded degree $(1,0)$, and generators
$x_{0j},x_{1j}$ of bigraded degree $(1,a_{j})$, respectively $(1,-a_{j})$ for
$j=1,\ldots,d$. There are also $d$ relations of bidegree $(2,0)$ which form a
regular sequence. The statement follows after factoring the numerator in the
resulting Hilbert series, and cancelling one factor of $(1-t)$. ∎

Notice that the only time the right hand side of 2.2 simplifies further is if
there exists $1\leq j\leq d$ with $a_{j}=0$. Apart from the trivial case ${\bf
a}={\mathbf{0}}$, the only time all negative terms in the numerator cancel is

$E^{(0,\ldots,0,1)}_{\Xi_{d}}(t,q)=\frac{(1+t)^{d}}{(1-qt)(1-q^{-1}t)(1-t)^{d-1}}.$

In this case, we do not get $q$-deformations of the numerator of (1.1).
Indeed, recent work of Tielker [11] uses a different, non-linear refined count
of lattice points inside dilations of the cross-polytope to get an interesting
$q$-refinement of its Ehrhart numerator. Our last example will behave very
differently in this regard.

###  2.4. The hypercube

Finally consider $\Sigma=C_{d}$, the unit $d$-hypercube. At the $r$-th
dilation, we have a hypercube of side length $r$, whose $(1+r)^{d}$ lattice
points are weighted by $a_{j}$ in direction $j$, so

$L^{\bf a}_{C_{d},r}(q)=\prod_{j=1}^{d}\left(\sum_{i=0}^{r}q^{ia_{j}}\right).$

Summing each of these geometric series, we deduce

(2.3) $E^{\bf
a}_{C_{d}}(t,q)=\sum_{r=0}^{\infty}\left(\prod_{j=1}^{d}[r+1]_{q^{a_{i}}}\\!\\!\right)t^{r},$

where for a positive integer $k$,

$[k]_{q}=\frac{1-q^{k}}{1-q}$

is the corresponding $q$-integer.

For symmetric weights ${\bf a}=(1,\ldots,1)$, we get

(2.4)
$E^{(1,\ldots,1)}_{C_{d}}(t,q)=\sum_{r=0}^{\infty}[r+1]_{z}^{d}\,t^{r}=\frac{\tilde{A}_{d}(t,q)}{\prod_{j=0}^{d}(1-q^{j}t)}.$

Here the second equality is a formula of MacMahon and Carlitz [6, 3, 5], and
involves the two-variable refined Eulerian polynomial

$\tilde{A}_{d}(t,q)=\sum_{\sigma\in
S_{d}}q^{\mathrm{maj}(\sigma)}t^{\mathrm{desc}(\sigma)}.$

This enumerates permutations $\sigma\in S_{d}$ according to descent and major
index, where the latter is defined as the sum of all the descent positions in
$\sigma$. The first few examples of this refined Eulerian polynomial are

$\tilde{A}_{2}(t,z)=1+qt,\ \tilde{A}_{3}(t,q)=1+2(1+q)qt+q^{3}t^{2},$

and

$\tilde{A}_{4}(t,q)=1+(3+5q+3q^{2})qt+(3+5q+3q^{2})q^{3}t^{2}+q^{6}t^{3}.$

At $q=1$, the refined Eulerian polynomial specialises to the original one,
while (2.4) specialises to (1.2).

To evaluate more cases of (2.3), consider first $d=2$, so work with
$\Sigma=C_{2}$, the unit square. Recall that in this case, the algebra
$R_{C_{2}}$ can be described as

$R_{C_{2}}\cong\mathbb{C}[x_{00},x_{01},x_{10},x_{11}]/\langle
x_{00}x_{11}-x_{01}x_{10}\rangle.$

For a refinement with weights ${\bf a}=(a,b)$, the bigraded algebra
$R_{C_{2}}$ is generated by four generators of bidegrees $(1,0)$, $(1,a)$,
$(1,b)$ and $(1,a+b)$ with a single relation in bidegree $(2,a+b)$. We deduce

######  Proposition 2.3.

For arbitrary integer weights ${\bf a}=(a,b)$, we have

$E^{(a,b)}_{C_{2}}(t,q)=\frac{1-q^{a+b}t^{2}}{(1-t)(1-q^{a}t)(1-q^{b}t)(1-q^{a+b}t)}.$

Note that the right hand side does not simplify for general ${\bf a}$.
However, for ${\bf a}=(1,1)$, there is an obvious cancellation of one term,
and we get

$E^{(1,1)}_{C_{2}}(t,q)=\frac{1+qt}{(1-t)(1-qt)(1-q^{2}t)},$

the expected $q$-refinement of

$E_{C_{2}}(t)=\frac{1+t}{(1-t)^{3}}.$

On the other hand, for arbitrary weights, we deduce the following identity, a
generalization of the classical MacMahon–Carlitz identity (2.4).

######  Corollary 2.4.

For arbitrary integers $(a,b)$, we have

$\sum_{r=0}^{\infty}[r+1]_{q^{a}}[r+1]_{q^{b}}t^{r}=\frac{1-q^{a+b}t^{2}}{(1-t)(1-q^{a}t)(1-q^{b}t)(1-q^{a+b}t)}.$

Next let $\Sigma=C_{3}$, the unit cube. Recall that here we have

$R_{C_{3}}\cong\mathbb{C}[x_{000},x_{001},x_{010},x_{011},x_{100},x_{101},x_{110},x_{111}]/I_{C_{3}},$

with $I_{C_{3}}$ a codimension $4$ ideal discussed in Example 1.5. For an
arbitrary refinement with weights ${\bf a}=(a,b,c)$, we get the following.

######  Proposition 2.5.

$E^{\bf a}_{C_{3}}(t,q)=\frac{N^{\bf
a}_{C_{3}}(t,q)}{(1-t)\left(1-q^{a}t\right)\left(1-q^{b}t\right)\left(1-q^{c}t\right)\left(1-q^{a+b}t\right)\left(1-q^{a+c}t\right)\left(1-q^{b+c}t\right)\left(1-q^{a+b+c}t\right)},$

where the numerator is

$\begin{array}[]{ll}N^{\bf a}_{C_{3}}(t,q)=&1\\\
&-\left(3q^{a+b+c}+q^{2a+b+c}+q^{a+2b+c}+q^{a+b+2c}+q^{a+b}+q^{a+c}+q^{b+c}\right)t^{2}\\\
&+\left(2q^{a+b+c}+2q^{2a+b+c}+2q^{a+2b+c}+2q^{a+b+2c}+2q^{2a+2b+c}+2q^{2a+b+2c}+2q^{a+2b+2c}+2q^{2(a+b+c)}\right)t^{3}\\\
&-\left(q^{2a+2b+c}+q^{2a+b+2c}+q^{a+2b+2c}+3q^{2(a+b+c)}+q^{3a+2b+2c}+q^{2a+3b+2c}+q^{2a+2b+3c}\right)t^{4}\\\
&+q^{3(a+b+c)}t^{6}.\end{array}$

###### Proof.

This can easily be obtained from the known (or machine computable) resolution
of the ideal
$I_{C_{3}}\lhd\mathbb{C}[x_{000},x_{001},x_{010},x_{011},x_{100},x_{101},x_{110},x_{111}]$,
which is compatible with the bigrading. We omit the details. ∎

From this formula, it is easy to derive the specialisation ${\bf a}=(1,1,1)$
corresponding to symmetric weights:

$E^{(1,1,1)}_{C_{3}}(t,q)=\frac{1+2q(1+q)t+q^{3}t^{2}}{(1-t)(1-qt)(1-q^{2}t)(1-q^{3}t)},$

recovering the refined Eulerian polynomial $\tilde{A}_{3}(t,q)$ in the
numerator. A slightly more general, but attractive, specialisation is ${\bf
a}=(a,1,1)$ which gives the following.

######  Corollary 2.6.

For any integer $a$, we have

$\sum_{r=0}^{\infty}[r+1]_{q^{a}}[r+1]_{q}^{2}\,t^{r}=\frac{1+q\left(q^{a}+1\right)t-2\left(q^{2}+q+1\right)q^{a+1}t^{2}+\left(q^{a}+1\right)q^{a+3}t^{3}+q^{2a+4}t^{4}}{(1-t)(1-qt)\left(1-q^{2}t\right)\left(1-q^{a}t\right)\left(1-q^{a+1}t\right)\left(1-q^{a+2}t\right)}.$

For $d\geq 4$, the graded algebra $R_{C_{d}}=S_{d}/I_{C_{d}}$ corresponding to
the embedding
$\nu_{d}\colon(\mathbb{P}^{1})^{d}\hookrightarrow\mathbb{P}^{2^{d}-1}$ is
getting more and more complicated. Already for $d=4$, one has a codimension
$11$ embedding, defined by $55$ quadrics; the formula corresponding to that of
Proposition 2.5, even if computable, would have thousands of terms.

Turning the ideas around, let us see whether we can learn anything about the
ideal $I_{C_{d}}$ of the embedding $\nu_{d}$ in general. For general $d$,
formula (1.3) becomes

$E_{C_{d}}(t)=\frac{A_{d}(t)}{(1-t)^{d+1}}=\frac{1-c_{d}t^{2}+O(t^{3})}{(1-t)^{2^{d}}},$

for a certain positive integer $c_{d}$. On the right hand side, we have the
Hilbert series of a variety embedded in $\mathbb{P}^{2^{d}-1}$, whose ideal is
cut out by $c_{d}$ quadrics. In the middle, remembering that all coefficients
of $A_{d}(t)$ are positive, we have the Hilbert series of a finite length
module over a free polynomial algebra in $d+1$ generators, corresponding to a
free subalgebra of $R_{C_{d}}$ formed of elements of a regular sequence; this
is basically a form of Noether normalization, remembering that the image
variety is projective of dimension $d$. See [9, Section 4] for more details on
these ideas in the present context.

The fact that the formula of Proposition 2.3 does not simplify for general
weights tells us that this procedure is not compatible with general bigradings
even for $d=2$. The fact that there is simplification for symmetric weights
suggests that it is compatible with that particular bigrading. On the other
hand, we have the general MacMahon–Carlitz formula

$E^{(1,\ldots,1)}_{C_{d}}(t,q)=\frac{\tilde{A}_{d}(t,q)}{\prod_{j=0}^{d}(1-q^{j}t)},$

with numerator with positive coefficients only. This suggests that for all
$d$, $R_{C_{d}}$ should be finite as a bigraded module over a bigraded
subalgebra generated by a regular sequence of elements of bidegrees
$(1,0),(1,1),\ldots,(1,d)$.

To find this subalgebra explicitly, for $0\leq k\leq d$ let

$y_{k}=\sum_{\begin{subarray}{c}{\bf v}\in C_{d}\cap\mathbb{Z}^{d}\\\ {\bf
v}\cdot(1,...,1)=k\end{subarray}}x_{{\bf v},1}\in S_{d}.$

These elements are bihomogeneous of degree $(1,k)$ for the symmetric
bigrading.

######  Proposition 2.7.

The images $[y_{0}],\ldots,[y_{d}]\in R_{C_{d}}$ form a regular sequence. The
quotient

$T_{C_{d}}=R_{C_{d}}/\langle[y_{0}],\ldots,[y_{d}]\rangle$

is a finite dimensional bigraded $\mathbb{C}$-vector space with Hilbert series
equal to the refined Eulerian polynomial $\tilde{A}_{d}(t,q)$.

###### Proof.

111We would like to thank Mateusz Michalek for suggesting this proof. The
graded ring $R_{C_{d}}$ is Cohen–Macaulay (in fact even Gorenstein). Thus in
order to show that its elements $[y_{0}],\ldots,[y_{d}]\in R_{C_{d}}$ form a
regular sequence, it is enough to show that the ideal they generate has
codimension $(d+1)$. The latter is equivalent to the fact that the projective
intersection of $\operatorname{im}\nu_{C_{d}}$ with the $d+1$ hyperplanes in
$\mathbb{P}^{2^{d}-1}$ given by the conditions $y_{k}=0$ for $0\leq k\leq d$
is empty. Assume otherwise: choose a point
$([a_{1}:b_{1}],[a_{2}:b_{2}],...,[a_{d}:b_{d}])\in(\mathbb{P}^{1})^{d}$ whose
image under $\nu_{C_{d}}$ is on the vanishing locus of all these hyperplanes.
The vanishing of the first linear form tells us that the product of all the
$a_{i}$ is zero; we may assume $a_{1}=0$, and hence $b_{1}\neq 0$. The
vanishing of the second linear form tells us that the sum of $d$-products of
all but one of the $a_{i}$ and the complimentary $b_{j}$ is 0. Given that
$a_{1}=0$, we deduce $b_{1}a_{2}a_{3}\cdots a_{d}=0$, so without loss of
generality $a_{2}=0$ and so $b_{2}\neq 0$. Continuing in this way, we get that
on the vanishing locus of all the hyperplanes, all the $a_{i}$ coordinates are
zero. But the vanishing of the last linear form tells us that the product of
all the $b_{i}$ is also zero, which is a contradiction. This proves the first
statement.

It follows in particular that $R_{C_{d}}/\langle[y_{0}],\ldots,[y_{d}]\rangle$
is a finite dimensional bigraded vector space. By the MacMahon–Carlitz
formula, its bigraded Hilbert series indeed equals the refined Eulerian
polynomial $\tilde{A}_{d}(t,q)$. ∎

## References

  * [1] P. Adeyemo, Grassmannians in the lattice points of dilations of standard simplex, manuscript 2022.
  * [2] M. Beck and S. Robins, Computing the Continuous Discretely. Integer-Point Enumeration in Polyhedra, UTM 154, Springer, New York, 2007.
  * [3] L. Carlitz, A combinatorial property of $q$-Eulerian numbers, Amer. Math. Monthly 82 (1975), 51–54.
  * [4] F. Chapoton, $q$-analogues of Ehrhart polynomials, Proc. Edinb. Math. Soc. (2) 59 (2016), 339–358.
  * [5] C.-O. Chow and I. M. Gessel, On the descent numbers and major indices for the hyperoctahedral group, Advances in Applied Mathematics 38 (2007), 275–301.
  * [6] P. A. MacMahon, Two applications of general theorems in combinatory analysis: (1) to the theory of inversion of permutations; (2) to the ascertainment of the numbers of terms in the development of a determinant which has amongst its elements an arbitrary number of zeros, Proc. London Math. Soc. (2) 15 (1916), 314–321.
  * [7] M. Michalek and B. Sturmfels, Invitation to non-linear algebra, Graduate Studies in Mathematics 172, AMS, Providence, R.I., 2021.
  * [8] R. P. Stanley, Ordered structures and partitions, Memoirs of the AMS 119, Providence, R.I., 1972.
  * [9] R. P. Stanley, Magic labelings of graphs, symmetric magic squares, systems of parameters, and Cohen-Macaulay rings, Duke M. J. 43 (1976), 511–531.
  * [10] B. Sturmfels, Grobner bases and convex polytopes, University Lecture Series 8, AMS, Providence, R.I., 1996.
  * [11] E. Tielker, Weighted Ehrhart series and a type-B analogue of a formula of MacMahon, arXiv:2203.15774.