# The generalized Ermakov conservative system: A discussion

Antonios Mitsopoulos1,a) and Michael Tsamparlis1,b)

1Faculty of Physics, Department of Astronomy-Astrophysics-Mechanics,
University of Athens, Panepistemiopolis, Athens 157 83, Greece
a)Email<EMAIL_ADDRESS>
b)Email<EMAIL_ADDRESS>
###### Abstract

Using older and recent results on the integrability of two-dimensional (2d)
dynamical systems, we prove that the results obtained in a recent publication
concerning the 2d generalized Ermakov system can be obtained as special cases
of a more general approach. This approach is geometric and can be used to
study efficiently similar dynamical systems.

##  1 Introduction

The two-dimensional (2d) generalized Ermakov system has attained attention in
90’s ([1], [2], [3], [4], [5], [6]) where most of its properties have been
revealed. A review of these studies can be found in [7]. However, the interest
in the topic is still alive and a recent article has appeared in this journal
[8] presenting new results. The purpose of the present discussion is to show
that these latter results can be obtained as special cases of older and more
recent results on the integrability of 2d dynamical systems.

##  2 The 2d generalized Ermakov system

The 2d generalized Ermakov system is defined by the equations

$\displaystyle\ddot{x}$ $\displaystyle=$
$\displaystyle-\omega^{2}(t)x+\frac{1}{x^{2}y}f\left(\frac{y}{x}\right)$ (1)
$\displaystyle\ddot{y}$ $\displaystyle=$
$\displaystyle-\omega^{2}(t)y+\frac{1}{xy^{2}}g\left(\frac{x}{y}\right)$ (2)

where $f,g$ are arbitrary functions. This system admits the Ermakov first
integral (FI)

$I_{0}=\frac{1}{2}(x\dot{y}-y\dot{x})^{2}+\int^{y/x}f(u)du+\int^{x/y}g(v)dv$
(3)

where $u=v^{-1}=\frac{y}{x}$. It is well-known that the 2d generalized Ermakov
system generalizes the one-dimensional (1d) time-dependent oscillator.

Introducing the functions $F,G$ by the relations

$f\left(\frac{y}{x}\right)=\frac{y}{x}F\left(\frac{y}{x}\right),\enskip
g\left(\frac{x}{y}\right)=\frac{x}{y}G\left(\frac{y}{x}\right)$

equations (1), (2) take the equivalent form

$\displaystyle\ddot{x}$ $\displaystyle=$
$\displaystyle-\omega^{2}(t)x+\frac{1}{x^{3}}F\left(\frac{y}{x}\right)$ (4)
$\displaystyle\ddot{y}$ $\displaystyle=$
$\displaystyle-\omega^{2}(t)y+\frac{1}{y^{3}}G\left(\frac{y}{x}\right)$ (5)

while the Ermakov FI (3) becomes

$I_{0}=\frac{1}{2}(x\dot{y}-y\dot{x})^{2}+\int^{y/x}\left[uF(u)-u^{-3}G(u)\right]du.$
(6)

If one introduces the variables [4]

$T=\int\rho^{-2}dt,\enskip X=\rho^{-1}x,\enskip Y=\rho^{-1}y$ (7)

where $\rho(t)$ is a solution of the 1d time-dependent oscillator

$\ddot{\rho}+\omega^{2}(t)\rho=0$ (8)

then equations (4), (5) become the autonomous system

$\displaystyle X^{\prime\prime}$ $\displaystyle=$
$\displaystyle\frac{1}{X^{3}}F\left(\frac{Y}{X}\right)$ (9) $\displaystyle
Y^{\prime\prime}$ $\displaystyle=$
$\displaystyle\frac{1}{Y^{3}}G\left(\frac{Y}{X}\right)$ (10)

and the Ermakov FI

$I_{0}=\frac{1}{2}\left(XY^{\prime}-YX^{\prime}\right)^{2}+\int^{Y/X}\left[uF(u)-u^{-3}G(u)\right]du$
(11)

where we use the notation $f^{\prime}\equiv\frac{df(T)}{dT}$ and
$\dot{f}\equiv\frac{df(t)}{dt}$. For general functions $F,G$ the autonomous
dynamical system (9) - (10) is not conservative. In the following we determine
the family of conservative Ermakov systems together with their FIs using
collineations of the metric defined by these equations.

##  3 Integrability of the 2d generalized Ermakov system

Since the system (9) - (10) is autonomous, the second FI will be the
Hamiltonian $H$. To find $H$ we do not have to do any new calculations,
because in [9], [10] all the integrable and superintegrable 2d autonomous
conservative systems have been determined. From these results, we find (see
e.g. section 8 case (1) in [10]) that the 2d integrable potential

$V_{21}=\frac{F_{1}(u)}{X^{2}+Y^{2}}+F_{2}(X^{2}+Y^{2})$ (12)

where $u=\frac{Y}{X}$ and $F_{1},F_{2}$ are arbitrary functions of their
arguments, admits the FI

$I_{11}=\frac{1}{2}\left(XY^{\prime}-YX^{\prime}\right)^{2}+F_{1}(u).$ (13)

If we consider $F_{1}(u)=(u^{2}+1)N(u)$ and $F_{2}=0$, we find that
$V_{21}=\frac{N(u)}{X^{2}}$ while the resulting equations are

$\displaystyle X^{\prime\prime}$ $\displaystyle=$
$\displaystyle\frac{2N+u\frac{dN}{du}}{X^{3}}$ (14) $\displaystyle
Y^{\prime\prime}$ $\displaystyle=$
$\displaystyle-\frac{1}{X^{3}}\frac{dN}{du}.$ (15)

Therefore, if we define the functions

$F(u)=2N+u\frac{dN}{du},\enskip G(u)=-u^{3}\frac{dN}{du}$ (16)

then equations (14), (15) become the Ermakov equations (9), (10) while
$I_{0}=I_{11}$.

We conclude that the family of the conservative 2d Ermakov systems is defined
by the potential $V=\frac{N(u)}{X^{2}}$ where $N(u)$ is an arbitrary function
while the Hamiltonian is given by the expression

$H=\frac{1}{2}\left(X^{\prime 2}+Y^{\prime 2}\right)+\frac{N(u)}{X^{2}}.$ (17)

In the original coordinates the system (14) - (15) becomes

$\displaystyle\ddot{x}$ $\displaystyle=$
$\displaystyle-\omega^{2}(t)x+\frac{2N(u)+u\frac{dN}{du}}{x^{3}}$ (18)
$\displaystyle\ddot{y}$ $\displaystyle=$
$\displaystyle-\omega^{2}(t)y-\frac{1}{x^{3}}\frac{dN}{du}$ (19)

where $u=\frac{y}{x}=\frac{Y}{X}$.

For $N(u)=\frac{u^{-2}}{2}$ we find respectively

$\displaystyle\ddot{x}$ $\displaystyle=$ $\displaystyle-\omega^{2}(t)x$ (20)
$\displaystyle\ddot{y}$ $\displaystyle=$
$\displaystyle-\omega^{2}(t)y+\frac{1}{y^{3}}$ (21)

while the Ermakov FI becomes the well-known Lewis invariant [11]

$I_{0}=\frac{1}{2}(x\dot{y}-y\dot{x})^{2}+\frac{1}{2}\left(\frac{x}{y}\right)^{2}.$
(22)

Equation (20) is the 1d time-dependent harmonic oscillator and (21) is the
auxiliary equation with which one determines the frequency $\omega(t)$ for a
given function $y(t)$. These justify the characterization of the Ermakov
system as a generalization of the harmonic oscillator.

##  4 The FIs of the conservative Ermakov system

There are two ways to find the FIs of the Ermakov system. One way is to use
the results of [9], [10] and read the FIs for the potential $V_{21}$ given in
(12) for $F_{1}(u)=(u^{2}+1)N(u)$, $F_{2}=0$. Here, we shall follow another
way which can be useful in many similar problems. We shall use Theorem 2 of
[12] where it is stated that the generators of the Noether point symmetries of
autonomous conservative systems are the elements of the homothetic algebra of
the metric defined by the kinetic energy (kinetic metric). In the Ermakov
case, this metric is the Euclidean 2d metric $\gamma_{ab}=diag(1,1)$.

For the convenience of the reader we state Theorem 2 of [12].

###### Theorem 1

Autonomous conservative dynamical systems of the form

$\ddot{q}^{a}=-\Gamma^{a}_{bc}\dot{q}^{b}\dot{q}^{c}-V^{,a}(q)$ (23)

where $\Gamma_{bc}^{a}$ are the Riemannian connection coefficients determined
from the kinetic metric $\gamma_{ab}$ (kinetic energy) and $V(q)$ the
potential of the system, admit the following point Noether symmetries.

Case 1. The point Noether symmetry

$\mathbf{A}_{1}=\partial_{t},\enskip f_{1}=const\equiv 0$ (24)

which produces the Noether FI (Hamiltonian)

$H=\frac{1}{2}\gamma_{ab}\dot{q}^{a}\dot{q}^{b}+V(q).$ (25)

Case 2. The point Noether symmetry

$\mathbf{A}_{2}=2\psi_{B}t\partial_{t}+B^{a}\partial_{q^{a}},\enskip
f_{2}=c_{1}t$ (26)

where $c_{1}$ is an arbitrary constant and $B^{a}$ is a KV ($\psi_{B}=0$) or
the HV ($\psi_{B}=1$) such that

$B_{a}V^{,a}+2\psi_{B}V+c_{1}=0.$ (27)

The associated Noether FI is

$I_{2}=2\psi_{B}tH-B_{a}\dot{q}^{a}+c_{1}t.$ (28)

Case 3. The point Noether symmetry

$\mathbf{A}_{3}=2\psi\int
C(t)dt\partial_{t}+C(t)\Phi^{,a}\partial_{q^{a}},\enskip
f_{3}=C_{,t}\Phi(q)+D(t)$ (29)

where $\Phi^{,a}(q)$ is a gradient KV ($\psi=0$) or a gradient HV ($\psi=1$)
such that ($c_{2},c_{3}$ are arbitrary constants)

