Let us detail a practical application of the formalism above. Consider the
product of $k$ isotropic three-dimensional oscillators given by the equations

$\frac{\mathrm{d}^{2}x^{i}_{\alpha}}{\mathrm{d}t^{2}}=-b_{\alpha}^{2}x^{i}_{\alpha}\,,\qquad\alpha=1,\dotsc,k\,,\qquad
i=1,2,3,$

where the $b_{\alpha}>0$, with $\alpha=1,\ldots,k$, are a series of constants.
The above system of second-order differential equations can be written as a
first-order system of differential equations

$\begin{dcases}\frac{\mathrm{d}x^{i}_{\alpha}}{\mathrm{d}t}=p^{i}_{\alpha}\,,\\\
\frac{\mathrm{d}p^{i}_{\alpha}}{\mathrm{d}t}=-b_{\alpha}^{2}x^{i}_{\alpha}\,,\end{dcases}\qquad\alpha=1,\dotsc,k\,,\qquad
i=1,2,3\,,$ (5.23)

on the product manifold $P=(\mathrm{T}^{\ast}\mathbb{R}^{3})^{k}$. The
$\alpha$-th factor $\mathrm{T}^{\ast}\mathbb{R}^{3}$ in $P$ is a symplectic
manifold equipped with the symplectic form

$\omega^{\alpha}=\sum_{i=1}^{3}\mathrm{d}x^{i}_{\alpha}\wedge\mathrm{d}p^{i}_{\alpha}\,,$

where we stress that there is no sum over the index $\alpha$. Then, $P$ is a
$k$-polysymplectic manifold when endowed with the $\mathbb{R}^{k}$-valued form
${\bm{\omega}}=\sum_{\alpha=1}^{k}\omega^{\alpha}\otimes e_{\alpha}$, where
$\omega^{1},\ldots,\omega^{k}$ are considered as pulled back to $P$ in the
natural way. Moreover, (5.23) describes the integral curves of the vector
field

$X_{\bm{h}}=\sum_{\alpha=1}^{k}\sum_{i=1}^{3}\left(p_{\alpha}^{i}\frac{\partial}{\partial
x^{i}_{\alpha}}-b_{\alpha}^{2}x_{\alpha}^{i}\frac{\partial}{\partial
p_{\alpha}^{i}}\right)\,,$

which is ${\bm{\omega}}$-Hamiltonian admitting an ${\bm{\omega}}$-Hamiltonian
function

${\bm{h}}=\frac{1}{2}\sum_{\alpha=1}^{k}\left(p_{\alpha}^{2}+b_{\alpha}^{2}x_{\alpha}^{2}\right)\otimes
e_{\alpha}\,,\qquad p_{\alpha}^{2}=\sum_{i=1}^{3}(p_{\alpha}^{i})^{2},\qquad
x_{\alpha}^{2}=\sum_{i=1}^{3}(x_{\alpha}^{i})^{2}.$ (5.24)

Let us consider a Lie group action
$\Phi^{\alpha}:\operatorname{SO}(3)\times(\mathrm{T}^{\ast}\mathbb{R}^{3})_{\alpha}\rightarrow(\mathrm{T}^{\ast}\mathbb{R}^{3})_{\alpha}$,
where each $\Phi^{\alpha}$ is the lift of the natural Lie group action
$\Psi:\operatorname{SO}(3)\times\mathbb{R}^{3}\rightarrow\mathbb{R}^{3}$
induced by rotations on $\mathbb{R}^{3}$ to the $\alpha$-th copy of
$\mathrm{T}^{\ast}\mathbb{R}^{3}$ in $P$. Then, the resulting Lie group action
$\Phi$ on $(\mathrm{T}^{\ast}\mathbb{R}^{3})^{k}$ given by (3.6) reads

$\Phi:\operatorname{SO}(3)^{k}\times(\mathrm{T}^{\ast}\mathbb{R}^{3})^{k}\longrightarrow(\mathrm{T}^{\ast}\mathbb{R}^{3})^{k}\,.$

The Lie algebra of fundamental vector fields of $\Phi$ is spanned by the basis
of vector fields on $P$ of the form

$\displaystyle\xi^{1}_{\alpha P}=\left(x_{\alpha}^{1}\frac{\partial}{\partial
x_{\alpha}^{2}}-x_{\alpha}^{2}\frac{\partial}{\partial
x_{\alpha}^{1}}+p_{\alpha}^{1}\frac{\partial}{\partial
p_{\alpha}^{2}}-p_{\alpha}^{2}\frac{\partial}{\partial
p_{\alpha}^{1}}\right),\qquad\xi^{2}_{\alpha
P}=\left(x_{\alpha}^{2}\frac{\partial}{\partial
x_{\alpha}^{3}}-x_{\alpha}^{3}\frac{\partial}{\partial
x_{\alpha}^{2}}+p_{\alpha}^{2}\frac{\partial}{\partial
p_{\alpha}^{3}}-p_{\alpha}^{3}\frac{\partial}{\partial
p_{\alpha}^{2}}\,\right),$ (5.25) $\displaystyle\xi^{3}_{\alpha
P}=\left(x_{\alpha}^{3}\frac{\partial}{\partial
x_{\alpha}^{1}}-x_{\alpha}^{1}\frac{\partial}{\partial
x_{\alpha}^{3}}+p_{\alpha}^{3}\frac{\partial}{\partial
p_{\alpha}^{1}}-p_{\alpha}^{1}\frac{\partial}{\partial
p_{\alpha}^{3}}\right)\,,$ (5.26)

with $\alpha=1,\ldots,k$. These vector fields are Lie symmetries of
$\bm{\omega}$ and $\bm{h}$. Moreover, a $k$-polysymplectic momentum map
associated with $\Phi$ is given by
$\mathbf{J}:(\mathrm{T}^{\ast}\mathbb{R}^{3})^{k}\rightarrow[(\mathfrak{so}_{3}^{k})^{*}]^{k}$
such that

$\mathbf{J}(\bm{q}_{1},\ldots,\bm{q}_{k})=(0,0,0;\ldots;J_{\alpha}^{1},J_{\alpha}^{2},J_{\alpha}^{3};\ldots;0,0,0)\otimes
e_{\alpha}$

where
$\bm{q}_{\alpha}=(x^{1}_{\alpha},x^{2}_{\alpha},x^{3}_{\alpha},p^{1}_{\alpha},p^{2}_{\alpha},p^{3}_{\alpha})\in\mathrm{T}^{*}\mathbb{R}^{3}$
for $\alpha=1,\ldots,k$, while

$(J^{1}_{\alpha},J^{2}_{\alpha},J^{3}_{\alpha})=(x_{\alpha}^{1}p_{\alpha}^{2}-x_{\alpha}^{2}p_{\alpha}^{1},x_{\alpha}^{2}p_{\alpha}^{3}-x_{\alpha}^{3}p_{\alpha}^{2},x_{\alpha}^{3}p_{\alpha}^{1}-x_{\alpha}^{1}p_{\alpha}^{3})\,,$

and $\alpha=1,\ldots,k$. Note that the elements of $\mathfrak{so}_{3}^{*}$ are
represented by the coordinates given in an appropriate basis. The function
$x_{\alpha}^{1}p_{\alpha}^{2}-x_{\alpha}^{2}p_{\alpha}^{1}$ is the angular
momentum, $p_{\alpha\varphi}$, of the $\alpha$-th particle in the
corresponding spherical coordinates
$\\{r_{\alpha},\theta_{\alpha},\varphi_{\alpha}\\}$. Meanwhile,
$L_{\alpha}^{2}=(J_{\alpha}^{1})^{2}+(J_{\alpha}^{2})^{2}+(J_{\alpha}^{3})^{2}$
is the square of the total angular momentum of the $\alpha$-th particle. Both
quantities are conserved by the evolution of $X_{\bm{h}}$.

The momentum map $\mathbf{J}$ is $(\operatorname{Ad}^{*})^{k}$-equivariant.
Recall that $\mathbf{J}=(0,\ldots,\mathbf{J}^{\alpha},\ldots,0)\otimes
e_{\alpha}$. Then,
$\bm{\mu}=(0,0,0;\ldots;J_{\alpha}^{1},J_{\alpha}^{2},J_{\alpha}^{3};\ldots;0,0,0)\otimes
e_{\alpha}$ is weakly regular value of $\mathbf{J}$ if and only if each triple
$\mu^{\alpha}=(J_{\alpha}^{1},J_{\alpha}^{2},J_{\alpha}^{3})\in\mathfrak{so}_{3}^{*}$
is a weakly regular value of $\mathbf{J}^{\Phi^{\alpha}}$, where
$\alpha=1,\ldots,k$. Let us fix some weakly regular $\bm{\mu}$. Then,

$\mathrm{T}_{\bm{q}}(\mathbf{J}^{-1}(\bm{\mu}))=\bigoplus_{\alpha=1}^{k}\mathrm{T}_{{\bm{q}}_{\alpha}}(\mathbf{J}_{\alpha}^{\Phi^{\alpha}-1}(\mu^{\alpha}))\,,\qquad\forall{\bf
q}=({\bf q}_{1},\ldots,{\bf q}_{k})\in P\,.$

Moreover,

$\xi_{\alpha
P}^{i}J^{j}_{\beta}=-\delta_{\alpha\beta}\epsilon_{ijk}J^{k}_{\beta}\,,\qquad
i,j=1,2,3,\qquad\alpha,\beta=1,\ldots,k.$

The isotropy subgroup of $\Phi$ at $\bm{\mu}$ is given by the Cartesian
product of all the isotropy subgroups corresponding to each $\mu^{\alpha}$
relative to $\Phi^{\alpha}$ and $\alpha=1,\ldots,k$. To obtain
$G_{\mu^{\alpha}}$, one may verify when $\sum_{i=1}^{3}\lambda_{i}(\xi_{\alpha
P}^{i})$ belongs to $\mathrm{T}_{{\bm{q}_{\alpha}}}({\bf
J}^{\Phi^{\alpha}-1}(\mu_{\alpha}))$, namely
$\sum_{i=1}^{3}\lambda_{i}(\xi_{P}^{i})_{\alpha}J^{j}_{\alpha}=0$ for
$j=1,2,3$ (with no summation over $\alpha$), which occurs if and only if

$\begin{pmatrix}0&-J_{\alpha}^{3}&J_{\alpha}^{2}\\\
J_{\alpha}^{3}&0&-J_{\alpha}^{1}\\\
-J_{\alpha}^{2}&J_{\alpha}^{1}&0\end{pmatrix}\begin{pmatrix}\lambda_{1}\\\
\lambda_{2}\\\ \lambda_{3}\end{pmatrix}=\begin{pmatrix}0\\\ 0\\\
0\end{pmatrix}\,.$

The matrix of coefficients has rank two for $L^{2}_{\alpha}\neq 0$. Moreover,
$\mathbf{J}^{\Phi^{\alpha}}$ has a regular value at $\mu_{\alpha}$ when
$L^{2}_{\alpha}\neq 0$ for every $\alpha=1,\ldots,k$. Let us restrict to that
case. Each isotropy subgroup $G_{\mu^{\alpha}}$ has always dimension one.
Hence, the reduced manifold ${\bf J}^{\Phi-1}({\bm{\mu}})/G_{\bm{\mu}}$ has
dimension $6k-3k-k=2k$. The conditions for the $k$-polysymplectic
Marsden–Weinstein reduction results, as already commented, from the ones for
the symplectic reduction on each component, which are satisfied. Hence, the
$k$-polysymplectic Marsden–Weinstein reduction exists.

Note that $\bm{h}=\frac{1}{2}\sum_{\alpha=1}^{k}(p_{\alpha
r}^{2}+p_{\alpha\varphi}^{2}/(r_{\alpha}^{2}\sin^{2}\theta_{\alpha})+p_{\alpha\theta}^{2}/r_{\alpha}^{2}+b_{\alpha}^{2}r_{\alpha}^{2})\otimes
e_{\alpha}$ and $\omega^{\alpha}=\mathrm{d}r_{\alpha}\wedge\mathrm{d}p_{\alpha
r}+\mathrm{d}\theta_{\alpha}\wedge\mathrm{d}p_{\alpha\theta}+\mathrm{d}\varphi_{\alpha}\wedge\mathrm{d}p_{\alpha\varphi}$,
for each $\alpha=1,\ldots,k$, in spherical coordinates for the $k$ component
manifolds of $(\mathrm{T}^{\ast}\mathbb{R}^{3})^{k}$. Then, the differential
equations for the integral curves of $X_{\bm{h}}$ read

$\begin{gathered}\frac{\mathrm{d}p_{\alpha
r}}{\mathrm{d}t}=\frac{p_{\alpha\varphi}^{2}}{r_{\alpha}^{3}\sin^{2}\theta_{\alpha}}+\frac{p_{\alpha\theta}^{2}}{r_{\alpha}^{3}}-b^{2}_{\alpha}r_{\alpha}\,,\qquad\frac{\mathrm{d}p_{\alpha\varphi}}{\mathrm{d}t}=0\,,\qquad\frac{\mathrm{d}p_{\alpha\theta}}{\mathrm{d}t}=\frac{p_{\alpha\varphi}^{2}\cos\theta_{\alpha}}{r_{\alpha}^{2}\sin^{3}\theta_{\alpha}}\,,\\\
\frac{\mathrm{d}r_{\alpha}}{\mathrm{d}t}=p_{\alpha
r}\,,\qquad\frac{\mathrm{d}\theta_{\alpha}}{\mathrm{d}t}=\frac{p_{\alpha\theta}}{r_{\alpha}^{2}}\,,\qquad\frac{\mathrm{d}\varphi_{\alpha}}{\mathrm{d}t}=\frac{p_{\alpha\varphi}}{r_{\alpha}^{2}\sin^{2}\theta_{\alpha}}\,.\end{gathered}$
(5.27)

