# Cloud on-demand emulation of quantum dynamics with tensor networks

Kemal Bidzhiev<EMAIL_ADDRESS>Aleksander Wennersteen
<EMAIL_ADDRESS>Mourad Beji  Mario Dagrada  Mauro
D’Arcangelo  Sebastián Grijalva  Anne-Claire Le Hénaff  Anton Quelle  Alvin
Sashala Naik  PASQAL, 7 rue Léonard de Vinci, 91300 Massy, France

(August 28, 2024)

###### Abstract

We introduce a tensor network based emulator, simulating a programmable analog
quantum processing unit (QPU). The software package is fully integrated in a
cloud platform providing a common interface for executing jobs on a HPC
cluster as well as dispatching them to a QPU device. We also present typical
emulation use cases in the context of Neutral Atom Quantum Processors, such as
evaluating the quality of a state preparation pulse sequence, and solving
Maximum Independent Set problems by applying a parallel sweep over a set of
input pulse parameter values, for systems composed of a large number of
qubits.

##  I Introduction

Quantum technology proposes leveraging the laws of quantum mechanics to
process information in a framework that can enable solving hard computational
problems more efficiently than classical alternatives. A new scientific age of
quantum information has been precipitated by this idea, and constructing the
first performant quantum processing units (QPUs) has become a key goal in
recent years. As QPUs gain in effectiveness and reliability, experimental
demonstrations have come forward solving specific computational tasks which
may show an advantage over classical methods, such as sampling problems [1, 2,
3, 4], and quantum simulation implementations [5, 6]. Platforms for QPUs have
been proposed based on diverse physical architectures, including neutral
atoms, superconducting circuits, trapped ions and photons.

As well as continuing the experimental development of QPUs, it is becoming
increasingly important to simulate the behavior of particular physical QPU
architectures, referred to as _emulation_. Emulators encode realistic physical
constraints including hardware-specific time dynamics and noise effects.
Furthermore, they can inform about the expected performance of the QPU and
benchmarking [7]. This allows the user to circumvent overhead by performing
test jobs e.g. on a High-Performance Computing (HPC) backend before forwarding
the full job on to the QPU. Moreover, integrating QPU and emulator backends
through cloud access, offering the choice of backend to the user, allows
efficient workflows as part of a full-stack quantum computer. The hybrid
workflow would send suitable routines to the quantum processor, while
executing larger computational tasks on the emulator backend. Systems
supporting such hybrid workflows with emulated quantum devices already provide
valuable services to researchers [8, 9, 10, 11, 12].

In most cases, the simulation of quantum systems has a complexity that grows
exponentially with the system size. A significant instrument which has emerged
to tackle this challenge are tensor networks [13], providing a powerful
structure to represent complex systems efficiently. The framework of tensor
networks has allowed the efficient simulation and study of quantum systems,
such as, most importantly for our purposes, the numerical study of time-
evolution of quantum states which are ground states of gapped local
Hamiltonians. Further, they have found uses in the identification of new
phases of matter [14, 15], as a structure for processing data [16, 17], or the
numerical benchmarking of state-of-the-art QPU experiments with large qubit
numbers [5, 18].

Figure 1: Role of emulation in a workflow for large quantum systems. Pulse-
level design tools are complemented with tensor network representations of the
corresponding quantum states. Algorithms performing time evolution take into
account Hamiltonian dynamics as well as realistic hardware conditions, which
are simulated on a HPC backend. These can be used to prepare optimized
instructions on a QPU. We present a cloud-based platform integrating these
tasks on-demand.

The paper is structured as follows. In Section II we review the elementary
setup for the numerical simulation of quantum dynamics and introduce relevant
tensor-network methods. Then, in Section III we describe the framework for
integrating HPC-supported numerical algorithms into the cloud. We demonstrate
applications of the emulator in IV, simulating a neutral atom quantum
processor in three examples: producing 1D $\mathbb{Z}_{n}$-ordered states,
preparing 2D antiferromagnetic states and performing a parallel sweep of pulse
parameters to search for Maximum Independent Sets on a particular graph. We
conclude with comments about extensions to the system and a discussion about
emulated environments in heterogeneous computing.

##  II Numerical simulation of Quantum Dynamics

###  II.1 Time-evolution of a Quantum System

_Qubits_ , or 2-level quantum systems, are the fundamental building blocks of
the quantum systems that we will study. A _quantum state_ for an $N$-qubit
system is an element of a $2^{N}$-dimensional Hilbert space. It can be
represented by a complex vector $|\psi\rangle=\sum_{x}c_{x}|x\rangle$ of size
$2^{N}$, where we choose $\\{|x\rangle\\}$ to be an orthonormal basis of the
space. We concentrate on an orthonormal basis that is built from the
eigenstates of the $\hat{\sigma}^{z}$ operator on each site, which is called
the _computational basis_. Its elements can be written as $|s_{1}s_{2}\cdots
s_{N}\rangle$, where each $s_{i}\in\\{0,1\\}$. The strings $s_{1}s_{2}\cdots
s_{N}$ are often called _bitstrings_. A general quantum state is thus written
as:

$|\psi\rangle=\sum_{s_{1}\cdots s_{N}}c_{s_{1}\cdots s_{N}}|s_{1}\cdots
s_{N}\rangle,$ (1)

where the $c_{s_{1}\cdots s_{N}}$ is a complex number corresponding to a
probability amplitude for a given set of indices $s_{1}s_{2}\cdots s_{N}$
representing the outcome of a measurement. The QPU takes an initial state and
evolves it through a series of gates and/or analog pulses. To obtain the final
state, it makes a projective measurement on it with respect to a certain
measurement basis, repeating the cycle until a certain condition is met. Of
these operations, the most computationally expensive calculation is to compute
the final time-evolved state, which we describe next.

The time evolution of a quantum state is governed by the _Schrödinger
Equation_

$i\frac{d|\psi\rangle}{dt}=\hat{H}(t)|\psi\rangle,$ (2)

where we consider a particular time-dependent _Hamiltonian_ operator
$\hat{H}(t)$ which is a sum of $J$ time-dependent _control_ terms and
$N(N-1)/2$ two-body _interaction_ terms $V_{ij}=V(\mathbf{r}_{ij})$ that
depend on the relative positions $\mathbf{r}_{ij}$ of each pair of qubits
$i,j$ in $\\{1,\ldots,N\\}$:

$\hat{H}(t)=\sum_{j=1}^{J}\hat{H}_{\text{ctrl}}^{(j)}(t)+\sum_{i>j}V_{ij}\hat{P}_{i}\hat{P}_{j}.$
(3)

with $\hat{P}_{i}$ a local operator constructed with Pauli matrices111The
Pauli matrices are $\hat{\sigma}^{x}=\begin{pmatrix}0&1\\\
1&0\end{pmatrix},\hat{\sigma}^{y}=\begin{pmatrix}0&-i\\\
i&0\end{pmatrix},\text{ and }\hat{\sigma}^{z}=\begin{pmatrix}1&0\\\
0&-1\end{pmatrix}$. We represent their action on the $i^{th}$ site of the
quantum system by a Kronecker product of matrices:
$\hat{P}_{i}=\mathbb{I}\otimes\cdots\otimes\hat{P}\otimes\cdots\mathbb{I}$,
for
$\hat{P}\in\\{\mathbb{I},\hat{\sigma}^{x},\hat{\sigma}^{y},\hat{\sigma}^{z}\\}$
occupying the $i$-th entry . Composite operators such as
$\hat{P}_{i}\hat{P}_{j}$ are then formed via the matrix product. and
$\hat{H}_{\text{ctrl}}^{(j)}(t)$ the $j$-th control term, a time-dependent
function representing the action of the experimental parameters affecting the
system, usually constructed using waveform generators and applied through
laser or microwave fields.