$\Phi_{,a}V^{,a}+2\psi V=c_{2}\Phi+c_{3}$ (30)

and the functions $C(t)$, $D(t)$ are determined by the relations ($C_{,t}\neq
0$)

$C_{,tt}=-c_{2}C,\enskip D_{,t}=-c_{3}C.$ (31)

The associated Noether FI is

$I_{3}=2\psi H\int C(t)dt-C(t)\Phi_{,a}\dot{q}^{a}+C_{,t}\Phi-c_{3}\int
C(t)dt.$ (32)

We apply Theorem 1 in the case of the autonomous integrable Ermakov system
(14) - (15) which has potential $V=\frac{N(u)}{X^{2}}$, where $u=Y/X$, and
kinetic metric $\gamma_{ab}=diag(1,1)$.

The homothetic algebra of $\gamma_{ab}$ consists of two gradient Killing
vectors (KVs) $\partial_{X},\partial_{Y}$, one non-gradient KV (rotation)
$Y\partial_{X}-X\partial_{Y}$ and the gradient homothetic vector (HV)
$X\partial_{X}+Y\partial_{Y}$.

For each case of Theorem 1 we have the following.

###  4.1 The vector $\partial_{T}$

Case 1. In this case the point Noether symmetry $\mathbf{A}_{1}=\partial_{T}$,
$f_{1}=0$ produces the Hamiltonian (as expected)

$H=\frac{1}{2}\left(X^{\prime 2}+Y^{\prime
2}\right)+V=\frac{1}{2}\left(X^{\prime 2}+Y^{\prime
2}\right)+\frac{N(u)}{X^{2}}.$ (33)

###  4.2 The gradient HV $X\partial_{X}+Y\partial_{Y}$

Case 2. Consider the gradient HV $B^{a}=(X,Y)$ with homothety factor
$\psi_{B}=1$.

Substituting in condition (27) we find that

$XV_{,X}+YV_{,Y}+2V+c_{1}=0\implies-\frac{2N+u\frac{dN}{du}}{X^{2}}+\frac{u}{X^{2}}\frac{dN}{du}+2\frac{N}{X^{2}}+c_{1}=0\implies
c_{1}=0.$

Therefore, the point Noether symmetry is

$\mathbf{A}_{2}=2T\partial_{T}+X\partial_{X}+Y\partial_{Y},\enskip f_{2}=0$
(34)

and the associated Noether FI

$I_{2}=2TH-(XX^{\prime}+YY^{\prime}).$ (35)

It can be shown that the three FIs $I_{0},H,I_{2}$ are independent; therefore,
the conservative generalized Ermakov system is superintegrable.

Although the remaining FIs will be expressible in terms of the
$I_{0},H,I_{2}$, we continue in order to show that we recover the results of
[5] which were obtained using Lie symmetries.

Case 3. The point Noether symmetry is

$\mathbf{A}_{3}=T^{2}\partial_{T}+TX\partial_{X}+TY\partial_{Y},\enskip
f_{3}=\frac{X^{2}+Y^{2}}{2}$ (36)

with associated Noether integral

$I_{3}=T^{2}H-T(XX^{\prime}+YY^{\prime})+\frac{X^{2}+Y^{2}}{2}=\frac{I_{2}^{2}+2I_{0}}{4H}.$
(37)

We observe that the Lie symmetries (2.9a), (2.9b), (2.9c) found in [5] are the
point Noether symmetries $\mathbf{A}_{1}$, (34), (36). Concerning the
remaining FIs of [5] we have: (4.11) $I^{\prime}=2I_{0}$, (4.12)
$J_{1}^{\prime}=2H$, (4.13) $J^{\prime}_{2}=I_{2}$ and (4.14)
$J^{\prime}_{3}=2I_{3}$. Using these relations, equation (4.18) is equivalent
to the expression (37).

###  4.3 The gradient KV $b_{1}\partial_{X}+b_{2}\partial_{Y}$

The potential becomes111This is a superintegrable potential of the form
$F(b_{1}Y-b_{2}X)$ (see section 7 in [10]).
$V_{1}=\frac{k}{(b_{1}Y-b_{2}X)^{2}}$ where $k,b_{1},b_{2}$ are arbitrary
constants.

Case 2. The Noether generator, the Noether function and the FI are

$\mathbf{A}_{21}=b_{1}\partial_{X}+b_{2}\partial_{Y},\enskip f_{21}=0,\enskip
I_{21}=b_{1}X^{\prime}+b_{2}Y^{\prime}.$

Case 3. The Noether generator, the Noether function and the FI are

$\mathbf{A}_{31}=Tb_{1}\partial_{X}+Tb_{2}\partial_{Y},\enskip
f_{31}=b_{1}X+b_{2}Y,\enskip
I_{31}=b_{1}(-TX^{\prime}+X)+b_{2}(-TY^{\prime}+Y).$

In order to compare these results with the ones of [8] we use polar
coordinates $X=r\cos\theta$ and $Y=r\sin\theta$. We find:

$V_{1}=\frac{k}{r^{2}(b_{1}\sin\theta-b_{2}\cos\theta)^{2}}$
$\mathbf{A}_{21}=\left(b_{1}\cos\theta+b_{2}\sin\theta\right)\partial_{r}+\frac{1}{r}\left(b_{2}\cos\theta-
b_{1}\sin\theta\right)\partial_{\theta},\enskip f_{21}=0$
$\mathbf{A}_{31}=\left(b_{1}T\cos\theta+b_{2}T\sin\theta\right)\partial_{r}+\frac{1}{r}\left(b_{2}T\cos\theta-
b_{1}T\sin\theta\right)\partial_{\theta},\enskip
f_{31}=r(b_{1}\cos\theta+b_{2}\sin\theta)$

and

$\displaystyle I_{21}$ $\displaystyle=$ $\displaystyle
b_{1}\left(r^{\prime}\cos\theta-r\theta^{\prime}\sin\theta\right)+b_{2}\left(r^{\prime}\sin\theta+r\theta^{\prime}\cos\theta\right)$
$\displaystyle=$ $\displaystyle
b_{1}\left(\bar{p}_{1}\cos\theta-\frac{\bar{p}_{2}\sin\theta}{r}\right)+b_{2}\left(\bar{p}_{1}\sin\theta+\frac{\bar{p}_{2}\cos\theta}{r}\right)$
$\displaystyle I_{31}$ $\displaystyle=$ $\displaystyle
b_{1}\left(-Tr^{\prime}\cos\theta+Tr\theta^{\prime}\sin\theta+r\cos\theta\right)+b_{2}\left(-Tr^{\prime}\sin\theta-
Tr\theta^{\prime}\cos\theta+r\sin\theta\right)$ $\displaystyle=$
$\displaystyle
b_{1}\left(-T\bar{p}_{1}\cos\theta+\frac{T\bar{p}_{2}\sin\theta}{r}+r\cos\theta\right)+b_{2}\left(-T\bar{p}_{1}\sin\theta-\frac{T\bar{p}_{2}\cos\theta}{r}+r\sin\theta\right)$

where $\bar{p}_{a}=\bar{\gamma}_{ab}\bar{q}^{b\prime}$ are the generalized
momenta. Replacing with $\bar{q}_{a}=(r,\theta)$ and
$\bar{\gamma}_{ab}=diag(1,r^{2})$ we find that $\bar{p}_{1}=r^{\prime}$ and
$\bar{p}_{2}=r^{2}\theta^{\prime}$.

It is straightforward to show that the point Noether symmetries
$\mathbf{A}_{21},\mathbf{A}_{31}$, the Noether functions $f_{21},f_{31}$ and
the FIs $I_{21},I_{31}$ are the symmetries $\mathbf{X}_{9},\mathbf{X}_{10}$,
the functions $B_{9},B_{10}$ and the FIs $I_{9},I_{10}$ respectively of [8];
while
\- for $b_{1}=0,b_{2}=1$ they reduce to the symmetries
$\mathbf{X}_{5},\mathbf{X}_{6}$, the functions $B_{5},B_{6}$ and the FIs
$I_{5},I_{6}$ respectively of [8] and
\- for $b_{1}=1$, $b_{2}=0$ they reduce to the symmetries
$\mathbf{X}_{7},\mathbf{X}_{8}$, the functions $B_{7},B_{8}$ and the FIs
$I_{7},I_{8}$ respectively of [8].

As expected, the non-gradient KV (rotation) $Y\partial_{X}-X\partial_{Y}$
leads to the linear first integral of angular momentum.

Finally, using the three FIs $H,I_{0},I_{2}$, we integrate the system (14) -
(15) and find that in polar coordinates the solution is

$\displaystyle r^{2}(T)$ $\displaystyle=$ $\displaystyle\frac{1}{2H}\left(2HT-
I_{2}\right)^{2}+\frac{I_{0}}{H}$ (38)
$\displaystyle\int\frac{d\theta}{\sqrt{I_{0}-\bar{F}(\theta)}}$
$\displaystyle=$
$\displaystyle\pm\int\frac{\sqrt{2}}{r^{2}(T)}dT=\pm\frac{1}{\sqrt{I_{0}}}\tan^{-1}\left[\frac{1}{\sqrt{2I_{0}}}(2HT-
I_{2})\right]$ (39)

where $\bar{F}(\theta)=(\tan^{2}\theta+1)N(\tan\theta)$. In Table 3 of [8] the
corresponding formula of (39) gives $\tanh^{-1}$ instead of
$\tan^{-1}=\arctan$ which is the correct result.

##  5 Conclusions

Using recent results on the integrability of 2d conservative dynamical
systems, we proved that the generalized Ermakov system is superintegrable and
determined all the quadratic FIs. We showed that the recent results of [8] can
be obtained from the more general method outlined in [10] and [12]. Obviously,
the methods discussed in the present work can be used by other authors in the
study of similar dynamical systems. As a final remark, we note that an
extension of this method, where one computes the time-dependent and autonomous
quadratic FIs without the use of Noether symmetries, has appeared recently in
[13].

## Data Availability

All data that supports the findings of this study are available within the
article.

## Conflict of interest