$k$-Polysymplectic relative equilibrium points are given by those points for
which the vector field $X_{\bm{h}}$ on $(\mathrm{T}^{\ast}\mathbb{R}^{3})^{k}$
corresponding to the dynamics is proportional to one of the fundamental vector
fields of $\Phi$. In particular, let us take $z_{e}\in P$ such that
$z_{e}=(r_{\alpha},\theta_{\alpha}\\!\\!=\\!\\!\frac{\pi}{2},\varphi_{\alpha},p_{\alpha
r}\\!\\!=\\!\\!0,p_{\alpha\theta}=0,p_{\alpha\varphi})$ and
$L_{\alpha}=b_{\alpha}r^{2}_{\alpha}=p_{\alpha\varphi}$ for every
$\alpha=1,\ldots,k$ on analysed points. Then, the ${\bm{\omega}}$-Hamiltonian
vector field $X_{{\bm{h}}}$ at such points is

$X_{\bm{h}}=\sum_{\alpha=1}^{k}\frac{p_{\alpha\varphi}}{r_{\alpha}^{2}}\frac{\partial}{\partial\varphi_{\alpha}}.$

This implies that $z_{e}\in P$ is a $k$-polysymplectic relative equilibrium
point of $X_{{\bm{h}}}$. Let us demonstrate this by applying Theorem 4.2. This
theorem ensures that $z_{e}$ is a $k$-polysymplectic relative equilibrium
point of $X_{{\bm{h}}}$ if and only if there exists
$\xi\in\mathfrak{so}_{3}^{k}$ such that
$\bm{h}_{\xi}={\bm{h}}-\langle\mathbf{J}^{\Phi}-\bm{\mu}_{e},\xi\rangle$ has a
critical point at $z_{e}$. Indeed, for

$\xi=(p_{1\varphi}/r^{2}_{1},0,0;\ldots;p_{k\varphi}/r^{2}_{k},0,0)\in\mathfrak{so}_{3}^{k},$

the $\mathbb{R}^{k}$-valued function

$\bm{h}_{\xi}=\bm{h}-\langle\mathbf{J}^{\Phi}-\bm{\mu}_{e},\xi\rangle=\sum_{\alpha=1}^{k}(h^{\alpha}-\langle(0,\ldots,{\bf
J}^{\Phi_{\alpha}}-(L_{\alpha},0,0),\ldots,0),\xi\rangle)\otimes
e_{\alpha}\,,$

has a critical point at $z_{e}$. Therefore, $z_{e}$ is a $k$-polysymplectic
relative equilibrium point.

By Theorem 3.2, the reduced manifolds is $(\mathrm{T}^{\ast}\mathbb{R})^{k}$
with coordinates $\\{r_{\alpha},p_{\alpha r}\\}$ for $\alpha=1,\ldots,k$. The
reduced $k$-polysymplectic form on the reduced manifold is given by

$\bm{\omega}_{\bm{\mu}}=\sum^{k}_{\alpha=1}\mathrm{d}r_{\alpha}\wedge\mathrm{d}p_{\alpha
r}\otimes e_{\alpha},$

and the reduced ${\bm{\omega}}_{\bm{\mu}}$-Hamiltonian reads

${\bm{f}}_{\bm{\mu}}=\frac{1}{2}\sum_{\alpha=1}^{k}\left(p_{\alpha
r}^{2}+\frac{L_{\alpha}^{2}}{r^{2}_{\alpha}}+b_{\alpha}^{2}r_{\alpha}^{2}\right)\otimes
e_{\alpha}\,.$

Furthermore, one has that

$\frac{\mathrm{d}p_{\alpha
r}}{\mathrm{d}t}=-b_{\alpha}^{2}r_{\alpha}+\frac{L_{\alpha}^{2}}{r^{3}_{\alpha}}\,,\qquad\frac{\mathrm{d}r_{\alpha}}{\mathrm{d}t}=p_{\alpha
r}\,,\qquad\alpha=1,\dotsc,k\,.$

Thus, the equilibrium points of $X_{{\bm{f}}_{\bm{\mu}}}$ have $p_{\alpha
r}=0$ and

$-b^{2}_{\alpha}r_{\alpha}+\frac{L_{\alpha}^{2}}{r^{3}_{\alpha}}=0\,$

for $\alpha=1,\ldots,k$. Note that this point is the projection of a
$k$-polysymplectic relative equilibrium point $z_{e}\in P$.

The Hessian of the functions $f_{\bm{\mu}}^{\alpha}$ is positive-definite in a
supplementary to the kernel of ${\bm{\omega}}_{\mu^{\alpha}}$ at the
equilibrium point. Indeed,

$\operatorname{Hess}(f_{\bm{\mu}}^{\alpha})=\begin{pmatrix}1&0\\\
0&4b_{\alpha}^{2}\end{pmatrix}\,.$

Moreover, the function $\sum_{\alpha=1}^{k}f^{\alpha}_{\bm{\mu}}$ has a
positive-definite Hessian and the equilibrium point becomes a strict minimum.
This means that the reduced $k$-polysymplectic relative equilibrium point is
stable. In the original manifold, the orbits around $k$-polysymplectic
relative equilibrium points remain in the anti-image in ${\bf
J}^{\Phi-1}(\bm{\mu})$ of an open neighbourhood of the projection of the
$k$-polysymplectic relative equilibrium points.

###  5.3 $k$-polysymplectic affine Lie systems

Let us apply our techniques to a family of affine inhomogeneous systems of
first-order differential equations. It is worth noting that all such systems
are Lie systems [13]. We will hereafter call such differential equations
affine Lie systems. Many such systems appear in control theory and other
relevant disciplines [16]. In particular, we are here concerned with affine
Lie systems admitting a Lie algebra of Hamiltonian vector fields relative to a
$k$-polysymplectic form. We call them $k$-polysymplectic affine Lie systems.

Although our techniques could be extended to other affine Lie systems, let us
restrict ourselves to the particular case

$\frac{\mathrm{d}}{\mathrm{d}t}\begin{pmatrix}x_{1}\\\ x_{2}\\\ x_{3}\\\
x_{4}\\\ x_{5}\end{pmatrix}=\begin{pmatrix}b_{1}(t)\\\ b_{2}(t)\\\ b_{3}(t)\\\
b_{4}(t)\\\ b_{5}(t)\end{pmatrix}+b_{6}(t)\begin{pmatrix}0&0&0&0&0\\\
0&0&0&0&0\\\ 0&0&0&0&1\\\ 0&0&0&0&0\\\
0&0&-1&0&0\end{pmatrix}\begin{pmatrix}x_{1}\\\ x_{2}\\\ x_{3}\\\ x_{4}\\\
x_{5}\\\ \end{pmatrix}\,,$ (5.28)

where $b_{1}(t),\ldots,b_{6}(t)$ are arbitrary $t$-dependent functions. The
above system is the system of differential equations describing the integral
curves of the $t$-dependent vector field

$X=\sum_{\alpha=1}^{6}b_{\alpha}(t)X_{\alpha}\,,$

where

$X_{1}=\frac{\partial}{\partial x_{1}}\,,\quad X_{2}=\frac{\partial}{\partial
x_{2}}\,,\quad X_{3}=\frac{\partial}{\partial x_{3}}\,,\quad
X_{4}=\frac{\partial}{\partial x_{4}}\,,\quad X_{5}=\frac{\partial}{\partial
x_{5}}\,,\quad X_{6}=x_{5}\frac{\partial}{\partial
x_{3}}-x_{3}\frac{\partial}{\partial x_{5}}\,.$

These vector fields span a six-dimensional Lie algebra of vector fields $V$,
with the following non-vanishing commutation relations

$[X_{3},X_{6}]=-X_{5}\,,\qquad[X_{5},X_{6}]=X_{3}\,.$

Consider the case where $b_{1}(t),\ldots,b_{6}(t)$ are constants, denoted as
$c_{1},\ldots,c_{6}\in\mathbb{R}$, respectively. Since the vector fields
$X_{1}\wedge\dotsb\wedge X_{6}=0$, the methods presented in Section 5.1 for
describing $k$-polysymplectic forms compatible with Lie systems do not apply
to (5.28). Nevertheless, there exists a two-polysymplectic form on
$\mathbb{R}^{5}$ given by

$\bm{\omega}=(\mathrm{d}x_{3}\wedge\mathrm{d}x_{5}+\mathrm{d}x_{4}\wedge\mathrm{d}x_{1})\otimes
e_{1}+(\mathrm{d}x_{3}\wedge\mathrm{d}x_{5}+\mathrm{d}x_{4}\wedge\mathrm{d}x_{2})\otimes
e_{2}$

turning all the vector fields of $V$ into ${\bm{\omega}}$-Hamiltonian vector
fields. Indeed, ${\bm{\omega}}$-Hamiltonian functions for $X_{1},\ldots,X_{6}$
have the form

$\displaystyle{\bm{h}}_{1}$ $\displaystyle=-x_{4}\otimes e_{1}\,,$
$\displaystyle{\bm{h}}_{2}$ $\displaystyle=-x_{4}\otimes e_{2}\,,$
$\displaystyle{\bm{h}}_{3}$ $\displaystyle=x_{5}\otimes e_{1}+x_{5}\otimes
e_{2}\,,$ (5.29) $\displaystyle{\bm{h}}_{4}$ $\displaystyle=x_{1}\otimes
e_{1}+x_{2}\otimes e_{2}\,,$ $\displaystyle{\bm{h}}_{5}$
$\displaystyle=-x_{3}\otimes e_{1}-x_{3}\otimes e_{2}\,,$
$\displaystyle{\bm{h}}_{6}$
$\displaystyle=\frac{1}{2}(x_{3}^{2}+x_{5}^{2})\otimes
e_{1}+\frac{1}{2}(x_{3}^{2}+x_{5}^{2})\otimes e_{2}\,.$ (5.30)

The flow of the vector field $X_{4}$ gives rise to a two-polysymplectic Lie
group action $\Phi:\mathbb{R}\times\mathbb{R}^{5}\rightarrow\mathbb{R}^{5}$.
Moreover, $X_{4}$, which spans the space of fundamental vector fields of
$\Phi$, is a Lie symmetry of system (5.28). Then, a two-polysymplectic
momentum map associated with $\Phi$ reads

${\bf
J}^{\Phi}:(x_{1},x_{2},x_{3},x_{4},x_{5})\in\mathbb{R}^{5}\mapsto(x_{1},x_{2})=\bm{\mu}\in\mathbb{R}^{*2}.$

Note that $\bm{\mu}\in\mathbb{R}^{*2}$ is a regular value of ${\bf J}^{\Phi}$,
and ${\bf J}^{\Phi}$ is $\operatorname{Ad}^{*2}$-equivariant two-
polysymplectic momentum map. It can be proved that the example satisfies the
conditions (3.1) and (3.2). Hence, Theorem 3.2 can be applied. The vector
field $X_{4}$ is tangent to $\mathbf{J}^{\Phi-1}(\bm{\mu})$ and
$\mathrm{T}_{x}(G_{\bm{\mu}}x)=\langle\frac{\partial}{\partial x_{4}}\rangle$
for $x\in\mathbb{R}^{5}$. Therefore,
$P_{\bm{\mu}}=\mathbf{J}^{\Phi-1}(\bm{\mu})/\mathbb{R}$ is a two-dimensional
manifold and the variables $\\{x_{3},x_{5}\\}$ can be considered in a natural
manner as variables on $P_{\bm{\mu}}$. The reduced two-polysymplectic form
reads

$\bm{\omega}_{\bm{\mu}}=\omega^{1}_{\bm{\mu}}\otimes
e_{1}+\omega^{2}_{\bm{\mu}}\otimes
e_{2}=\mathrm{d}x_{3}\wedge\mathrm{d}x_{5}\otimes
e_{1}+\mathrm{d}x_{3}\wedge\mathrm{d}x_{5}\otimes e_{2}\,.$

To apply Theorem 3.3, the affine Lie system must be tangent to ${\bf
J}^{\Phi-1}(\bm{\mu})$, which can be ensured by assuming that its associated
$\bm{\omega}$-Hamiltonian function has to be invariant relative to $X_{4}$.
These conditions are satisfied by imposing $c_{1}=c_{2}=0$. The resulting
vector field, $X_{\bm{\omega}}=c_{3}X_{3}+c_{4}X_{4}+c_{5}X_{5}+c_{6}X_{6}$,
projects onto $P_{\bm{\mu}}={\bf J}^{\Phi-1}(\bm{\mu})/\mathbb{R}$ giving rise
to an $\bm{\omega}_{\bm{\mu}}$-Hamiltonian vector field of the form

$X_{{\bm{\mu}}}=c_{6}\left(x_{5}\frac{\partial}{\partial
x_{3}}-x_{3}\frac{\partial}{\partial
x_{5}}\right)+c_{3}\frac{\partial}{\partial
x_{3}}+c_{5}\frac{\partial}{\partial x_{5}}\,.$

The $\bm{\omega}_{\bm{\mu}}$-Hamiltonian function of $X_{\bm{\mu}}$ reads

$\bm{f}_{\bm{\mu}}=\left(c_{3}x_{5}-c_{5}x_{3}+c_{6}\left(\frac{x_{3}^{2}}{2}+\frac{x_{5}^{2}}{2}\right)\right)\otimes
e_{1}+\left(c_{3}x_{5}-c_{5}x_{3}+c_{6}\left(\frac{x_{3}^{2}}{2}+\frac{x_{5}^{2}}{2}\right)\right)\otimes
e_{2}\,.$

Next, the methods introduced in Section 4 will be employed to find the two-
polysymplectic relative equilibrium points of the $\bm{\omega}$-Hamiltonian
vector field

$Y=X_{4}+X_{6}$

and study their stability. According to Theorem 4.2, a two-polysymplectic
relative equilibrium point $z_{e}\in P$ is a point for which there exists
$\xi\in\mathfrak{g}\simeq\mathbb{R}$ such that $z_{e}$ is a critical point of
each component of the $\mathbb{R}^{2}$-valued function

${\bm{h}}_{\xi}=\left(x_{1}-\xi(x_{1}-\mu^{1})+\frac{1}{2}(x_{3}^{2}+x_{5}^{2})\right)\otimes
e_{1}+\left(x_{2}-\xi(x_{2}-\mu^{2})+\frac{1}{2}(x_{3}^{2}+x_{5}^{2})\right)\otimes
e_{1}\,.$

This happens for $\xi=1$ and
$z_{e}=(x_{1},x_{2},x_{3}=0,x_{4},x_{5}=0)\in\mathbb{R}^{5}$, where
$x_{1},x_{2},x_{4}$ are arbitrary.