The solution of the differential equation (2) can be obtained by a variety of
methods [19]: one can approximate the time-evolution operator by
$\mathcal{U}(t)\approx\prod\hat{U}(\delta t)$, with $\hat{U}(\delta
t)=\exp(-i\hat{H}\delta t)$ for small time intervals $\delta t$, so that
$\hat{H}$ can be considered time independent during $\delta t$. However,
computing the exponential of an operator can be prohibitive for large systems
[20]. One can also expand $U(\delta t)$ into a product of simpler terms up to
a desired degree of accuracy, depending on the form of $\hat{H}$ [21]. In
addition, one can try to numerically solve the full differential equation,
such as by using Runge-Kutta methods [22]. Another family of methods consists
of implementing an algorithm simulating the time evolution in a quantum
computer, for example, by decomposing the exponential into evolution blocks
using the Trotter product formulas [23, 24, 25], or by using tools from
quantum machine learning [26, 27]. One could also try to _learn_ the final
state vector using tools from machine learning [28], or apply Monte Carlo
techniques [29, 30, 31].

Finally, a very successful approach is given by tensor network techniques,
which is the approach we choose to implement as part of our emulation
framework. For readers unfamiliar with tensor networks, we briefly review them
below. An in-depth detailed description can be found for example in [16, 32,
33, 13].

###  II.2 Representing quantum states and their evolution using tensor
networks

Simualting time evolution of non-equilibrium interacting systems is
challenging due to the growth of the entanglement entropy with time [34].
Tensor network operations are a matter of linear algebra, and the size of the
matrices involved scales with the entanglement of the system under
consideration. Different types of tensor networks exhibit different scaling,
depending on the geometry and dimension of the physical system [35, 36, 37,
38], however, for all of these methods, the matrices involved eventually
become prohibitively large. Furthermore, only some existing time evolution
algorithms can efficiently encode the dynamics of long-range Hamiltonians [39,
40], which are necessary to simulate two-body interactions (3) with
Hamiltonians contains terms of the form $V_{ij}\sim
1/|\mathbf{r}_{ij}|^{\alpha}$.

In order to simulate a quantum many-body system (1), we encode the wave
function into a _matrix-product state_ (MPS) [13, 41]. MPS allow efficient
representations of the $N$-order tensor $c_{s_{1}\cdots s_{N}}$ as a product
of $N$ smaller tensors $A[\sigma]^{s_{\sigma}}_{i_{\sigma-1}i_{\sigma}}$:

$|\psi\rangle=\sum_{\\{s\\}}\sum_{\\{i\\}}A[1]^{s_{1}}_{i_{1}}A[2]^{s_{2}}_{i_{1}i_{2}}\cdots
A[N]^{s_{N}}_{i_{N-1}}|s_{1}\cdots s_{N}\rangle,$ (4)

where each tensor index $i_{\sigma}\in{i}$ runs from $1$ to $\chi_{\sigma}$.
The bond dimension $\chi_{\sigma}$ controls the size of each $A[\sigma]$,
which determines the computation complexity of the simulation [34].

One of the most successful algorithms able to treat long-range interactions,
while maintaining a sufficiently small bond dimension, is the Time-Dependent
Variational Principle (TDVP) [42]. TDVP constrains the time evolution to the
MPS manifold $\mathcal{M}(\psi[A])$ with a given $\chi$ by projecting the
right hand side of the Schrödinger equation onto the tangent space
$T\mathcal{M}(\psi[A])$ :

$\frac{d|\psi[A]\rangle}{dt}=-i\hat{\Pi}_{T\mathcal{M}(\psi[A])}\hat{H}|\psi\rangle.$
(5)

The resulting equations (5) are nonlinear and couple all degrees of freedom in
the MPS. Approaches to solving (5) have different ways of controlling accuracy
and convergence. In the emulator presented here, we implement 2-site TDVP to
deal with the core numerical simulations. Details of the implementation are
given in Appendix B.

##  III Platform architecture

In this section we discuss the full architecture of our platform. We describe
how the HPC cluster is integrated with cloud services to provide our quantum
device emulation. The emulator, which we will refer to as EMU-TN, includes the
constraints of a particular QPU as detailed in Section II.1, implementing the
tensor network algorithms of Section II.2 as the main numerical backend.
Below, we describe the input to the emulator, the pre-processing, as well as
post-processing before an output is returned. We finally discuss the cloud
infrastructure, including orchestration, scheduling and dispatching.

###  III.1 Encoding the dynamics information

Our platform takes as input an encoded abstract representation of the control
parameters $\hat{H}_{\text{ctrl}}(t)$. It then performs the required quantum
dynamic evolution, applies the measurements to the final state, and returns
the readout data to the design tool.

EMU-TN includes a JSON-parser, which is important to uniformize the
information sent to the cluster. We take as initial expression the Hamiltonian
in eq. (3), with a set of $J$ control fields and the positions
$\\{\mathbf{r}_{ij}\\}$ of the qubits in the register. Each control field
$\hat{H}_{\text{ctrl}}^{(j)}(t)$ acts on a subset $I^{(j)}$ of qubits in the
register,

$\hat{H}_{\text{ctrl}}^{(j)}(t)=\Omega^{(j)}(t)\sum_{i\in{I^{(j)}}}\hat{P}^{(j)}_{i}.$
(6)

For example, if the $j$-th control is a global field of magnitude
$\Omega^{(j)}(t)$ in the $z$-axis, then $\hat{P}^{(j)}=\hat{\sigma}_{i}^{z}$
and $\hat{H}_{\text{ctrl}}^{(j)}=\Omega^{(j)}(t)\sum_{i\in\\{1\ldots
N\\}}\hat{\sigma}^{z}_{i}$. To parse this information, a design tool will have
to specify $\Omega^{(j)}(t)$ as either an array of numbers, the initial and
final values plus a total duration, or a generating function. Additionally,
locality has to be specified with the indices $I^{(j)}$.

Finally one needs to inform about the initial state preparation and
measurement of the state. The initial state is usually experimentally
convenient to produce, such as a _product state_ indicated typically by a
bitstring. Platforms such as neutral atom quantum processors further allow
choosing the positions of the qubits, determining the interaction strength
between them. More generically, initial state conditions can be provided as
labels produced by the design tool. The measurement protocol can consist of a
basis (or a set of bases) where the measurements take place and a number of
repetitions of the experiment, in order to generate statistics to estimate
observables. This can be specified through a number of _runs_ of the job
submission, together with a label for the desired basis, as long as it is
compatible with the emulated QPU device.

There is an increasing list of software packages that offer tools for the
design and customization of pulse sequences for different types of qubit
architectures and quantum processing tasks [43, 44, 45, 46, 47, 48]. These
have been used in the context of the design of quantum algorithms, of quantum
optimal control as well as the research of the dynamical properties of many-
body quantum systems. A great majority of them use directly or indirectly a
state-vector representation of the quantum state, and are available as open
source tools for personal use. EMU-TN can be included as an alternative of the
solver, thereby extending the types of tasks that can be handled by these
design tools.

###  III.2 Cloud infrastructure

We present how the cloud platform dispatches jobs to the classical emulation
and quantum computing backends.

####  III.2.1 Orchestration and scheduling

The infrastructure, illustrated in Fig. 2, is split between public cloud
services hosted on OVHCloud and orchestrated through Kubernetes and an HPC
cluster hosted in a private data center running a Slurm scheduler. We describe
below their integration for providing quantum computing as a service.

Figure 2: Architecture of the cloud platform providing quantum computing as a
service. The cloud services are deployed on a Kubernetes cluster and use a
microservice architecture, each microservice exposing a public REST API to end
users. The HPC cluster where EMU-TN runs, consists of 10 DGX A100 systems from
NVIDIA (See Appendix A for more details) and uses Slurm, an HPC cluster
management and job scheduling system, for resource management and job
submission