The authors declare no conflict of interest.

## References

  * [1] V.P. Ermakov, Appl. Anal. Discrete Math. 2, 123 (2008).
  * [2] J.R. Ray and J.L. Reid, Phys. Lett. A 71(4), 317 (1979).
  * [3] C. Athorne, Phys. Lett. A 151(8), 407 (1990).
  * [4] P.G.L. Leach, Phys. Lett. A 158(3,4), 102 (1991).
  * [5] K.S. Govinder and P.G.L. Leach, Phys. Lett. A 186, 391 (1994).
  * [6] J. Goedert and F. Haas, Phys. Lett. A 239, 348 (1998).
  * [7] P.G.L. Leach and K. Andriopoulos, Appl. Anal. Discrete Math. 2, 146 (2008).
  * [8] R. Naz, Eur. Phys. J. Plus 135, 641 (2020).
  * [9] J. Hietarinta, Phys. Rep. 147(2), 87 (1987).
  * [10] A. Mitsopoulos, M. Tsamparlis and A. Paliathanasis, Symmetry 12, 1655 (2020).
  * [11] H.R. Lewis, J. Math. Phys. 9(11), 1976 (1968).
  * [12] M. Tsamparlis and A. Paliathanasis, J. Phys. A: Math. Theor. 44, 175202 (2011).
  * [13] A. Mitsopoulos and M. Tsamparlis, Mathematics 9(13), 1503 (2021).


#  pySYD: Automated measurements of global asteroseismic parameters

Ashley Chontos  Daniel Huber  Maryum Sayeed  Pavadol Yamsiri