To examine the relative stability of $z_{e}\in\mathbb{R}^{5}$, note that the
supplementary spaces to
$\mathrm{T}_{z_{e}}(G_{{\bm{\mu}}_{e}}z_{e})+\ker\omega^{1}_{z_{e}}$ and
$\mathrm{T}_{z_{e}}(G_{{\bm{\mu}}_{e}}z_{e})+\ker\omega^{2}_{z_{e}}$ in
$\mathrm{T}_{z_{e}}({\bf J}^{\Phi-1}(\bm{\mu}_{e}))$ have the form

$\mathcal{S}^{1}=\left\langle\frac{\partial}{\partial
x_{3}},\frac{\partial}{\partial
x_{5}}\right\rangle\,,\qquad\mathcal{S}^{2}=\left\langle\frac{\partial}{\partial
x_{3}},\frac{\partial}{\partial x_{5}}\right\rangle\,,$

respectively. Then,
$\mathcal{S}^{1}+\mathcal{S}^{2}+\mathrm{T}_{z_{e}}(G_{\bm{\mu}_{e}}z_{e})=\mathrm{T}_{z_{e}}({\bf
J}^{\Phi-1}(\bm{\mu}_{e}))$ and the Hessian of
$(\delta^{2}h^{1}_{\xi})_{z_{e}}$ at $z_{e}$ in the subspace $\mathcal{S}^{1}$
and the Hessian of $(\delta^{2}h^{2}_{\xi})_{z_{e}}$ in the subspace
$\mathcal{S}^{2}$ are definite-positive. Therefore, by our criterion, a two-
polysymplectic relative equilibrium point $z_{e}\in\mathbb{R}^{5}$ is
relatively stable, namely its projection to $P_{\bm{\mu}_{e}}$ is stable. More
specifically, the reduced system has an
${\bm{\omega}}_{\bm{\mu}_{e}}$-Hamiltonian function whose components,
$f^{1}_{\bm{\mu}_{e}},f^{2}_{\bm{\mu}_{e}}$ are such that their Hessians at
equilibrium points $\pi_{\bm{\mu}_{e}}(z_{e})=(x_{3}=0,x_{5}=0)$ are definite-
positive in the directions of $\ker(\omega^{1}_{\bm{\mu}_{e}})_{z_{e}}$ and
$\ker(\omega^{2}_{\bm{\mu}_{e}})_{z_{e}}$, respectively. Indeed, the reduced
${\bm{\omega}}_{\bm{\mu}_{e}}$-Hamiltonian function reads

$\bm{f}_{\bm{\mu}_{e}}=\frac{1}{2}(x_{3}^{2}+x_{5}^{2})\otimes
e_{1}+\frac{1}{2}(x_{3}^{2}+x_{5}^{2})\otimes e_{2}$

and the function

$f_{\bm{\mu}_{e}}^{1}+f_{\bm{\mu}_{e}}^{2}=x_{3}^{2}+x_{5}^{2}\,,$

is invariant under the dynamics of $Y_{\bm{\mu}_{e}}$ and has a strict minimum
at $\pi_{\bm{\mu}_{e}}(z_{e})=(x_{3}=0,x_{5}=0)$. Hence, the reduced two-
polysymplectic Hamiltonian system is stable.

###  5.4 Quantum quadratic Hamiltonian operators

Next, let us analyse an example based upon the Wei–Norman equations for the
automorphic Lie system related to quantum mechanical systems described by
quadratic Hamiltonian operators, which describe as particular cases quantum
harmonic oscillators with/without dissipation [13, 56]. In this case, the
system of differential equations under study is the one determining the
integral curves of the time-dependent vector field

$X=\sum_{\alpha=1}^{6}b_{\alpha}(t)X_{\alpha}^{R},$ (5.31)

for certain $t$-dependent functions $b_{1}(t),\ldots,b_{6}(t)$ and the vector
fields

$\displaystyle X^{R}_{1}=\frac{\partial}{\partial
v_{1}}+v_{5}\frac{\partial}{\partial
v_{4}}-\frac{1}{2}v_{5}^{2}\frac{\partial}{\partial v_{6}}\,,$ $\displaystyle
X^{R}_{2}$ $\displaystyle=v_{1}\frac{\partial}{\partial
v_{1}}+\frac{\partial}{\partial
v_{2}}+\frac{1}{2}v_{4}\frac{\partial}{\partial
v_{4}}-\frac{1}{2}v_{5}\frac{\partial}{\partial v_{5}}\,,$ (5.32)
$\displaystyle X^{R}_{3}=v_{1}^{2}\frac{\partial}{\partial
v_{1}}+2v_{1}\frac{\partial}{\partial v_{2}}+e^{v_{2}}\frac{\partial}{\partial
v_{3}}-v_{4}\frac{\partial}{\partial
v_{5}}+\frac{1}{2}v_{4}^{2}\frac{\partial}{\partial v_{6}}\,,$ $\displaystyle
X^{R}_{4}$ $\displaystyle=\frac{\partial}{\partial v_{4}}\,,$ $\displaystyle
X^{R}_{5}=\frac{\partial}{\partial v_{5}}-v_{4}\frac{\partial}{\partial
v_{6}}\,,$ $\displaystyle X^{R}_{6}$ $\displaystyle=\frac{\partial}{\partial
v_{6}}\,.$

The commutation relations between the above vector fields read

$\displaystyle[X^{R}_{1},X^{R}_{2}]=X^{R}_{1}\,,$ (5.33)
$\displaystyle[X^{R}_{1},X^{R}_{3}]=2\,X^{R}_{2}\,,$
$\displaystyle[X^{R}_{2},X^{R}_{3}]$ $\displaystyle=X^{R}_{3}\,,$
$\displaystyle[X^{R}_{1},X^{R}_{4}]=0\,,$ $\displaystyle[X^{R}_{2},X^{R}_{4}]$
$\displaystyle=-\frac{1}{2}\,X^{R}_{4}\,,$
$\displaystyle[X^{R}_{3},X^{R}_{4}]$ $\displaystyle=X^{R}_{5}\,,$
$\displaystyle[X^{R}_{1},X^{R}_{5}]=-X^{R}_{4}\,,$
$\displaystyle[X^{R}_{2},X^{R}_{5}]$ $\displaystyle=\frac{1}{2}\,X^{R}_{5}\,,$
$\displaystyle[X^{R}_{3},X^{R}_{5}]$ $\displaystyle=0\,,$
$\displaystyle[X^{R}_{4},X^{R}_{5}]$ $\displaystyle=-X^{R}_{6}\,,$
$\displaystyle[X^{R}_{1},X^{R}_{6}]=0\,,$ $\displaystyle[X^{R}_{2},X^{R}_{6}]$
$\displaystyle=0\,,$ $\displaystyle[X^{R}_{3},X^{R}_{6}]$ $\displaystyle=0\,,$
$\displaystyle[X^{R}_{4},X^{R}_{6}]$ $\displaystyle=0\,,$
$\displaystyle[X^{R}_{5},X^{R}_{6}]$ $\displaystyle=0\,.$

It is known that the Lie algebra of Lie symmetries of $\langle
X^{R}_{1},\ldots,X_{6}^{R}\rangle$ is spanned by

$\begin{gathered}X^{L}_{1}=e^{v_{2}}\frac{\partial}{\partial
v_{1}}+2v_{3}\frac{\partial}{\partial v_{2}}+v_{3}^{2}\frac{\partial}{\partial
v_{3}}\,,\qquad X^{L}_{2}=\frac{\partial}{\partial
v_{2}}+v_{3}\frac{\partial}{\partial v_{3}}\,,\qquad
X^{L}_{3}=\frac{\partial}{\partial v_{3}}\,,\\\
X^{L}_{4}=e^{-v_{2}/2}(e^{v_{2}}-v_{1}v_{3})\frac{\partial}{\partial
v_{4}}-e^{-v_{2}/2}v_{3}\frac{\partial}{\partial
v_{5}}-e^{-v_{2}/2}(e^{v_{2}}-v_{1}v_{3})v_{5}\frac{\partial}{\partial
v_{6}}\,,\\\ X^{L}_{5}=v_{1}e^{-v_{2}/2}\frac{\partial}{\partial
v_{4}}+e^{-v_{2}/2}\frac{\partial}{\partial
v_{5}}-v_{1}v_{5}e^{-v_{2}/2}\frac{\partial}{\partial v_{6}}\,,\quad
X^{L}_{6}=\frac{\partial}{\partial v_{6}}\,.\end{gathered}$ (5.34)

In particular, let us focus on systems (5.31) with constant coefficients and,
in particular,

$X^{R}_{5}=\frac{\partial}{\partial v_{5}}-v_{4}\frac{\partial}{\partial
v_{6}}\,.$

Then, a Lie symmetry of our system is given by

$Y=\frac{\partial}{\partial v_{5}}\,.$

A two-polysymplectic form on $\mathbb{R}^{6}$ can be defined in the following
way

$\displaystyle\bm{\omega}$ $\displaystyle=\omega^{1}\otimes
e_{1}+\omega^{2}\otimes e_{2}$
$\displaystyle=\left(\mathrm{d}v_{1}\wedge\mathrm{d}v_{3}+\mathrm{d}v_{2}\wedge\mathrm{d}v_{4}+\mathrm{d}v_{5}\wedge\mathrm{d}v_{1}+\mathrm{d}v_{4}\wedge\mathrm{d}v_{6}\right)\otimes
e_{1}+\left(\mathrm{d}v_{4}\wedge\mathrm{d}v_{6}-\mathrm{d}v_{3}\wedge\mathrm{d}v_{5}\right)\otimes
e_{2}\,.$

Note that

$\ker\omega^{1}=\left\langle\frac{\partial}{\partial
v_{3}}+\frac{\partial}{\partial v_{5}},\frac{\partial}{\partial
v_{2}}+\frac{\partial}{\partial
v_{6}}\right\rangle\,,\qquad\ker\omega^{2}=\left\langle\frac{\partial}{\partial
v_{1}},\frac{\partial}{\partial v_{2}}\right\rangle,$

and $\ker\omega^{1}\cap\ker\omega^{2}=0$ and $(\mathbb{R}^{6},\bm{\omega})$
becomes a two-polysymplectic manifold. The vector field $Y$ is a Lie symmetry
of the two-polysymplectic form, i.e. $\mathscr{L}_{Y}{\bm{\omega}}=0$. Then,

$\iota_{Y_{3}}\omega^{1}=\mathrm{d}v_{1}\,,\qquad\iota_{Y_{3}}\omega^{2}=\mathrm{d}v_{3}\,,$

and a two-polysymplectic momentum map $\mathbf{J}^{\Phi}$ associated with the
Lie group action given by the flow of $Y$ reads

$\mathbf{J}^{\Phi}:\mathbb{R}^{6}\ni
x\longmapsto(v_{1},v_{3})=\bm{\mu}\in(\mathfrak{g}^{*})^{2}\simeq\mathbb{R}^{*2}\,.$

Note that every $\bm{\mu}=(\mu^{1},\mu^{2})\in\mathbb{R}^{*2}$ is a weakly
regular value of $\mathbf{J}^{\Phi}$, which is ${\rm Ad}^{*2}$-equivariant.
The isotropy group for every $\bm{\mu}\in\mathbb{R}^{*2}$ reads
$G_{\bm{\mu}}=\mathbb{R}$. Hence, $\mathbf{J}^{\Phi-1}(\bm{\mu})$ is a
submanifold, as well as $\mathbf{J}_{1}^{\Phi-1}(\mu^{1})$ and
$\mathbf{J}_{2}^{\Phi-1}(\mu^{2})$. Since $Y_{5}$ is tangent to
$\mathbf{J}^{\Phi-1}(\bm{\mu})$, then
$P_{\bm{\mu}}=\mathbf{J}^{\Phi-1}(\bm{\mu})/G_{\bm{\mu}}$ can be locally
coordinated by the functions $\\{v_{2},v_{4},v_{6}\\}$.

The vector field $X_{5}^{R}$ is $\bm{\omega}$-Hamiltonian with

$\iota_{X_{5}^{R}}\bm{\omega}=\iota_{X_{5}^{R}}\omega^{1}\otimes
e_{1}+\iota_{X_{5}^{R}}\omega^{2}\otimes
e_{2}=\mathrm{d}\left(v_{1}+\frac{v_{4}^{2}}{2}\right)\otimes
e_{1}+\mathrm{d}\left(v_{3}+\frac{v_{4}^{2}}{2}\right)\otimes
e_{2}=\mathrm{d}{\bm{h}}^{R}_{5}\,.$ (5.35)

Then, the reduced two-forms read

$\omega^{1}_{\bm{\mu}}=\mathrm{d}v_{2}\wedge\mathrm{d}v_{4}+\mathrm{d}v_{4}\wedge\mathrm{d}v_{6}\,,\qquad\omega^{2}_{\bm{\mu}}=\mathrm{d}v_{4}\wedge\mathrm{d}v_{6}\,.$

Furthermore, one has

$\ker\omega^{1}_{\bm{\mu}}=\left\langle\frac{\partial}{\partial
v_{2}}+\frac{\partial}{\partial
v_{6}}\right\rangle\,,\qquad\ker\omega^{2}_{\bm{\mu}}=\left\langle\frac{\partial}{\partial
v_{2}}\right\rangle\,,$

and $\omega_{\bm{\mu}}^{1},\omega_{\bm{\mu}}^{2}$ define a two-polysymplectic
form on $P_{\bm{\mu}}$. Moreover, the ${\bm{\omega}}$-Hamiltonian function of
$X^{R}_{5}$ is invariant relative to $Y$. Then, Theorem 3.3 ensures that the
projection of $X^{R}_{5}$ onto $P_{\bm{\mu}}$ exists and is given by

$X^{R}=-v_{4}\frac{\partial}{\partial v_{6}},$

which is the $\bm{\omega}_{\bm{\mu}}$-Hamiltonian vector field of the
$\bm{\omega}_{\bm{\mu}}$-Hamiltonian function

${\bm{f}}_{\bm{\mu}}=\left(\mu^{1}+\frac{v_{4}^{2}}{2}\right)\otimes
e_{1}+\left(\mu^{2}+\frac{v_{4}^{2}}{2}\right)\otimes e_{2}\,,$