Kubernetes [49] is a container orchestration system centered around cloud
computing, public, private, or hybrid. It simplifies the deployment and
operation of containerised applications on a set of resources, as well as the
scaling of said resources. Slurm [50], on the other hand, was created with HPC
workloads in mind. It provides prioritization aware scheduling, job queues,
topology aware job allocations, reservations and backfill policies. These
features led to software like Slurm often being chosen to manage HPC clusters
over newer solutions bespoke for other compute paradigms. Slurm submits a
general Bash script to run its job, which can also include a sequence of
containerised applications.

It is possible to run the same workloads on a Kubernetes pod, either as a
long-running service responsible for recieving and executing a job or
scheduled using Kubernetes compatible schedulers. However, a downside of these
approaches is that it would be challenging to integrate the on-premise HPC
cluster with the Kubernetes scheduler. Not scheduling the jobs on the cluster,
which would mean just running a job execution service in a Kubernetes pod,
would put limitations both on how effectively the resources in the cluster are
used, and ultimately on the performance of the HPC applications.

As cloud computing and HPC become increasingly mixed, extensions between
Kubernetes and Slurm have recently seen a lot of progress, see for example
[51, 52, 53]. However, these approaches have generally centered around
submitting Kuberenetes pods, i.e. containers, to Slurm systems. This involves
relying on Kubernetes to schedule jobs. The Kubernetes scheduler does not
provide out-of-the-box support to work with often changing and highly
heterogeneous environments consisting of several QPU types, CPUs and GPUs. In
this work we use a custom scheduler for our pulse sequences. Thus, rather than
relying on mechanisms similar to those proposed in [51, 52, 53] a service that
directly submits regular Slurm jobs was chosen as this integrates better with
our existing hardware and is more similar to the way jobs are scheduled on the
QPU.

Moreover, we think that for the foreseeable future hybrid quantum-classical
workloads will constitute a large part of quantum jobs. The classical part of
these can often be GPU-accelerated, motivating the need to be able to request
such compute resources when needed. For example, NVIDIA recently introduced a
platform for hybrid quantum-classical computing, known as Quantum Optimized
Device Architecture (QODA) [54].

####  III.2.2 Job dispatching

To connect our two clusters, we used the recently introduced Slurm REST API
[55] to submit jobs to Slurm from the regular Backend service. A “Slurm
connector” service has been developed to act as a bridge between the existing
PASQAL cloud Services and the HPC cluster. The service creates a control
script on-demand for the job which is submitted to Slurm as the batch script.
This script then takes care of executing all the higher-level logic of the
job. The service is also responsible for setting up the job with appropriate
resources and inserting all necessary information for future communication
with the cloud service.

Each pulse sequence is encoded in a JSON format and sent to PASQAL cloud
Services by making an HTTP POST request to the job submission endpoint. The
pulse sequence format is the same for QPUs and emulators. We extend the body
of the HTTP request with a Configuration argument that allows the user to
control the parameters of the emulator, see Code Sample 1. This design allows
executing the same quantum program seamlessly on both QPUs and emulators,
while also allowing control of the numerical simulation parameters such as the
maximum bond dimension and the number of cores.

The internal logic of the job submission endpoint validates the pulse sequence
and the emulation configuration. The validation step checks that the request
and the underlying data such as the sequence and emulator configuration are
correctly formatted, that the requested resources are accessible to the user
and, for a real device, that the sequence is physically valid. It then
dispatches the request to the appropriate scheduler. In Code Sample 1 below,
we show an example of sending a serialized pulse sequence to PASQAL Cloud
Services for execution. To build the sequences we have used Pulser [46], an
open source package for designing neutral atom QPU sequences. To communicate
with cloud services we have used the PASQAL Cloud Services Python SDK [56],
which provides helper functions to the end user to request the APIs.

Code Sample 1: Running jobs through the cloud service.

⬇

1from pulser import Pulse, Sequence, Register

2from pulser.devices import MockDevice

3from sdk import Configuration, Endpoints, SDK

4import json

5

6# 1. Define Qubit Register

7N = 6 # Qubits per side

8reg = Register.square(side=N, spacing=7, prefix="q")

9

10# 2. Define Pulse Sequence and encode to JSON

11seq = Sequence(reg, MockDevice)

12seq.declare_channel(

13 name="ch0",

14 channel_id="rydberg_global"

15)

16seq.add(

17 Pulse.ConstantPulse(

18 duration=1000, # ns

19 amplitude=2, # rad/us

20 detuning=-6, # rad/us

21 phase=0),

22 channel="ch0"

23)

24encoded_seq = seq.to_abstract_repr() # JSON file

25

26# 3. Send encoded sequence to PASQAL Cloud Services

27# Get API keys at

28# https://portal.pasqal.cloud/api-keys

29

30endpoints = Endpoints(

31 core="https://apis.pasqal.cloud/core",

32 account="https://apis.pasqal.cloud/account"

33)

34

35sdk = SDK(

36 client_id="my_client_id",

37 client_secret="my_client_secret",

38 endpoints=endpoints,

39)

40

41# Configure Job

42my_job = {"runs": 1000} # Number of shots

43config = Configuration( # Configuration of emulator

44 dt=10.0,

45 precision="normal",

46 extra_config={"max-bond-dim": 100}

47)

48

49# Create batch on the Cloud service, wait until the execution completes, and
get results

50batch = sdk.create_batch(

51 encoded_seq,

52 device_type="EMU_TN",

53 jobs=[my_job],

54 configuration=config,

55 fetch_results=True,

56)

57res = [job.result for job in batch.jobs.values()]

In this context, EMU-TN is built as a Singularity container and executed as a
Slurm job. The Slurm job consists of the numerical simulation and is preceded
and followed by custom-made scripts to communicate with the cloud platform.
Since the pulse sequence has been previously validated, the pre-processing
script only translates the user configuration into the format required by EMU-
TN. After the completion of the simulation, the post-processing script
collects the results and communicates with the cloud service to ensure the
results are correctly stored and accessible by the frontend. This architecture
allows the emulators to be stand-alone and also to be easily executable
outside the overall platform. It also allows for easy extension for more
frequent two-way communication between the platform and the emulator, should
it be desirable in the future.

####  III.2.3 Other considerations

Storing bistring results or heavily compressed tensor-states in a columnar
database works very well as they can be serialised to short text strings. For
larger datasets this can cause performance degradation in the database and
eventually the dataset will be too large for the specific database
implementation. Thus, object storage may be useful in such cases, such as for
the tensor network representations produced by the emulator. Since obtaining
such a state is computationally costly, the platform makes this state
available as a serialised object that can be de-serialised for further
processing as required by the numerical algorithm.

A second consideration is that not every user request constitutes a feasible
computational objective. Some pulse sequences (e.g. those that introduce steep
variations) usually require a larger matrix size to simulate. It could also be
that the algorithm does not converge to the desired accuracy. There is an
unavoidable iterative process where a candidate pulse sequence and quantum
system is studied by the user on physical grounds and successive verification
with exact diagonalization methods for smaller systems is crucial. Once a
promising job is selected, the data is deserialized in the cloud core service
where a validation script checks that it describes a valid quantum program,
for a given emulated device based on existing QPUs. This ensures that the
pulse sequence is physically implementable on the hardware. The Slurm-
connector also does some simple validation of Slurm configuration.

##  IV Applications

In this section we discuss some applications of the EMU-TN framework,
motivated by typical tasks that can be performed with a Neutral Atoms Quantum
Processor [57, 58, 59]. The native Hamiltonian will be taken as the following
particular form of (3):

