# Compressed Video Quality Assessment for Super-Resolution: a Benchmark and a
Quality Metric††thanks: This study was supported by Russian Science Foundation
under grant 22-21-00478, https://rscf.ru/en/project/22-21-00478/

Evgeney Bogatyrev1 Ivan Molodetskikh12 Dmitriy Vatolin12  Lomonosov Moscow
State University1
MSU Institute for Artificial Intelligence2
{evgeney.zimin, ivan.molodetskikh<EMAIL_ADDRESS>
###### Abstract

We developed a super-resolution (SR) benchmark to analyze SR’s capacity to
upscale compressed videos. Our dataset employed video codecs based on five
compression standards: H.264, H.265, H.266, AV1, and AVS3. We assessed 17
state-of-the-art SR models using our benchmark and evaluated their ability to
preserve scene context and their susceptibility to compression artifacts. To
get an accurate perceptual ranking of SR models, we conducted a crowd-sourced
side-by-side comparison of their outputs. The benchmark is publicly available
at https://videoprocessing.ai/benchmarks/super-resolution-for-video-
compression.html.

We also analyzed benchmark results and developed an objective-quality-
assessment metric based on the current best-performing objective metrics. Our
metric outperforms others, according to Spearman correlation with subjective
scores for compressed video upscaling. It is publicly available at
https://github.com/EvgeneyBogatyrev/super-resolution-metric.

###### Index Terms:

video processing, super-resolution, dataset, video quality metric

##  I Introduction

Super-resolution (SR) has garnered extensive research in recent years, with
new articles appearing monthly. Some state-of-the-art SR methods can restore
details that are unclear in the original (lower-resolution) clip while working
in real time [30]. Neighboring frames can help fill gaps when upscaling,
because small movements caused by camera tremor may provide enough information
to accurately increase the resolution, as demonstrated using a Google Pixel 3
camera [29].

Video compression that can reduce bandwidth consumption with minimal changes
to visual quality is more critical than ever. Some recent codec [7, 12]
downscale a video before compression to cut the bitrate and then upscale it to
its original resolution using SR. This approach decreases bandwidth
consumption and preserves the video’s perceptual quality. Not all SR methods
are suitable for downsample-based video compression, however, since few real-
time SR models can generate acceptable-quality video.

To analyze which SR models work better with each compression standard, and to
help researchers find the best model for their codec, we present our Super-
Resolution for Video Compression benchmark. To develop it we selected 17 SR
models with different architectures and assessed their compressed-video-
restoration capabilities on our dataset, which includes videos compressed
using five codecs. Our effort employed objective metrics and subjective
evaluation to assess model quality. In addition, analysis of the correlation
between objective-quality metrics and subjective scores aided us in developing
our objective metric for analyzing SR’s capacity to upscale compressed videos.

Figure 1: Example videos from the dataset. The dataset includes real-world
sequences, animation, and clips from games.

Our main contributions are as follows:

  1. 1.

We developed a comprehensive SR benchmark to test the ability of SR models to
upscale and restore videos compressed by video codecs of different standards.
The benchmark is publicly available at
https://videoprocessing.ai/benchmarks/super-resolution-for-video-
compression.html.

  2. 2.

We provided an analysis of 6 video quality metrics by their correlation with
subjective scores on our dataset.

  3. 3.

We developed a new objective quality metric for assessing compressed video
quality that has the highest correlation with subjective scores among other
metrics. The metric is publicly available at
https://github.com/EvgeneyBogatyrev/super-resolution-metric.

Figure 2: Distribution of Google SI/TI [25] features for videos we considered
when creating our dataset. Chosen videos appear in orange, others in blue.

##  II Related Work

In this section we provide an overview of existing SR methods and video-
quality metrics.

###  II-A Super-Resolution Methods

SR has received extensive attention since neural networks were first applied
to this area, resulting in many approaches.

Some SR algorithms rely on the temporal redundancy of video frames, allowing
them to restore a single high-resolution frame from a series of low-resolution
ones. RBPN [8] integrates spatial and temporal contexts from a video using a
recurrent encoder-decoder module. COMISR [16] upscales compressed videos; it
employs bidirectional recurrent warping for detail-preserving flow estimation,
and it applies Laplacian enhancement. BasicVSR++ [6] also adopts bidirectional
propagation and spatial alignment. VRT [17] extracts video features, upscales
them, and then reconstructs HQ frames on the basis of these features using a
transformer network.

Generative adversarial networks (GANs) serve widely in deep learning and
especially in SR. ESRGAN [24] modifies the SRGAN architecture by adding
residual-in-residual dense block as well as improved adversarial loss and
perceptual loss. Real-ESRGAN [23] enhances this approach by incorporating
high-order degradation modeling to simulate real-world degradation.

Given the limited number of SR models designed to work with compressed video,
assessing the performance of existing SR models on compressed video remains a
critical task.

###  II-B Super-Resolution Quality Metrics

Peak signal-to-noise ratio (PSNR) uses mean squared error and maximum pixel
values. Structural similarity index measure (SSIM)[26] calculates average,
variance, and covariance pixel values for image windows. Both still commonly
appear in super-resolution papers despite their simple techniques and poor
correlation with subjective scores [35].

Deep-learning approaches to video-quality metrics are growing in popularity.
VMAF [2] is an advanced, objective full-reference video quality assessment
metric developed by Netflix that utilizes SVM-based regression to accurately
and consistently predict the perceived video quality score based on various
video features. LPIPS [32] is a metric that uses deep features of various
neural networks to compare images. It additionally serves as a loss function
for some SR models [3].

The ERQA metric [13] assesses images in terms of detail restoration. To do so,
it detaches object edges and matches them with their counterparts in the
reference image.

##  III Benchmark

In this section we present our Super-Resolution for Video Compression
benchmark, including a description of its method and an overview of our
evaluation results for SR models.

Our benchmark is publicly available at
https://videoprocessing.ai/benchmarks/super-resolution-for-video-
compression.html. We welcome SR researchers to contribute to it by submitting
SR models.

###  III-A Dataset Preparation

To ensure the benchmark dataset is diverse enough to test various aspects of
SR models, we collected videos from multiple sources:

  * •

Vimeo: We gathered 50 FullHD sequences, including both real world and
animation. We split them into scenes using the Scene Change Detector plugin
for
VQMT111https://www.compression.ru/video/quality_measure/video_measurement_tool.html.

  * •

Camera: We shot several videos using a Canon EOS 7D. The settings aimed to
minimize blur and achieve the appropriate brightness — the ISO was 4000 and
the aperture 400. Those settings provided clear ground-truth (GT) videos
without blur or noise. We shot 20 indoor videos and 30 outdoor videos. The
indoor ones consist of synthetically crafted scenes containing objects from
everyday life. Each scene includes either moving objects or parallel camera
motion.

  * •

Games: We recorded 20 clips from various 2D and 3D videogames.

We then obtained the following features for each video: Google SI/TI features
[25], frames per second (FPS), colorfulness [9], and maximal number of
faces222https://github.com/ageitgey/face_recognition throughout the video. On
the basis of these features we separated all videos into nine clusters using
the K-Means clustering and selected one video from each cluster. We refer to
these nine selections as source videos. A preview of them appears in Fig. 1.

###  III-B Super-Resolution Models

To select SR models for evaluation we used SR benchmarks that target two
tasks: detail restoration333https://videoprocessing.ai/benchmarks/video-super-
resolution.html and perceptual-quality
improvement444https://videoprocessing.ai/benchmarks/video-upscalers.html. We
excluded similar SR methods, focusing on similar metric values.

We selected the following 17 methods: BasicVSR++ [6], COMISR [16], DBVSR [34],
EGVSR [5], LGFN [21], RBPN [8], Real-ESRGAN [23], RealSR [11], RSDN [10], SOF-
VSR-BD [22], SOF-VSR-BI [22], SwinIR [18], TMNet [31], VRT [17], ahq-11 [1],
amq-12 [1], and bicubic interpolation.

TABLE I: Average BSQ-rate [33] over subjective scores for each SR+codec pair.
Lower BSQ-rate is better. The best result appears in bold, the second best is
underlined, and the third best is in italics.

SR codec | x264 | x265 | aomenc | vvenc | uavs3e
---|---|---|---|---|---
No SR | 1.000 | 1.000 | 1.000 | 1.000 | _1.000_
RealSR[11] | 0.196 | 0.502 | _1.513_ | 4.470 | 0.639
ahq-11[1] | 0.271 | 0.724 | 2.806 | _1.665_ | 1.750
SwinIR[18] | _0.304_ | 0.346 | 1.505 | 2.502 | 0.640
Real-ESRGAN[23] | 0.335 | _0.640_ | 2.411 | 4.368 | 2.430
COMISR[16] | 0.367 | 0.741 | 2.799 | 0.701 | 1.603
VRT[17] | 1.245 | 3.175 | 4.157 | 4.185 | 3.663
BasicVSR++[6] | 1.971 | 3.390 | 4.199 | 4.185 | 4.250
RBPN[8] | 1.979 | 3.434 | 4.226 | 4.246 | 4.470

###  III-C Objective-Quality Estimation

First, we downscaled the source video to $480\times 270$ resolution using
FFmpeg with the flags=bicubic option. We then compressed the low-resolution
video using each of five video codecs with six target bitrates: 0.1, 0.3, 0.6,
1.0, 2.0, and 4.0 Mbps. We chose these bitrates to be relatively low and to
form a logarithmic curve. Our selection of video codecs included x264, x265,
aomenc, vvenc [28], and uavs3e555https://github.com/uavs3/uavs3e. All employed
the medium preset during compression. Compressed videos underwent transcoding
to PNG sequences using FFmpeg and were inputs to an SR model. We applied image
SR models to each frame individually; video SR models received the path to the
directory containing frames in the correct order. We tested a 4× upscale using
our benchmark, but some SR models can only handle 2×. In this latter case, we
applied the model twice.

TABLE II: Comparison of SR+codec pairs by subjective score and objective metrics for “Restaurant.” The best result appears in bold, the second best is underlined, and the third best appears in italics. Rank | SR + codec | Subj. score $\uparrow$ | ERQA $\uparrow$ | LPIPS $\downarrow$ | VMAF $\uparrow$ | PSNR $\uparrow$ | MS-SSIM $\uparrow$
---|---|---|---|---|---|---|---
1 | SwinIR[18] \+ x264 | 5.855 | 0.601 | 0.237 | 73.921 | 24.961 | 0.932
2 | RealSR[11] \+ x264 | 5.838 | 0.565 | 0.268 | 69.906 | 25.449 | 0.933
3 | Real-ESRGAN[23] \+ x264 | 5.142 | 0.560 | 0.238 | 71.277 | 25.083 | 0.934
4 | ahq-11[1] \+ x264 | 5.049 | 0.579 | 0.217 | 74.954 | 26.209 | 0.939
5 | COMISR[16] \+ x264 | 4.966 | 0.550 | 0.256 | 42.892 | 24.417 | 0.909
6 | SwinIR[18] \+ x265 | 4.801 | 0.585 | 0.231 | 75.784 | 25.034 | 0.936
7 | RealSR[11] \+ x265 | 4.738 | 0.584 | 0.260 | 71.991 | 25.519 | 0.936
8 | Real-ESRGAN[23] \+ x265 | 4.312 | 0.576 | 0.232 | 73.013 | 25.113 | 0.937
9 | SwinIR[18] \+ uavs3e | 4.206 | 0.597 | 0.228 | 80.147 | 24.954 | 0.933
10 | SwinIR[18] \+ aomenc | 3.843 | 0.598 | 0.198 | 89.183 | 25.245 | 0.952

After super-resolving videos, we calculated the following objective-video-
quality metrics on the results: PSNR, MS-SSIM[26], VMAF[2], LPIPS [32], and
ERQA [13]. We ranked the SR methods by calculating BSQ-rate (bitrate-for-the-
same-quality rate) [33] for each SR+codec pair relative to base codec
performance, where the base codec is the one we used to compress low-
resolution video.

###  III-D Subjective Comparison

To subjectively rank SR models, we conducted a crowdsourced comparison through
Subjectify.us service. Because detail loss and compression artifacts can be
difficult to notice in a full frame, the subjective evaluation employed crops.
First, we generated saliency maps for each source video using a method
proposed in [14]. Second, we averaged the saliency maps over all frames and
applied a Gaussian-blur kernel to the result in order to determine the video’s
most salient region. Third, we took distorted videos from the benchmark
dataset, which were compressed to three target bitrates (0.6, 1.0, and 2.0
Mbps), and cut one $480\times 270$ crop from each one, with the most salient
area at the center of the crop. We chose these three bitrates because they
represent the logarithmic structure of an RD curve. Our crop resolution is
guaranteed to fit the screens of subjective-comparison participants. We
evaluated objective metrics on these crops to determine the correlation with
the subjective scores.

We split our comparison into five sections by codec, using only the 10 best SR
models as determined by the LPIPS value.

During the experiment we showed each participant a pair of videos from two
random SR models and asked them to choose the video that looks more realistic
and has fewer compression artifacts (“indistinguishable” is also an option).
Every video pair was viewed by 15 participants. Each participant compared 25
pairs total.

Among the 25 questions were three verification ones, which had obvious
predefined answers. We excluded the results from any participant who failed to
correctly answer one or more of the verification questions.

A total of 5,662 people participated in our subjective evaluation. We excluded
the results from 265 of them because they failed to correctly answer
verification questions. Our calculation of the final subjective scores, using
the Bradley-Terry model [4], employed the remaining 120,316 responses.

TABLE III: Mean Pearson (PLCC) and Spearman (SRCC) correlation coefficients between metrics and subjective-comparison results. Metric | PLCC | SRCC
---|---|---
MS-SSIM [27] | 0.146 | 0.151
PSNR | 0.187 | 0.285
VMAF [2] | 0.344 | 0.448
LPIPS [32] | 0.414 | 0.431
ERQA [13] | 0.582 | 0.624
MDTVSFA [15] | 0.634 | 0.644
ERQAxMDTVSFA | 0.770 | 0.801
Proposed | 0.821 | 0.838

###  III-E Benchmark Results

The results of our subjective evaluation of each SR+codec pair appears in Tab.
I. As the table shows, the best SR model differs by codec, proving that no
single SR model can handle distortion from all compression standards with
equal effectiveness.

The ”No SR” method applies the video codec to the source video without
downscaling or super-resolving. This method serves as a reference during the
BSQ-rate calculation. We see that it exhibits the best results for the aomenc
codec, because AV1-based codecs have a mode that encodes frames at low
resolution and applies an upsampler when decoding. This mode is normally used
at low bitrate, so it is nearly identical to our evaluation method.

Tab. II lists the ranking of SR+codec pairs for ”Restuarant” video. It shows
the top 10 on the basis of subjective scores. The Pearson (PLCC) and Spearman
(SRCC) correlations of objective metrics with subjective score appear in Tab.
III.

##  IV Quality Metric

Judging by the subjective benchmark comparison (Tab. III), existing video-
quality metrics correlate poorly with subjective scores and are therefore
unsuitable for assessing the results of downscale-based video coding.

###  IV-A Dataset

We expanded our benchmark dataset so we could train the metric. When selecting
videos we used the same method described in Sec. III-A, but we collected 20
new ones. Fig. 2 shows the distribution of SI/TI features in the dataset.

We selected four SR methods (RealSR, Real-ESRGAN, COMISR, and BasicVSR++) as
well as three video codecs (x264, x265, and aomenc) to generate distorted
videos. We then downscaled the source videos and compressed them using each
codec at two target bitrates: 0.5 Mbps and 2.0 Mbps. The next step was to
apply each SR method to each video, thereby producing the dataset of distorted
videos. We also compressed the source videos without downscaling and added the
results to our dataset.

We conducted a subjective comparison like the one in Sec. III-D. Because
subjective scores differ in range from one video to the next, we rescaled them
to range from 0 to 1 in each video.

###  IV-B Proposed Approach

Although individual metrics perform poorly on this task, MDTVSFA, ERQA, and
LPIPS have the highest correlation among them. We decided to combine these
metrics to boost the correlation. We calculated MDTVSFA, ERQA, and LPIPS
values over our distorted dataset and then added ERQAxMDTVSFA and ERQAxLPIPS
features, which result from multiplying the ERQA value by the MDTVSFA and
LPIPS values, respectively. We tried other metrics such as NIQE [36] and
SSIM[26] in combination with ERQA, but they showed worse results. To
generalize the information about the video, we also calculated other features:
Google SI/TI [25], colorfulness [9], and bitrate.

Our approach used min-max normalization to normalize the data before training.
We then divided all videos in the distorted dataset using threefold cross-
validation. To train the metric we selected the SVR model with a linear kernel
from Scikit-Learn 1.2.1. Our choice of this model was because it is easy to
train and has produced better results than other regression models.

To optimize the metric’s run time we reduced the number of features. We
iterated through each pair of features and trained the metric without them.
Our approach used the same three splits, and for each pair it involved
calculating the worst SRCC value for those splits. The metric showed the
highest correlation without MDTVSFA and bitrate features. Removing other
features caused much higher correlation losses.

The metric’s final version uses the features ERQA, LPIPS, ERQAxLPIPS,
ERQAxMDTVSFA, SI, TI, and colorfulness.

TABLE IV: Comparison of VQA-algorithm performance by SROCC. The best result
appears in bold, the second best is underlined, and the third best appears in
italics. Most of the data is from [19].

Prediction Model | Wireless | H.264 | MPEG-2 | All Data
---|---|---|---|---
PSNR | 0.4334 | 0.4296 | 0.3588 | 0.3684
SSIM | 0.5233 | 0.6514 | 0.5545 | 0.5257
Speed SSIM | 0.5630 | 0.7086 | 0.6185 | 0.5849
VSNR | 0.7019 | 0.6460 | 0.5915 | 0.6755
V-VIF | 0.5507 | 0.6807 | 0.6116 | 0.5710
Spatial MOVIE | 0.7928 | 0.7066 | 0.6911 | 0.7270
Temporal MOVIE | 0.8114 | 0.7797 | 0.8170 | 0.8055
MOVIE | 0.8109 | 0.7664 | 0.7733 | 0.7890
Proposed | 0.7934 | 0.8465 | 0.7554 | 0.7946

###  IV-C Experiments

We tested the proposed metric on our benchmark, where it outperformed the
other metrics in both Pearson and Spearman correlation, as Tab. III shows.

We also evaluated the metric on the Live Video Quality Database [19, 20]. Tab.
IV shows the Spearman correlation of our metric with DMOS values for each of
the four test cases, as well as the metric’s overall Spearman correlation.

The results demonstrate that our metric exhibits the highest correlation on
H.264 compression and the second-highest correlation overall. Correlation in
the case of MPEG-2 compression is lower, probably because the training set
lacked any videos encoded using that standard.

##  V Conclusion

In this paper we proposed a new benchmark for super-resolution (SR)
compression restoration. Our work assessed 17 SR models applied to five codecs
using both objective-quality metrics and subjective evaluation. We proposed a
new quality metric that has better correlation with subjective scores than do
other metrics on our benchmark’s dataset, and its results more closely
resemble ranking based on human judgment.

Our research shows that SR models, such as RealSR [11] and SwinIR [18], can
serve in downscale-based codecs on the decoder side to enhance the
subjectively perceived quality of videos with low bitrates.

## References

  * [1] Topaz Video Enhance AI.  https://www.topazlabs.com/topaz-video-ai.  Online; accessed 2022-11-08.
  * [2] VMAF - Video Multi-Method Assessment Fusion.  https://github.com/Netflix/vmaf.  Online; accessed 2022-08-06.
  * [3] Andreas Aakerberg, Kamal Nasrollahi, and Thomas B Moeslund.  Real-world super-resolution of face-images from surveillance cameras.  IET Image Processing, 16(2):442–452, 2022.
  * [4] Ralph Allan Bradley and Milton E Terry.  Rank analysis of incomplete block designs: I. The method of paired comparisons.  Biometrika, 39(3/4):324–345, 1952.
  * [5] Yanpeng Cao, Chengcheng Wang, Changjun Song, Yongming Tang, and He Li.  Real-time super-resolution system of 4k-video based on deep learning.  In 2021 IEEE 32nd International Conference on Application-specific Systems, Architectures and Processors (ASAP), pages 69–76. IEEE, 2021.
  * [6] Kelvin CK Chan, Shangchen Zhou, Xiangyu Xu, and Chen Change Loy.  BasicVSR++: Improving video super-resolution with enhanced propagation and alignment.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 5972–5981, 2022.
  * [7] Ying Chen, Qing Li, Aoyang Zhang, Longhao Zou, Yong Jiang, Zhimin Xu, Junlin Li, and Zhenhui Yuan.  Higher quality live streaming under lower uplink bandwidth: an approach of super-resolution based video coding.  In Proceedings of the 31st ACM Workshop on Network and Operating Systems Support for Digital Audio and Video, pages 74–81, 2021.
  * [8] Muhammad Haris, Gregory Shakhnarovich, and Norimichi Ukita.  Recurrent back-projection network for video super-resolution.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 3897–3906, 2019.
  * [9] David Hasler and Sabine E Suesstrunk.  Measuring colorfulness in natural images.  In Human vision and electronic imaging VIII, volume 5007, pages 87–95. SPIE, 2003.
  * [10] Takashi Isobe, Xu Jia, Shuhang Gu, Songjiang Li, Shengjin Wang, and Qi Tian.  Video super-resolution with recurrent structure-detail network.  In European conference on computer vision, pages 645–660. Springer, 2020.
  * [11] Xiaozhong Ji, Yun Cao, Ying Tai, Chengjie Wang, Jilin Li, and Feiyue Huang.  Real-world super-resolution via kernel estimation and noise injection.  In proceedings of the IEEE/CVF conference on computer vision and pattern recognition workshops, pages 466–467, 2020.
  * [12] Mehrdad Khani, Vibhaalakshmi Sivaraman, and Mohammad Alizadeh.  Efficient Video Compression via Content-Adaptive Super-Resolution.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 4521–4530, 2021.
  * [13] Anastasia Kirillova, Eugene Lyapustin, Anastasia Antsiferova, and Dmitry Vatolin.  ERQA: Edge-restoration quality assessment for video super-resolution.  arXiv preprint arXiv:2110.09992, 2021.
  * [14] Alexander Kroner, Mario Senden, Kurt Driessens, and Rainer Goebel.  Contextual encoder–decoder network for visual saliency prediction.  Neural Networks, 129:261–270, 2020.
  * [15] Dingquan Li, Tingting Jiang, and Ming Jiang.  Unified quality assessment of in-the-wild videos with mixed datasets training.  International Journal of Computer Vision, 129(4):1238–1257, 2021\.
  * [16] Yinxiao Li, Pengchong Jin, Feng Yang, Ce Liu, Ming-Hsuan Yang, and Peyman Milanfar.  Comisr: Compression-informed video super-resolution.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 2543–2552, 2021.
  * [17] Jingyun Liang, Jiezhang Cao, Yuchen Fan, Kai Zhang, Rakesh Ranjan, Yawei Li, Radu Timofte, and Luc Van Gool.  Vrt: A video restoration transformer.  arXiv preprint arXiv:2201.12288, 2022.
  * [18] Jingyun Liang, Jiezhang Cao, Guolei Sun, Kai Zhang, Luc Van Gool, and Radu Timofte.  SwinIR: Image restoration using swin transformer.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 1833–1844, 2021.
  * [19] Kalpana Seshadrinathan, Rajiv Soundararajan, Alan C Bovik, and Lawrence K Cormack.  A subjective study to evaluate video quality assessment algorithms.  In Human Vision and Electronic Imaging XV, volume 7527, pages 128–137. SPIE, 2010.
  * [20] Kalpana Seshadrinathan, Rajiv Soundararajan, Alan Conrad Bovik, and Lawrence K Cormack.  Study of subjective and objective quality assessment of video.  IEEE transactions on Image Processing, 19(6):1427–1441, 2010.
  * [21] Dewei Su, Hua Wang, Longcun Jin, Xianfang Sun, and Xinyi Peng.  Local-global fusion network for video super-resolution.  IEEE Access, 8:172443–172456, 2020.
  * [22] Longguang Wang, Yulan Guo, Zaiping Lin, Xinpu Deng, and Wei An.  Learning for video super-resolution through HR optical flow estimation.  In Asian Conference on Computer Vision, pages 514–529. Springer, 2018.
  * [23] Xintao Wang, Liangbin Xie, Chao Dong, and Ying Shan.  Real-ESRGAN: Training Real-World Blind Super-Resolution with Pure Synthetic Data.  In International Conference on Computer Vision Workshops (ICCVW), 2021.
  * [24] Xintao Wang, Ke Yu, Shixiang Wu, Jinjin Gu, Yihao Liu, Chao Dong, Yu Qiao, and Chen Change Loy.  ESRGAN: Enhanced super-resolution generative adversarial networks.  In The European Conference on Computer Vision Workshops (ECCVW), September 2018.
  * [25] Yilin Wang, Sasi Inguva, and Balu Adsumilli.  YouTube UGC dataset for video compression research.  In 2019 IEEE 21st International Workshop on Multimedia Signal Processing (MMSP), pages 1–5. IEEE, 2019.
  * [26] Zhou Wang, Alan C Bovik, Hamid R Sheikh, and Eero P Simoncelli.  Image quality assessment: from error visibility to structural similarity.  IEEE transactions on image processing, 13(4):600–612, 2004.
  * [27] Zhou Wang, Eero P Simoncelli, and Alan C Bovik.  Multiscale structural similarity for image quality assessment.  In The Thrity-Seventh Asilomar Conference on Signals, Systems & Computers, 2003, volume 2, pages 1398–1402. Ieee, 2003.
  * [28] Adam Wieckowski, Jens Brandenburg, Tobias Hinz, Christian Bartnik, Valeri George, Gabriel Hege, Christian Helmrich, Anastasia Henkel, Christian Lehmann, Christian Stoffers, Ivan Zupancic, Benjamin Bross, and Detlev Marpe.  VVenC: An Open And Optimized VVC Encoder Implementation.  In Proc. IEEE International Conference on Multimedia Expo Workshops (ICMEW), pages 1–2.
  * [29] Bartlomiej Wronski, Ignacio Garcia-Dorado, Manfred Ernst, Damien Kelly, Michael Krainin, Chia-Kai Liang, Marc Levoy, and Peyman Milanfar.  Handheld Multi-Frame Super-Resolution.  ACM Trans. Graph., 38(4), July 2019.
  * [30] Yanze Wu, Xintao Wang, Gen Li, and Ying Shan.  AnimeSR: Learning Real-World Super-Resolution Models for Animation Videos.  arXiv preprint arXiv:2206.07038, 2022.
  * [31] Gang Xu, Jun Xu, Zhen Li, Liang Wang, Xing Sun, and Ming-Ming Cheng.  Temporal modulation network for controllable space-time video super-resolution.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 6388–6397, 2021.
  * [32] Richard Zhang, Phillip Isola, Alexei A Efros, Eli Shechtman, and Oliver Wang.  The unreasonable effectiveness of deep features as a perceptual metric.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 586–595, 2018.
  * [33] Anastasia V Zvezdakova, Dmitriy L Kulikov, Sergey V Zvezdakov, and Dmitriy S Vatolin.  BSQ-rate: a new approach for video-codec performance comparison and drawbacks of current solutions.  Programming and computer software, 46(3):183–194, 2020.
  * [34] Jinshan Pan, Haoran Bai, Jiangxin Dong, Jiawei Zhang, and Jinhui Tang.  Deep blind video super-resolution.  In Proceedings of the IEEE/CVF International Conference on Computer Vision, pages 4811–4820, 2021.
  * [35] Quan Huynh-Thu and Mohammed Ghanbari.  The accuracy of PSNR in predicting video quality for different video scenes and frame rates.  Telecommunication Systems, 49(1):35–48, 2012.
  * [36] Anish Mittal, Rajiv Soundararajan, and Alan C. Bovik.  Making a “Completely Blind” Image Quality Analyzer.  IEEE Signal Processing Letters, 20(3):209–212, 2013.


# Dynamical Photon Condensation into Wannier-Stark States

Arkadiusz Kosior  0000-0002-5039-1789 Karol Gietka  0000-0001-7700-3208 Farokh
Mivehvar  0000-0003-4776-1352 Helmut Ritsch  0000-0001-7013-5208 Institut für
Theoretische Physik, Universität Innsbruck, A-6020 Innsbruck, Austria

###### Abstract

Strongly coupled light-matter systems can exhibit nonequilibrium collective
phenomena due to loss and gain processes on the one hand and effective photon-
photon interactions on the other hand. Here we study a photonic lattice system
composed of a linear array of driven-dissipative coupled cavities (or cavity
modes) with linearly increasing resonance frequencies across the lattice. The
model amounts to a driven-dissipative Bose-Hubbard model in a tilted potential
_without_ the particle-conservation constraint. We predict a diverse range of
stationary and non-stationary states resulted from the interplay of the tilt,
tunneling, on-site interactions, and the loss and gain processes. Our key
finding is that, under weak on-site interactions, photons mostly Bose condense
into a selected, single-particle Wannier-Stark state, instead of exhibiting
expected Bloch oscillations. As the strength of the photon-photon interactions
increase, a non-stationary regime emerges which is marked surprisingly by
periodic Bloch-type oscillations. These intriguing, nontrivial effects are a
direct consequence of the driven-dissipative nature of the system.

_Introduction._ —In classical physics, matter and light appear as distinct
entities. Matter is composed of elementary particles, while light is described
as an electromagnetic wave. In quantum mechanics this sharp division between
matter and light becomes less stringent owing to the particle-wave duality
[1]. On the one hand, matter particles are assigned a de Broglie wavelength
and can exhibit wave-like phenomena such as diffraction and interference. On
the other hand, light is considered to be composed of discrete
quanta—photons—which can exhibit particle-like behavior on detection. While
direct photon-photon interaction is typically negligible, when sufficiently
enhanced it can lead to intriguing collective phenomena in an analogy to
correlated matter.

Nevertheless there is a major difference between massive many-particle and
many-photon systems at low energy scales: while the matter-particle number is
strictly conserved, photons can appear and disappear due to absorption,
spontaneous or stimulated emission, and external photon sources. This causes
open photonic systems to be inherently out of equilibrium [2] and even
stationary sates are typically not determined simply by temperature and
entropy but rather by the dynamical balance of gain and loss. Intriguing
nonequilibrium phenomena thus appear in composite light-matter systems [3, 4,
5, 6] and in particular in quantum fluids of light [7]. The most notable
example is the observation of the quasi-equilibrium Bose-Einstein condensate
(BEC) of exciton polaritons—bosonic quasiparticles composed of a mixture of an
exciton (an electron-hole pair) and a cavity photon—in a semiconductor
microcavity [8, 9, 10, 11] and the BEC of photons interacting via molecules in
a multimode optical microcavity [12]. Despite the driven-dissipative nature of
these systems, they still exhibit an effective thermalization process to which
one can attribute an effective temperature. This stands in a sharp contrast to
a typical laser operation, where the thermalization is completely ineffective
and the photon gas is far out-of-equilibrium.

Substantial interest has been focused on the investigation of strongly
correlated many-body effects with photons [13, 14]. In earlier investigations,
the focus was on establishing connections between driven-dissipative steady
states and equilibrium many-body phases. These included the prediction of a
phase transition from a superfluid to a Mott-insulator state for photons via
the photon-blockade effect in coupled cavities [15, 16, 17, 18, 19, 20, 21].
More recently the focus was shifted towards the intriguing realm of the
driven-dissipative regime, where nonequilibrium steady-state phases exhibit
distinct properties in comparison to thermally equilibrium cases [22, 23, 24].
For example, the boundary between monostable and bistable phases in a driven-
dissipative model resembles characteristic Mott insulator lobes but the mean
photon density is not constant within these regions [22]. In a wider context
various ideas and schemes have also been put forward to simulate geometric
phases and gauge potentials for photons, opening the possibility for realizing
nonequilibrium topological photonic states [25, 26].

In this Letter, we investigate a driven-dissipative array of coupled Kerr
cavities (or cavity modes) with linearly increasing resonant frequencies
mimicking a 1D bosonic optical lattice subjected to an external force. In our
generic model photons (bosons) are continuously injected into the system via a
coherent pump to a selected cavity. Photons are then redistributed via nearest
neighbor mode couplings until they eventually leak out through imperfect
cavity mirrors. The dynamics of the system can be effectively captured by a
driven-dissipative Bose-Hubbard model in a tilted potential _without_ particle
conservation. Despite the conceptual simplicity, the model, as we demonstrate,
exhibits a variety of intriguing stationary and non-stationary, nonequilibrium
phenomena controllable by the system parameters.

A central result of our study is that for sufficiently weak on-site
interactions, photons dynamically condense due to the explicit $U(1)$ symmetry
breaking into a selected, spatially-localized state (Fig. 1), instead of
exhibiting expected Bloch oscillations [27, 28]. Specifically, we show that
the condensate wavefunction is often close to a single Wannier-Stark (WS)
state with only small contributions from a few neighboring WS states; see Fig.
2. Intriguingly, by increasing the strength of local on-site interactions the
system enters a non-stationary regime (Fig. 3) where the photon density
undergoes periodic Bloch-type oscillations over time as shown in Fig. 4 [29,
30, 31]. This finding strongly contrasts with the irreversible decay of Bloch
oscillations of interacting atoms in a 1D tilted lattice [32, 33, 34, 35, 36].
Let us emphasize here that both of these regimes appear independent of initial
conditions and the choice of the pumped cavity, and are a direct consequence
of the driven-dissipative nature of the system. As our model is quite generic,
a corresponding experimental realization should be implementable in various
platforms including superconducting circuits [37, 38, 39], photonic crystal
structures [40], waveguide-coupled optical cavities [41], and atom-filled
transverse multi-mode cavities [42].

_Model and its Hamiltonian._ –Consider a linear array of standing-wave coupled
cavities (labeled by $j\in\mathbb{Z}$) with linearly increasing resonant
frequencies $\omega_{j}\propto j$, each containing a Kerr-like non-linear
medium. Coherent pumps with the frequency $\omega_{p}$ continuously inject
photons into the cavity modes. Each cavity is coupled to two adjacent cavities
due to the photon leakage through the cavity mirrors, which leads to a
coherent hopping of photons in the cavity lattice. The Hamiltonian of the
system in the rotating frame of the coherent pumps is given by [43],

$\displaystyle\hat{H}$
$\displaystyle=\sum_{j}\left[\hbar\Delta_{j}\hat{a}_{j}^{\dagger}\hat{a}_{j}-J\big{(}\hat{a}_{j}^{\dagger}\hat{a}_{j+1}+\text{H.c.}\big{)}+\chi\hat{a}_{j}^{\dagger
2}\hat{a}_{j}^{2}\right]$
$\displaystyle+\hbar\sum_{j}\eta_{j}\left(\hat{a}_{j}+\hat{a}_{j}^{\dagger}\right),$
(1)

with $\hat{a}_{j}$, $\hat{a}_{j}^{\dagger}$ being bosonic operators
annihilating and creating a photon in the $j$-th cavity, respectively. Here we
have defined

$\Delta_{j}=\omega_{j}-\omega_{p}\equiv\Delta\omega\left(j-j_{0}\right),$ (2)

as the cavity-pump detuning. Moreover, $J$ is the nearest-neighbor photon
tunneling-amplitude rate, $\chi$ the on-site photon-photon interaction
strength due to the effective Kerr non-linearity, and $\eta_{j}$ the pumping
rate of the $j$-th cavity. Photon losses, assumed to be the same for all the
cavities, $\kappa$ are taken into account via the quantum Heisenberg equations
of motion,

$\frac{d\hat{a}_{j}}{dt}=\frac{i}{\hbar}[\hat{H},\hat{a}_{j}]-\kappa\hat{a}_{j}.$
(3)

The first line of the Hamiltonian (Dynamical Photon Condensation into Wannier-
Stark States) describes the familiar equilibrium Bose-Hubbard model in a
tilted lattice. While in the non-interacting limit the equilibrium model
features well-known Bloch oscillations [27, 28], it has been recently shown
that strong interactions can lead to disorder-free many-body localization [44,
45] (for related experiments see Ref. [46]). The second line of the
Hamiltonian (Dynamical Photon Condensation into Wannier-Stark States)
introduces a coherent photon pumping, which along with the photon decay into
the environment $\kappa$ [see Eq. (3)] explicitly breaks the $U(1)$ symmetry
of the system associated with the particle number conservation [43]. In the
following we will show that this lack of the particle conservation and
explicitly broken $U(1)$ symmetry due to the loss and gain processes has
fundamental consequences in both static and dynamics of the system.

Figure 1: Nonequilibrium dynamics of the system for a weak photon-photon
interaction strength $\chi=10^{-2}$ and the lattice tilt $\Delta\omega=0.5$.
The system reaches a spatially localized steady state in long time. The
expectation values of the cavity particle number operators
$\langle\hat{n}_{j}(t)\rangle$ in the course of time evolution in the first-
order cumulant expansion (i.e., mean field) for (a) $j_{0}=0$ and (b)
$j_{0}=5$. (c) The time evolution of the expectation value of the total photon
number operator $\langle\hat{N}\rangle$ for $j_{0}=0$ in both first- and
second-order cumulant expansion. (d) The distribution of
$\langle\hat{n}_{j}(t_{f})\rangle$ over lattice sites in the stationary state
for $j_{0}=0$ (cf. panel a). To a very good approximation, the distribution is
proportional to the probability density of a single WS state (see also Fig.
2), signalling a nonequilibrium photon Bose-condensation into a WS state. Mean
field is quite accurate in this weakly interacting regime.

_Consequences of explicit U(1) symmetry breaking: State selection_.—In order
to gain some physical intuition, let us start with the the non-interacting
limit, $\chi=0$. In the WS basis, the Heisenberg equations of motions read
[43],

$i\frac{d\hat{b}_{n}}{dt}=\left(\Delta_{n}-i\kappa\right)\hat{b}_{n}+\tilde{\eta}_{n},$
(4)

where $\hat{b}_{n}=\sum_{j}\beta_{n,j}\hat{a}_{n}$,
$\tilde{\eta}_{n}=\sum_{j}\beta_{n,j}\eta_{j}$, and
$\beta_{n,j}=\mathcal{J}_{j-n}\left({2J}/{\Delta\omega}\right)$, with
$\mathcal{J}_{k}$ being the Bessel function of the first kind of order $k$.
The equations of motion (4) readily yield

$\hat{b}_{n}(t)=e^{-it(\Delta_{n}-i\kappa)}\,\hat{b}_{n}(0)+\frac{e^{-it(\Delta_{n}-i\kappa)}-1}{\Delta_{n}-i\kappa}\,\tilde{\eta}_{n}.$
(5)

If $\forall_{n}\tilde{\eta}_{n}=0$, the solutions correspond to damped Bloch
oscillations [43]. If additionally $\kappa=0$, one then recovers the common
Bloch oscillations and the choice of $j_{0}$ in Eq. (2) becomes arbitrary as
it only adds an irrelevant phase factor.

In the long-time limit $t\gg\kappa^{-1}\gg J^{-1}$, regardless of the choice
of initial conditions $\hat{b}_{n}(0)$, the Bloch oscillations are completely
damped out and the system reaches a steady state with the WS mode occupations,

$\langle\hat{n}_{n}\rangle=\frac{\tilde{\eta}_{n}^{2}}{(\Delta\omega)^{2}(n-j_{0})^{2}+\kappa^{2}}.$
(6)

As can be seen from Eq. (6), now $j_{0}$ plays an essential role. In
particular, by properly choosing $j_{0}\in\mathbb{Z}$ one can select
dynamically a single WS mode ($n=j_{0}$) to be microscopically occupied.
However, we note this approach does not work when $\tilde{\eta}_{j_{0}}=0$,
which happens for some specific ratios of $J/\Delta\omega$. In these cases,
one encounters a sequence of pumping anti-resonances, which lead to the
occupation a few adjacent WS states instead. Furthermore, these anti-
resonances are responsible for non-trivial phase boundaries between stationary
and non-stationary states as shown in Fig. 3, which we will delve into it in
the last section.

The above simple analysis is valid qualitatively also for sufficiently small
photon-photon interactions, which we confirm numerically in the following. In
numerical simulations we set $J=\hbar=1$ (as the unit of energy),
$\kappa=10^{-2}$, $\forall_{j}\hat{a}_{j}(0)=0$ , and calculate the
expectation values of the relevant operators using the cumulant expansion in
both first (i.e., mean-field) and second order [47, 43]. Moreover, for the
sake of simplicity and without loss of generality, we consider only a single
cavity pumping $\eta_{j}=\eta\delta_{j,0}$, but different choices do not
change our main conclusions [43]. Finally, with the exception of Fig. 1, we
also fix $j_{0}=0$.

_Nonequilibrium photon condensation in the weakly interacting regime._ —Now,
we turn our attention to the weakly interacting regime, $\chi\ll 1$. As in the
non-interacting case, in the weakly interacting regime the many-photon system
still occupies macroscopically a single or a few one-particle WS states with a
high degree of coherence. In Fig. 1, we show the expectation value of cavity
photon-number operators $\hat{n}_{j}=\hat{a}_{j}^{\dagger}\hat{a}_{j}$ as well
as the total photon-number operator $\hat{N}=\sum_{j}\hat{n}_{j}$. As can be
seen, in each case the system reaches a spatially-localised stationary steady
state which, for the chosen parameters, is proportional to a single WS state
with $n=j_{0}$.

To further quantify this observation, we calculate the fidelity
$P_{n}(t)=|\langle\Psi_{n}|\psi(t)\rangle|^{2}$ between the mean-field
wavefunction $|\psi(t)\rangle$ and the WS basis states $|\Psi_{n}\rangle$. The
distribution of $P_{n}$ is illustrated in Fig. 2(a) in steady states, showing
it is centered around $n=j_{0}=0$ [cf. Eq. (6)]. In general two distinct
scenarios are possible: (i) The mean-field wavefunction consists predominantly
of a single WS state, or (ii) it has contributions from a few different WS
states. As Fig. 2(a) shows, for $\Delta\omega=0.5$ [as in Fig. 1(a) and (d)]
the wavefunction is close to the central WS state, while for
$\Delta\omega=0.35$ two additional modes $n=\pm 1$ are also significantly
occupied at the macroscopic level. This is because of the aforementioned
pumping anti-resonances, where the population of $n=\pm 1$ ($n=0$) mode is
completely suppressed at $\Delta\omega\approx 0.522$ ($\Delta\omega\approx
0.362$) due to the hitting a zero of the $\mathcal{J}_{1}$ ($\mathcal{J}_{0}$)
Bessel function.

Figure 2:  (a) Steady-sate fidelity $P_{n}(t_{f})$ between the mean-field
wavefunction $|\psi(t_{f})\rangle$ and the WS basis states $|\Psi_{n}\rangle$.
To a very good approximation, the condensate wavefunction is either
proportional to only one WS state or is a superposition of a few WS states as
in an anti-resonant case (see the discussion in the main text). (b) The
dominant eigenvalue $N_{0}$ of the single-particle density matrix as a
function of $\Delta\omega$ remains close to the total number of photons $N$,
indicating a high condensate fraction.  Figure 3: Mean-field nonequilibrium
many-body phase diagram of the system in the parameter plane of the force
$\Delta\omega$ vs. the photon-photon interaction strength $\chi$. (a) Relative
change of the total photon number over a long-time evolution $\Delta n$ [Eq.
(7)] reveals three regimes: stationary steady-state phase (deep blue),
dynamically unstable chaotic regime (light blue), and non-stationary regular
oscillatory states (warm colors). (b) Time-averaged maximal fidelity between
the mean-field wavefunction $|\psi\rangle$ and the WS basis states
$|\Psi_{n}\rangle$ is complementary to panel (a) and reveals particularly a
series of narrow bands that can be explained on a single-particle level as
pumping anti-resonances; see the discussion in the main text. The dashed
(dotted) lines correspond to zeros of the Bessel function $\mathcal{J}_{0}$
($\mathcal{J}_{1}$). Note that interactions slightly shift the positions of
the anti-resonances.

In order to confirm that photons are Bose condensed, we consider the single-
particle reduced density matrix
$\hat{\rho}_{1}(x,x^{\prime})=\langle\hat{\psi}^{\dagger}(x)\hat{\psi}(x^{\prime})\rangle$,
whose eigenvalues determine the occupation probabilities of the natural
orbitals [48]. The highest occupation number quantifies the level of coherence
in the system. Although we can expand the field operators in any orthogonal
basis, we choose the Wannier basis $\hat{\psi}(x)=\sum_{j}w_{j}(x)\hat{a}_{j}$
and calculate the eigenvalues of
$\langle\hat{a}^{\dagger}_{j}\hat{a}_{l}\rangle$ in a steady state of the
system in the second-order cumulant expansion [43]. Indeed, as can be seen
from Fig. 2(b) the highest eigenvalue $N_{0}$ is very close to the total
number of photons $N=\langle\hat{N}\rangle$, supporting the interpretation of
a photon BEC. In contrast to the distribution of $P_{n}$ as shown in Fig.
2(a), the highest eigenvalue of the the reduced density matrix is only weakly
affected by the lattice tilt $\Delta\omega$. Hence, we infer that the
condensate wavefunction is either close to a single WS states or is a
superposition of a few WS states.

Figure 4: Mean-field dynamics of the expectation values of the cavity photon
number operators $\langle\hat{n}_{j}(t)\rangle$ for different interaction
strengths: (a) $\chi=0.13$ in the oscillatory regime, and (b) $\chi=0.135$ in
the chaotic regime. Stationary and different non-stationary solutions are
distinguished by monitoring (c) the relative change of the total photon number
over a long-time evolution and (d) the time-averaged maximal fidelity. The
lattice tilt is set to $\Delta\omega=0.5$ for all panels.

_Non-stationary phases in the strongly interacting regime._ —Above we showed
that in the non-interacting and weakly interacting regimes photons can
condensate into a one or a few selected WS states and reach a steady state.
However, with increasing the interaction strength $\chi$, the stationary
states lose their stability and intriguing non-stationary solutions appear.
Figure 3(a) depicts the nonequilibrium many-body phase diagram of the system
in the parameter plane of $\\{\chi,\Delta\omega\\}$ and contains three main
regimes: stationary steady states (deep blue), dynamically unstable chaotic
regime (light blue), and non-stationary regular oscillatory states (warm
colors). The three phases are distinguished by the relative change of the
total photon number over a long-time evolution, defined as

$\Delta
n=\frac{\max_{t}\langle\hat{N}(t)\rangle-\min_{t}\langle\hat{N}(t)\rangle}{\text{avg}_{t}\langle\hat{N}(t)\rangle},$
(7)

with $\max_{t}$, $\min_{t}$, and $\text{avg}_{t}$ denoting, respectively, the
maximal, minimal, and average value of a quantum-averaged observable during a
long-time evolution (where the initial transient dynamics of the system has
been neglected). Although the phase diagram is dominated by the stationary and
dynamically unstable chaotic regimes, regions of regular oscillatory dynamics
appear mostly on boundaries between the two regimes; see also Fig. 3(b)
showing the time-averaged maximal fidelity and Supplemental Material [43].

Interestingly, photons in the regular oscillatory regime tunnel between a few
neighboring cavities performing spatially-confined oscillations reminiscent of
standard Bloch oscillations; see Fig. 4. Unlike standard Bloch oscillations
that tend to decay in interacting systems [32, 33, 34, 35, 36], in our driven-
dissipative model oscillatory solutions are in fact induced by the on-site
photon-photon interactions and persist over long time dynamics. These
solutions can be conceptualized as multimode limit cycles (as shown in
Supplemental Material [43]) within driven-dissipative systems, where the
steady-state solutions of equations of motion exhibit dynamic instability [49,
50, 51, 52, 53, 54, 55, 56, 57].

_Summary, conclusions, and perspectives._ —In summary, we have studied an
array of coupled Kerr cavities with linearly increasing resonant cavity
frequencies, which is equivalent to a $U(1)$-symmetry broken driven-
dissipative Bose-Hubbard model with a tilted potential. The model reveals a
range of both stationary and non-stationary nonequilibrium phenomena. Notably,
photons condense into selected, stationary WS states under sufficiently weak
on-site interactions. As the strength of photon-photon interactions increases,
we observe a transition to a non-stationary dynamical regime marked by
periodic Bloch-like oscillations over time. Unlike the equilibrium
counterpart, these regular oscillations are induced by the photon-photon
interactions and do not decay over time.

Our research sheds light on the intricate behavior of driven-dissipative
coupled cavity systems, emphasizing the role of explicit $U(1)$ symmetry
breaking and interactions in shaping their dynamics. These findings hold
promise for applications such as coherent light storage [58, 59, 60], light
confinement [61, 62, 63], generation of non-classical many-photon states [64,
65, 66], and distributed quantum sensing [67, 68, 69]. While our focus in this
Letter has been on weakly interacting systems, we underscore the importance of
delving deeper into dynamics beyond perturbation regimes. Specifically, it is
interesting to explore the non-stationary states in many-body regimes,
particularly in the context of ergodicity breaking [70, 71] and Stark many-
body localization [72]. Another intriguing scenario is to extend our driven-
dissipative model to topological setting, where exotic nonequilibrium
topological effects are expected to appear [73, 74, 75, 76].

The data presented in this article is available from [77].

###### Acknowledgements.

We acknowledge stimulating discussions with Darrick Chang and Thomas Pohl.
This research was funded in whole or in part by the Austrian Science Fund
(FWF) [grant DOIs: 10.55776/ESP171, 10.55776/M3304 and 10.55776/P35891]. For
open access purposes, the authors have applied a CC BY public copyright
license to any author accepted manuscript version arising from this
submission. F. M. additionally acknowledges the support of the ESQ-Discovery
Grant of the Austrian Academy of Sciences (ÖAW).

## References

  * Bohm [1989] D. Bohm,  _Quantum Theory_ , Dover Books on Physics (Dover Publications, Mineola, NY, 1989).
  * Walls and Milburn [2008] D. F. Walls and G. F. Milburn,  _Quantum Optics_ (Springer Berlin Heidelberg, 2008).
  * Deng  _et al._ [2010] H. Deng, H. Haug, and Y. Yamamoto, Exciton-polariton Bose-Einstein condensation, Rev. Mod. Phys. 82, 1489 (2010).
  * Keeling and Berloff [2011] J. Keeling and N. G. Berloff, Exciton–polariton condensation, Contemp. Phys. 52, 131–151 (2011).
  * Byrnes  _et al._ [2014] T. Byrnes, N. Y. Kim, and Y. Yamamoto, Exciton–polariton condensates, Nat. Phys 10, 803–813 (2014).
  * Mivehvar  _et al._ [2021] F. Mivehvar, F. Piazza, T. Donner, and H. Ritsch, Cavity QED with quantum gases: new paradigms in many-body physics, Adv. Phys. 70, 1 (2021).
  * Carusotto and Ciuti [2013] I. Carusotto and C. Ciuti, Quantum fluids of light, Rev. Mod. Phys. 85, 299 (2013).
  * Deng  _et al._ [2002] H. Deng, G. Weihs, C. Santori, J. Bloch, and Y. Yamamoto, Condensation of semiconductor microcavity exciton polaritons, Science 298, 199–202 (2002).
  * Kasprzak  _et al._ [2006] J. Kasprzak, M. Richard, S. Kundermann, A. Baas, P. Jeambrun, J. M. J. Keeling, F. M. Marchetti, M. H. Szymańska, R. André, J. L. Staehli, V. Savona, P. B. Littlewood, B. Deveaud, and L. S. Dang, Bose–Einstein condensation of exciton polaritons, Nature 443, 409–414 (2006).
  * Balili  _et al._ [2007] R. Balili, V. Hartwell, D. Snoke, L. Pfeiffer, and K. West, Bose-Einstein Condensation of Microcavity Polaritons in a Trap, Science 316, 1007–1010 (2007).
  * Amo  _et al._ [2009] A. Amo, J. Lefrère, S. Pigeon, C. Adrados, C. Ciuti, I. Carusotto, R. Houdré, E. Giacobino, and A. Bramati, Superfluidity of polaritons in semiconductor microcavities, Nat. Phys 5, 805–810 (2009).
  * Klaers  _et al._ [2010] J. Klaers, J. Schmitt, F. Vewinger, and M. Weitz, Bose–Einstein condensation of photons in an optical microcavity, Nature 468, 545–548 (2010).
  * Noh and Angelakis [2016] C. Noh and D. G. Angelakis, Quantum simulations and many-body physics with light, Rep. Prog. Phys. 80, 016401 (2016).
  * Hartmann [2016] M. J. Hartmann, Quantum simulation with interacting photons, J. Opt. 18, 104005 (2016).
  * Hartmann  _et al._ [2006] M. J. Hartmann, F. G. S. L. Brandão, and M. B. Plenio, Strongly interacting polaritons in coupled arrays of cavities, Nat. Phys 2, 849–855 (2006).
  * Greentree  _et al._ [2006] A. D. Greentree, C. Tahan, J. H. Cole, and L. C. L. Hollenberg, Quantum phase transitions of light, Nat. Phys 2, 856–861 (2006).
  * Hartmann  _et al._ [2008] M. Hartmann, F. Brandão, and M. Plenio, Quantum many‐body phenomena in coupled cavity arrays, Laser & Photon. Rev. 2, 527–556 (2008).
  * Schmidt and Blatter [2009] S. Schmidt and G. Blatter, Strong Coupling Theory for the Jaynes-Cummings-Hubbard Model, Phys. Rev. Lett. 103, 086403 (2009).
  * Tomadin  _et al._ [2010] A. Tomadin, V. Giovannetti, R. Fazio, D. Gerace, I. Carusotto, H. E. Türeci, and A. Imamoglu, Signatures of the superfluid-insulator phase transition in laser-driven dissipative nonlinear cavity arrays, Phys. Rev. A 81, 061801 (2010).
  * Wu  _et al._ [2011] C.-W. Wu, M. Gao, Z.-J. Deng, H.-Y. Dai, P.-X. Chen, and C.-Z. Li, Quantum phase transition of light in a one-dimensional photon-hopping-controllable resonator array, Phys. Rev. A 84, 043827 (2011).
  * foo [a] Although these models share many similarities, such as bosons hopping coherently between nearest-neighbor sites with local nonlinearities, we note that quantitative differences may arise on different energy scales.
  * Le Boité  _et al._ [2013] A. Le Boité, G. Orso, and C. Ciuti, Steady-State Phases and Tunneling-Induced Instabilities in the Driven Dissipative Bose-Hubbard Model, Phys. Rev. Lett. 110, 233601 (2013).
  * Le Boité  _et al._ [2014] A. Le Boité, G. Orso, and C. Ciuti, Bose-Hubbard model: Relation between driven-dissipative steady states and equilibrium quantum phases, Phys. Rev. A 90, 063821 (2014).
  * Graß [2019] T. Graß, Excitations and correlations in the driven-dissipative Bose-Hubbard model, Phys. Rev. A 99, 043607 (2019).
  * Lu  _et al._ [2014] L. Lu, J. D. Joannopoulos, and M. Soljačić, Topological photonics, Nat. Photon 8, 821–829 (2014).
  * Ozawa  _et al._ [2019] T. Ozawa, H. M. Price, A. Amo, N. Goldman, M. Hafezi, L. Lu, M. C. Rechtsman, D. Schuster, J. Simon, O. Zilberberg, and I. Carusotto, Topological photonics, Rev. Mod. Phys. 91, 015006 (2019).
  * Wannier [1960] G. H. Wannier, Wave Functions and Effective Hamiltonian for Bloch Electrons in an Electric Field, Phys. Rev. 117, 432 (1960).
  * Hartmann  _et al._ [2004] T. Hartmann, F. Keck, H. J. Korsch, and S. Mossmann, Dynamics of Bloch oscillations, New J. Phys. 6, 2–2 (2004).
  * Peschel  _et al._ [1998] U. Peschel, T. Pertsch, and F. Lederer, Optical Bloch oscillations in waveguide arrays, Opt. Lett. 23, 1701 (1998).
  * Sapienza  _et al._ [2003] R. Sapienza, P. Costantino, D. Wiersma, M. Ghulinyan, C. J. Oton, and L. Pavesi, Optical Analogue of Electronic Bloch Oscillations, Phys. Rev. Lett. 91, 263902 (2003).
  * Yuan and Fan [2016] L. Yuan and S. Fan, Bloch oscillation and unidirectional translation of frequency in a dynamically modulated ring resonator, Optica 3, 1014 (2016).
  * Buchleitner and Kolovsky [2003] A. Buchleitner and A. R. Kolovsky, Interaction-Induced Decoherence of Atomic Bloch Oscillations, Phys. Rev. Lett. 91, 253002 (2003).
  * Gustavsson  _et al._ [2008] M. Gustavsson, E. Haller, M. J. Mark, J. G. Danzl, G. Rojas-Kopeinig, and H.-C. Nägerl, Control of Interaction-Induced Dephasing of Bloch Oscillations, Phys. Rev. Lett. 100, 080404 (2008).
  * Krimer  _et al._ [2009] D. O. Krimer, R. Khomeriki, and S. Flach, Delocalization and spreading in a nonlinear Stark ladder, Phys. Rev. E 80, 036201 (2009).
  * Eckstein and Werner [2011] M. Eckstein and P. Werner, Damping of Bloch Oscillations in the Hubbard Model, Phys. Rev. Lett. 107, 186406 (2011).
  * foo [b] However, note that bound states of a few strongly interacting particles can experience fractional Bloch oscillations at frequencies that are multiples of single-particle Bloch oscillations [78, 79, 80, 81] (for experiments see [82, 83]). Also, the decay of Blosch oscillations and the emergence of the directional motion in the excitation-phonon has been recently reported in [84].
  * Houck  _et al._ [2012] A. A. Houck, H. E. Türeci, and J. Koch, On-chip quantum simulation with superconducting circuits, Nat. Phys 8, 292–299 (2012).
  * Jin  _et al._ [2013] J. Jin, D. Rossini, R. Fazio, M. Leib, and M. J. Hartmann, Photon Solid Phases in Driven Arrays of Nonlinearly Coupled Cavities, Phys. Rev. Lett. 110, 163605 (2013).
  * Jin  _et al._ [2014] J. Jin, D. Rossini, M. Leib, M. J. Hartmann, and R. Fazio, Steady-state phase diagram of a driven QED-cavity array with cross-Kerr nonlinearities, Phys. Rev. A 90, 023827 (2014).
  * Majumdar  _et al._ [2012] A. Majumdar, A. Rundquist, M. Bajcsy, V. D. Dasika, S. R. Bank, and J. Vučković, Design and analysis of photonic crystal coupled cavity arrays for quantum simulation, Phys. Rev. B 86, 195312 (2012).
  * Lepert  _et al._ [2011] G. Lepert, M. Trupke, M. J. Hartmann, M. B. Plenio, and E. A. Hinds, Arrays of waveguide-coupled optical cavities that interact strongly with atoms, New J. Phys. 13, 113002 (2011).
  * Kollár  _et al._ [2015] A. J. Kollár, A. T. Papageorge, K. Baumann, M. A. Armen, and B. L. Lev, An adjustable-length cavity and bose–einstein condensate apparatus for multimode cavity qed, New J. Phys. 17, 043012 (2015).
  * [43] See Supplemental Materials for the derivation of the rotating frame Hamiltonian, discussion of Wannier-Stark states and the U(1) symmetry breaking terms, details of the cumulant expansion, and more examples of interaction-induced Bloch oscillations.
  * Yao and Zakrzewski [2020] R. Yao and J. Zakrzewski, Many-body localization of bosons in an optical lattice: Dynamics in disorder-free potentials, Phys. Rev. B 102, 104203 (2020).
  * Taylor  _et al._ [2020] S. R. Taylor, M. Schulz, F. Pollmann, and R. Moessner, Experimental probes of Stark many-body localization, Phys. Rev. B 102, 054206 (2020).
  * Morong  _et al._ [2021] W. Morong, F. Liu, P. Becker, K. S. Collins, L. Feng, A. Kyprianidis, G. Pagano, T. You, A. V. Gorshkov, and C. Monroe, Observation of Stark many-body localization without disorder, Nature 599, 393–398 (2021).
  * Plankensteiner  _et al._ [2022] D. Plankensteiner, C. Hotter, and H. Ritsch, QuantumCumulants.jl: A Julia framework for generalized mean-field equations in open quantum systems, Quantum 6, 617 (2022).
  * Pethick and Smith [2008] C. J. Pethick and H. Smith,  _Bose–Einstein Condensation in Dilute Gases_ (Cambridge University Press, 2008).
  * Piazza and Ritsch [2015] F. Piazza and H. Ritsch, Self-Ordered Limit Cycles, Chaos, and Phase Slippage with a Superfluid inside an Optical Resonator, Phys. Rev. Lett. 115, 163601 (2015).
  * Chan  _et al._ [2015] C.-K. Chan, T. E. Lee, and S. Gopalakrishnan, Limit-cycle phase in driven-dissipative spin systems, Phys. Rev. A 91, 051601 (2015).
  * Owen  _et al._ [2018] E. T. Owen, J. Jin, D. Rossini, R. Fazio, and M. J. Hartmann, Quantum correlations and limit cycles in the driven-dissipative Heisenberg lattice, New J. Phys. 20, 045004 (2018).
  * Keßler  _et al._ [2019] H. Keßler, J. G. Cosme, M. Hemmerling, L. Mathey, and A. Hemmerich, Emergent limit cycles and time crystal dynamics in an atom-cavity system, Phys. Rev. A 99, 053605 (2019).
  * Colella  _et al._ [2022] E. Colella, A. Kosior, F. Mivehvar, and H. Ritsch, Open Quantum System Simulation of Faraday’s Induction Law via Dynamical Instabilities, Phys. Rev. Lett. 128, 070603 (2022).
  * Gao  _et al._ [2023] P. Gao, Z.-W. Zhou, G.-C. Guo, and X.-W. Luo, Self-organized limit cycles in red-detuned atom-cavity systems, Phys. Rev. A 107, 023311 (2023).
  * Kosior  _et al._ [2023a] A. Kosior, H. Ritsch, and F. Mivehvar, Nonequilibrium phases of ultracold bosons with cavity-induced dynamic gauge fields, SciPost Phys. 15, 046 (2023a).
  * Mivehvar [2024] F. Mivehvar, Conventional and Unconventional Dicke Models: Multistabilities and Nonequilibrium Dynamics, Phys. Rev. Lett. 132, 073602 (2024).
  * Skulte  _et al._ [2024] J. Skulte, P. Kongkhambut, H. Keßler, A. Hemmerich, L. Mathey, and J. G. Cosme, Realizing limit cycles in dissipative bosonic systems (2024), arXiv:2401.05332 [cond-mat.quant-gas] .
  * Phillips  _et al._ [2001] D. F. Phillips, A. Fleischhauer, A. Mair, R. L. Walsworth, and M. D. Lukin, Storage of Light in Atomic Vapor, Phys. Rev. Lett. 86, 783 (2001).
  * Dudin  _et al._ [2013] Y. O. Dudin, L. Li, and A. Kuzmich, Light storage on the time scale of a minute, Phys. Rev. A 87, 031801 (2013).
  * Katz and Firstenberg [2018] O. Katz and O. Firstenberg, Light storage for one second in room-temperature alkali vapor, Nat. Commun 9, 2074 (2018).
  * Cao  _et al._ [2000] H. Cao, J. Y. Xu, D. Z. Zhang, S.-H. Chang, S. T. Ho, E. W. Seelig, X. Liu, and R. P. H. Chang, Spatial confinement of laser light in active random media, Phys. Rev. Lett. 84, 5584 (2000).
  * Notomi [2011] M. Notomi, Strong light confinement with periodicity, Proc. IEEE 99, 1768 (2011).
  * Riboli  _et al._ [2014] F. Riboli, N. Caselli, S. Vignolini, F. Intonti, K. Vynck, P. Barthelemy, A. Gerardino, L. Balet, L. H. Li, A. Fiore, M. Gurioli, and D. S. Wiersma, Engineering of light confinement in strongly scattering disordered media, Nat. Mater 13, 720 (2014).
  * Cooper  _et al._ [2013] M. Cooper, L. J. Wright, C. Söller, and B. J. Smith, Experimental generation of multi-photon Fock states, Opt. Express 21, 5309 (2013).
  * Chang  _et al._ [2014] D. E. Chang, V. Vuletić, and M. D. Lukin, Quantum nonlinear optics—photon by photon, Nat. Photon 8, 685 (2014).
  * Pizzi  _et al._ [2023] A. Pizzi, A. Gorlach, N. Rivera, A. Nunnenkamp, and I. Kaminer, Light emission from strongly driven many-body systems, Nat. Phys 19, 551 (2023).
  * Guo  _et al._ [2020] X. Guo, C. R. Breum, J. Borregaard, S. Izumi, M. V. Larsen, T. Gehring, M. Christandl, J. S. Neergaard-Nielsen, and U. L. Andersen, Distributed quantum sensing in a continuous-variable entangled network, Nature Physics 16, 281 (2020).
  * Zhang and Zhuang [2021] Z. Zhang and Q. Zhuang, Distributed quantum sensing, Quantum Sci. Technol. 6, 043001 (2021).
  * Pelayo  _et al._ [2023] J. C. Pelayo, K. Gietka, and T. Busch, Distributed quantum sensing with optical lattices, Phys. Rev. A 107, 033318 (2023).
  * Luca and Scardicchio [2013] A. D. Luca and A. Scardicchio, Ergodicity breaking in a model showing many-body localization, Europhys Lett. 101, 37003 (2013).
  * Gietka  _et al._ [2019] K. Gietka, J. Chwedeńczuk, T. Wasak, and F. Piazza, Multipartite entanglement dynamics in a regular-to-ergodic transition: Quantum fisher information approach, Phys. Rev. B 99, 064303 (2019).
  * Schulz  _et al._ [2019] M. Schulz, C. A. Hooley, R. Moessner, and F. Pollmann, Stark many-body localization, Phys. Rev. Lett. 122, 040606 (2019).
  * Hafezi  _et al._ [2013] M. Hafezi, S. Mittal, J. Fan, A. Migdall, and J. M. Taylor, Imaging topological edge states in silicon photonics, Nat. Photon 7, 1001–1005 (2013).
  * Mittal  _et al._ [2014] S. Mittal, J. Fan, S. Faez, A. Migdall, J. M. Taylor, and M. Hafezi, Topologically robust transport of photons in a synthetic gauge field, Phys. Rev. Lett. 113, 087403 (2014).
  * Gong  _et al._ [2018] Z. Gong, Y. Ashida, K. Kawabata, K. Takasan, S. Higashikawa, and M. Ueda, Topological phases of non-hermitian systems, Phys. Rev. X 8, 031079 (2018).
  * Colella  _et al._ [2019] E. Colella, F. Mivehvar, F. Piazza, and H. Ritsch, Hofstadter butterfly in a cavity-induced dynamic synthetic magnetic field, Phys. Rev. B 100, 224306 (2019).
  * [77] A. Kosior, K. Gietka, F. Mivehvar, and H. Ritsch, Dynamical photon condensation into wannier-stark states.
  * Claro  _et al._ [2003] F. Claro, J. F. Weisz, and S. Curilef, Interaction-induced oscillations in correlated electron transport, Phys. Rev. B 67, 193101 (2003).
  * Dias  _et al._ [2007] W. S. Dias, E. M. Nascimento, M. L. Lyra, and F. A. B. F. de Moura, Frequency doubling of Bloch oscillations for interacting electrons in a static electric field, Phys. Rev. B 76, 155124 (2007).
  * Khomeriki  _et al._ [2010] R. Khomeriki, D. O. Krimer, M. Haque, and S. Flach, Interaction-induced fractional Bloch and tunneling oscillations, Phys. Rev. A 81, 065601 (2010).
  * Wiater  _et al._ [2017] D. Wiater, T. Sowiński, and J. Zakrzewski, Two bosonic quantum walkers in one-dimensional optical lattices, Phys. Rev. A 96, 043629 (2017).
  * Corrielli  _et al._ [2013] G. Corrielli, A. Crespi, G. Della Valle, S. Longhi, and R. Osellame, Fractional Bloch oscillations in photonic lattices, Nat. Commun 4, 1555 (2013).
  * Preiss  _et al._ [2015] P. M. Preiss, R. Ma, M. E. Tai, A. Lukin, M. Rispoli, P. Zupancic, Y. Lahini, R. Islam, and M. Greiner, Strongly correlated quantum walks in optical lattices, Science 347, 1229–1233 (2015).
  * Kosior  _et al._ [2023b] A. Kosior, S. Kokkelmans, M. Lewenstein, J. Zakrzewski, and M. Płodzień, Phonon-assisted coherent transport of excitations in Rydberg-dressed atom arrays, Phys. Rev. A 108, 043308 (2023b).

##  I Supplemental Material

###  I.1 Derivation of the rotation frame Hamiltonian

Consider an array of optical Kerr cavities labeled by $j=-L/2,\ldots,L/2$.
Each cavity in the array is assumed to have a distinct primary mode with a
unique frequency denoted as $\omega_{j}$. The neighboring cavities are coupled
through photon exchange, described as a nearest-neighbor tunneling process,
where a photon leaks out of one cavity and enters the adjacent one. The full
Hamiltonian in the lab frame is given by

$\displaystyle\hat{H}_{\text{Lab}}=-J\sum_{j}\left(\hat{a}_{j}^{\dagger}\hat{a}_{j+1}+\text{H.c.}\right)+\sum_{j}\omega_{j}\hat{a}_{j}^{\dagger}\hat{a}_{j}+\chi\sum_{j}\hat{a}_{j}^{\dagger
2}\hat{a}_{j}^{2}+\sum_{j}\eta_{j}\left(\hat{a}_{j}e^{i\omega_{p}t}+\hat{a}_{j}^{\dagger}e^{-i\omega_{p}t}\right),$
(S1)

where $\hat{a}_{j}$, $\hat{a}_{j}^{\dagger}$ are bosonic operators for
annihilating and creating photons in the $j$-th cavity with frequency
$\omega_{j}=\omega_{0}+j\Delta\omega$. Here, $J$ is the tunneling amplitude,
$\chi$ is the Kerr non-linearity leading to on-site photon-photon interaction,
and $\eta_{j}$ is the pumping rate of the $j$th cavity with a pumping laser of
frequency $\omega_{p}$.

To switch to a rotating frame with slowly varying variables, a unitary
transformation with a time-dependent operator is applied

$U(t)=\prod_{j}\exp\left(i\omega_{p}t\hat{a}_{j}^{\dagger}\hat{a}_{j}\right),$
(S2)

transforming the bosonic operators to
$\hat{a}^{\prime}_{j}=U\hat{a}_{j}U^{\dagger}=\exp(-i\omega_{p}t)\hat{a}_{j}$
and the Hamiltonian to
$\hat{H}=U\hat{H}_{\text{Lab}}U^{\dagger}+i\left(\partial_{t}U\right)U^{\dagger}$,
i.e.,

$\displaystyle\hat{H}=-J\sum_{j}\left(\hat{a}_{j}^{\dagger}\hat{a}_{j+1}+\text{H.c.}\right)+\sum_{j}\Delta_{j}\hat{a}_{j}^{\dagger}\hat{a}_{j}+\chi\sum_{j}\hat{a}_{j}^{\dagger
2}\hat{a}_{j}^{2}+\sum_{j}\eta_{j}\left(\hat{a}_{j}+\hat{a}_{j}^{\dagger}\right),$
(S3)

where $\Delta_{j}=\Delta\omega\left(j-j_{0}\right)$ with
$j_{0}=(\omega_{p}-\omega_{0})/\Delta\omega$.

###  I.2 Bloch oscillations and consequences of U(1) symmetry

Here, for pedagogical reasons, we first consider a non-interacting bosonic
Hamiltonian in a tilted potential and show the system exhibits Bloch
oscillations. Then, we introduce a pumping term
$\sum_{j}\eta_{j}\left(\hat{a}_{j}+\hat{a}_{j}^{\dagger}\right)$ and discuss
the consequences of U(1) symmetry breaking in weakly interacting systems.

####  I.2.1 Bloch oscillations

The non-interacting Hamiltonian of our interest reads

$\hat{H}_{0}=-J\sum_{j}\left(\hat{a}_{j}^{\dagger}\hat{a}_{j+1}+\text{H.c.}\right)+\sum_{j}\Delta_{j}\hat{a}_{j}^{\dagger}\hat{a}_{j}$
(S4)

with $\Delta_{j}=\Delta\omega\left(j-j_{0}\right)$, where the choice of
$j_{0}$ is irrelevant due to the conservation of the particle number. The
Hamiltonian is quadratic in $\hat{a}_{j}^{\dagger}$, $\hat{a}_{j}$ and can be
readily diagonalised by a unitary transformation

$\hat{b}_{n}=\sum_{j=-\infty}^{\infty}\beta_{n,j}\hat{a}_{j},\quad\beta_{n,j}=\mathcal{J}_{j-n}\left(\gamma=\frac{2J}{\Delta\omega}\right)$
(S5)

with $\mathcal{J}_{m}$ being the Bessel function of the first kind of order
$m$. Since the transformation is unitary, then $\beta_{n,j}$ must fulfill

$\sum_{j}\beta_{n,j}\beta_{m,j}=\delta_{n,m}.$ (S6)

Indeed, using the Bessel function addition theorem

$\sum_{m}\mathcal{J}_{m}(x)\mathcal{J}_{n-m}(y)=\mathcal{J}_{n}(x+y)$ (S7)

as well as the property $\mathcal{J}_{n}(x)=\mathcal{J}_{-n}(-x)$, we obtain

$\sum_{j}\mathcal{J}_{j-n}(\gamma)\mathcal{J}_{j-m}(\gamma)=\mathcal{J}_{m-n}(0)=\delta_{n,m}.$
(S8)

The inverse transformation then reads

$\hat{a}_{j}=\sum_{n}\beta_{n,j}\hat{b}_{n}.$ (S9)

Figure S1:  (a)-(c) The expectation values of particle number operators
$\langle\hat{n}_{j}(t)\rangle$ in early times of evolution in a very weakly
interacting case with a different pumping $\eta_{j}=\delta_{j,j^{\prime}}$:
(a) $j^{\prime}=0$, (b) $j^{\prime}=7$, (c) $j^{\prime}=20$. (d)-(e) The
respective site occupation distributions after long time evolution in
stationary states. In all panels $j_{0}=0$, $\Delta\omega=0.3$ and
$\chi=10^{-2}$.

The eigenstates $|\Psi_{n}\rangle=\hat{b}_{n}^{\dagger}|0\rangle$ are called
the Wannier-Stark (WS) states. In practise, the infinite summation can be
truncated, as the Bessel function $\mathcal{J}_{j-n}(\gamma)$ is mainly
localised in the interval $|j-n|<|\gamma|$. Although the formula holds for an
infinite lattice, the WS states have finite spatial size, and therefore it can
be used with a very good accuracy to describe bulk states of a finite-size
lattice. Using the WS basis, the Hamiltonian (S4) can be recast into a
diagonal form,

$\hat{H}_{0}=\Delta\omega\sum_{n}n\,\hat{b}_{n}^{\dagger}\hat{b}_{n}.$ (S10)

The time evolution of the Wannier state,
$|j\rangle=\hat{a}_{j}^{\dagger}|0\rangle$, reads

$|j(t)\rangle=\exp(i\hat{H}_{0}\,t)|j\rangle=\sum_{n}\beta_{n,j}e^{in\Delta\omega
t}|\Psi_{n}\rangle=\sum_{n,k}\beta_{n,j}\beta_{n,k}e^{in\Delta\omega
t}|k\rangle,$ (S11)

which is a $T$-periodic function with period $T=2\pi/\Delta\omega$.
Consequently, a state which is initially maximally localised on a tilted
lattice performs Bloch oscillations around the initial position.

####  I.2.2 U(1) symmetry

Because the total number of particles in the Hamiltonian $\hat{H}_{0}$ of Eq.
(S4) is conserved, we were able to neglect the term
$\hat{K}=-j_{0}\sum_{j}\hat{a}_{j}^{\dagger}\hat{a}_{j}$, that commutes with
the Hamiltonian. Equivalently, $\hat{K}$ can be gauged out by applying the
time-dependent transformation

$V(t)=\Pi_{j}\exp\left(i\chi(t)\hat{a}_{j}^{\dagger}\hat{a}_{j}\right),\quad\chi(t)=\int^{t}j_{0}(\tau)\mbox{d}\tau,$
(S12)

where we assumed that $j_{0}$ can be explicitly time dependent. The
transformation only modifies the global phase of the bosonic operators
$\hat{a}^{\prime}_{j}=V\hat{a}_{j}\tilde{V}^{\dagger}=\exp[-i\chi(t)]\hat{a}_{j}$
and does not contribute to the dynamics. It is a mere manifestation of the
$U(1)$ symmetry of the model.

The full Hamiltonian of our system (S3), which we study in the main text,
consist of two additional terms,

$\hat{H}=\hat{H}_{0}+\chi\sum_{j}\hat{a}_{j}^{\dagger
2}\hat{a}_{j}^{2}+\sum_{j}\eta_{j}\left(\hat{a}_{j}+\hat{a}_{j}^{\dagger}\right),$
(S13)

where the first term describes contact boson-boson interactions, and the
second accounts for coherent injection of bosons into the lattice. Although
the interacting term is invariant under the $U(1)$ symmetry of Eq. (S12), in
general it is responsible for the decay of Bloch oscillations [32, 33, 34,
35]. The second term explicitly breaks the $U(1)$ symmetry of the model, and
therefore $\hat{K}=-j_{0}\sum_{j}\hat{a}_{j}^{\dagger}\hat{a}_{j}$ does not
commute with $\hat{H}$ and cannot be neglected anymore. As described in the
main text, the presence of $U(1)$ symmetry breaking terms in the Hamiltonian
has enormous consequences to the system’s dynamics and is responsible for the
dynamical photon condensation.

Additionally, in Fig. S1(a)-(c) we show the early-time evolution of
expectation values of particle number operators
$\langle\hat{n}_{j}\rangle=\langle\hat{a}_{j}^{\dagger}\hat{a}_{j}\rangle$ in
each lattice site with different pumping schemes in the non-interacting case,
while in Fig. S1(d)-(e) we show the corresponding long-time occupation numbers
of photons over lattice sites. The final state is proportional to a central WS
state irrespective of $j^{\prime}$, as long as $j^{\prime}$ lies within the WS
localization range, i.e., approximately, $|j^{\prime}|<|\gamma|$.

###  I.3 Cumulant expansion

In this section, we briefly introduce the cumulant expansion method. For an
open system, the time evolution of an average of an observable
$\hat{\mathcal{O}}$ can be calculated according to

$\frac{d\langle\hat{\mathcal{O}}\rangle}{dt}=\frac{i}{\hbar}\langle[\hat{H},\hat{\mathcal{O}}]\rangle+\sum_{j}\kappa_{j}\langle
2\hat{c}_{j}\hat{\mathcal{O}}\hat{c}_{j}^{\dagger}-\hat{c}_{j}^{\dagger}\hat{c}_{j}\hat{\mathcal{O}}-\hat{\mathcal{O}}\hat{c}_{j}^{\dagger}\hat{c}_{j}\rangle,$
(S14)

where $\kappa_{j}$ characterizes the $j$th decay channel and $c_{j}$ the
corresponding jump operator. Writing down the equations explicitly, we
typically find out an infinite hierarchy of equations for products of
operators, which can be attributed to the non-commutativity of various
operators. A well-established approach to deal with the infinite set of
equations is to neglected quantum correlations. In such a case, an average
value of operator product can be rewritten as product of average values of
operators,
$\langle\hat{O}_{1}\hat{O}_{2}\rangle=\langle\hat{O}_{1}\rangle\langle\hat{O}_{2}\rangle$,
which is typically referred to as the mean field approach. More generally, it
is possible to neglect quantum correlations of an arbitrary order. A
systematic approach to such a truncation is realized by a so-called cumulant
expansion, which relies on decomposing average values of arbitrary operator
products into products of average values of operators of a given lower order,
leading to a closed set of equations.

###  I.4 Interaction-induced Bloch oscillations

In this section we focus on interaction-induced Bloch-like oscillations, which
arise in our system whenever stationary solutions are not dynamically stable,
in a similarity to the limit cycle solutions in driven dissipative systems
[add citations]. Indeed, these oscillations can be interpreted as multi-mode
limit cycles as shown in Fig. S2(a). [On the other hand, stronger on-site
interactions result in chaotic behavior, see Fig. S2(b).]

In the main text of this manuscript, we argue that the nonequilibrium phases
of the system can be identified by two observables:

  * •

$\Delta n$—relative change of the total photon number over a long-time
evolution, defined as

$\Delta
n=\frac{\max_{t}\langle\hat{N}(t)\rangle-\min_{t}\langle\hat{N}(t)\rangle}{\text{avg}_{t}\langle\hat{N}(t)\rangle},$
(S15)

with $\max_{t}$, $\min_{t}$ and $\text{avg}_{t}$ denoting maximal, minimal and
the average value of an observable during a long-time evolution (where the
initial transient dynamics of the system has been neglected), and

  * •

$\langle\max[P_{n}(t)]\rangle_{av}$—time-averaged maximal fidelity (i.e.
square modulus of the overlap) between the mean-field wavefunction
$|\psi\rangle$ and the WS basis
$|\Psi_{n}\rangle=\hat{b}_{n}^{\dagger}|0\rangle=\sum_{j}\beta_{n,j}\hat{a}_{j}^{\dagger}|0\rangle$.

In order to find the stable islands of interaction-induced Bloch oscillation
it is convenient to define a third observable, as a simple product of the
previous two observables,

$\Theta=\Delta n\cdot\mathrm{avg}_{t}\langle\max[P_{n}(t)]\rangle$ (S16)

see Fig. S3, where we show that $\Theta$ is close to zero both in the
stationary and chaotic phases, and therefore it is a good observable to
identify non-stationary regular dynamics in the system.

Finally, in Fig. S4 we illustrate the expectation values of particle number
operators $\langle\hat{n}_{j}(t)\rangle$ in the course of time evolution for
different parameters of the system, showing more examples of interaction-
induced oscillations.

Figure S2:  Evolution of the photonic field amplitudes
$\langle\hat{a}_{j}\rangle$ on the complex plane for $\Delta\omega=0.5$ and
(a) $\chi=0.13$, (b) $\chi=0.135$ [as in Fig. 4 in the main text]. (a) In the
regime corresponding to the regular oscillatory behavior, the photonic fields
evolve along closed trajectories. (b) In the unstable regime, photonic fields
evolve chaotically. For the clarity of the presentation only selected modes
are shown.  Figure S3: The density plot of $\Theta$ which identifies stability
islands of interaction-induced regular oscillations (bright colors). Panels
(b) and (c) depict cuts through the density plot at the positions of red
dashed lines. (b) $\Delta\omega=0.525$, (c) $\Delta\omega=0.49$.  Figure S4:
Panels present the expectation values of particle number operators
$\langle\hat{n}_{j}(t)\rangle$ in the course of time evolution showing more
examples of interaction-induced oscillations for different parameters of the
system. (a) $\chi=0.11$, $\Delta\omega=0.525$ (b) $\chi=0.05$,
$\Delta\omega=0.2$, (c) $\chi=0.03$, $\Delta\omega=0.45$, (d) $\chi=0.09$,
$\Delta\omega=0.525$. Other parameters as in the main text of this Letter.