which has a critical point at every point $(v_{4}=0,v_{6})$, where $v_{6}$ is
arbitrary. Such points are not stable equilibrium points. In particular, this
$\bm{\omega}_{\bm{\mu}}$-Hamiltonian function does not satisfy that
$f^{1}_{\bm{\mu}}+f^{2}_{\bm{\mu}}$ has a strict minimum at the equilibrium
point: it has only a minimum. Note that the points in ${\bf
J}^{\Phi-1}(\bm{\mu})$ that project onto the above-mentioned equilibrium
points are two-polysymplectic relative equilibrium points. The analysis of
(5.35) with our two-polysymplectic energy-momentum methods at the mentioned
two-polysymplectic relative equilibrium points suggests the same results.

###  5.5 Equilibrium points and vector fields with polynomial coefficients

Let us illustrate certain aspects of our $k$-polysymplectic energy-momentum
method by studying vector fields with a polynomial behaviour. Moreover, our
example will illustrate some features of weakly regular points of
$k$-polysymplectic momentum maps and the character of their associated
$k$-polysymplectic Marsden–Weinstein reductions.

Consider coordinates $\\{x_{1},x_{2},x_{3},x_{4},x_{5},x_{6},x_{7},x_{8}\\}$
on $\mathbb{R}^{8}$ and the vector field $X$ on $\mathbb{R}^{8}$ given by

$X=x_{6}^{a}\frac{\partial}{\partial x_{2}}+x_{4}^{b}\frac{\partial}{\partial
x_{3}}-x_{3}^{c}\frac{\partial}{\partial
x_{4}}+x_{8}^{d}\frac{\partial}{\partial
x_{7}}-x_{7}^{e}\frac{\partial}{\partial x_{8}}\,,$

where $a,b,c,d,e\in\mathbb{N}$. Define the two-polysymplectic for
${\bm{\omega}}$ on $\mathbb{R}^{8}$ of the form

$\bm{\omega}=\omega^{1}\otimes e_{1}+\omega^{2}\otimes
e_{2}=(\mathrm{d}x_{3}\wedge\mathrm{d}x_{4}+\mathrm{d}x_{1}\wedge\mathrm{d}x_{5})\otimes
e_{1}+(\mathrm{d}x_{2}\wedge\mathrm{d}x_{6}+\mathrm{d}x_{7}\wedge\mathrm{d}x_{8})\otimes
e_{2}\,.$

Then,

$\ker\omega_{x}^{1}=\left\langle\frac{\partial}{\partial
x_{2}},\frac{\partial}{\partial x_{6}},\frac{\partial}{\partial
x_{7}},\frac{\partial}{\partial
x_{8}}\right\rangle,\qquad\ker\omega^{2}_{x}=\left\langle\frac{\partial}{\partial
x_{1}},\frac{\partial}{\partial x_{3}},\frac{\partial}{\partial
x_{4}},\frac{\partial}{\partial
x_{5}}\right\rangle,\qquad\ker\omega^{1}_{x}\cap\ker\omega^{2}_{x}=0$

for any $x\in\mathbb{R}^{8}$, and thus $\bm{\omega}$ becomes a two-
polysymplectic form.

The vector field $X$ admits the Lie symmetries,
$Y_{1}=\frac{\partial}{\partial x_{2}}\,,\ Y_{2}=\frac{\partial}{\partial
x_{1}}\,,\ Y_{3}=\frac{\partial}{\partial x_{5}}$, which span a three-
dimensional abelian Lie algebra of vector fields. These Lie symmetries are the
infinitesimal generators of the translations along the $x_{2}$, $x_{1}$, and
$x_{5}$ coordinates, and they also leave the two-polysymplectic form
invariant, i.e. $\mathscr{L}_{Y_{i}}\omega^{\alpha}=0$ for $i=1,2,3$ and
$\alpha=1,2$. They give rise to a Lie group action
$\Phi:\mathbb{R}^{3}\times\mathbb{R}^{8}\rightarrow\mathbb{R}^{8}$ that leaves
invariant $\bm{\omega}$.

Since

$\displaystyle\iota_{Y_{1}}\omega^{1}$ $\displaystyle=0\,,$
$\displaystyle\iota_{Y_{2}}\omega^{1}$ $\displaystyle=\mathrm{d}x_{5}\,,$
$\displaystyle\iota_{Y_{3}}\omega^{1}$ $\displaystyle=-\mathrm{d}x_{1}\,,$
(5.36) $\displaystyle\iota_{Y_{1}}\omega^{2}$
$\displaystyle=\mathrm{d}x_{6}\,,$ $\displaystyle\iota_{Y_{2}}\omega^{2}$
$\displaystyle=0\,,$ $\displaystyle\iota_{Y_{3}}\omega^{2}$
$\displaystyle=0\,,$ (5.37)

a two-polysymplectic momentum map $\mathbf{J}^{\Phi}$ can be defined by
setting

$\mathbf{J}^{\Phi}:\mathbb{R}^{8}\ni
x\longmapsto\mathbf{J}^{\Phi}(x)=(0,x_{5},-x_{1};x_{6},0,0)\in(\mathbb{R}^{3*})^{2}\simeq(\mathbb{R}^{3})^{2}\,.$

Then,
$\mathrm{T}_{x}\mathbf{J}^{\Phi-1}(\bm{\mu})=\left\langle\frac{\partial}{\partial
x_{2}},\frac{\partial}{\partial x_{3}},\frac{\partial}{\partial
x_{4}},\frac{\partial}{\partial x_{7}},\frac{\partial}{\partial
x_{8}}\right\rangle$ for each $x\in\mathbf{J}^{\Phi-1}(\bm{\mu})$ and
$\bm{\mu}\in(\mathbb{R}^{3})^{2}$. The two-polysymplectic momentum map
$\mathbf{J}^{\Phi}$ is $\operatorname{Ad}^{*2}$-equivariant and every
$\bm{\mu}\in(\mathbb{R}^{3})^{2}$ is a weakly regular value of
$\mathbf{J}^{\Phi}$ since each $\mathbf{J}^{\Phi-1}(\bm{\mu})$ is a five-
dimensional submanifold of $\mathbb{R}^{8}$ and its tangent space at each
point coincides with the kernel of ${\bf J}^{\Phi}$ at that point. Moreover,
${\bf J}^{\Phi}$ has no regular points.

Note that $Y_{2}$ and $Y_{3}$ do not take values at $x$ in
$\mathrm{T}_{x}(G_{\bm{\mu}}x)$ but $Y_{1}$ does. The assumptions of Theorem
3.2 are satisfied, and the quotient space
$\mathrm{T}_{x}\mathbf{J}^{\Phi-1}(\bm{\mu})/\mathrm{T}_{x}\left(G_{\bm{\mu}}x\right)$
is a two-dimensional subspace, where

$\mathrm{T}_{x}\mathbf{J}^{\Phi-1}(\bm{\mu})/\mathrm{T}_{x}\left(G_{\bm{\mu}}x\right)=\left\langle\frac{\partial}{\partial
x_{3}},\frac{\partial}{\partial x_{4}},\frac{\partial}{\partial
x_{7}},\frac{\partial}{\partial x_{8}}\right\rangle\,,\qquad\forall x\in{\bf
J}^{\Phi-1}(\bm{\mu})$

and

$\bm{\omega}_{\bm{\mu}}=\omega^{1}_{\bm{\mu}}\otimes
e_{1}+\omega^{2}_{\bm{\mu}}\otimes
e_{2}=(\mathrm{d}x_{3}\wedge\mathrm{d}x_{4})\otimes
e_{1}+(\mathrm{d}x_{7}\wedge\mathrm{d}x_{8})\otimes e_{2}\,.$

The vector field $X$ is $\bm{\omega}$-Hamiltonian relative to

$\displaystyle\mathrm{d}\bm{h}$
$\displaystyle=\iota_{X}{\bm{\omega}}=\iota_{X}\omega^{1}\otimes
e_{1}+\iota_{X}\omega^{2}\otimes e_{2}$
$\displaystyle=\mathrm{d}\left(\frac{1}{1+b}x^{b+1}_{4}+\frac{1}{c+1}x_{3}^{c+1}\right)\otimes
e_{1}+\mathrm{d}\left(\frac{1}{1+a}x^{a+1}_{6}+\frac{1}{d+1}x_{8}^{d+1}+\frac{1}{1+e}x_{7}^{e+1}\right)\otimes
e_{2}\,.$

Moreover, $\bm{h}$ is invariant relative to the Lie symmetries
$Y_{1},Y_{2},Y_{3}$. By Theorem 3.3, the vector field $X$ projects onto the
quotient manifold and its projection $X_{\bm{\mu}}$ is given by

$X_{\bm{\mu}}=x_{4}^{b}\frac{\partial}{\partial
x_{3}}-x_{3}^{c}\frac{\partial}{\partial
x_{4}}+x_{8}^{d}\frac{\partial}{\partial
x_{7}}-x_{7}^{e}\frac{\partial}{\partial x_{8}}\,,$

which is an $\bm{\omega}_{\bm{\mu}}$-Hamiltonian vector field since

$\mathrm{d}\bm{f}_{\bm{\mu}}=\iota_{X_{\bm{\mu}}}\bm{\omega}_{\bm{\mu}}=\mathrm{d}\left(\frac{1}{1+b}x^{b+1}_{4}+\frac{1}{c+1}x_{3}^{c+1}\right)\otimes
e_{1}+\mathrm{d}\left(\frac{1}{d+1}x_{8}^{d+1}+\frac{1}{1+e}x_{7}^{e+1}\right)\otimes
e_{2}\,.$

According to Theorem 4.2, a point $z_{e}$ is a two-polysymplectic relative
equilibrium point if it is a critical point of $\bm{h}_{\xi}$ for some
$\xi=(\xi_{1},\xi_{2},\xi_{3})\in\mathfrak{g}\simeq\mathbb{R}^{3}$. Then, one
has that

$\displaystyle\mathrm{d}\bm{h}_{\xi}$
$\displaystyle=\mathrm{d}f^{1}_{\xi}\otimes e_{1}+\mathrm{d}f^{2}_{\xi}\otimes
e_{2}$
$\displaystyle=\left(x_{4}^{b}\mathrm{d}x_{4}+x^{c}_{3}\mathrm{d}x_{3}-\xi_{2}\mathrm{d}x_{5}+\xi_{3}\mathrm{d}x_{1}\right)\otimes
e_{1}+\left((x^{a}_{6}-\xi_{1})\mathrm{d}x_{6}+x_{8}^{d}\mathrm{d}x_{8}+x_{7}^{e}\mathrm{d}x_{7}\right)\otimes
e_{2}\,.$

It follows that $\xi_{2}=\xi_{3}=0$ and we have two-polysymplectic relative
equilibrium points of $X$ of the form
$z_{e}=(x_{1},x_{2},0,0,x_{5},x_{6},0,0)$ for $x^{a}_{6}=\xi_{1}$. In fact,
$x_{1},x_{2},x_{5},x_{6}$ are arbitrary. Indeed,
$(X_{\bm{\mu}_{e}})_{[z_{e}]}=0$ for $\bm{\mu}_{e}={\bf J}^{\Phi}(z_{e})$.

To analyse the stability of the above-mentioned two-polysymplectic relative
equilibrium points, let us analyse the second derivatives of $\bm{h}_{\xi}$ at
such points $z_{e}$. Then,

$\displaystyle(\delta^{2}\bm{h}_{\xi})_{z_{e}}$
$\displaystyle=(\delta^{2}h^{1}_{\xi})_{z_{e}}\otimes
e_{1}+(\delta^{2}h^{2}_{\xi})_{z_{e}}\otimes e_{2}$
$\displaystyle=\left(cx_{3}^{c-1}\mathrm{d}x_{3}\otimes\mathrm{d}x_{3}+bx_{4}^{b-1}\mathrm{d}x_{4}\otimes\mathrm{d}x_{4}\right)\otimes
e_{1}+\left(ex_{7}^{e-1}\mathrm{d}x_{7}\otimes\mathrm{d}x_{7}+dx_{8}^{d-1}\mathrm{d}x_{8}\otimes\mathrm{d}x_{8}\right)\otimes
e_{2}\,.$

Taking into account that the supplementary spaces to
$\mathrm{T}_{z_{e}}G_{\bm{\mu}_{e}}+\ker\omega^{1}_{z_{e}}\cap\mathrm{T}_{z_{e}}({\bf
J}^{\Phi-1}(\bm{\mu}_{e}))$ and
$\mathrm{T}_{z_{e}}G_{\bm{\mu}_{e}}+\ker\omega^{2}_{z_{e}}\cap\mathrm{T}_{z_{e}}({\bf
J}^{\Phi-1}(\bm{\mu}_{e}))$ can be chosen
$\mathcal{S}_{z_{e}}^{1}=\langle\frac{\partial}{\partial
x_{3}},\frac{\partial}{\partial x_{4}}\rangle$ and
$\mathcal{S}^{2}_{z_{e}}=\langle\frac{\partial}{\partial
x_{7}},\frac{\partial}{\partial x_{8}}\rangle$, respectively, Definition 4.6
and the posterior explanation give that the $z_{e}$ are stable two-
polysymplectic relative equilibrium points if

$(\delta^{2}h^{1}_{\xi})_{z_{e}}(v_{z_{e}},v_{z_{e}})>0\,,\qquad\forall
v_{z_{e}}\in\mathcal{S}^{1}_{z_{e}}\setminus\\{0\\}\,,$

and

$(\delta^{2}h^{2}_{\xi})_{z_{e}}(v_{z_{e}},v_{z_{e}})>0\,,\qquad\forall
v_{z_{e}}\in\mathcal{S}^{2}_{z_{e}}\setminus\\{0\\}\,.$

These inequalities hold if and only if $b,c,d,e=1$. Hence, the $z_{e}$ are
formally stable two-polysymplectic relative equilibrium points of $X$ for
$b,c,d,e=1$. Indeed, it is immediate that previous conditions of stability
imply that in the reduced space close to the equilibrium point, the
coordinates $x_{3},x_{4},x_{7},x_{8}$ are bounded for every motion close
enough to the equilibrium point, which ensures real stability.

##  6 Conclusions and Outlook