$\hat{H}(t)=\frac{\Omega(t)}{2}\sum_{i=1}^{N}\hat{\sigma}^{x}_{i}-\delta(t)\sum_{i=1}^{N}\hat{n}_{i}+\sum_{i>j}\frac{C}{|\mathbf{r}_{ij}|^{6}}\hat{n}_{i}\hat{n}_{j}$
(7)

where $\hat{n}=|1\rangle\langle 1|=(\mathbb{I}-\hat{\sigma}^{z})/2$ is the
projector into the excited state and $C$ is a constant that depends on the
properties of the targeted excited state. After the quantum evolution and
readout, one obtains a bitstring of 0’s and 1’s representing excited and
ground-state qubits, respectively. To calculate the expectation value of a
given observable, many cycles need to be performed in order to generate enough
statistics. The output of the emulator provides bistrings, which simulates the
actual output of the QPU.

We remark that EMU-TN can be applied to other Hamiltonians, for example long-
range “XY” Hamiltonians that can be constructed with neutral atom [60] and
trapped ion [61] arrays.

###  IV.1 1D Regular Lattice: Rydberg Crystals

One notable property of the Hamiltonian (7) is that, in certain regimes, the
simultaneous excitation of two neighbouring atoms is suppressed. This
phenomenon is called Rydberg blockade [62, 63], and it takes place when the
interaction term $C/r^{6}$ exceeds the Rabi frequency $\Omega$ (it is
therefore distance dependent). The Rydberg blockade mechanism is the main
source of entanglement in neutral atom systems, since the maximally entangled
state $|01\rangle+|10\rangle$ becomes energetically favored over the product
state $|11\rangle$.

The Rydberg blockade has been used in [64] to create highly regular excitation
patterns in 1D chains of neutral atoms. These structures are known as Rydberg
crystals and, for experimentally reasonable values of atomic spacing and Rabi
frequency, they can be prepared in such a way as to display $\mathbb{Z}_{n}$
order for $n=2,3,4$. The situation is depicted schematically in Fig. 3.

Figure 3: Examples of 1D Rydberg crystals for a chain of 13 atoms. A solid
circle represents an excited atom, while an empty one represents an atom in
the ground state. The $\mathbb{Z}_{2}$, $\mathbb{Z}_{3}$ and $\mathbb{Z}_{4}$
orders are characterized by excitations separated by respectively 1, 2, or 3
atoms in the ground state.

We first present the results of a test simulation on a chain of 16 atoms. For
the pulse sequence shown in the top part of Fig. 4 with
$\Omega_{\text{max}}=6.3\ \text{rad}/\mu s$ and a spacing of 4 $\mu m$, we
expect the final ground state of the system to be in the $\mathbb{Z}_{3}$
order. The excitation probability for each qubit in its final state is also
reported for both EMU-TN and an exact Schrödinger equation solver in Fig. 4,
together with a heatmap representing the time evolution of the excitation
probability for each qubit given by the exact solver and EMU-TN respectively.

Figure 4: Comparison between an exact Schrödinger equation solver (using QuTiP
[65]) and the emulator (EMU-TN) for $\mathbb{Z}_{3}$ Rydberg crystal
preparation on a 1D chain of 16 atoms placed 4 $\mu m$ apart. The driving
sequence represents a pulse of $t_{tot}=10\mu s$ with
$\Omega_{\text{max}}=6.3\ \text{rad}/\mu s$. The histogram represents the
excitation probability for each qubit in the final state, while the heatmaps
represent the evolution over time of the excitation probability obtained with
the exact solver and EMU-TN respectively.

Next, we present the results for the same simulation but with a chain of 37
atoms, which is intractable by an exact solver. The pulse in this case is
chosen to be twice as long in order to allow correlations to spread from the
borders all the way to the center of the chain and observe a clearer
alternation of the excitations. The excitation probability in the final state
and its time evolution are reported in Fig. 5.

Figure 5: Same simulation as Fig. 4, but with a chain of 37 qubits and a pulse
of 20 $\mu s$.

###  IV.2 2D Regular Lattice: Antiferromagnetic State Preparation

One of the most common tasks one can try to simulate is the preparation of
particular states with a regular register. This has been used in milestone
implementations in programmable neutral-atom arrays of hundreds of qubits [5,
6].

A typical pulse sequence (Fig. 6, above) represents a path through the phase
diagram in the thermodynamic limit that ends in a point of the expected
antiferromagnetic phase, which has been analytically studied before. We
present the results from the sampling of the evolved MPS as well as from a
straightforward implementation of exact diagonalization solved numerically on
a local computer. Using exact diagonalization it is possible to run
simulations just above the 20-qubit range, but it soon becomes impractical
once the number of qubits increases. On the other hand, adjusting the bond
dimension of the tensor network algorithm, one can aim to explore the behavior
of sequences in the range of 20 to 60 qubits in a comparably short time.
Moreover, the flexibility of the cluster approach allows for parallelization
of tasks which provide information about parameter sweeps without the user
having to spend time adapting their code. In Figure 7, we include information
about the elapsed time of the performed simulations, for different array sizes
and bond dimensions, with and without access to GPUs.

Figure 6: Comparison of sampling results ($10^{3}$ samples) from the prepared
quantum state between exact diagonalization (ED) and the emulator (EMU-TN) for
increasing array sizes. The driving sequence (same structure as the one in
Fig. 4, above) represents a pulse of $t_{\text{tot}}=2.7\mu s$ moving in the
phase diagram of the Hamiltonian towards the antiferromagnetic phase. The
$5\times 7$ and $5\times 9$ arrays required $10^{4}$ samples to resolve better
the peaks (1% and 0.4% respectively), since the short pulse breaks
adiabaticity for large arrays. Figure 7: Computation time for increasing array
sizes. We compare only CPU use and then an added GPU backend, for different
maximal bond dimension, $\chi_{\text{max}}$.

###  IV.3 Maximum Independent Set on Unit Disk graphs

We present a typical use case of neutral atom devices for two-dimensional
registers. The possibility of placing atoms in arbitrary 2D configurations
allows to solve certain classes of hard graph problems on the QPU. An
interesting application is solving the Maximum Independent Set (MIS) problem
on Unit Disk (UD) graphs [66].

A graph is an object with nodes and connections between nodes. A UD graph is a
graph where two nodes are connected if an only if they are closer than a
certain minimal distance. By representing the nodes of a graph with neutral
atoms and matching the Rydberg blockade radius with the UD graph minimal
distance, one can establish a direct mapping where a connection exists between
atoms that are within a blockade distance of each other. The quantum evolution
of such a system is naturally restricted to those sectors of the Hilbert space
where excitations of connected atoms are forbidden. These configurations
translate to independent sets of a graph, i.e. subsets of nodes that are not
directly connected to each other. Driving the quantum evolution in such a way
as to produce as many excitations as possible, one would obtain with high
probability as the outcome of a measurement an independent set of high
cardinality, representing a good candidate solution to the MIS problem.

The best pulse sequence to find the MIS of a graph is based on the same
adiabatic pulses already shown in Figs. 4 and 6, partly since Rydberg crystals
and antiferromagnetic states can be seen as the MIS of regular 1D and 2D
lattices. As shown in [67], however, for arbitrary graphs it is often
necessary to optimize the pulse shape in order to find the MIS with a high
enough probability. We present in Fig. 8 a one-parameter family of pulses
where the parameter $t_{c}$ controls at which time the detuning goes from
negative to positive. Rather than sweeping over the free parameter
sequentially, one can send multiple simultaneous jobs for each value to be
tested, and finally select the parameter giving the best result. Parallelizing
this task is a direct advantage of the cluster infrastructure, and the cloud
service provides a dedicated interface in which all jobs can be organized and
retrieved, saving large amounts of time. We performed a sweep over 9
parameters for a pulse of 5 $\mu s$ on a graph with 30 nodes. The results are
shown in Fig. 8, below, where the optimal zero-crossing for the detuning is
found to be at around 1.6 $\mu s$.