††margin:  DOI: https://doi.org/DOI unavailablehttps://doi.orgDOI unavailable
Software • N/Ahttps://doi.orgReview • NO_REPOSITORYhttps://doi.orgRepository •
DOI unavailablehttps://doi.orgArchive Editor:
https://example.comhttps://doi.orgPending Editor
Reviewers: • https://github.com/Pending Reviewershttps://doi.org@Pending
Reviewers Submitted: N/A
Published: N/A License
Authors of papers retain copyright and release the work under a Creative
Commons Attribution 4.0 International License
(http://creativecommons.org/licenses/by/4.0/https://doi.orgCC BY 4.0).

## Summary

Asteroseismology, the study of stellar oscillations, is a powerful tool for
studying the interiors of stars and determining their fundamental properties
(Aerts, 2021). For stars with temperatures that are similar to the Sun,
turbulent near-surface convection excites sound waves that propagate within
the stellar cavity (Bedding, 2014). These waves penetrate into different
depths within the star and therefore provide powerful constraints on stellar
interiors that would otherwise be inaccessible. Asteroseismology is well-
established in astronomy as the gold standard for characterizing fundamental
properties like masses, radii, densities, and ages for single stars, which has
broad impacts on several fields in astronomy. For example, ages of stars are
important to reconstruct the formation history of the Milky Way (so-called
galactic archaeology). For exoplanets that are discovered indirectly through
changes in stellar observables, precise and accurate stellar masses and radii
are critical for learning about the planets that orbit them.

## Statement of Need

The NASA space telescopes _Kepler_ , K2 and TESS have recently provided very
large databases of high-precision light curves of stars. By detecting
brightness variations due to stellar oscillations, these light curves allow
the application of asteroseismology to large numbers of stars, which requires
automated software tools to efficiently extract observables. Several tools
have been developed for asteroseismic analyses (e.g., A2Z, Mathur et al.,
2010; COR, Mosser & Appourchaux, 2009; OCT, Hekker et al., 2010; SYD, Huber et
al., 2009), but many of them are closed-source and therefore inaccessible to
the general astronomy community. Some open-source tools exist (e.g., DIAMONDS
and FAMED, Corsaro & De Ridder, 2014; PBjam, Nielsen et al., 2021; lightkurve,
Lightkurve Collaboration et al., 2018), but they are either optimized for
smaller samples of stars or have not yet been extensively tested against
closed-source tools.

pySYD is adapted from the framework of the IDL-based SYD pipeline (Huber et
al., 2009; hereafter referred to as SYD), which has been used frequently to
measure asteroseismic parameters for _Kepler_ stars and has been extensively
tested against other closed-source tools on _Kepler_ data (Hekker et al.,
2011; Verner et al., 2011). Papers based on asteroseismic parameters measured
using the SYD pipeline include Huber et al. (2011), Bastien et al. (2013),
Chaplin et al. (2014), Serenelli et al. (2017), and Yu et al. (2018). pySYD
was developed using the same well-tested methodology, but has improved
functionality including automated background model selection and parallel
processing as well as improved flexibility through a user-friendly interface,
while still maintaining its speed and efficiency. Well-documented, open-source
asteroseismology software that has been benchmarked against closed-source
tools are critical to ensure the reproducibility of legacy results from the
_Kepler_ mission (Borucki et al., 2010). The combination of well-tested
methodology, improved flexibility and parallel processing capabilities will
make pySYD a promising tool for the broader community to analyze current and
forthcoming data from the NASA TESS mission (Ricker et al., 2015).

## The pySYD library

The excitation mechanism for solar-like oscillations is stochastic and modes
are observed over a range of frequencies. Oscillation modes are separated by
the so-called large frequency spacing ($\Delta\nu$) with an approximately
Gaussian-shaped power excess centered on the frequency of maximum power
($\rm\nu_{max}$). The observables $\rm\nu_{max}$ and $\Delta\nu$ are directly
related to fundamental properties such as surface gravity, density, mass and
radius (Kjeldsen & Bedding, 1995).

pySYD is a Python package for detecting solar-like oscillations and measuring
global asteroseismic parameters. Derived parameters include $\rm\nu_{max}$ and
$\Delta\nu$, as well as characteristic amplitudes and time scales of various
granulation processes.

A pySYD pipeline Target class object has two main methods:

  * •

The first module searches for signatures of solar-like oscillations by
implementing a frequency-resolved, collapsed autocorrelation (ACF) method. The
output from this routine provides an estimate for $\rm\nu_{max}$, which is
used as an initial guess for the main pipeline function (i.e. second module).
However if $\rm\nu_{max}$ is already known, the user can provide an estimate
and hence bypass this module.

  * •

The second routine begins by masking out the region in the power spectrum (PS)
with the power excess in order to characterize the stellar background. pySYD
optimizes and selects the stellar background model that minimizes the BIC
(Bayesian Information Criterion; Schwarz, 1978). The best-fit background model
is then subtracted from the PS, where the peak of the smoothed, background-
corrected PS is $\rm\nu_{max}$. An ACF is computed using the region in the
power spectrum centered on $\rm\nu_{max}$ and the peak in the ACF closest to
the expected value for the large frequency separation is adopted as
$\Delta\nu$.

The pySYD software depends on a number of powerful libraries, including
Astropy (Astropy Collaboration et al., 2018, 2013), Matplotlib (Hunter, 2007),
Numpy (Harris et al., 2020), and SciPy (Virtanen et al., 2020). pySYD has been
tested against SYD using results from the _Kepler_ sample for $\sim 100$ stars
(Figure 1). The comparisons show no significant systematic differences, with a
median offset and scatter of $0.2\%$ and $0.4\%$ for $\rm\nu_{max}$ as well as
$0.002\%$ and $0.09\%$ for $\Delta\nu$, which is smaller or comparable to the
typical random uncertainties (Huber et al., 2011).

Figure 1: Comparison of global parameters $\rm\nu_{max}$ (left) and
$\Delta\nu$ (right) measured by pySYD and IDL SYD for $\sim 100$ _Kepler_
stars. The bottom panels show the fractional residuals.

## Documentation & Examples

The main documentation for the pySYD software is hosted at
https://pysyd.readthedocs.iohttps://doi.orgpysyd.readthedocs.io. pySYD
provides a convenient setup feature that will download data for three example
stars and automatically create the recommended data structure for an easy
quickstart. The features of the pySYD output results are described in detail
in the documentation.

## Acknowledgements

We thank Dennis Stello, Jie Yu, Marc Hon, Yifan Chen, Yaguang Li, and other
pySYD users for discussion and suggestions which helped with the development
of this code.

We also acknowledge support from:

  * •

The National Science Foundation (DGE-1842402, AST-1717000)

  * •

The National Aeronautics and Space Administration (80NSSC19K0597)

  * •

The Alfred P. Sloan Foundation

## References

Aerts, C. (2021). Probing the interior physics of stars through
asteroseismology. _Reviews of Modern Physics_ , _93_(1), 015001\.
https://doi.org/10.1103/RevModPhys.93.015001

Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., Günther, H. M.,
Lim, P. L., Crawford, S. M., Conseil, S., Shupe, D. L., Craig, M. W.,
Dencheva, N., Ginsburg, A., VanderPlas, J. T., Bradley, L. D., Pérez-Suárez,
D., de Val-Borro, M., Aldcroft, T. L., Cruz, K. L., Robitaille, T. P.,
Tollerud, E. J., … Astropy Contributors. (2018). The Astropy Project: Building
an Open-science Project and Status of the v2.0 Core Package. _156_(3), 123.
https://doi.org/10.3847/1538-3881/aabc4f

Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., Greenfield, P.,
Droettboom, M., Bray, E., Aldcroft, T., Davis, M., Ginsburg, A., Price-Whelan,
A. M., Kerzendorf, W. E., Conley, A., Crighton, N., Barbary, K., Muna, D.,
Ferguson, H., Grollier, F., Parikh, M. M., Nair, P. H., … Streicher, O.
(2013). Astropy: A community Python package for astronomy. _558_ , A33.
https://doi.org/10.1051/0004-6361/201322068

Bastien, F. A., Stassun, K. G., Basri, G., & Pepper, J. (2013). An
observational correlation between stellar brightness variations and surface
gravity. _500_(7463), 427–430. https://doi.org/10.1038/nature12419

Bedding, T. R. (2014). Solar-like oscillations: An observational perspective.
In P. L. Pallé & C. Esteban (Eds.), _Asteroseismology_ (p. 60).

Borucki, W. J., Koch, D., Basri, G., Batalha, N., Brown, T., Caldwell, D.,
Caldwell, J., Christensen-Dalsgaard, J., Cochran, W. D., DeVore, E., Dunham,
E. W., Dupree, A. K., Gautier, T. N., Geary, J. C., Gilliland, R., Gould, A.,
Howell, S. B., Jenkins, J. M., Kondo, Y., … Prsa, A. (2010). Kepler Planet-
Detection Mission: Introduction and First Results. _Science_ , _327_ , 977.
https://doi.org/10.1126/science.1185402

Chaplin, W. J., Basu, S., Huber, D., Serenelli, A., Casagrande, L., Silva
Aguirre, V., Ball, W. H., Creevey, O. L., Gizon, L., Handberg, R., Karoff, C.,
Lutz, R., Marques, J. P., Miglio, A., Stello, D., Suran, M. D., Pricopi, D.,
Metcalfe, T. S., Monteiro, M. J. P. F. G., … Salabert, D. (2014).
Asteroseismic Fundamental Properties of Solar-type Stars Observed by the NASA
Kepler Mission. _210_ , 1. https://doi.org/10.1088/0067-0049/210/1/1

Corsaro, E., & De Ridder, J. (2014). DIAMONDS: A new Bayesian nested sampling
tool. Application to peak bagging of solar-like oscillations. _571_ , A71.
https://doi.org/10.1051/0004-6361/201424181

Harris, C. R., Millman, K. J., Walt, S. J. van der, Gommers, R., Virtanen, P.,
Cournapeau, D., Wieser, E., Taylor, J., Berg, S., Smith, N. J., Kern, R.,
Picus, M., Hoyer, S., Kerkwijk, M. H. van, Brett, M., Haldane, A., Río, J. F.
del, Wiebe, M., Peterson, P., … Oliphant, T. E. (2020). Array programming with
NumPy. _Nature_ , _585_(7825), 357–362.
https://doi.org/10.1038/s41586-020-2649-2

Hekker, S., Broomhall, A.-M., Chaplin, W. J., Elsworth, Y. P., Fletcher, S.
T., New, R., Arentoft, T., Quirion, P.-O., & Kjeldsen, H. (2010). The Octave
(Birmingham-Sheffield Hallam) automated pipeline for extracting oscillation
parameters of solar-like main-sequence stars. _402_(3), 2049–2059.
https://doi.org/10.1111/j.1365-2966.2009.16030.x

Hekker, S., Elsworth, Y., De Ridder, J., Mosser, B., García, R. A., Kallinger,
T., Mathur, S., Huber, D., Buzasi, D. L., Preston, H. L., Hale, S. J., Ballot,
J., Chaplin, W. J., Régulo, C., Bedding, T. R., Stello, D., Borucki, W. J.,
Koch, D. G., Jenkins, J., … Christensen-Dalsgaard, J. (2011). Solar-like
oscillations in red giants observed with Kepler: comparison of global
oscillation parameters from different methods. _525_ , A131.
https://doi.org/10.1051/0004-6361/201015185

Huber, D., Bedding, T. R., Stello, D., Hekker, S., Mathur, S., Mosser, B.,
Verner, G. A., Bonanno, A., Buzasi, D. L., Campante, T. L., Elsworth, Y. P.,
Hale, S. J., Kallinger, T., Silva Aguirre, V., Chaplin, W. J., De Ridder, J.,
García, R. A., Appourchaux, T., Frandsen, S., … Smith, J. C. (2011). Testing
Scaling Relations for Solar-like Oscillations from the Main Sequence to Red
Giants Using Kepler Data. _743_ , 143.
https://doi.org/10.1088/0004-637X/743/2/143

Huber, D., Stello, D., Bedding, T. R., Chaplin, W. J., Arentoft, T., Quirion,
P.-O., & Kjeldsen, H. (2009). Automated extraction of oscillation parameters
for Kepler observations of solar-type stars. _Communications in
Asteroseismology_ , _160_ , 74. http://arxiv.org/abs/0910.2764

Hunter, J. D. (2007). Matplotlib: A 2D graphics environment. _Computing in
Science & Engineering_, _9_(3), 90–95. https://doi.org/10.1109/MCSE.2007.55

Kjeldsen, H., & Bedding, T. R. (1995). Amplitudes of stellar oscillations: the
implications for asteroseismology. _293_ , 87–106. http://arxiv.org/abs/astro-
ph/9403015

Lightkurve Collaboration, Cardoso, J. V. de M., Hedges, C., Gully-Santiago,
M., Saunders, N., Cody, A. M., Barclay, T., Hall, O., Sagear, S., Turtelboom,
E., Zhang, J., Tzanidakis, A., Mighell, K., Coughlin, J., Bell, K., Berta-
Thompson, Z., Williams, P., Dotson, J., & Barentsen, G. (2018). _Lightkurve:
Kepler and TESS time series analysis in Python_ (p. ascl:1812.013).

Mathur, S., García, R. A., Régulo, C., Creevey, O. L., Ballot, J., Salabert,
D., Arentoft, T., Quirion, P.-O., Chaplin, W. J., & Kjeldsen, H. (2010).
Determining global parameters of the oscillations of solar-like stars. _511_ ,
A46. https://doi.org/10.1051/0004-6361/200913266

Mosser, B., & Appourchaux, T. (2009). On detecting the large separation in the
autocorrelation of stellar oscillation times series. _508_(2), 877–887.
https://doi.org/10.1051/0004-6361/200912944

Nielsen, M. B., Davies, G. R., Ball, W. H., Lyttle, A. J., Li, T., Hall, O.
J., Chaplin, W. J., Gaulme, P., Carboneau, L., Ong, J. M. J., García, R. A.,
Mosser, B., Roxburgh, I. W., Corsaro, E., Benomar, O., Moya, A., & Lund, M. N.
(2021). PBjam: A Python Package for Automating Asteroseismology of Solar-like
Oscillators. _161_(2), 62. https://doi.org/10.3847/1538-3881/abcd39

Ricker, G. R., Winn, J. N., Vanderspek, R., Latham, D. W., Bakos, G. Á., Bean,
J. L., Berta-Thompson, Z. K., Brown, T. M., Buchhave, L., Butler, N. R.,
Butler, R. P., Chaplin, W. J., Charbonneau, D., Christensen-Dalsgaard, J.,
Clampin, M., Deming, D., Doty, J., De Lee, N., Dressing, C., … Villasenor, J.
(2015). Transiting Exoplanet Survey Satellite (TESS). _Journal of Astronomical
Telescopes, Instruments, and Systems_ , _1_(1), 014003.
https://doi.org/10.1117/1.JATIS.1.1.014003

Schwarz, G. (1978). Estimating the Dimension of a Model. _The Annals of
Statistics_ , _6_(2), 461–464. https://doi.org/10.1214/aos/1176344136

Serenelli, A., Johnson, J., Huber, D., Pinsonneault, M., Ball, W. H., Tayar,
J., Silva Aguirre, V., Basu, S., Troup, N., Hekker, S., Kallinger, T., Stello,
D., Davies, G. R., Lund, M. N., Mathur, S., Mosser, B., Stassun, K. G.,
Chaplin, W. J., Elsworth, Y., … Zamora, O. (2017). The First APOKASC Catalog
of Kepler Dwarf and Subgiant Stars. _233_(2), 23.
https://doi.org/10.3847/1538-4365/aa97df

Verner, G. A., Elsworth, Y., Chaplin, W. J., Campante, T. L., Corsaro, E.,
Gaulme, P., Hekker, S., Huber, D., Karoff, C., Mathur, S., Mosser, B.,
Appourchaux, T., Ballot, J., Bedding, T. R., Bonanno, A., Broomhall, A.-M.,
García, R. A., Handberg, R., New, R., … Fanelli, M. N. (2011). Global
asteroseismic properties of solar-like oscillations observed by Kepler: a
comparison of complementary analysis methods. _415_(4), 3539–3551.
https://doi.org/10.1111/j.1365-2966.2011.18968.x

Virtanen, P., Gommers, R., Oliphant, T. E., Haberland, M., Reddy, T.,
Cournapeau, D., Burovski, E., Peterson, P., Weckesser, W., Bright, J., van der
Walt, S. J., Brett, M., Wilson, J., Millman, K. J., Mayorov, N., Nelson, A. R.
J., Jones, E., Kern, R., Larson, E., … SciPy 1.0 Contributors. (2020). SciPy
1.0: Fundamental Algorithms for Scientific Computing in Python. _Nature
Methods_ , _17_ , 261–272. https://doi.org/10.1038/s41592-019-0686-2

Yu, J., Huber, D., Bedding, T. R., & Stello, D. (2018). Predicting radial-
velocity jitter induced by stellar oscillations based on Kepler data. _480_ ,
L48–L53. https://doi.org/10.1093/mnrasl/sly123


# Experimental Evaluation of Multi-operator RIS-assisted Links in Indoor
Environment

Mir Lodro1, Jean Baptiste Gros2, Steve Greedy1, Geoffroy Lerosey2, Anas Al
Rawi3, Gabriele Gradoni1,

Mir Lodro, Steve Greedy and Gabriele Gradoni are with George Green Institute
for Electromagnetics Research-GGIEMR, the University of Nottingham, UK. Jean
Baptiste Gros and Geoffroy Lerosey are with Greenerwave, Paris 1University of
Nottingham, NG7 2RD, United Kingdom  2Greenerwave, 75005 Paris, France
3OFCOM, London, SE1 9HA, United Kingdom

###### Abstract

In this work, we present reconfigurable intelligent surface (RIS)-assisted
optimization of the multiple links in the same indoor environment. Multiple
RISs from different operators can co-exists and handle independent robust
communication links in the same indoor environment. We investigated the key
performance metrics with the help of two simultaneously operating RIS-
empowered robust communication links at different center frequencies in the
same indoor environment. We found with the help of bit error rate (BER) and
error vector magnitude (EVM) measurements that two operators can co-exist in
the same RF environment without seriously impacting quality of service of
users.

###### Index Terms:

Co-channel interference, Reconfigurable intelligent Surface, SDR, USRP.

##  I Introduction

Reconfigurable intelligent surface is a promising future technology that
consists of a large number of sub-wavelength and low-cost passive elements
whose phase and amplitude can be optimized to increase signal coverage,
enhance spectral-efficiency [1], perform accurate RF sensing [2] and increase
PHY layer security [3]. Each element in the RIS can introduce phase shift to
impinging electromagnetic waves and the collective phase shifts introduced by
all the elements can perform passive beamforming to the desired user. Coverage
extension [4] and localization are the two important benefits of RIS that
shall pave the way for sixth generation (6G) wireless communication [5]. In
order to demonstrate RIS-assisted communication system different experimental
setup and prototypes have been created at sub-6GHz and at mmWaves. For
example, authors in [6] have proposed a 256 2-bit element RIS and presented
experimental results at 2.3 GHz and 28.5 GHz. They have performed RIS-assisted
measurements using single Tx-Rx link by transmitting OFDM modulated waveforms
over-the-air in an indoor environment. Several other works have provided their
experimental setups designed for a single-carrier frequency at sub-6 GHz and
at mmWaves [7],[8][9][10][11]. Authors in [10] have performed single-carrier
RIS-assisted over-the-air test in both indoor and outdoor environment at
operating frequency of 5.8 GHz. Similarly, authors in [11] have performed RIS-
assisted measurements in an indoor environment and performed KPI optimization.
Recent experimental research on RIS-assisted communication systems remains
focused on passive beamforming and KPI optimization of single Tx-Rx links.
However, the question of understanding the interference control among multiple
RIS-assisted links in the same environment remains outstanding. We use an
experimental approach to address the issue of interference among co-existing
multiple operator communication links in the same indoor environment. In the
same RF environment both the RISs from two operators are exposed to
transmitted EM waves from both transmitters at the same time. In this setup
where an independent phase optimization of one link may degrade the quality of
another communication link. It’s highly important to experimentally understand
the effect of multiple RIS-assisted links in the same RF environment.

We have divided our work into six sections. After the introduction in section
I, we have discussed indoor RF environment in section II, where RIS-assisted
measurements are conducted. Section III discusses communication system KPIs
that have been optimized. Section IV gives detailed discussion of the
experimental setup and the measurement procedure. Section V is dedicated to
results and discussion section where the details about the two configurations
of RIS are presented. Section VI is about the lessons learned from the
measurement campaign, the conclusion and future aspects of the research.

##  II Indoor Environment

The RIS-empowered two link optimization was performed in an indoor office
environment. The floor plan of the RF environment is shown in Fig. 1. The
indoor environment is an office environment in a multi-story building and the
measurement setup is created in an office room. The office room consists of
wooden furniture including desks, chairs, small cabinets and IT equipment such
as monitors and the host PCs. The indoor environment for sub-6GHz is a rich-
scattering environment where the two Tx-Rx links are assisted using RISs. Each
communication link is assisted independently by a separate RIS. There exist
direct link between transmitter and receiver in addition to cascaded channel
Tx-RIS-Rx. The direct link is weak enough to be considered negligible and
hence the KPI optimization is performed only using indirect i.e. Tx-RIS-Rx
channel.

Figure 1: Indoor environment floor plan.

##  III EVM and BER Measurement

EVM and BER are two key performance indicators (KPIs) in digital wireless
communication systems allowing greater understanding of modulation quality and
the end-to-end communication system performance metric. In baseband model, the
EVM can be measured at different points in the digital receiver before frame
synchronization and data recovery stage. The modulation techniques with
constant magnitude such as QPSK, the EVMs are normalized to constellation
maximum and those with multiple magnitude levels the EVM is normalized to RMS
level of I/Q reference symbols [12]. Bit error rate (BER) is an end-to-end
communication system performance metric. The BER is measured after all the RF
impairments have been corrected. Additionally, the BER is measured after the
frame equalization, phase ambiguity and data recovery steps have been
correctly applied. BER is an ultimate performance indicator which guarantees
all the transmitter steps have been reversed at the receiver including RF
impairment compensation.

$BER=\frac{\textrm{number of received bits in error}}{\textrm{total number of
transmitted bits}}$ (1)

##  IV Experimental Setup

In this work, we used two full-duplex X310 USRPs from Ettus Research. Each of
the X310 USRP was equipped with two UBX-160 RF cards that can support 2x2 MIMO
links. We used separate USRPs for the transmitter and the receiver in order to
avoid inter-board cross-talk. The two RISs and the Tx USRPs and the Rx USRPs
were connected to the host PC. We transmitted QPSK modulated waveforms from
each channel of the Tx USRP independently at different center frequencies in
the indoor environment. The two links operate at adjacent center frequencies
and are a single-input single-output (SISO) links. The complex samples were
captured by another Rx USRP whose Rx channels were tuned to the corresponding
channels of the Tx USRPs. The transmitted waveforms and the received waveforms
are generated and processed in a separate MATLAB session. The complex samples
are transferred to and from the host PC at a sample rate of 400 kS/s. Hence,
the maximum data-transfer rate is 12.8 Mbps for the sample size of 16-bit for
I and Q samples. The captured QPSK waveforms are down-converted and processed
in the baseband model where the baseband waveforms are compensated for RF
impairments such as coarse frequency compensation, fine-frequency compensation
and the timing synchronization. The measured waveforms and the EVMs of two
links are recorded after all the RF impairments have been rectified. We
recorded our previous experience with current baseband model in [13]. For
customizing the radio channel we used two 5.2 GHz 1-bit RISs patented by
Greenerwave, Paris that resonates at 5.2 GHz. The RISs have been used in
previous experiments on link-optimization and changing the boundary conditions
in a chaotic metal-enclosure exhibiting rich-scattering [14][15]. Fig. 2 shows
experimental setup of two separated RIS-assisted links measurement setup. Each
metasurface consists of 76 binary pixels whose reflection co-efficient can be
changed by controlling the biase voltage of PIN diodes. Binary phase shift is
applied to RIS elements using RIS controller. The RIS controller connected to
the host PC using USB2.0 cable and it receives PIN diode mapping instructions
from a MATLAB program. The two RIS-assisted links use monopoles for the Tx and
Rx links such that the two RISs are exposed to transmitted signals from the
cross-links.

Figure 2: Photograph of the experimental setup of the two separated RISs in
relation to Tx1 and Tx2.

##  V Discussion and Results

We conducted EVM and BER optimization by optimizing the configurations of two
RISs where each RIS supporting separate links at designated center frequency.
The RISs are designed to operate at center frequency of 5.2 GHz. Multi-
operator link optimization is performed by optimizing configurations of each
RIS supporting independent link at adjacent frequencies. Two scenarios at
different center frequencies are investigated. Co-located RISs when the two
RISs are placed side by side and they are at the same distance from the
transmitters. In separated RISs scenario the two RISs are placed along the
axes that are orthogonal to each other.

###  V-A Co-located RISs

We considered two co-located RISs supporting links at the center frequencies
of 5.1 GHz and 5.2 GHz respectively. In co-located scenario we placed two RIS
side by side and established two independent communication links. The Tx side
of both the links were placed closely to co-located RIS. The receivers of both
the links were also co-located and the receivers were placed at a distance of
4 meters from the RISs. Each RISs assisted an independent link operating at
different carrier frequencies with a minimum gap of 50 MHz. For co-located
scenario, we performed KPI optimization at 5.1 GHz and 5.2 GHz respectively.
We started with poor channel condition for both the links. It was visually
observed from the constellation diagrams at the starting position the received
QPSK symbols were centered around zero on the IQ plane as shown in Fig.3a.
After RIS optimization the QPSK symbols were received in their corresponding
quadrant on the IQ plane (see Fig. 3b).

(a)

(b)

(c)

(d)

Figure 3: Typical QPSK constellations (a) typical initial QPSK constellation
at 5.1 GHz (b) typical optimized QPSK constellation at 5.1 GHz (c) typical
initial QPSK constellation at 5.2 GHz (d) typical optimized QPSK constellation
at 5.2 GHz.

Fig. 4 shows the instantaneous RMS EVM and the optimized RMS EVM at 5.1 GHz
and 5.2 GHz links respectively in a co-located RISs scenario. Instantaneous
EVM of the two links decreases and the RMS EVM for 5.2 GHz link is
consistently less than the 5.1 GHz link. Fig. 5 shows instantaneous BER and
the optimized BER of two RIS-assisted links in the co-located RISs scenario.
Fig. 5a shows the instantaneous BER of 5.2 GHz link is consistently lower than
the 5.1 GHz link because of the RIS peak response.

(a)

(b)

Figure 4: Co-located RISs scenario: QPSK EVM at 5.1 GHz and 5.2 GHz
respectively. The center frequency difference is 100 MHz (a) instantaneous EVM
(b) optimized EVM.

(a)

(b)

Figure 5: Co-located RISs scenario: QPSK BER at 5.1 GHz and 5.2 GHz
respectively. The center frequency difference is 100 MHz (a) instantaneous BER
(b) optimized BER.

###  V-B Separated RISs

We considered separated RISs placed along x-axis and y-axis respectively. We
repeated the measurements for several frequencies centered around RISs
operating frequencies. We observed RISs optimization was more effective at
frequencies closer to its operating frequencies. For the fixed USRP transmit
gains we observed low-edge frequencies with better RF propagation specially
when the two links were separately operating with separation of 100 MHz i.e.
at 5.1 GHz and 5.2 GHz respectively. We further reduced the gap of two RF
frequencies to 50 MHz. The two links were centered around 5.2 GHz and we
observed almost same RF effects. We visually observed this from the QPSK
constellation diagram. The EVM optimization of two RISs links closely followed
each other. Fig. 7 shows instantaneous BER and the optimized BER of two RIS-
assisted links in the separated RISs scenario. Similar to Fig. 5a, Fig. 7a
shows the instantaneous BER of 5.2 GHz link is consistently lower than the 5.3
GHz link because of RIS peak response.

(a)

(b)

Figure 6: Separated RISs scenario: QPSK EVM at 5.2 GHz and 5.3 GHz
respectively. The center frequency difference is 100 MHz (a) instantaneous EVM
(b) optimized EVM.

(a)

(b)

Figure 7: Separated RISs scenario: QPSK BER at 5.2 GHz and 5.3 GHz
respectively. The center frequency difference is 100 MHz (a) instantaneous BER
(b) optimized BER.

We also measured QPSK EVM and BER in separated RISs scenario by changing
center frequency of the SDR to 5.175 GHz and 5.225 GHz. The center frequency
difference of 50 MHz. The separated RISs scenario optimization of QPSK EVM and
QPSK BER was performed. We observed the similar pattern as in previous cases.
We observed BER only at initial RISs iterations then no BERs were observed as
the link started with favourable propagation but the EVM optimization followed
same repetitive patterns. The measured QPSK EVM and the optimized QPSK EVM in
separated RIS scenario is shown in Fig.8. Both links have been optimized.

(a)

(b)

Figure 8: Separated RISs QPSK EVM at 5.175 GHz and 5.225 GHz respectively. The
center frequency difference is 50 MHz (a) instantaneous EVM (b) optimized EVM.

From Fig. 4a, Fig. 6a and Fig. 8a we understand that the RIS effectively
optimizes the links at 5.2 GHz or closer to 5.2 GHz. The instantaneous EVM is
higher for center frequency centered on or closer to 5.2 GHz independent of
whether the two links started with poor or favourable propagation conditions.

##  VI Conclusion

An experimental study was carried out to understand multiple RIS-assisted
links in the same indoor environment. We conducted RIS-assisted multiple link
optimization in different scenarios using SDR based testbed. In all the
experiments RIS configurations were optimized to decrease the EVM and the BER
of each communication link. From our measurements experience, we understand
that multi-operator RIS-assisted links can be sustained in the same
environment without affecting or degrading the QoS of served users.

## Acknowledgment

This work has been supported by the European Commission through the H2020
RISE-6G Project under Grant 101017011.

## References

  * [1] L. You, J. Xiong, D. W. K. Ng, C. Yuen, W. Wang, and X. Gao, “Energy efficiency and spectral efficiency tradeoff in ris-aided multiuser mimo uplink transmission,” _IEEE Transactions on Signal Processing_ , vol. 69, pp. 1407–1421, 2020.
  * [2] J. Hu, H. Zhang, B. Di, L. Li, K. Bian, L. Song, Y. Li, Z. Han, and H. V. Poor, “Reconfigurable intelligent surface based rf sensing: Design, optimization, and implementation,” _IEEE Journal on Selected Areas in Communications_ , vol. 38, no. 11, pp. 2700–2716, 2020.
  * [3] J. Chen, Y.-C. Liang, Y. Pei, and H. Guo, “Intelligent reflecting surface: A programmable wireless environment for physical layer security,” _IEEE Access_ , vol. 7, pp. 82 599–82 612, 2019.
  * [4] V. Popov, M. Odit, J.-B. Gros, V. Lenets, A. Kumagai, M. Fink, K. Enomoto, and G. Lerosey, “Experimental demonstration of a mmwave passive access point extender based on a binary reconfigurable intelligent surface,” _arXiv preprint arXiv:2107.02087_ , 2021.
  * [5] E. C. Strinati, G. C. Alexandropoulos, V. Sciancalepore, M. Di Renzo, H. Wymeersch, D.-T. Phan-Huy, M. Crozzoli, R. D’Errico, E. De Carvalho, P. Popovski, P. Di Lorenzo, L. Bastianelli, M. Belouar, J. E. Mascolo, G. Gradoni, S. Phang, G. Lerosey, and B. Denis, “Wireless environment as a service enabled by reconfigurable intelligent surfaces: The rise-6g perspective,” in _2021 Joint European Conference on Networks and Communications & 6G Summit (EuCNC/6G Summit)_, 2021, pp. 562–567.
  * [6] L. Dai, B. Wang, M. Wang, X. Yang, J. Tan, S. Bi, S. Xu, F. Yang, Z. Chen, M. Di Renzo _et al._ , “Reconfigurable intelligent surface-based wireless communications: Antenna design, prototyping, and experimental results,” _IEEE Access_ , vol. 8, pp. 45 913–45 923, 2020.
  * [7] J.-B. Gros, V. Popov, M. A. Odit, V. Lenets, and G. Lerosey, “A reconfigurable intelligent surface at mmwave based on a binary phase tunable metasurface,” _IEEE Open Journal of the Communications Society_ , vol. 2, pp. 1055–1064, 2021.
  * [8] G. C. Alexandropoulos, N. Shlezinger, and P. Del Hougne, “Reconfigurable intelligent surfaces for rich scattering wireless communications: Recent experiments, challenges, and opportunities,” _IEEE Communications Magazine_ , vol. 59, no. 6, pp. 28–34, 2021.
  * [9] R. Fara, D.-T. Phan-Huy, P. Ratajczak, A. Ourir, M. Di Renzo, and J. De Rosny, “Reconfigurable intelligent surface-assisted ambient backscatter communications–experimental assessment,” in _2021 IEEE International Conference on Communications Workshops (ICC Workshops)_. IEEE, 2021, pp. 1–7.
  * [10] X. Pei, H. Yin, L. Tan, L. Cao, Z. Li, K. Wang, K. Zhang, and E. Björnson, “Ris-aided wireless communications: Prototyping, adaptive beamforming, and indoor/outdoor field trials,” _IEEE Transactions on Communications_ , vol. 69, no. 12, pp. 8627–8640, 2021.
  * [11] M. M. Amri, N. M. Tran, and K. W. Choi, “Reconfigurable intelligent surface-aided wireless communications: Adaptive beamforming and experimental validations,” _IEEE Access_ , vol. 9, pp. 147 442–147 457, 2021.
  * [12] H. Chen, W. Lin, S. Wang, W. Che, and Q. Xue, “A calibrated over-the-air measurement method for error vector magnitude characterization,” _IEEE Transactions on Instrumentation and Measurement_ , vol. 71, pp. 1–8, 2022.
  * [13] M. Lodro, C. Smart, G. Gradoni, A. Vukovic, D. Thomas, and S. Greedy, “Near-field ber and evm measurement at 5.8 ghz in mode-stirred metal enclosure,” _The Applied Computational Electromagnetics Society Journal (ACES)_ , pp. 1080–1088, 2020.
  * [14] M. Lodro, G. Gradoni, J.-B. Gros, S. Greedy, and G. Lerosey, “Reconfigurable intelligent surface-assisted bluetooth low energy link in metal enclosure,” _Frontiers in Communications and Networks_ , p. 44, 2021.
  * [15] J.-B. Gros, P. del Hougne, and G. Lerosey, “Tuning a regular cavity to wave chaos with metasurface-reconfigurable walls,” _Physical Review A_ , vol. 101, no. 6, p. 061801, 2020.


# Damping law of photocount distribution in a dissipative channel††thanks:
This work was supported by the National Natural Science Foundation of China
(Grant Nos. 11175113 and 11264018), and the Young Talents Foundation of
Jiangxi Normal University.

Hong-Yi Fan1, Sen-Yue Lou1 and Li-Yun Hu2† 1Department of Physics, Ningbo
University, Ningbo 315211, P. R. China
2Department of physics, Jiangxi Normal University, Nanchang 330022, P. R.
China
${\dagger}$Email<EMAIL_ADDRESS>
###### Abstract

For a dissipative channel governed by the master equation of density operator
$d\rho/dt=\kappa\left(2a\rho a^{\dagger}-a^{\dagger}a\rho-\rho
a^{\dagger}a\right),$ we find that photocount distribution formula at time
$t,$ $p\left(n,t\right)=Tr\left\\{\rho\left(t\right)\mathbf{\colon}\left(\xi
a^{\dagger}a\right)^{n}e^{-\xi a^{\dagger}a}/n!\colon\right\\},$ becomes
$p\left(n,t\right)=Tr\left[\rho\left(0\right)\mathbf{\colon}\left(\xi
e^{-2\kappa t}a^{\dagger}a\right)^{n}e^{-\xi e^{-2\kappa
t}a^{\dagger}a}/n!\colon\right],$ as if the quantum efficiency $\xi$ of the
detector becomes $\xi e^{-2\kappa t}$. This law greatly simplifies the
theoretical study of photocount distribution for quantum optical field.

In nature, syetems we concerned usually are surrounded by thermo reservoir.
Decoherence is an important topic in the fields of quantum optics, quantum
computing and quantum information processing. Decoherence of a quantum optical
field can be judged by photocounting its final state comparing with its
initial state. In this paper we discuss time-dissipation of photocount
distribution in a dissipative reservoir. We hope to know photon counting law
for $\rho(t)$, when the initial density operator $\rho_{0}=\rho(t=0)$ evolves
with time in a photon-loss channel. The corresponding master equation is r1

$\frac{d\rho}{dt}=\kappa\left(2a\rho a^{\dagger}-a^{\dagger}a\rho-\rho
a^{\dagger}a\right),$ (1)

where $\kappa$ is the rate of dissipation, $\left[a,a^{\dagger}\right]=1$. We
shall deduce a new damping formula of photon counting after an optical quantum
field passing through a photon-loss channel.

_Converting_ $\rho$_’s master equations into pure state’s evolution
equations._ In order to solve equation Eq.(1), instead of converting the
operator equation into its $c$-number equation by the Wigner function approach
or by P-representation approach, we introduce thermo entangled state
representation r2

$\left|\tau\right\rangle=\exp\left(-\frac{1}{2}|\tau|^{2}+\tau
a^{\dagger}-\tau^{\ast}\tilde{a}^{\dagger}+a^{\dagger}\tilde{a}^{\dagger}\right)\left|0\tilde{0}\right\rangle,$
(2)

which is complete, to make up the pure state $\left|\rho\right\rangle$,

$\left|\rho\right\rangle\equiv\rho\left|I\right\rangle,\text{
}\left|I\right\rangle\equiv\left|\tau=0\right\rangle=e^{a^{\dagger}\tilde{a}^{\dagger}}\left|0\tilde{0}\right\rangle=\sum_{n=0}^{\infty}\left|n\tilde{n}\right\rangle,$
(3)

where $\tilde{a}^{\dagger}$ is a fictitious mode accompanying the real photon
creation operator $a^{\dagger},$ $\left|\tilde{0}\right\rangle$ is annihilated
by $\tilde{a},$ $\left[\tilde{a},\tilde{a}^{\dagger}\right]=1,$
$\left|n\tilde{n}\right\rangle=\left(a^{\dagger}\tilde{a}^{\dagger}\right)^{n}\left|0\tilde{0}\right\rangle/n!,\
$ $\left[a,\tilde{a}^{\dagger}\right]=0$. From Eq.(3) it is easy to see that

$\displaystyle a\left|I\right\rangle$ $\displaystyle=$
$\displaystyle\tilde{a}^{\dagger}\left|I\right\rangle,$ $\displaystyle
a^{\dagger}\left|I\right\rangle$ $\displaystyle=$
$\displaystyle\tilde{a}\left|I\right\rangle,$
$\displaystyle\left(a^{\dagger}a\right)^{n}\left|I\right\rangle$
$\displaystyle=$
$\displaystyle\left(\tilde{a}^{\dagger}\tilde{a}\right)^{n}\left|I\right\rangle,$
(4)

which indicates that the transformation relation between the real mode and the
fictitious mode. Acting the both sides of Eq.(1) on the state
$\left|I\right\rangle,$ we have

$\frac{d}{dt}\left|\rho\right\rangle=\kappa\left(2a\rho
a^{\dagger}-a^{\dagger}a\rho-\rho a^{\dagger}a\right)\left|I\right\rangle.$
(5)

Using Eq.(4) and noting that $\rho$ is defined in the real space which is
commutative with the fictitious operators, we can convert Eq.(5) into

$\frac{d}{dt}\left|\rho\right\rangle=\kappa\left(2a\tilde{a}-a^{\dagger}a-\tilde{a}^{\dagger}\tilde{a}\right)\left|\rho\right\rangle.$
(6)

Its formal solution is given by

$\left|\rho\right\rangle=\exp\\{\kappa
t\left(2a\tilde{a}-a^{\dagger}a-\tilde{a}^{\dagger}\tilde{a}\right)\\}\left|\rho_{0}\right\rangle,$
(7)

where $\left|\rho_{0}\right\rangle\equiv\rho_{0}\left|I\right\rangle$. Thus by
introducing a fictitious mode $\tilde{a}^{\dagger}$ to be a counterpart of the
system mode, one can convert density operators’ master equations into
evolution equations of pure states, which leads to the formal solution (7).

On the other hand, noticing that the commutative relation
$\left[(a^{\dagger}a+\tilde{a}^{\dagger}\tilde{a})/2,a\tilde{a}\right]=-a\tilde{a},$
thus we can use the operator identity

$e^{\lambda\left(A+\sigma B\right)}=e^{\lambda
A}\exp\left[\sigma\left(1-e^{-\lambda\tau}\right)B/\tau\right],$ (8)

which is valid for $\left[A,B\right]=\tau B$, to decompose the exponential
operator in Eq.(7) as r3

$e^{-2\kappa
t(\frac{a^{\dagger}a+\tilde{a}^{\dagger}\tilde{a}}{2}-a\tilde{a})}=e^{-\kappa
t\left(a^{\dagger}a+\tilde{a}^{\dagger}\tilde{a}\right)}e^{Ta\tilde{a}},$ (9)

where $T=1-e^{-2\kappa t}.$ On substituting Eq.(9) into Eq.(7) yields

$\displaystyle\left|\rho\left(t\right)\right\rangle$ $\displaystyle=$
$\displaystyle\exp\left[-\kappa
t\left(a^{\dagger}a+\tilde{a}^{\dagger}\tilde{a}\right)\right]\exp\left[Ta\tilde{a}\right]\left|\rho_{0}\right\rangle$
(10) $\displaystyle=$ $\displaystyle\exp\left[-\kappa
t\left(a^{\dagger}a+\tilde{a}^{\dagger}\tilde{a}\right)\right]\sum_{n=0}^{\infty}\frac{T^{n}}{n!}a^{n}\tilde{a}^{n}\rho_{0}\left|I\right\rangle.$

Noticing that the system mode is commutative with the fictitious operators,
and the transition relation shown in Eq.(4), we see that

$\tilde{a}^{n}\rho_{0}\left|I\right\rangle=\rho_{0}\tilde{a}^{n}\left|I\right\rangle=\rho_{0}a^{\dagger
n}\left|I\right\rangle,$ (11)

and

$\displaystyle\exp\left[-\kappa
t\tilde{a}^{\dagger}\tilde{a}\right]a^{n}\rho_{0}a^{\dagger
n}\left|I\right\rangle$ (12) $\displaystyle=$ $\displaystyle
a^{n}\rho_{0}a^{\dagger n}\exp\left[-\kappa
t\tilde{a}^{\dagger}\tilde{a}\right]\left|I\right\rangle$ $\displaystyle=$
$\displaystyle a^{n}\rho_{0}a^{\dagger n}\exp\left[-\kappa
ta^{\dagger}a\right]\left|I\right\rangle,$

thus we have

$\displaystyle\rho\left(t\right)\left|I\right\rangle$ $\displaystyle=$
$\displaystyle e^{-\kappa
ta^{\dagger}a}\sum_{n=0}^{\infty}\frac{T^{n}}{n!}a^{n}\rho_{0}a^{\dagger
n}e^{-\kappa t\tilde{a}^{\dagger}\tilde{a}}\left|I\right\rangle$ (13)
$\displaystyle=$ $\displaystyle\sum_{n=0}^{\infty}\frac{T^{n}}{n!}e^{-\kappa
ta^{\dagger}a}a^{n}\rho_{0}a^{\dagger n}e^{-\kappa
ta^{\dagger}a}\left|I\right\rangle.$

Comparing with the standard form of the operator-sum representation r4 ,

$\rho\left(t\right)=\sum_{m=0}^{\infty}M_{m}\rho_{0}M_{m}^{\dagger}$ (14)

we see that the Kraus operator for this damping process is given by

$M_{m}=\sqrt{\frac{T^{m}}{m!}}e^{-\kappa ta^{\dagger}a}a^{m},$ (15)

and we can prove the unitrity $\sum_{m=0}^{\infty}M_{m}^{\dagger}M_{m}=1$.

_Damping formula of Photocount distribution in a dissipative channel._ Now we
discuss how does the photocount distribution formula of optical fields evolve
in the photon-loss channel?

The quantum mechanical photon counting formula was first derived by Kelley and
Kleiner r5 . As shown in Ref. r5 ; r6 ; r7 ; r8 ; r9  the probability
distribution $\mathfrak{p}\left(n\right)$ of registering $n$ photoelectrons in
time interval $\Delta\tau$ is given by

$\mathfrak{p}\left(n\right)=\mathtt{Tr}\left\\{\rho\mathbf{\colon}\frac{\left(\xi
a^{\dagger}a\right)^{n}}{n!}e^{-\xi a^{\dagger}a}\colon\right\\},$ (16)

where $\mathbf{\colon}\colon$ denotes normal ordering, $\xi\propto\Delta\tau$
is called the quantum efficiency (a measure) of the detector,
$\mathbf{\rho}\left(0\right)$ is the single-mode density operator of the light
field concerned at initial time. According to the operator-sum representation
at time $t$ in Eq.(14), we have

$\displaystyle\mathfrak{p}\left(n,t\right)$ $\displaystyle=$
$\displaystyle\mathtt{Tr}\left\\{\rho\left(t\right)\mathbf{\colon}\frac{\left(\xi
a^{\dagger}a\right)^{n}}{n!}e^{-\xi a^{\dagger}a}\colon\right\\}$ (17)
$\displaystyle=$
$\displaystyle\mathtt{Tr}\left\\{\rho_{0}\hat{O}\left(t\right)\right\\},$

where we have set

$\hat{O}\left(t\right)=\sum_{m=0}^{\infty}M_{m}^{\dagger}\mathbf{\colon}\frac{\left(\xi
a^{\dagger}a\right)^{n}}{n!}e^{-\xi a^{\dagger}a}\colon M_{m}.$ (18)

Substituting Eq.(15) into Eq.(18) and using the following operator identities:
$\mathbf{\colon}e^{-\xi
a^{\dagger}a}\colon=\exp\left[a^{\dagger}a\ln\left(1-\xi\right)\right]=\left(1-\xi\right)^{a^{\dagger}a}$
and $e^{\lambda a^{\dagger}a}a^{\dagger}e^{-\lambda
a^{\dagger}a}=a^{\dagger}e^{\lambda}$, $e^{\lambda a^{\dagger}a}ae^{-\lambda
a^{\dagger}a}=ae^{-\lambda}$, we can put Eq.(18) into the following form

$\displaystyle\hat{O}\left(t\right)$ (19) $\displaystyle=$
$\displaystyle\sum_{m=0}^{\infty}\frac{T^{m}}{m!}a^{\dagger m}e^{-\kappa
ta^{\dagger}a}\mathbf{\colon}\frac{\left(\xi
a^{\dagger}a\right)^{n}}{n!}e^{-\xi a^{\dagger}a}\mathbf{\colon}e^{-\kappa
ta^{\dagger}a}a^{m}$ $\displaystyle=$
$\displaystyle\sum_{m=0}^{\infty}\frac{T^{m}}{m!}a^{\dagger m}e^{-\kappa
ta^{\dagger}a}\frac{\xi^{n}a^{\dagger
n}}{n!}\left(1-\xi\right)^{a^{\dagger}a}a^{n}e^{-\kappa ta^{\dagger}a}a^{m}$
$\displaystyle=$ $\displaystyle\frac{\xi^{n}}{n!}e^{-2n\kappa
t}\sum_{m=0}^{\infty}\frac{T^{m}}{m!}a^{\dagger n+m}e^{-\kappa
ta^{\dagger}a}e^{a^{\dagger}a\ln\left(1-\xi\right)}e^{-\kappa
ta^{\dagger}a}a^{n+m}$ $\displaystyle=$
$\displaystyle\frac{\xi^{n}}{n!}e^{-2\kappa
tn}\sum_{m=0}^{\infty}\frac{T^{m}}{m!}a^{\dagger
n+m}\mathbf{\colon}e^{\left[e^{\left[\ln\left(1-\xi\right)-2\kappa
t\right]}-1\right]a^{\dagger}a}\mathbf{\colon}a^{n+m}.$

Further using the property that $a$ and $a^{\dagger}$ are commute within the
normal ordering symbol $:$ $:$, we then remember $T=1-e^{-2\kappa t}$ to have

$\displaystyle\hat{O}\left(t\right)$ (20) $\displaystyle=$
$\displaystyle\frac{\xi^{n}}{n!}e^{-2\kappa tn}a^{\dagger
n}\colon\sum_{m=0}^{\infty}\frac{\left(Ta^{\dagger}a\right)^{m}}{m!}e^{\left[e^{\left[\ln\left(1-\xi\right)-2\kappa
t\right]}-1\right]a^{\dagger}a}\colon a^{n}$ $\displaystyle=$
$\displaystyle\frac{\left(\xi e^{-2\kappa t}\right)^{n}}{n!}a^{\dagger
n}\colon\exp\\{-\xi e^{-2\kappa t}a^{\dagger}a\\}\colon a^{n}$
$\displaystyle=$ $\displaystyle\colon\frac{\left(\xi e^{-2\kappa
t}a^{\dagger}a\right)^{n}}{n!}e^{-\xi e^{-2\kappa t}a^{\dagger}a}\colon.$

Substituting Eq.(20) into (17) yields

$\mathfrak{p}\left(n,t\right)=\mathtt{Tr}\left\\{\rho_{0}\mathbf{\colon}\frac{\left(\xi
e^{-2\kappa t}a^{\dagger}a\right)^{n}}{n!}e^{-\xi e^{-2\kappa
t}a^{\dagger}a}\colon\right\\}.$ (21)

Comparing Eq.(21) with Eq.(16) we obtain the law of damping formula of
photocount distribution in a dissipative channel, that is

$\displaystyle\mathfrak{p}\left(n,0\right)$ $\displaystyle=$
$\displaystyle\mathtt{Tr}\left\\{\rho\left(0\right)\mathbf{\colon}\frac{\left(\xi
a^{\dagger}a\right)^{n}}{n!}e^{-\xi a^{\dagger}a}\colon\right\\}$
$\displaystyle\left.\rightarrow\right.\mathfrak{p}\left(n,t\right)$
$\displaystyle=$
$\displaystyle\mathtt{Tr}\left\\{\rho_{0}\mathbf{\colon}\frac{\left(\xi
e^{-2\kappa t}a^{\dagger}a\right)^{n}}{n!}e^{-\xi e^{-2\kappa
t}a^{\dagger}a}\colon\right\\},$ (22)

namely, this formula exhibits the law as if the quantum efficiency $\xi$ of
the detector becomes $\xi e^{-2\kappa t}$. So from Eq.(21) one can see that
the photocount distribution is related to the decay of quantum efficiency
$\xi^{\prime}=\xi e^{-2\kappa t}$. This is our main result.

_Example._ Now we check the formula through an example. When the initial state
is a pure number state, $\rho\left(0\right)=\left|m\right\rangle\left\langle
m\right|,$ using the following operator identities $a^{\dagger
m}a^{m}=N\left(N-1\right)\cdots\left(N-m+1\right)$, where $N=a^{\dagger}a$,
then the photocount distribution is

$\displaystyle\mathfrak{p}\left(n,0\right)$ (23) $\displaystyle=$
$\displaystyle\left\langle m\right|\mathbf{\colon}\frac{\left(\xi
a^{\dagger}a\right)^{n}}{n!}e^{-\xi
a^{\dagger}a}\mathbf{\colon}\left|m\right\rangle$ $\displaystyle=$
$\displaystyle\frac{\xi^{n}}{n!}\sum_{l=0}^{\infty}\frac{\left(-\xi\right)^{l}}{l!}\left\langle
m\right|a^{\dagger n+l}a^{n+l}\left|m\right\rangle$ $\displaystyle=$
$\displaystyle\frac{\xi^{n}}{n!}\sum_{l=0}^{\infty}\frac{\left(-\xi\right)^{l}}{l!}\left\langle
m\right|N\left(N-1\right)\cdots\left(N-n-l+1\right)\left|m\right\rangle$
$\displaystyle=$
$\displaystyle\sum_{l=0}^{m-n}\frac{\left(-\xi\right)^{l}}{l!}\frac{\xi^{n}}{n!}\frac{m!}{\left(m-n-l\right)!}$
$\displaystyle=$
$\displaystyle\frac{m!}{n!\left(m-n\right)!}\xi^{n}\left(1-\xi\right)^{m-n},$

which indicates that the probability of a photon being counted during the
period $\Delta\tau$ is the quantum efficiency $\xi$, thus the probability of
counting $n$ out of $m$ photons is proportional to the probability of counting
$n$ photons $\xi^{n}$ times the probability of counting $m-n$ photons
$\left(1-\xi\right)^{m-n}$. Then according to the new formula Eq.(22) we
immediately have the probability distribution of optical filed
$\left|m\right\rangle$ after damping at time $t,$ i.e.,

$\displaystyle\mathfrak{p}\left(n,t\right)$ $\displaystyle=$
$\displaystyle\mathtt{Tr}\left\\{\rho_{0}\mathbf{\colon}\frac{\left(\xi
e^{-2\kappa t}a^{\dagger}a\right)^{n}}{n!}e^{-\xi e^{-2\kappa
t}a^{\dagger}a}\mathbf{\colon}\right\\}$ (24) $\displaystyle=$
$\displaystyle\frac{m!\left(\xi e^{-2\kappa
t}\right)^{n}}{n!\left(m-n\right)!}\left(1-\xi e^{-2\kappa t}\right)^{m-n}.$

Now we check it through another approach. We see that using (14)-(15) the
initial state $\left|m\right\rangle\left\langle m\right|$ evolves into a
binomial state,

$\left|m\right\rangle\left\langle m\right|\rightarrow e^{-2m\kappa
t}\sum_{l=0}^{m}\frac{m!\left(e^{2\kappa
t}-1\right)^{l}}{l!\left(m-l\right)!}\left|m-l\right\rangle\left\langle
m-l\right|,$ (25)

then using Eq.(16) and Eq.(23) we have

$\displaystyle\mathfrak{p}\left(n,t\right)$ (26) $\displaystyle=$
$\displaystyle\sum_{l=0}^{m}\frac{m!\left(e^{2\kappa
t}-1\right)^{l}}{l!\left(m-l\right)!e^{2m\kappa t}}\left\langle
m-l\right|\mathbf{\colon}\frac{\left(\xi a^{\dagger}a\right)^{n}}{n!}e^{-\xi
a^{\dagger}a}\colon\left|m-l\right\rangle$ $\displaystyle=$ $\displaystyle
e^{-2m\kappa t}\sum_{l=0}^{m-n}\frac{m!\left(e^{2\kappa
t}-1\right)^{l}}{l!n!\left(m-n-l\right)!}\xi^{n}\left(1-\xi\right)^{m-l-n}$
$\displaystyle=$ $\displaystyle\frac{m!\xi^{n}e^{-2m\kappa
t}}{n!\left(m-n\right)!}\sum_{l=0}^{m}\binom{m-n}{m-n-l}\left(e^{2\kappa
t}-1\right)^{l}\left(1-\xi\right)^{m-n-l}$ $\displaystyle=$
$\displaystyle\frac{m!\left(\xi e^{-2\kappa
t}\right)^{n}}{n!\left(m-n\right)!}\left(1-\xi e^{-2\kappa t}\right)^{m-n},$

which is the same as Eq.(24). Thus our new formula is correct.

In summary, we conclude that for a dissipative system characteristic of
amplitude damping, the photocount distribution at any time can be considered
as the one at initial time with the time evolution quantum efficiency $\xi
e^{-2\kappa t},$ shown in Eq.(22). This law greatly simplifies the theoretical
study of photocount distribution for quantum optical field. In addition, in
Ref.r10 , we have derived two new quantum-mechanical photocount formulas which
build the relation from Wigner function and Q-function to photocount
distribution. Thus our result can be applied directly to the case.

Acknowledgements: This work was supported by the National Natural Science
Foundation of China (Grant Nos. 11175113 and 11264018), and the Young Talents
Foundation of Jiangxi Normal University.

## References

  * (1) C. Gardiner and P. Zoller, Quantum Noise (Springer Berlin, 2000).
  * (2) Hong-yi Fan and Yue Fan, New representation of thermal states in thermal field dynamics, Phys. Lett. A 246, 242 (1998).
  * (3) Hong-yi Fan, Hai-liang Lu, Yue Fan, Newton–Leibniz integration for ket–bra operators in quantum mechanics and derivation of entangled state representations, Ann. Phys. 321, 480 (2006).
  * (4) Preskill J. Lecture Notes for Physics 229: Quantum Information and Computation, California Institution of Technology, 1998.
  * (5) P. L. Kelley and W. H. Kleiner, Theory of Electromagnetic Field Measurement and Photoelectron Counting, Phys. Rev. 136, 316 (1964).
  * (6) M. O. Scully and W. E. Lamb, Jr. Quantum Theory of an Optical Maser. III. Theory of Photoelectron Counting Statistics, Phys. Rev. 179, 368 (1969).
  * (7) B. R. Mollow, Quantum Theory of Field Attenuation, Phys. Rev. 168, 1896 (1968).
  * (8) M. Orszag, Quantum Optics, Springer-Verlag, Berlin, 2000.
  * (9) R. Loudon, The Quantum Theory of Light, Oxford University Press, 2nd, 1983.
  * (10) H.Y. Fan and L.Y. Hu, Two quantum-mechanical photocount formulas, Opt. Lett. 33, 443 (2008).