In this work, we have devised a new energy-momentum method for systems of
ordinary differential equations given by Hamiltonian vector fields with
respect to a $k$-polysymplectic form. This led to defining and characterising
$k$-polysymplectic relative equilibrium points and introducing new techniques
to study stability through $k$-polysymplectic geometry. In this respect, we
have also reviewed several aspects and mistakes in previous Marsden–Weinstein
reductions for $k$-polysymplectic forms and systems [7, 27, 46, 37], which is
a relevant part in energy-momentum methods. To illustrate our new energy-
momentum method and theoretical results, we have studied several relevant
examples in detail: complex Schwarz equations, the product of several
symplectic manifolds, with a family of particles subjected to the effect of an
isotropic potential for each of them, affine homogeneous differential
equations (with potential applications to control theory and Lie systems
theory), a quantum harmonic oscillator, integrable symplectic systems, and
some dynamical systems with polynomial coefficients.

A non-autonomous analogue of the methods devised in this paper could be
accomplished by using the Lyapunov theory depicted in [25]. Note that the
stability with respect to $k$-polysymplectic forms is a topic that requires
further development. The criteria here used are enough for the family of
examples to be studied, but a deeper study with an analysis of all
possibilities is in order. Recall also that, even in the one-polysymplectic
case, the criterion for the stability of the energy-momentum method, which we
here recover as a particular case and is classical [2], is not a necessary
condition for the stability of the reduced system.

The study of complex Schwarz equations and the Schwarzian derivative could be
more appropriately studied through a complex Lie system formalism. We aim to
study this possibility in further work.

## Acknowledgements

J. de Lucas, X. Rivas and B.M. Zawora acknowledge partial financial support
from the Nowe Idee 2B-POB II project PSP: 501-D111-20-2004310 funded by the
“Inicjatywa Doskonałości - Uczelnia Badawcza” (IDUB) program. J. de Lucas
acknowledges a CRM-Simons professorship funded by the Simons Foundation and
the Centre de Recherches Mathématiques (CRM) of the Université Montréal. X.
Rivas acknowledges funding from J. de Lucas’s CRM-Simons professorship to
accomplish a stay at the CRM which allowed us to obtain several relevant
findings of the present work. X. Rivas also acknowledges funding from the
Spanish Ministry of Science and Innovation, grants PID2021-125515NB-C21, and
RED2022-134301-T of AEI, and Ministry of Research and Universities of the
Catalan Government, project 2021 SGR 00603. B.M. Zawora acknowledges funding
from the IDUB mikrogrant program to accomplish a research stay at the CRM. L.
Colombo acknowledges financial support from Grant PID2022-137909NB-C21 funded
by MCIN/AEI/ 10.13039/501100011033. B.M. Zawora and X. Rivas acknowledge
partial funding from the Simons Foundation to take part in a conference in the
Simons Center for Geometry and Physics of Stony Brook University, where the
last parts of this paper were accomplished and several relevant results
concerning the stability of reduced $k$-polysymplectic Hamiltonian systems and
their applications were devised.

## References

  * [1] H. Abarbanel and D. Holm.  Nonlinear stability analysis of inviscid flows in three dimensions: incompressible fluids and barotropic fluids.  Phys. Fluids, 30(11):3369–3382, 1987.  10.1063/1.866469.
  * [2] R. Abraham and J. E. Marsden.  Foundations of mechanics, volume 364 of AMS Chelsea publishing.  Benjamin/Cummings Pub. Co., New York, 2nd edition, 1978.  10.1090/chel/364.
  * [3] C. Albert.  Le théorème de réduction de Marsden–Weinstein en géométrie cosymplectique et de contact.  J. Geom. Phys., 6(4):627–649, 1989.  10.1016/0393-0440(89)90029-6.
  * [4] A. Awane.  $k$-symplectic structures.  J. Math. Phys., 33(12):4046, 1992.  10.1063/1.529855.
  * [5] A. Awane and M. Goze.  Pfaffian systems, $k$-symplectic systems.  Springer, Dordrecht, 1st edition, 2000.  10.1007/978-94-015-9526-1.
  * [6] C. Bai and H. Zhang.  A weak energy-momentum method for stochastic instability induced by dissipation and random excitations.  Probabilistic Eng. Mech., 37:35–40, 2014.  10.1016/j.probengmech.2014.03.009.
  * [7] C. Blacker.  Polysymplectic reduction and the moduli space of flat connections.  J. Phys. A: Math. Theor., 52(33):335201, 2019.  10.1088/1751-8121/ab2eed.
  * [8] A. M. Blaga.  Connections on $k$-symplectic manifolds.  Balk. J. Geom. Appl., 14(2):28–33, 2009.  https://www.emis.de/journals/BJGA/v14n2/B14-2-bl.pdf.
  * [9] G. Blankenstein and T. S. Ratiu.  Singular reduction of implicit Hamiltonian systems.  Rep. Math. Phys., 53(2):211–260, 2004.  10.1016/S0034-4877(04)90013-4.
  * [10] A. M. Bloch.  Nonholonomic mechanics and control, volume 24 of Interdisciplinary Applied Mathematics.  Springer, New York, 2003.  10.1007/978-1-4939-3017-3.
  * [11] Y. D. Bozhkov and P. R. da Conceição.  On the Generalizations of the Kummer–Schwarz Equation.  Nonlinear Anal., 192:111691, 2020.  10.1016/j.na.2019.111691.
  * [12] L. Búa, I. Bucataru, M. de León, M. Salgado, and S. Vilariño.  Symmetries in Lagrangian field theory.  Rep. Math. Phys., 75(3):333–357, 2015.  10.1016/S0034-4877(15)30010-0.
  * [13] J. F. Cariñena and J. de Lucas.  Lie systems: theory, generalisations, and applications.  Dissertationes Math., 479:1–162, 2011.  10.4064/dm479-0-1.
  * [14] J. F. Cariñena, J. Grabowski, and G. Marmo.  Lie–Scheffers systems: a geometric approach.  Napoli Series on Physics and Astrophysics. Bibliopolis, Naples, 2000.
  * [15] J. F. Cariñena, J. Grabowski, and G. Marmo.  Superposition rules, lie theorem, and partial differential equations.  Rep. Math. Phys., 60(2):237–258, 2007.  10.1016/S0034-4877(07)80137-6.
  * [16] J. F. Cariñena and A. Ramos.  Applications of Lie systems in quantum mechanics and control theory.  Banach Center Publ., 59:143–162, 2003.  10.4064/bc59-0-7.
  * [17] M. de León, E. Merino, J. A. Oubiña, and M. Salgado.  Stable almost cotangent structures.  Bolletino Unione Mat. Ital. B (7), 11(3):509–529, 1997.
  * [18] M. de León, I. Méndez, and M. Salgado.  $p$-almost tangent structures.  Rend. Circ. Mat. Palermo, 37(2):282–294, 1988.  10.1007/BF02844526.
  * [19] M. de León, I. Méndez, and M. Salgado.  Regular $p$-almost cotangent structures.  J. Korean Math. Soc., 25(2):273–287, 1988.  https://jkms.kms.or.kr/journal/view.html?spage=273&volume=25&number=2.
  * [20] M. de León, M. Salgado, and S. Vilariño.  Methods of Differential Geometry in Classical Field Theories.  World Scientific, 2015.  10.1142/9693.
  * [21] J. de Lucas, A. Maskalaniec, and B. M. Zawora.  A cosymplectic energy-momentum method with applications.  2302.05827, 2023.
  * [22] J. de Lucas, X. Rivas, S. Vilariño, and B. M. Zawora.  On $k$-polycosymplectic Marsden–Weinstein reductions.  J. Geom. Phys., 191:104899, 2023.  10.1016/j.geomphys.2023.104899.
  * [23] J. de Lucas and C. Sardón.  A Guide to Lie Systems with Compatible Geometric Structures.  World Scientific Publishing Co. Pte. Ltd., Singapore, 2020.  10.1142/q0208.
  * [24] J. de Lucas and S. Vilariño.  $k$-symplectic Lie systems: theory and applications.  J. Differ. Equ., 258(6):2221–2255, 2015.  10.1016/j.jde.2014.12.005.
  * [25] J. de Lucas and B. M. Zawora.  A time-dependent energy-momentum method.  J. Geom. Phys., 170:104364, 2021.  10.1016/j.geomphys.2021.104364.
  * [26] A. Echeverría-Enríquez, M. C. Muñoz-Lecanda, and N. Román-Roy.  Geometry of Lagrangian first-order classical field theories.  Fortschr. Phys., 44(3):235–280, 1996.  10.1002/prop.2190440304.
  * [27] E. García-Toraño Andrés and T. Mestdag.  Conditions for symmetry reduction of polysymplectic and polycosymplectic structures.  J. Phys. A: Math. Theor., 56(33):335202, 2023.  10.1088/1751-8121/ace74c.
  * [28] G. Giachetta, L. Mangiarotti, and G. A. Sardanashvily.  New Lagrangian and Hamiltonian Methods in Field Theory.  World Scientific, River Edge, 1997.  10.1142/2199.
  * [29] X. Gràcia, J. de Lucas, M. C. Muñoz-Lecanda, and S. Vilariño.  Multisymplectic structures and invariant tensors for Lie systems.  J. Phys. A: Math. Theor., 52(21):215201, 2019.  10.1088/1751-8121/ab15f2.
  * [30] L. Guieu and C. Roger.  L’algèbre et le Groupe de Virasoro. Aspects géométriques et algébriques, généralisations.  Les Publications CRM, Montréal, 2007.
  * [31] C. Günther.  The polysymplectic Hamiltonian formalism in field theory and calculus of variations I: The local case.  J. Differ. Geom., 25(1):23–53, 1987.  10.4310/jdg/1214440723.
  * [32] E. Hille.  Ordinary differential equations in the complex domain.  Dover Publications, Mineola, NY, 1997.
  * [33] I. V. Kanatchikov.  Canonical structure of classical field theory in the polymomentum phase space.  Rep. Math. Phys., 41(1):49–90, 1998.  10.1016/S0034-4877(98)80182-1.
  * [34] J. M. Lee.  Manifolds and differential geometry, volume 107 of Graduate Studies in Mathematics.  American Mathematical Society, Providence, RI, 2009.
  * [35] J. M. Lee.  Introduction to Smooth Manifolds, volume 218 of Graduate Texts in Mathematics.  Springer New York Heidelberg Dordrecht London, 2nd edition, 2012.  10.1007/978-1-4419-9982-5.
  * [36] O. Lehto.  Remarks on Nehari’s theorem about the Schwarzian derivative and schlicht functions.  J. Anal. Math., 36:184–190, 1979.  10.1007/BF02798778.
  * [37] J. C. Marrero, N. Román-Roy, M. Salgado, and S. Vilariño.  Reduction of polysymplectic manifolds.  J. Phys. A: Math. Theor., 48(5):055206, 2015.  10.1088/1751-8113/48/5/055206.
  * [38] J. C. Marrero, N. Román-Roy, M. Salgado, and S. Vilariño.  Symmetries, Noether’s theorem and reduction in k-cosymplectic field theories.  AIP Conference Proceedings, 1260(1):173–179, 2010.  10.1063/1.3479319.
  * [39] J. E. Marsden, T. A. Posbergh, and J. C. Simo.  Stability of coupled rigid body and geometrically exact rods: block diagonalization and the energy–momentum method.  Phys. Rep., 193:280–360, 1990.  10.1016/0370-1573(90)90125-L.
  * [40] J. E. Marsden and T. Ratiu.  Introduction to Mechanics and Symmetry. A basic exposition of classical mechanical systems, volume 17 of Texts in Applied Mathematics.  Springer-Verlag, New York, 1999.  10.1007/978-0-387-21792-5.
  * [41] J. E. Marsden and J. C. Simo.  The energy momentum method.  Act. Acad. Sci. Tau., 1(124):245–268, 1988.  https://resolver.caltech.edu/CaltechAUTHORS:20101019-093814651.
  * [42] J. E. Marsden, J. C. Simo, D. Lewis, and T. Posbergh.  A block diagonalization theorem in the energy-momentum method.  In Dynamics and control of multibody systems, volume 97 of Contemp. Math., page 297–313. Amer. Math. Soc., Providence, RI, 1989.  10.1090/conm/097/1021043.
  * [43] N. Martinez.  Poly-symplectic groupoids and poly-Poisson structures.  Lett. Math. Phys., 105(5):693–721, 2015.
  * [44] M. McLean and L. K. Norris.  Covariant field theory on frame bundles of fibered manifolds.  J. Math. Phys., 41(10):6808–6823, 2000.  10.1063/1.1288797.
  * [45] E. E. Merino.  Geometría $k$-simpléctica y $k$-cosimpléctica. Aplicaciones a las teorías clásicas de campos.  PhD thesis, Universidade de Santiago de Compostela, 1997.
  * [46] F. Munteanu, A. M. Rey, and M. Salgado.  The Günther’s formalism in classical field theory: momentum map and reduction.  J. Math. Phys., 45(5):1730–1751, 2004.  10.1063/1.1688433.
  * [47] L. K. Norris.  Generalized symplectic geometry on the frame bundle of a manifold.  In Proc. Symp. Pure Math., volume 54.2, page 435–465. Amer. Math. Soc., Providence RI, 1993.  10.1090/pspum/054.2.
  * [48] J. Ortega, V. Planas-Bielsa, and T. Ratiu.  Asymptotic and Lyapunov stability of constrained and Poisson equilibria.  J. Differ. Equ., 214(1):92–127, 2005.  10.1016/j.jde.2004.09.016.
  * [49] J. P. Ortega and T. S. Ratiu.  Momentum maps and Hamiltonian reduction, volume 222 of Progress in Mathematics.  Birkhäuser Boston, Inc., Boston, 2004.  10.1007/978-1-4757-3811-7.
  * [50] N. Román-Roy, M. Salgado, and S. Vilariño.  Symmetries and conservation laws in the Günther $k$-symplectic formalism of field theory.  Rev. Math. Phys., 19(10):1117–1147, 2007.  10.1142/S0129055X07003188.
  * [51] G. A. Sardanashvily.  Generalized Hamiltonian formalism for field theory. Constraint systems.  World Scientific Publishing Company, River Edge, NJ, 1995.  10.1142/2550.
  * [52] J. Simo, D. Lewis, and J. Marsden.  Stability of relative equilibria. Part I: The reduced energy-momentum method.  Arch. Rat. Mech. Anal., 115:15–59, 1991.  10.1007/BF01881678.
  * [53] J. Simo and N. Tarnow.  The discrete energy-momentum method. Conserving algorithms for nonlinear elastodynamics.  Z. Angew. Math. Phys., 43:757–792, 1992.  10.1007/BF00913408.
  * [54] M. Vidyasagar.  Nonlinear systems analysis, volume 42 of Classics in Applied Mathematics.  Society for Industrial and Applied Mathematics (SIAM), Philadelphia, PA, 2002.  10.1137/1.9780898719185.
  * [55] L. Wang and P. Krishnaprasad.  Gyroscopic control and stabilization.  J. Nonlinear Sci., 2:367–415, 1992.  10.1007/BF01209527.
  * [56] J. Wei and E. Norman.  Lie Algebraic Solution of Linear Differential Equations.  J. Math. Phys., 4(4):575–581, 1963.  10.1063/1.1703993.
  * [57] P. Winternitz.  Lie groups and solutions of nonlinear differential equations.  In K. B. Wolf, editor, Nonlinear Phenomena, volume 189 of Lecture Notes in Physics, page 263–305. Springer, Berlin, Heidelberg, 1983.  10.1007/3-540-12730-5_12.
  * [58] B. M. Zawora.  A time-dependent energy-momentum method.  Master’s thesis, University of Warsaw, Faculty of Physics, 2021.
  * [59] D. V. Zenkov, A. M. Bloch, and J. E. Marsden.  The energy-momentum method for the stability of non-holonomic systems.  Dyn. Stab. Syst., 13(2):123–165, 1998.  10.1080/02681119808806257.