Figure 8: Parameter sweep. We arrange 30 qubits in 2D (above), where the
dashed circles show half the blockade radius. The pulse sequence (middle,
$t_{\text{tot}}=5\mu s$) contains a tunable parameter $t_{c}$ controlling the
time when the detuning goes from negative to positive. We sweep in parallel
for 9 values of $t_{c}$. The results (below) show the average size of the MIS
found for each value of $t_{c}$. The entire sweep, took about 2 hours of cloud
compute time.

While the example shown here is centered on pulse parameters, other tasks
suited for parallelization are also implementable, e.g. averages over noise
and disorder realizations. Being able to quickly set up and run studies of
this kind for systems of many qubits opens the door for creative algorithm
design that is at the same time compatible with a QPU. We remark that in 1D
systems, the TDVP algorithm has been used to study dynamics for
$\mathcal{O}(10^{2})$ spins [68, 69].

##  V Discussion

Due to the rising interest in developing and using quantum processing units
integrated in a modern computational infrastructure, we envision the emergence
of workflows requiring intensive use of emulators. They are not just a
placeholder for the actual quantum device, but serve for diagnostics, resource
estimation and management. Given the current price and availability of QPUs,
emulation is crucial in the initial phases of an application, providing a
faster feedback loop during the R&D phase, as well as enforcing a realistic
level of constraints based on a given quantum computing architecture. This is
also relevant for designing upgrades to a quantum device, as the relevant
figures of merit for expected behavior can be compared quickly.

By providing a seamless interface between emulators running locally, emulators
of large systems in a HPC environment and the QPUs, the researcher can
distribute the most convenient compute option. By further centralising this in
a workflow and allowing the scheduler to manage the optimal software resource
one can both reduce the total time and cost to solution. Workflows will also
have to take into account whether the computational task should be executed in
the shortest walltime, the most energy efficient way, the most cost-effective
way or some weighted combination. The computational community has taken great
strides in workflow management, and the ever increasing use of machine
learning has created new fields like MLOps (Machine Learning Operations). It
is expected that the quantum computing community strives for such state of the
art workflows, and to be able to integrate into existing ones, to achieve real
industrial advantage. This will include not only execution of algorithms, but
also providing applicable training, and testing.

In this paper we described the components at the basis of such workflows. By
leveraging tensor network algorithms for quantum dynamics that can exploit the
resources provided by an HPC cluster, we discussed how a cloud on-demand
service can make use of a cluster backend to simplify access to the required
computational operations. The numerics can be easily adapted to different
architectures and devices, separating the role of design tools from the
computation of quantum dynamics. Other techniques for the study of dynamical
quantum systems can be adapted as well: future objectives include considering
the low qubit ($<20$) range, where state-vector evolution methods allow
implementing a microscopic model of hardware noise effects and other state-
preparation and readout errors [70]. This includes implementing
parallelization routines on the algorithm management side, controlled
exploitation of the cluster resources for heavy computational loads and
efficient solvers for other dynamical equations (quantum master equations and
stochastic Schrödinger equation). A suite of benchmarks for common tasks in
the low qubit range is also desirable for understanding speed bottlenecks, and
will be the subject of a follow up work. In the case of larger arrays, the use
of tensor network structures like PEPS or Tree Tensor Networks, provide other
convenient representations for two-dimensional systems and their application
to a time evolution algorithms is under active study [71, 72, 73, 74, 75]

Developing workflows is much easier with a fully fledged cloud backend to
orchestrate the tasks. This enables creating workflows that automatically
assess the validity of the submitted pulse sequence for the QPU. Simulating
the dynamics from a particular pulse sequence requires attention to tuning and
understanding how to set the algorithm parameters. The cloud platform can
provide the user suggestions for future hyperparameters based on live
monitoring of the simulation. Finally, the worker nodes on the cluster have
large amounts of local storage in a so-called “scratch” directory. This allows
storing large amounts of data cheaply on a temporary basis. A future extension
of the platform could allow subsequent job-steps to use previous results,
recycling previous computational work for new tasks.

We believe the platform presented here can greatly improve reproducibility and
benchmarking of research in the developing quantum industry, as well as
allowing a broader scope of researchers, engineers and developers to
understand, explore and design the behavior of quantum systems and quantum
devices.

###### Acknowledgements.

We thank Loïc Henriet, Caroline de Groot and Matthieu Moreau for valuable
conversations, suggestions for this manuscript and collaboration in related
work.