# Point Spread Function Estimation for Wide Field Small Aperture Telescopes
with Deep Neural Networks and Calibration Data

Peng Jia1,4, Xuebo Wu1, Zhengyang Li 2 , Bo Li2, Weihua Wang1, Qiang Liu1,
robinmartin20<EMAIL_ADDRESS>Adam Popowicz3 and Dongmei Cai1
1College of Physics and Optoelectronics, Taiyuan University of Technology,
Taiyuan, 030024, China
2Nanjing Institute of Astronomical Optics and Technology CAS, Nanjing,
Jiangsu, 210042, China
3Department of Electronics, Electrical Engineering and Microelectronics,
Silesian University of Technology, Akademicka 16, 44-100 Gliwice, Poland
4Department of Physics, Durham University, DH1 3LE, UK
<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

The point spread function (PSF) reflects states of a telescope and plays an
important role in development of data processing methods, such as PSF based
astrometry, photometry and image restoration. However, for wide field small
aperture telescopes (WFSATs), estimating PSF in any position of the whole
field of view is hard, because aberrations induced by the optical system are
quite complex and the signal to noise ratio of star images is often too low
for PSF estimation. In this paper, we further develop our deep neural network
(DNN) based PSF modelling method and show its applications in PSF estimation.
During the telescope alignment and testing stage, our method collects system
calibration data through modification of optical elements within engineering
tolerances (tilting and decentering). Then we use these data to train a DNN
(Tel–Net). After training, the Tel–Net can estimate PSF in any field of view
from several discretely sampled star images. We use both simulated and
experimental data to test performance of our method. The results show that the
Tel–Net can successfully reconstruct PSFs of WFSATs of any states and in any
positions of the FoV. Its results are significantly more precise than results
obtained by the compared classic method - Inverse Distance Weight (IDW)
interpolation. Our method provides foundations for developing of deep neural
network based data processing methods for WFSATs, which require strong prior
information of PSFs.

###### keywords:

telescopes – methods: numerical – techniques: image processing

††pubyear: 2021††pagerange: Point Spread Function Estimation for Wide Field
Small Aperture Telescopes with Deep Neural Networks and Calibration Data–Point
Spread Function Estimation for Wide Field Small Aperture Telescopes with Deep
Neural Networks and Calibration Data

##  1 Introduction

Telescopes with wide FoV (normally with a diameter greater than 1 deg) and
small apertures (less than 2 metre) are called WFSATs. They are low-cost,
light-weighted and commonly used for wide field sky surveys in time domain
astronomy (Burd et al., 2005; Ma et al., 2007; Yuan et al., 2008; Cui et al.,
2008; Pablo et al., 2016; Ping & Zhang, 2017; Ratzloff et al., 2019; Sun & Yu,
2019). Because apertures of WFSATs are small, they are normally used to
observe bright targets. To further increase detection and observation
abilities of WFSATs, new data processing methods with strong prior conditions
have recently been developed and proven their superiority over existing image
processing techniques (Jia et al., 2020b; Glazier et al., 2020).

The PSF is the impulse response of an optical system and it reflects states of
the whole telescope. For optical observation data processing methods, PSF is a
commonly used prior condition. Prior conditions of PSFs mainly include shape
and variations of PSFs within the whole FoV. For general purpose sky survey
telescopes, different data processing methods have been recently tested on
images with various PSFs and have successfully proved its robustness (González
et al., 2018a; Duev et al., 2019; Burke et al., 2019a; He et al., 2020).
However, PSFs of WFSATs are quite different from that of general purpose sky
survey telescopes, because WFSATs are often working remotely (lack of
maintenance) and images obtained by them have both low sampling rate and
fluctuating quality (e.g. due to the high sensitivity of some optical elements
to temperature and/or humidity). PSFs of WFSATs have significant deformation
in different regions of FoV (Piotrowski & Żarnecki, 2010) and they can evolve
with time (Jia et al., 2017). It is very hard to design a fully deterministic
image processing method that can model spatially and temporally PSFs, so more
sophisticated PSF estimation methods are required.

The statistical PSF modelling method is commonly used for WFSATs nowadays. The
statistical methods require large amount of star images as samples of PSFs and
obtain effective components from these images through either principal
component analysis (Jee et al., 2007; Jia et al., 2017; Popowicz, 2018; Sun et
al., 2020) or denoising autoencoders (Jia et al., 2020c). Because large amount
of star images with adequate signal to noise ratio (SNR) are required (Wang et
al., 2018), application of such methods is frequently unsuitable for WFSATs
that are used for fast sky survey.

The forward PSF modelling method, which was firstly applied in space based
telescopes (Krist et al., 2011; Perrin et al., 2012) and later in telescopes
with adaptive optics systems (Martin et al., 2016; Fétick et al., 2019;
Beltramo-Martin et al., 2019; Fusco et al., 2020; Jia et al., 2021b), shows
good performance in modelling PSFs generated by the atmospheric turbulence
with the help of DNN (Jia et al., 2020d; Jia et al., 2020e). The forward DNN
based PSF modelling method (PSF–NET) builds the Monte Carlo model of the
imaging process and generates huge amount of PSFs to train the DNN. After
training, the PSF–NET learns a function that can reflect the repose of a
telescope (PSF) to complex external environments (such as the shape of
wavefront distortion due to the atmospheric turbulence).

For WFSATs, because all lenses and mirrors are manufactured and assembled with
predefined tolerances, PSFs of WFSATs can be interpreted as the response of an
imaging system to the actual state of optical misalignment together with the
actual quality of fabricated optical elements (Li et al., 2015). The
distribution of PSFs with different shapes also reflects the impact of
external conditions during imaging, like temperature, humidity or wind
strength, which can lead to very complex deformations of PSFs (Jia et al.,
2020a). Therefore, we decided to design the PSF estimation technique using
current capabilities of artificial intelligence since such algorithms have
already proven to be successful in realization of complex tasks in plenty of
various fields.

In this paper, we propose the Tel–Net as a novel PSF modelling method for
WFSATs. With tolerances of each optical element in a WFSAT as prior
information, the Tel-Net is trained on various possible realizations of the
optical alignment. After training, the Tel–Net can output PSFs of a WFSAT in
any state and in any position within the FoV. We introduce the concept and
structure of the Tel–Net in Section 2. The performance of the Tel–Net is
tested with simulated data in Section 3 and with real data in Section 4. We
give the conclusions and propose the future work in Section 5.

##  2 The Tel-Net

The philosophy of the Tel–Net is quite simple. PSFs in the whole FoV are
representations of the telescope state (Li et al., 2015) which can be modelled
by a DNN (Jia et al., 2020a). After we train the Tel–Net with some PSFs as
inputs and all PSFs within the whole FoV as outputs, the Tel–Net will be able
to output all PSFs of a WFSAT in different states within the whole FoV with
only some PSFs given at its input. The structure of the Tel–Net is adopted
from the PSF–NET proposed in Jia et al. (2020a) and the DAE–NET proposed in
Jia et al. (2020c). As shown in figure 1, the Tel–Net has an encoding and
decoding structure with DNN blocks adopted to increase its representation
ability.

Figure 1: The structure of the Tel–Net. The Tel–Net consists of convolutional
layers (Conv2d in blue), Bias layers (Bias in brown), Activation layers (Relu
in purple), Conv-transpose layers (ConvTranspose2d in red), multiply layers
(multiplier in grey) and Fixup blocks (Fixup block in green). The structure of
Fixup block is shown in the right part of this figure.

Basic components in the encoding part of Tel–Net is convolutional neural
network (CNN). The CNN uses trainable shared-weight architecture to extract
features from two multi-dimensional data, which makes it suitable for multi-
dimensional data processing, particularly for images. The CNN is widely used
for tasks such as: astronomical target detection (González et al., 2018b;
Burke et al., 2019b; Jia et al., 2020b; Hausen & Robertson, 2020),
classification (Jia et al., 2019a; Turpin et al., 2020; He et al., 2020) and
image restoration (Schawinski et al., 2017; Jia et al., 2019b, 2021a). In the
Tel–Net, the CNN is used to extract features of PSFs between pixels in the
same PSF and pixels of different PSFs in different positions of FoV. A bias
layer and a Relu layer are added after each convolutional layer. The bias
layer adds additional values to outputs of the previous layer to make the
activation function of the Relu layer effective. The Relu layer would applies
function $f(x)=max(0,x)$ to inputs and it would add non–linear properties to
the DNN (where $x$ stands for the input of Relu layer and $max$ stands for
function to return maximal value between $x$ and 0).

In the decoding part of the Tel–Net, transposed convolutional layer
($ConvTranspose2d$), bias layer and Relu layer are arranged to form basic
structure. The transpose convolutional layer also has shared-weight
architecture and it generates outputs with dimensions that are larger than
dimensions of inputs. With the encoding and decoding structure, the Tel–Net
can estimate PSFs in different positions of the FoV from input PSFs.

The deeper the neural networks are, the better it would be able to model
complex functions (Goodfellow et al., 2016). To increase the representation
ability of the Tel–Net, we introduce Fixup blocks in the Tel–Net as DNN
blocks. The structure of Fixup blocks is shown in right panel of figure 1,
which uses fixed–update initialization, also known as Fixup method proposed by
Zhang et al. (2019). The Fixup method is a regularization method that rescales
weights of all layers in DNNs properly at the beginning of training. Comparing
with other ordinary regularization methods, such as the normalization layers
(Ba et al., 2016), the Fixup method can accelerate convergence and improve
generalization without changing greyscale values of input images. Introduction
of the Fixup method would be benefit in applications of astronomical image
processing, because grey values of astronomical images have a wider range than
that in typical 8-bit images and original grey values instead of their
distributions are important for scientific data processing.

For the Tel–Net, PSFs of WFSATs in a given state are represented by a data
cube (PSF-Cube) with dimensions of $N\times N\times M$, where $N$ stands for
width and length of star images and $M$ stands for the number of PSFs that are
used to represent a given FoV. In this paper, we divide the FoV of WFSATs into
$M$ small sections with equal size (called patches) and we assume that PSFs
are uniform in each patch (La Camera et al., 2015). These PSFs are arranged
from left to right and top to bottom according to positions of patches in the
FoV. Each PSF in a PSF-Cube is aligned with its brightest pixel in the centre.
Besides, to reduce effects brought by sub-pixel uncertainties and to extend
the training set, we randomly shift PSFs by a fraction of a pixel (within -0.5
to 0.5 pixels with uniform probability distribution) when we obtain PSF–Cubes.

The input and output of the Tel–Net are both PSF-Cubes of the same size.
However, there are only several PSFs in input PSF-Cubes and others are zero
matrices. For output PSF–Cubes, all matrices are PSFs. Because the Tel–Net
should be able to predict all PSFs from several known PSFs, the input and
output PSF–Cubes that correspond to the same state of the WFSAT are used as a
pair of training data for the Tel–Net. During the training stage, the batch
size is set to 4 (due to limited available memory in the GPU – GTX 1080 in
this paper) and we choose the Adam Optimizer (Kingma & Ba, 2014) as the
optimization algorithm. The learning rate is set as 0.01 and the loss function
– $L_{1}$ norm is defined in equation 1, where $n$ equals to the batch size,
$TelNet$ $(PSFCube_{in})$ is the output PSF-Cube obtained from the Tel–Net and
the $PSFCube_{out}$ is expected output PSF–Cube, the norm $||$ utilizes the
absolute difference, $n$ stands for the total number of PSFs in the PSF-Cube
and $k$ stands for index of these PSFs.

$L=\frac{1}{n}\sum_{k=1}^{n}|TelNet(PSFCube_{in})-PSFCube_{out}|\\\ $ (1)

There are two possible ways to obtain PSFs to train the Tel–Net:
computation–based PSFs or data–driven PSFs. We can obtain PSF of any position
in the whole FoV for a WFSAT in any state through physical computation (Jia et
al., 2020a). However, because there are differences between PSFs obtained by
physical computation and PSFs in real observed images, the former can only be
used as pre–trained data. PSFs that are obtained from real observations
(data–driven PSFs) are required for real applications. Due to hardware
limitations or time consumptions in obtaining PSFs in any position of the FoV
during real observation, data–driven PSFs can only be obtained for a WFSAT
with finite states in several predefined positions in the FoV and with finite
accuracy. These problems will be further discussed in Section 3 and Section 4.

##  3 Training and testing the Tel-Net with simulated data

In this section, we simulate a classic WFSAT (Baker super-Schmidt telescope)
with ZEMAX to test the performance of the Tel–Net. Parameters of this WFSAT
are shown in figure 2. This WFSAT has a field of view of $10^{\circ}$ and we
divide it equally into $61\times 61$ patches. In each patch, variations of
PSFs are small thus a single PSF is representative. There are four optical
elements in this telescope and we assume that the second lens has alignment
error. There are six free dimensions for the second lens: decenter and tilt
along x, y and z directions. To make sure that star images in different places
within the FoV can be obtained by the WFSAT, we set the tolerance of the
second lens as: -3 millimetres to 3 millimetres with step of 1.5 millimetres
for x, y and z decenter, -0.2 degree to 0.2 degree with step of 0.1 degree for
x and y tilt. In total, there are 2000 different states for this element and
we could obtain 2000 PSF-Cubes. In each PSF-Cube, there are $61\times 61$ PSFs
and we obtain 744200 PSFs in total as the training set of the Tel–Net. The
size of PSF window is 15.33 um with $64\times 64$ pixels and pixelscale of
0.27 arcsec/pixel. Meanwhile, we generate 100 PSF-Cubes with the same WFSAT
and this time we set the second mirror in random states within its
misalignment configurations as the test set.

Surface | Surface Type | Radius (mm) | Separation (mm) | Glass | Semi-Diameter (mm)
---|---|---|---|---|---
1 | Spherical | 225.0 | 42.8 | BK7 | 164.5
2 | Spherical | 191.3 | 197.5 | Air | 143.8
3 | Even Aspherical | Infinity | 14.0 | F2 | 98.6
| (0, $-5.6\times 10^{-10}$, $-2.4\times 10^{-14}$) |  |  |  |
4 | Spherical | -1355.0 | 13.0 | SK2 | 101.1
5 | Plane | Infinity | 239.2 | Air | 105.1
6 | Spherical | -293.1 | 42.8 | BK7 | 177.5
7 | Spherical | -324.2 | 251.0 | Air | 196.3
8 | Spherical | -550.0 | -251.0 | Mirror | 290.8
9 | Spherical | -324.2 | -42.8 | BK7 | 119.9
10 | Spherical | -293.1 | -14.5 | Air | 95.1
Focal Plane | Spherical | -241.8 | - | - | 82.6

Figure 2: Parameters of a simulated WFSAT with Baker super-Schmidt design. The
third surface is an even aspherical surface and the aspherical coefficients
are also shown in this table

It should be noted that in essence PSF-Cubes in the training set are sparse
samples of all PSFs that correspond to the WFSAT within tolerances of its
secondary mirror. Because we sampled tolerances with relatively big steps,
there are still large number of PSFs that can be obtained. Assuming all
elements have the same number of free states as the second element, there
would be tens of thousands of PSF–Cubes to be obtained, which is impossible to
be processed in real applications. However, for PSF estimation task, many of
these states are redundant (although the telescope belongs to two different
states, PSFs in a given point within FoV are almost the same). This is why we
decided to make a trade–off between complexity and accuracy, and so we use
sparsely sampled PSF–Cubes to train the Tel–Net. In the future, it would be
necessary to further investigate possible methods of the optimal reduction
(minimization) of the training set to its most informative subset.

###  3.1 Training and testing the Tel-Net with simulated data of small size

First of all, we test the robustness of the Tel–Net using a small data set,
which includes only 200 PSF–Cubes reflecting states of the WFSAT along tilt of
x and y and decenter of x. We select 150 from 200 PSF–Cubes to train the
Tel–Net. For the input PSF–Cubes, we randomly select 60 PSFs as known PSFs and
set all other PSFs to zero matrices. For the output PSF–Cubes, we use the same
PSF–Cube with all PSFs as valid PSFs. After training, we select the remaining
50 PSF-Cubes to test the Tel–Net.

Meanwhile we use direct interpolation method for comparison. We select the
Inverse Distance Weight (IDW) interpolation method (Lu & Wong, 2008) as
exemplary classic interpolation method. The IDW interpolation method is a
commonly used weighted average interpolation method in which the weights
depend on the reciprocal of the distance between interpolated points. We use
the mean square error (MSE) between predicted PSFs and original PSFs as the
evaluation function.

Since the Tel–Net is very deep and the number of PSF-Cubes is small, there
would be some risks that the Tel–Net would over–fit these PSFs. Therefore, we
plot the MSE of PSFs (both in the training set and the test set) as function
of epochs as shown in figure 3. In order to save storage space, we save
weights for every 10 epochs. We have trained the Tel–Net for 480 epochs in
total and calculate the mean MSE for every 4 epochs. From this figure, we can
conclude that the Tel–Net could obtain effective representation of PSFs
without over–fitting.

Figure 3: MSE of PSFs in the training set and the test set as function of epochs, when the training set has smaller number of PSFs. Table 1: MSE between predicted PSFs (estimated by the Tel–Net and the IDW) and original PSFs. The Tel–Net is trained by 150 training PSF–Cubes and tested by other 50 test PSF–Cubes. Method | MSE mean | MSE variance
---|---|---
Tel-Net | $3.77\times 10^{-8}$ | $2.76\times 10^{-16}$
IDW | $8.01\times 10^{-7}$ | $2.09\times 10^{-13}$

We use the trained Tel–Net to reconstruct PSFs for WFSATs. The statical
results of the MSE are shown in table 1 and histograms of the MSE are exposed
in figure 4. We can find that with only a limited number of PSFs, the Tel–Net
could still be able to predict PSFs with relatively high accuracy, while the
existing IDW method presents a magnitude higher MSE. We further plot the heat
map of MSE between PSFs reconstructed by different methods and original PSFs
in the whole FoV as shown in figure 5. We can find that the MSE of
reconstructed PSFs is smaller in the centre of the FoV and PSFs in corner of
the FoV are not well reconstructed. This trend is the same for the IDW method,
because in the corner there would be less reference PSFs. But the MSE of PSFs
reconstructed by the Tel–Net is smaller than PSFs reconstructed by IDW method.
One of original PSFs and PSFs reconstructed by different methods are shown in
figure 6. We can conclude the following results:

  * 1)

PSFs of WFSATs tested in this part have irregular shapes;

  * 2)

PSFs obtained by the IDW method can not retain complex shapes, which is
probably the reason why our method has much better prediction results than
that of the IDW method;

  * 3)

Tel-Net has both low mean MSE and low variance of MSE which indicates high
repeatability of good results of this estimator.

Figure 4: Histogram of MSE between PSFs predicted by Tel–Net and original
PSFs. The blue histogram stands for MSE of PSFs predicted by the Tel–Net and
orange histogram stands for MSE of PSFs predicted by the IDW method. The
Tel–Net is trained by 150 training PSF–Cubes and tested by the other 50 test
PSF–Cubes.

(a) The heat map of mean MSE of PSFs reconstructed by the Tel–Net

(b) The heat map of mean MSE of PSFs reconstructed by the IDW

Figure 5: The heat map of mean MSE between original PSFs and PSFs
reconstructed by different methods.

(a) Original PSF

(b) Tel-Net predicted PSF

(c) IDW predicted PSF

Figure 6: Several examples of original PSFs and corresponding PSFs
reconstructed by different methods with their MSE values, when the number of
PSFs in the training set is small. (a) are original PSFs, (b) are PSFs
predicted by the Tel–Net and (c) are PSFs predicted by the IDW with with four
nearby PSFs.

###  3.2 Training and testing the Tel-Net with simulated data of big size

In the next step, we test the performance of the Tel–Net with a bigger
dataset. This time we use the entire 2000 PSF-Cubes that are defined in the
beginning of this section. We use 100 PSF–Cubes with random misalignments of
the secondary mirror as the test set. The Tel–Net is trained through 10 epochs
which is enough to achieve reasonable abilities. The MSE of the training set
and the test set of different epochs is shown in figure 7. The MSE mean values
are also calculted for every 4 epochs. As shown in this figure, the Tel–Net
converges after 40 epochs. The statistical results of the MSE are shown in
table 2 and the histograms of the MSE are shown in figure 8. The results show
that the Tel–Net has better ability to model PSFs in fine details than the
interpolation method. The residual error is more than one order smaller than
the IDW method. It indicates us that with adequate sampled PSFs as training
set, the Tel–Net could almost predict any PSFs of WFSATs in any states. In
comparison to results obtained on the smaller data set, we can find that the
mean MSE has been reduced by half. The heat map of MSE in different positions
of the FoV is shown in figure 9. We can find that with more PSFs as training
set, the Tel–Net can better predict unknown PSFs. We show original and
reconstructed PSFs in figure 10. As shown in this figure, we can find that
with more PSFs as training set, the Tel–Net could better reconstruct PSFs.

Figure 7: MSE of PSFs in the training set and the test set as function of
epochs, when the training set has bigger number of PSFs. Figure 8: Histogram
of MSE between PSFs predicted by Tel–Net and original PSFs. The blue histogram
stands for MSE of PSFs predicted by the Tel–Net and orange histogram stands
for MSE of PSFs predicted by the IDW method. The Tel–Net is trained by 2000
training PSF–Cubes and tested by 100 test PSF–Cubes. Figure 9: The heat map of
mean MSE between original PSFs and PSFs reconstructed by the Tel–Net.

(a) Original PSF

(b) Tel-Net predicted PSF

(c) IDW predicted PSF

Figure 10:  Several examples of original PSFs and corresponding PSFs reconstructed by different methods with their MSE values, when the number of PSFs in the training set is large. (a) are original PSFs, (b) are PSFs predicted by the Tel–Net and (c) are PSFs predicted by the IDW with with four nearby PSFs. Table 2: MSE between predicted PSFs (estimated by the Tel–Net and the IDW) and original PSFs. The Tel–Net is trained by 2000 training PSF–Cubes and tested by 100 test PSF–Cubes. Method | MSE mean | MSE variance
---|---|---
Tel-Net | $1.49\times 10^{-8}$ | $5.81\times 10^{-17}$
IDW | $7.79\times 10^{-7}$ | $2.04\times 10^{-13}$

###  3.3 Training and testing the Tel-Net with PSFs of different size

To cover a large FoV in a cost effective way, WFSATs used in fast sky survey
would obtain images with low spatial sampling rate. In these images, PSFs have
very limited number of pixels. Meanwhile, for photometry sky surveys, PSFs of
WFSATs may be well sampled (normally each PSF has tens to hundreds pixels). To
test robustness of the Tel–Net in processing images with different number of
pixels, we generate PSFs of different sampling rates ($32\times 32$ pixels,
$16\times 16$ pixels and $8\times 8$ pixels). The number and the strategy we
used in this subsection is the same as those used in subsection 3.2. The
results of MSE for reconstructed PSFs of different sampling rates is shown in
table 3. We can find that PSFs with different sampling rates can be well
reconstructed by a trained Tel–Net. Meanwhile, as the pixel number increases,
the MSE would decrease.

Table 3: MSE between predicted PSFs (estimated by the Tel–Net) and original PSFs for PSFs with different sizes. The Tel–Net is trained by 1900 training PSF–Cubes and tested by other 100 test PSF–Cubes. Number of pixels (pixelscale) | MSE mean | MSE variance
---|---|---
$32\times 32$ (0.54 arcsec/pixel) | $1.49\times 10^{-8}$ | $5.81\times 10^{-17}$
$16\times 16$ (1.08 arcsec/pixel) | $3.03\times 10^{-8}$ | $1.03\times 10^{-16}$
$8\times 8$ (2.16 arcsec/pixel) | $8.52\times 10^{-8}$ | $5.72\times 10^{-16}$

###  3.4 Training and testing the Tel-Net with star images of different SNR

Since the aperture size of WFSATs is small, images obtained by WFSATs may have
relatively low SNR. Stars with low SNR would introduce additional error to the
Tel–Net. To test the robustness of Tel–Net, we add different levels of Poisson
noise to simulated PSFs ($32\times 32$ pixels) as observed star images. Then,
these simulated star images are used as PSFs to train and test the Tel–Net.
The results are shown in table 4. Because PSFs are all normalized to 1 for the
Tel–Net, we multiply 10 or 80 with PSFs to generate star images with different
levels of signal and then add the same level of Poisson noise as background
noise. Therefore, 80 stand for star images with photon counts that are 8 times
larger than stars with 10. As shown in this table, we can find that as noise
levels increases, the performance of the Tel–Net would drop down. But the
performance is still better than the IDW method.

Table 4: MSE between predicted PSFs (estimated by the Tel–Net and the IDW) and original PSFs for stars with different SNRs. Counts of photons | Methods | MSE mean | MSE variance
---|---|---|---
$80$ | Tel-Net | $2.77\times 10^{-8}$ | $3.12\times 10^{-16}$
$80$ | IDW | $7.04\times 10^{-3}$ | $3.84\times 10^{-5}$
$10$ | Tel-Net | $4.18\times 10^{-8}$ | $7.69\times 10^{-16}$
$10$ | IDW | $9.62\times 10^{-3}$ | $4.32\times 10^{-5}$

##  4 Training and testing Tel–Net with experimental data

To further test the performance of the Tel–Net, we set up an optical
experimental bench in the laboratory according to the schema given in figure
11 and with the experiment setup presented in figure 12. It includes a
parabolic mirror (M3), an Ritcher-Chretien (RC) telescope with wavefront error
(RMS value) less than $1/30\lambda$ (M1 and M2), a fibre laser as a point-like
source of light and a camera (CCD). We tilt the parabolic mirror (M3) to
introduce additional aberrations and then move the camera together with the
telescope to obtain PSFs in different positions of the FoV. The experiment is
carried out in the optical laboratory located in Nanjing Institute of
Astronomical Optics & Technology, Chinese Academy of Sciences. In the
laboratory, the temperature variation is smaller than 5 degrees and and the
variation of humidity is below $10\%$. It took around 6 days to obtain all
PSFs and in this paper, we assume variations of temperature and humidity would
not introduce signifiant variations of PSFs. In total we generate 40 different
levels of M3 tilts. For each levels of tilt, we obtain 25 PSFs equally
distributed in $5\times 5$ patches. An exemplary PSF with aberrations obtained
at the edge of the FoV is shown in figure 13.

Figure 11: The concept design of the optical path in our experiment. Figure
12: The optical path utilized in our experiment.