## References

  * Arute et al. [2019] Frank Arute, Kunal Arya, Ryan Babbush, Dave Bacon, Joseph C Bardin, Rami Barends, Rupak Biswas, Sergio Boixo, Fernando GSL Brandao, David A Buell, et al.  Quantum supremacy using a programmable superconducting processor.  _Nature_ , 574(7779):505–510, 2019.
  * Zhong et al. [2020] Han-Sen Zhong, Hui Wang, Yu-Hao Deng, Ming-Cheng Chen, Li-Chao Peng, Yi-Han Luo, Jian Qin, Dian Wu, Xing Ding, Yi Hu, et al.  Quantum computational advantage using photons.  _Science_ , 370(6523):1460–1463, 2020.
  * Wu et al. [2021] Yulin Wu, Wan-Su Bao, Sirui Cao, Fusheng Chen, Ming-Cheng Chen, Xiawei Chen, Tung-Hsun Chung, Hui Deng, Yajie Du, Daojin Fan, et al.  Strong quantum computational advantage using a superconducting quantum processor.  _Physical review letters_ , 127(18):180501, 2021\.
  * Madsen et al. [2022] Lars S Madsen, Fabian Laudenbach, Mohsen Falamarzi Askarani, Fabien Rortais, Trevor Vincent, Jacob FF Bulmer, Filippo M Miatto, Leonhard Neuhaus, Lukas G Helt, Matthew J Collins, et al.  Quantum computational advantage with a programmable photonic processor.  _Nature_ , 606(7912):75–81, 2022.
  * Scholl et al. [2021] Pascal Scholl, Michael Schuler, Hannah J Williams, Alexander A Eberharter, Daniel Barredo, Kai-Niklas Schymik, Vincent Lienhard, Louis-Paul Henry, Thomas C Lang, Thierry Lahaye, et al.  Quantum simulation of 2d antiferromagnets with hundreds of rydberg atoms.  _Nature_ , 595(7866):233–238, 2021.
  * Ebadi et al. [2021] Sepehr Ebadi, Tout T Wang, Harry Levine, Alexander Keesling, Giulia Semeghini, Ahmed Omran, Dolev Bluvstein, Rhine Samajdar, Hannes Pichler, Wen Wei Ho, et al.  Quantum phases of matter on a 256-atom programmable quantum simulator.  _Nature_ , 595(7866):227–232, 2021.
  * Jaschke et al. [2022] Daniel Jaschke, Alice Pagano, Sebastian Weber, and Simone Montangero.  Ab-initio two-dimensional digital twin for quantum computer benchmarking.  _arXiv preprint arXiv:2210.03763_ , 2022.
  * Guerreschi et al. [2020] Gian Giacomo Guerreschi, Justin Hogaboam, Fabio Baruffa, and Nicolas PD Sawaya.  Intel quantum simulator: A cloud-ready high-performance simulator of quantum circuits.  _Quantum Science and Technology_ , 5(3):034007, 2020.
  * Humble et al. [2021] Travis S Humble, Alexander McCaskey, Dmitry I Lyakh, Meenambika Gowrishankar, Albert Frisch, and Thomas Monz.  Quantum computers for high-performance computing.  _IEEE Micro_ , 41(5):15–23, 2021.
  * Ravi et al. [2021] Gokul Subramanian Ravi, Kaitlin N Smith, Pranav Gokhale, and Frederic T Chong.  Quantum computing in the cloud: Analyzing job and machine characteristics.  In _2021 IEEE International Symposium on Workload Characterization (IISWC)_ , pages 39–50. IEEE, 2021.
  * Mandrà et al. [2021] Salvatore Mandrà, Jeffrey Marshall, Eleanor G Rieffel, and Rupak Biswas.  Hybridq: A hybrid simulator for quantum circuits.  In _2021 IEEE/ACM Second International Workshop on Quantum Computing Software (QCS)_ , pages 99–109. IEEE, 2021.
  * McCaskey et al. [2018] Alexander J McCaskey, Eugene F Dumitrescu, Dmitry Liakh, Mengsu Chen, Wu-chun Feng, and Travis S Humble.  A language and hardware independent approach to quantum–classical computing.  _SoftwareX_ , 7:245–254, 2018.
  * Schollwöck [2011] Ulrich Schollwöck.  The density-matrix renormalization group in the age of matrix product states.  _Annals of Physics_ , 326(1):96–192, January 2011\.
  * Iregui et al. [2014] Juan Osorio Iregui, Philippe Corboz, and Matthias Troyer.  Probing the stability of the spin-liquid phases in the kitaev-heisenberg model using tensor network algorithms.  _Physical Review B_ , 90(19):195102, 2014.
  * Jiang and Ran [2017] Shenghan Jiang and Ying Ran.  Anyon condensation and a generic tensor-network construction for symmetry-protected topological phases.  _Physical Review B_ , 95(12):125107, 2017.
  * Cichocki et al. [2016] Andrzej Cichocki, Namgil Lee, Ivan Oseledets, Anh-Huy Phan, Qibin Zhao, and Danilo P. Mandic.  Tensor networks for dimensionality reduction and large-scale optimization: Part 1 low-rank tensor decompositions.  _Foundations and Trends® in Machine Learning_ , 9(4-5):249–429, 2016.  doi: 10.1561/2200000059.
  * Stoudenmire and Schwab [2016] Edwin Stoudenmire and David J Schwab.  Supervised learning with tensor networks.  _Advances in Neural Information Processing Systems_ , 29, 2016.
  * Pan et al. [2022] Feng Pan, Keyang Chen, and Pan Zhang.  Solving the sampling problem of the sycamore quantum circuits.  _Physical Review Letters_ , 129(9):090502, 2022\.
  * Higham [2008] Nicholas J Higham.  _Functions of matrices: theory and computation_.  SIAM, 2008.
  * Moler and Van Loan [2003] Cleve Moler and Charles Van Loan.  Nineteen dubious ways to compute the exponential of a matrix, twenty-five years later.  _SIAM review_ , 45(1):3–49, 2003.
  * Barthel and Zhang [2020] Thomas Barthel and Yikang Zhang.  Optimized lie–trotter–suzuki decompositions for two and three non-commuting terms.  _Annals of Physics_ , 418:168165, July 2020.  doi: 10.1016/j.aop.2020.168165.  URL https://doi.org/10.1016/j.aop.2020.168165.
  * Hairer et al. [2010] Ernst Hairer, Syvert Paul Nørsett, and Gerhard Wanner.  _Solving ordinary differential equations. I: Nonstiff problems. 3rd corrected printing_ , volume 8.  2010\.
  * Cirstoiu et al. [2020] Cristina Cirstoiu, Zoe Holmes, Joseph Iosue, Lukasz Cincio, Patrick J Coles, and Andrew Sornborger.  Variational fast forwarding for quantum simulation beyond the coherence time.  _npj Quantum Information_ , 6(1):1–10, 2020.
  * Pastori et al. [2022] Lorenzo Pastori, Tobias Olsacher, Christian Kokail, and Peter Zoller.  Characterization and verification of trotterized digital quantum simulation via hamiltonian and liouvillian learning.  _arXiv preprint arXiv:2203.15846_ , 2022.
  * Childs et al. [2021] Andrew M Childs, Yuan Su, Minh C Tran, Nathan Wiebe, and Shuchen Zhu.  Theory of trotter error with commutator scaling.  _Physical Review X_ , 11(1):011020, 2021.
  * Caro et al. [2022] Matthias C Caro, Hsin-Yuan Huang, Nicholas Ezzell, Joe Gibbs, Andrew T Sornborger, Lukasz Cincio, Patrick J Coles, and Zoë Holmes.  Out-of-distribution generalization for learning quantum dynamics.  _arXiv preprint arXiv:2204.10268_ , 2022.
  * Gibbs et al. [2022] Joe Gibbs, Zoë Holmes, Matthias C Caro, Nicholas Ezzell, Hsin-Yuan Huang, Lukasz Cincio, Andrew T Sornborger, and Patrick J Coles.  Dynamical simulation via quantum machine learning with provable generalization.  _arXiv preprint arXiv:2204.10269_ , 2022.
  * Vicentini et al. [2022] Filippo Vicentini, Damian Hofmann, Attila Szabó, Dian Wu, Christopher Roth, Clemens Giuliani, Gabriel Pescia, Jannes Nys, Vladimir Vargas-Calderón, Nikita Astrakhantsev, et al.  Netket 3: Machine learning toolbox for many-body quantum systems.  _SciPost Physics Codebases_ , page 007, 2022.
  * Aoki et al. [2014] Hideo Aoki, Naoto Tsuji, Martin Eckstein, Marcus Kollar, Takashi Oka, and Philipp Werner.  Nonequilibrium dynamical mean-field theory and its applications.  _Reviews of Modern Physics_ , 86(2):779, 2014\.
  * De Vega and Alonso [2017] Inés De Vega and Daniel Alonso.  Dynamics of non-markovian open quantum systems.  _Reviews of Modern Physics_ , 89(1):015001, 2017\.
  * Plenio and Knight [1998] Martin B Plenio and Peter L Knight.  The quantum-jump approach to dissipative dynamics in quantum optics.  _Reviews of Modern Physics_ , 70(1):101, 1998\.
  * Cichocki et al. [2017] Andrzej Cichocki, Namgil Lee, Ivan Oseledets, Anh-Huy Phan, Qibin Zhao, Masashi Sugiyama, and Danilo P. Mandic.  Tensor networks for dimensionality reduction and large-scale optimization: Part 2 applications and future perspectives.  _Foundations and Trends® in Machine Learning_ , 9(6):249–429, 2017.  doi: 10.1561/2200000067.  URL https://doi.org/10.1561/2200000067.
  * Bridgeman and Chubb [2017] Jacob C Bridgeman and Christopher T Chubb.  Hand-waving and interpretive dance: an introductory course on tensor networks.  _Journal of physics A: Mathematical and theoretical_ , 50(22):223001, 2017.
  * Latorre [2007] José I Latorre.  Entanglement entropy and the simulation of quantum mechanics.  _Journal of Physics A: Mathematical and Theoretical_ , 40(25):6689, 2007.
  * Verstraete and Cirac [2004] F. Verstraete and J. I. Cirac.  Renormalization algorithms for quantum-many body systems in two and higher dimensions, 2004.  URL https://arxiv.org/abs/cond-mat/0407066.
  * Shi et al. [2006] Y.-Y. Shi, L.-M. Duan, and G. Vidal.  Classical simulation of quantum many-body systems with a tree tensor network.  _Physical Review A_ , 74(2), August 2006.  doi: 10.1103/physreva.74.022320.  URL https://doi.org/10.1103/physreva.74.022320.
  * Tagliacozzo et al. [2009] L. Tagliacozzo, G. Evenbly, and G. Vidal.  Simulation of two-dimensional quantum systems using a tree tensor network that exploits the entropic area law.  _Physical Review B_ , 80(23), December 2009.  doi: 10.1103/physrevb.80.235127.  URL https://doi.org/10.1103/physrevb.80.235127.
  * Orús [2014] Román Orús.  A practical introduction to tensor networks: Matrix product states and projected entangled pair states.  _Annals of Physics_ , 349:117–158, 2014.  ISSN 0003-4916.
  * Zaletel et al. [2015a] Michael P. Zaletel, Roger S. K. Mong, Christoph Karrasch, Joel E. Moore, and Frank Pollmann.  Time-evolving a matrix product state with long-ranged interactions.  _Physical Review B_ , 91(16), April 2015a.  doi: 10.1103/physrevb.91.165112.  URL https://doi.org/10.1103/physrevb.91.165112.
  * Paeckel et al. [2019] Sebastian Paeckel, Thomas Köhler, Andreas Swoboda, Salvatore R. Manmana, Ulrich Schollwöck, and Claudius Hubig.  Time-evolution methods for matrix-product states.  _Annals of Physics_ , 411:167998, December 2019.
  * Oseledets [2011] I. V. Oseledets.  Tensor-train decomposition.  _SIAM Journal on Scientific Computing_ , 33(5):2295–2317, January 2011.
  * Haegeman et al. [2016] Jutho Haegeman, Christian Lubich, Ivan Oseledets, Bart Vandereycken, and Frank Verstraete.  Unifying time evolution and optimization with matrix product states.  _Physical Review B_ , 94(16):165116, 2016.
  * Teske et al. [2022] Julian David Teske, Pascal Cerfontaine, and Hendrik Bluhm.  qopt: An experiment-oriented software package for qubit simulation and quantum optimal control.  _Physical Review Applied_ , 17(3):034036, 2022\.
  * Li et al. [2022a] Boxi Li, Shahnawaz Ahmed, Sidhant Saraogi, Neill Lambert, Franco Nori, Alexander Pitchford, and Nathan Shammah.  Pulse-level noisy quantum circuits with qutip.  _Quantum_ , 6:630, 2022a.
  * Alexander et al. [2020] Thomas Alexander, Naoki Kanazawa, Daniel J Egger, Lauren Capelluto, Christopher J Wood, Ali Javadi-Abhari, and David C McKay.  Qiskit pulse: Programming quantum computers through the cloud with pulses.  _Quantum Science and Technology_ , 5(4):044006, 2020.
  * Silvério et al. [2022] Henrique Silvério, Sebastián Grijalva, Constantin Dalyac, Lucas Leclerc, Peter J Karalekas, Nathan Shammah, Mourad Beji, Louis-Paul Henry, and Loïc Henriet.  Pulser: An open-source package for the design of pulse sequences in programmable neutral-atom arrays.  _Quantum_ , 6:629, 2022.
  * Rossignolo et al. [2022] Marco Rossignolo, Thomas Reisser, Alastair Marshall, Phila Rembold, Alice Pagano, Philipp J Vetter, Ressa S Said, Matthias M Müller, Felix Motzoi, Tommaso Calarco, et al.  QuOCS: The quantum optimal control suite.  _arXiv preprint arXiv:2212.11144_ , 2022.
  * Goerz et al. [2022] Michael H Goerz, Sebastián C Carrasco, and Vladimir S Malinovsky.  Quantum optimal control via semi-automatic differentiation.  _arXiv preprint arXiv:2205.15044_ , 2022.
  * [49] Cloud Native Computing Foundation.  Kubernetes.  URL https://kubernetes.io/.
  * Yoo et al. [2003] Andy B. Yoo, Morris A. Jette, and Mark Grondona.  Slurm: Simple linux utility for resource management.  In Dror Feitelson, Larry Rudolph, and Uwe Schwiegelshohn, editors, _Job Scheduling Strategies for Parallel Processing_ , pages 44–60, Berlin, Heidelberg, 2003. Springer Berlin Heidelberg.  ISBN 978-3-540-39727-4.
  * Wickberg [2022] Tim Wickberg.  Slurm and/or/vs kubernetes.  SuperComputing, 2022.  URL https://slurm.schedmd.com/SC22/Slurm-and-or-vs-Kubernetes.pdf.
  * Zhou et al. [2020] Naweiluo Zhou, Yiannis Georgiou, Li Zhong, Huan Zhou, and Marcin Pospieszny.  Container orchestration on hpc systems, 2020.  URL https://arxiv.org/abs/2012.08866.
  * Lublinsky et al. [2022] Boris Lublinsky, Elise Jennings, and Viktória Spišaková.  A kubernetes ’bridge’ operator between cloud and external resources, 2022\.  URL https://arxiv.org/abs/2207.02531.
  * [54] Tim Costa.  NVIDIA QODA: The platform for hybrid quantum-classical computing.  URL https://developer.nvidia.com/qoda/.
  * Rini [2020] Nate Rini.  Rest api.  Slurm User Group Meeting, 2020.  URL https://slurm.schedmd.com/SLUG20/REST_API.pdf.
  * Services [2021] PASQAL Cloud Services.  Pasqal cloud sdk.  https://github.com/pasqal-io/cloud-sdk, 2021.
  * Henriet et al. [2020] Loïc Henriet, Lucas Beguin, Adrien Signoles, Thierry Lahaye, Antoine Browaeys, Georges-Olivier Reymond, and Christophe Jurczak.  Quantum computing with neutral atoms.  _Quantum_ , 4:327, 2020.
  * Saffman et al. [2010] Mark Saffman, Thad G Walker, and Klaus Mølmer.  Quantum information with rydberg atoms.  _Reviews of modern physics_ , 82(3):2313, 2010\.
  * Morgado and Whitlock [2021] M Morgado and S Whitlock.  Quantum simulation and computing with rydberg-interacting qubits.  _AVS Quantum Science_ , 3(2):023501, 2021.
  * de Léséleuc et al. [2017] Sylvain de Léséleuc, Daniel Barredo, Vincent Lienhard, Antoine Browaeys, and Thierry Lahaye.  Optical control of the resonant dipole-dipole interaction between rydberg atoms.  _Physical Review Letters_ , 119(5):053202, 2017\.
  * Jurcevic et al. [2014] Petar Jurcevic, Ben P Lanyon, Philipp Hauke, Cornelius Hempel, Peter Zoller, Rainer Blatt, and Christian F Roos.  Quasiparticle engineering and entanglement propagation in a quantum many-body system.  _Nature_ , 511(7508):202–205, 2014.
  * Jaksch et al. [2000] D. Jaksch, J. I. Cirac, P. Zoller, S. L. Rolston, R. Côté, and M. D. Lukin.  Fast quantum gates for neutral atoms.  _Physical Review Letters_ , 85(10):2208–2211, September 2000.  doi: 10.1103/physrevlett.85.2208.  URL https://doi.org/10.1103/physrevlett.85.2208.
  * Browaeys and Lahaye [2020] Antoine Browaeys and Thierry Lahaye.  Many-body physics with individually controlled rydberg atoms.  _Nature Physics_ , 16(2):132–142, January 2020\.  doi: 10.1038/s41567-019-0733-z.  URL https://doi.org/10.1038/s41567-019-0733-z.
  * Bernien et al. [2017] Hannes Bernien, Sylvain Schwartz, Alexander Keesling, Harry Levine, Ahmed Omran, Hannes Pichler, Soonwon Choi, Alexander S. Zibrov, Manuel Endres, Markus Greiner, Vladan Vuletić, and Mikhail D. Lukin.  Probing many-body dynamics on a 51-atom quantum simulator.  _Nature_ , 551(7682):579–584, Nov 2017.  ISSN 1476-4687.
  * Johansson et al. [2012] J Robert Johansson, Paul D Nation, and Franco Nori.  Qutip: An open-source python framework for the dynamics of open quantum systems.  _Computer Physics Communications_ , 183(8):1760–1772, 2012.
  * Pichler et al. [2018] Hannes Pichler, Sheng-Tao Wang, Leo Zhou, Soonwon Choi, and Mikhail D. Lukin.  Quantum optimization for maximum independent set using rydberg atom arrays.  2018\.  doi: 10.48550/ARXIV.1808.10816.  URL https://arxiv.org/abs/1808.10816.
  * Ebadi et al. [2022] S. Ebadi, A. Keesling, M. Cain, T. T. Wang, H. Levine, D. Bluvstein, G. Semeghini, A. Omran, J.-G. Liu, R. Samajdar, X.-Z. Luo, B. Nash, X. Gao, B. Barak, E. Farhi, S. Sachdev, N. Gemelke, L. Zhou, S. Choi, H. Pichler, S.-T. Wang, M. Greiner, V. Vuletić, and M. D. Lukin.  Quantum optimization of maximum independent set using rydberg atom arrays.  _Science_ , 376(6598):1209–1215, 2022.  doi: 10.1126/science.abo6587.  URL https://www.science.org/doi/abs/10.1126/science.abo6587.
  * Yang and White [2020] Mingru Yang and Steven R White.  Time-dependent variational principle with ancillary krylov subspace.  _Physical Review B_ , 102(9):094315, 2020.
  * Li et al. [2022b] Jheng-Wei Li, Andreas Gleis, and Jan Von Delft.  Time-dependent variational principle with controlled bond expansion for matrix product states.  _arXiv preprint arXiv:2208.10972_ , 2022b.
  * De Léséleuc et al. [2018] Sylvain De Léséleuc, Daniel Barredo, Vincent Lienhard, Antoine Browaeys, and Thierry Lahaye.  Analysis of imperfections in the coherent optical excitation of single atoms to rydberg states.  _Physical Review A_ , 97(5):053803, 2018.
  * Zaletel et al. [2015b] Michael P Zaletel, Roger SK Mong, Christoph Karrasch, Joel E Moore, and Frank Pollmann.  Time-evolving a matrix product state with long-ranged interactions.  _Physical Review B_ , 91(16):165112, 2015b.
  * Bauernfeind and Aichhorn [2020] Daniel Bauernfeind and Markus Aichhorn.  Time dependent variational principle for tree tensor networks.  _SciPost Physics_ , 8(2):024, 2020.
  * Rams and Zwolak [2020] Marek M Rams and Michael Zwolak.  Breaking the entanglement barrier: Tensor network simulation of quantum transport.  _Physical review letters_ , 124(13):137701, 2020\.
  * Kloss et al. [2020] Benedikt Kloss, David Reichman, and Yevgeny Bar Lev.  Studying dynamics in two-dimensional quantum lattices using tree tensor network states.  _SciPost Physics_ , 9(5):070, 2020.
  * Vanhecke et al. [2021] Bram Vanhecke, Laurens Vanderstraeten, and Frank Verstraete.  Symmetric cluster expansions with tensor networks.  _Physical Review A_ , 103(2):L020402, 2021.
  * Fishman et al. [2022] Matthew Fishman, Steven R. White, and E. Miles Stoudenmire.  The ITensor Software Library for Tensor Network Calculations.  _SciPost Phys. Codebases_ , page 4, 2022.  doi: 10.21468/SciPostPhysCodeb.4.  URL https://scipost.org/10.21468/SciPostPhysCodeb.4.
  * Fishman and Stoudenmire [2022] Matthew Fishman and Miles Stoudenmire.  itensors tdvp.  https://github.com/ITensor/ITensorTDVP.jl, 2022.
  * Bezanson et al. [2017] Jeff Bezanson, Alan Edelman, Stefan Karpinski, and Viral B Shah.  Julia: A fresh approach to numerical computing.  _SIAM review_ , 59(1):65–98, 2017.  URL https://doi.org/10.1137/141000671.