Because automatic alignment devices such as electric precise lifting platforms
are not used in this experiment, tilt of M3 and capture of PSFs are all done
manually. Therefore, although we have spent a long time to collect PSFs, only
very limited data are obtained, comparing with the amount of data generated
previously by ZEMAX. Besides, due to limitations of human beings, these data
have obvious imperfections: different levels of tilts are not strictly equally
distributed and PSFs of the same patch corresponding to different tilts have
position uncertainties up to several pixels. These problems impose the need
for the use of automatic alignment instruments and related control methods
while calibrating the smart data processing algorithms.

Besides, limited PSF data would introduce over–fitting problem. The relation
between MSE and epochs is shown in figure 13. As shown in this figure, the
Tel–Net converges after 250 epochs but we could find that the Tel–Net is
over–fitting the training set. It is caused by limited number of PSF-Cubes,
which would require automatic alignment devices to obtain more data in our
future experiments. Although these PSF–Cubes have some issues, we are still
able to train the Tel–Net with these data to test its robustness. In total, we
obtain 40 PSF–Cubes and we randomly select 36 PSF–Cubes to train the Tel–Net.
8 PSFs are randomly selected as known PSFs in input PSF–Cubes. After training,
we use the Tel–Net to predict PSFs with PSF–Cubes from the test set. One of
PSFs predicted by the Tel–Net and that predicted by the IDW method are shown
in figure 14. The Tel–Net can predict PSFs in fine details, while the IDW
method can only predict main structure of the original PSF. Kernel density
estimation and statistical results of reconstruction error are shown in figure
LABEL:figure015 and table 5. We can find that the Tel–Net has better
performance than that of the IDW method, even with very limited data. It
should be noted that because the number of PSFs in the training data is small,
the results obtained in real applications are not as good as results obtained
with the simulated data. This problem is caused by limited number of training
data and could be solved through increasing the volumn of training data with
applying automatic alignment devices when we obtain training data.

Figure 13: MSE of PSFs in the training set and the test set as function of
epochs for real observation PSFs. As we can find from this figure, the Tel–Net
is over–fitting, because there is only very limited number of PSF-Cubes.

(a) Original PSF

(b) Tel-Net predicted PSF

(c) IDW predicted PSF

Figure 14:  (a) is the original PSF, (b) is PSF predicted by the Tel–Net and (c) is PSF predicted by the IDW with four PSFs as prior information. Table 5: MSE between predicted PSFs (estimated by the Tel–Net and the IDW) and original PSFs. The Tel–Net is trained by 36 training PSF–Cubes and tested by 4 PSF–Cubes. Method | MSE mean | MSE variance
---|---|---
Tel-Net | $1.29\times 10^{-6}$ | $3.07\times 10^{-12}$
IDW | $2.83\times 10^{-6}$ | $5.60\times 10^{-12}$
Figure 15: Histogram of MSE between PSFs predicted by Tel–Net and original
PSFs. To make the results easier to compare, we use the kernel density
estimation method to fit the Histogram. The blue curve stands for MSE of PSFs
predicted by the Tel–Net and orange curve stands for MSE of PSFs predicted by
the IDW method. The Tel–Net has better performance in PSF estimation.

##  5 Conclusions and future work

In this paper, we propose the Tel-Net as a new PSF modelling framework for
WFSATs. The Tel-Net assumes that PSFs in WFSATs can be understood as a
response of a complex system for a single state of an imaging instrument.
States of the whole optical system are sampled by PSF–Cubes during system test
procedure. Then we train the Tel-Net with PSF–Cubes and use it to predict
PSFs. We test the Tel-Net with simulated and real data and find that the Tel-
Net can obtain PSFs in any positions of the FoV with finite number of star
images. An example of classic PSF estimation method, the interpolation method
utilizing the Inverse Distance Weighting has significantly worse performance
and is not able to predict complex aberrations seen in both real and simulated
PSFs.

The Tel-Net could be added as a necessary part to the calibration method for
WFSATs, particularly for WFSATs that are lack of maintenance during
observations, such as space telescopes. The Tel–Net indicates the importance
of alignment data of a WFSAT pointing out the need for automatic
characterization of the optical instrument during commissioning phase. In the
future, we will use the Tel–Net as PSF modelling method to modify astronomical
target detection and photometry neural networks (Jia et al., 2020b) to
increase the observation accuracy of these telescopes when there are PSF
variations in data obtained by these WFSATs.

## Acknowledgements

Authors would like to thank the reviewer, Dr. Morgan, Robert for his kindly
suggestions that improve this paper a lot. Peng Jia would like to thank Dr.
Alastair Basden from Durham University, Professor Rongyu Sun from Purple
Mountain Observatory who provide very helpful suggestions for this paper. This
work is supported by National Natural Science Foundation of China
(NSFC)(11503018), the Joint Research Fund in Astronomy (U1631133, U1931207)
under cooperative agreement between the NSFC and Chinese Academy of Sciences
(CAS). Adam Popowicz was supported by grants: 02/140/RGJ21/0012 and
BK-225/RAu-11/2021. Authors acknowledge the French National Research Agency
(ANR) to support this work through the ANR APPLY (grant ANR-19-CE31-0011)
coordinated by B. Neichel. This work is also supported by Shanxi Province
Science Foundation for Youths (201901D211081), Research and Development
Program of Shanxi (201903D121161), Research Project Supported by Shanxi
Scholarship Council of China, the Scientific and Technological Innovation
Programs of Higher Education Institutions in Shanxi (2019L0225).
Data Availability Statements: the code in this paper can be downloaded from
aojp.lamost.org and after acceptance the code will be released in PaperData
Repository powered by China-VO with a DOI number.

## References

  * Ba et al. (2016) Ba J. L., Kiros J. R., Hinton G. E., 2016, arXiv preprint arXiv:1607.06450
  * Beltramo-Martin et al. (2019) Beltramo-Martin O., Correia C., Ragland S., Jolissaint L., Neichel B., Fusco T., Wizinowich P., 2019, Monthly Notices of the Royal Astronomical Society, 487, 5450
  * Burd et al. (2005) Burd A., et al., 2005, in Photonics Applications in Industry and Research IV. p. 59481H
  * Burke et al. (2019a) Burke C. J., Aleo P. D., Chen Y.-C., Liu X., Peterson J. R., Sembroski G. H., Lin J. Y.-Y., 2019a, MNRAS, 490, 3952
  * Burke et al. (2019b) Burke C. J., Aleo P. D., Chen Y.-C., Liu X., Peterson J. R., Sembroski G. H., Lin J. Y.-Y., 2019b, Monthly Notices of the Royal Astronomical Society, 490, 3952
  * Cui et al. (2008) Cui X., Yuan X., Gong X., 2008, in Ground-based and Airborne Telescopes II. p. 70122D
  * Duev et al. (2019) Duev D. A., et al., 2019, MNRAS, 486, 4158
  * Fétick et al. (2019) Fétick R., et al., 2019, Astronomy & Astrophysics, 628, A99
  * Fusco et al. (2020) Fusco T., et al., 2020, Astronomy & Astrophysics, 635, A208
  * Glazier et al. (2020) Glazier A. L., Howard W. S., Corbett H., Law N. M., Ratzloff J. K., Fors O., del Ser D., 2020, arXiv e-prints, p. arXiv:2006.14712
  * González et al. (2018a) González R. E., Muñoz R. P., Hernández C. A., 2018a, Astronomy and Computing, 25, 103
  * González et al. (2018b) González R. E., Munoz R. P., Hernández C. A., 2018b, Astronomy and computing, 25, 103
  * Goodfellow et al. (2016) Goodfellow I., Bengio Y., Courville A., Bengio Y., 2016, Deep learning. Vol. 1, MIT press Cambridge
  * Hausen & Robertson (2020) Hausen R., Robertson B. E., 2020, The Astrophysical Journal Supplement Series, 248, 20
  * He et al. (2020) He Z., et al., 2020, MNRAS,
  * Jee et al. (2007) Jee M. J., Blakeslee J. P., Sirianni M., Martel A. R., White R. L., Ford H. C., 2007, PASP, 119, 1403
  * Jia et al. (2017) Jia P., Sun R., Wang W., Cai D., Liu H., 2017, MNRAS, 470, 1950
  * Jia et al. (2019a) Jia P., Zhao Y., Xue G., Cai D., 2019a, The Astronomical Journal, 157, 250
  * Jia et al. (2019b) Jia P., Huang Y., Cai B., Cai D., 2019b, The Astrophysical Journal Letters, 881, L30
  * Jia et al. (2020a) Jia P., Wu X., Yi H., Cai B., Cai D., 2020a, AJ, 159, 183
  * Jia et al. (2020b) Jia P., Liu Q., Sun Y., 2020b, AJ, 159, 212
  * Jia et al. (2020c) Jia P., Li X., Li Z., Wang W., Cai D., 2020c, MNRAS, 493, 651
  * Jia et al. (2020d) Jia P., Wu X., Yang X., Huang Y., Cai B., Cai D., 2020d, in Advances in Optical Astronomical Instrumentation 2019. p. 112030Q, doi:10.1117/12.2541083
  * Jia et al. (2020e) Jia P., Ma M., Zhang M., Wu X., Wang W., Zhan Y., 2020e, in Adaptive Optics Systems VII. p. 114484A
  * Jia et al. (2021a) Jia P., Ning R., Sun R., Yang X., Cai D., 2021a, Monthly Notices of the Royal Astronomical Society, 501, 291
  * Jia et al. (2021b) Jia P., Ma M., Cai D., Wang W., Li J., Li C., 2021b, MNRAS, 503, 3194
  * Kingma & Ba (2014) Kingma D. P., Ba J., 2014, arXiv preprint arXiv:1412.6980
  * Krist et al. (2011) Krist J. E., Hook R. N., Stoehr F., 2011, in Proc. SPIE. p. 81270J, doi:10.1117/12.892762
  * La Camera et al. (2015) La Camera A., Schreiber L., Diolaiti E., Boccacci P., Bertero M., Bellazzini M., Ciliegi P., 2015, Astronomy & Astrophysics, 579, A1
  * Li et al. (2015) Li Z., Yuan X., Cui X., 2015, Monthly Notices of the Royal Astronomical Society, 449, 425
  * Lu & Wong (2008) Lu G. Y., Wong D. W., 2008, Computers & geosciences, 34, 1044
  * Ma et al. (2007) Ma Y., Zhao H., Yao D., 2007, in Valsecchi G. B., Vokrouhlický D., Milani A., eds, IAU Symposium Vol. 236, Near Earth Objects, our Celestial Neighbors: Opportunity and Risk. pp 381–384, doi:10.1017/S1743921307003468
  * Martin et al. (2016) Martin O., et al., 2016, in Adaptive Optics Systems V. p. 99091Q
  * Pablo et al. (2016) Pablo H., et al., 2016, PASP, 128, 125001
  * Perrin et al. (2012) Perrin M. D., Soummer R., Elliott E. M., Lallo M. D., Sivaramakrishnan A., 2012, in Proc. SPIE. p. 84423D, doi:10.1117/12.925230
  * Ping & Zhang (2017) Ping Y., Zhang C., 2017, Advances in Space Research, 60, 907
  * Piotrowski & Żarnecki (2010) Piotrowski L. W., Żarnecki A. F., 2010, in Leroy C., Rancoita P.-G., Barone M., Gaddi A., Price L., Ruchti R., eds, Astroparticle, Particle and Space Physics, Detectors and Medical Physics Applications. pp 313–318, doi:10.1142/9789814307529_0051
  * Popowicz (2018) Popowicz A., 2018, in Proc. SPIE. p. 1069820, doi:10.1117/12.2299498
  * Ratzloff et al. (2019) Ratzloff J. K., Law N. M., Fors O., Corbett H. T., Howard W. S., Ser D. D., Haislip J. B., 2019, Publications of the Astronomical Society of the Pacific, 131, 075001
  * Schawinski et al. (2017) Schawinski K., Zhang C., Zhang H., Fowler L., Santhanam G. K., 2017, Monthly Notices of the Royal Astronomical Society: Letters, 467, L110
  * Sun & Yu (2019) Sun R., Yu S., 2019, Astrophysics and Space Science, 364, 39
  * Sun et al. (2020) Sun R., Yu S., Jia P., Zhao C., 2020, Monthly Notices of the Royal Astronomical Society, 497, 4000
  * Turpin et al. (2020) Turpin D., et al., 2020, Monthly Notices of the Royal Astronomical Society, 497, 2641
  * Wang et al. (2018) Wang W., Jia P., Cai D., Liu H., 2018, MNRAS, 478, 5671
  * Yuan et al. (2008) Yuan X., et al., 2008, in Ground-based and Airborne Telescopes II. p. 70124G
  * Zhang et al. (2019) Zhang H., Dauphin Y. N., Ma T., 2019, arXiv preprint arXiv:1901.09321


estimated factors are smooth over time due to the kernel,
$\mathbf{K}_{\text{time}}$.  Figure C.6: Inferred time weights for California.
Figure C.7: Inferred unit weights for California.

###  C.2 Additional impact estimates

Figure C.8: Estimates from MTGP for Normal model Figure C.9: Estimates from
MTGP for Normal model with a linear adjustment for background covariates.
Figure C.10: Estimates from MTGP with Poisson link Figure C.11: Estimates from
MTGP with Poisson link and linear covariate adjustment. Figure C.12: Estimates
from SCM with and without an intercept shift, ridge ASCM, gsynth, and matrix
completion, plotted over the posterior distribution of the effect of APPS on
homicide rates with the rank 5 Poisson model. Figure C.13: The posterior
distribution for the rank 5 Poisson model when leaving out the homicide rate
series for each of the 49 comparison states. The posterior with all states is
shown in orange. Figure C.14: Posterior predictive correlation between gun and
non-gun related homicides in California, pre-treatment.

###  C.3 Select MCMC Convergence Diagnostics

Figure C.15: Predicted counterfactual homicides for California in 2010 for the
Rank 5 Poisson model. For this model, $\hat{R}\approx 1$ and the bulk
effective sample size is approximately 3500, indicating excellent mixing.
Figure C.16: Predicted counterfactual homicides for California in 2010 for the
Rank 5 Gaussian model. For this model, $\hat{R}\approx 1$ and the bulk
effective sample size is approximately 1000, indicating excellent mixing.