## Appendix

###  V.1 Computation Cluster

As an example of a modern HPC cluster infrastructure in which the algorithms
mentioned in the main text are expected to be executed, we describe PASQAL’s
computation cluster. Note that the amount of computing power of different
providers varies widely, as does their availability and cost. This is a
crucial element to consider when intensive computational tasks are needed.

Both CPUs and GPUs can be used to compute complex scientific algorithms,
although they are often used for different types of tasks. CPUs are generally
more versatile and can handle a wide range of tasks, including running the
operating system and other software, managing memory and input/output, and
performing complex calculations. They are well suited for tasks that require
sequential processing, such as executing instructions in a specific order or
following a set of rules. GPUs, on the other hand, are specialized processors
that are designed to handle the large number of calculations needed to render
images and video. They are particularly well suited for tasks that can be
parallelized, meaning that they can be broken down into smaller pieces
consisting of identical operations on potentially different data.

_GPU.-_ The NVIDIA A100 GPU delivers subsequent acceleration for a wide range
of applications in AI, data analytics, and high-performance computing (HPC). A
particularly parallel algorithm can efficiently scale to thousands of GPUs or,
with NVIDIA Multi-Instance GPU (MIG) technology, be partitioned into seven GPU
instances to accelerate workloads of all sizes. Furthermore, third-generation
Tensor Cores accelerate every precision for diverse workloads. At PASQAL, we
currently have 80 A100 NVIDIA GPUs distributed into 10 DGX A100 NVIDIA
systems. These servers are linked to each other through Mellanox’s
proprietatry 200Gbps Ethernet fibre optics connections. Each DGX A100 system
is itself comprised of 2 Dual 64-core AMD EPYC 7742 CPUs and 8 NVIDIA A100
SXM4 GPUs. Each GPU is connected through NVLINK to each other, with offers a
GPU memory pool of 320 GB of shared GPU memory for one node, in addition to
the 1 TB of RAM memory for the AMD CPUs.

_CPU.-_ With two 64-core EPYC CPUs and 1TB of system memory, the DGX A100
boasts high performance even before the GPUs are considered. The architecture
of the AMD EPYC “Rome” CPUs is outside the scope of this article, but offers
an elegant design of its own. Each CPU provides 64 processor cores (supporting
up to 128 threads), 256MB L3 cache, and eight channels of DDR4-3200 memory
(which provides the highest memory throughput of any mainstream x86 CPU). The
15TB high-speed NVMe SSD Scratch Space allows for very fast data read/write
speeds which definitely impacts total processing speeds when large data needs
to be uploaded from disk to RAM memory at different steps of the computation.
In addition to the on-node scratch space the cluster has a 15TB distributed
file system providing the users home directories as well as software and
containers. Using Slurm’s partitions we separate internal, external and
production usage coming from the cloud service herein.

###  V.2 TDVP implementation details

For the applications shown in the Section IV, Tensor network calculations are
implemented by using the ITensor.jl [76] and ITensorTDVP.jl [77] packages
written in Julia language [78]. As a time evolution algorithm we chose 2-site
TDVP with the time step $\delta t=10.0$ ns. To deal with the entanglement
growth during the time evolution, we use an adaptive bond dimension
$\chi_{\sigma}$ with a maximal default value of $\chi^{\rm max}$ = 400 that
may be changed according to the requirements of the given problem. The actual
value of $\chi_{\sigma}$ at each $A[\sigma]$ is controlled by the uniform
truncation parameter $\epsilon$, such that discarded singular values of the
SDV decomposition $\sum_{i\in\text{discarded}}\lambda_{i}^{2}<\epsilon$ . Our
setup allows users to chose between three compression levels
$\epsilon\in\\{10^{-8},10^{-10},10^{-12}\\}$.