$\tau_{\pm 1}=\mp\frac{\tau_{1}\pm i\tau_{2}}{\sqrt{2}},\ \tau_{0}=\tau_{3}.$
(144)

We obtain the $j$-representation of $\tau_{\alpha}$ in terms of the 3$j$
symbol, according to the Wigner–Eckart theorem, as

$D^{{}^{\prime}j}_{mm^{\prime}}(\tau_{\alpha})=iw_{j}\epsilon^{j}_{nm}\left(\begin{array}[]{ccc}j&j&1\\\
n&m^{\prime}&\alpha\end{array}\right)$ (145)

where $w_{j}=\sqrt{j(j+1)(2j+1)}$ and
$\epsilon^{j}_{nm}=(-1)^{j+m}\delta(n,-m)$ is the 2-$j$ symbol. The 2-$j$
symbol is also SU(2) invariant. By 3-$j$ symbol and $\epsilon^{j}_{nm}$, any
SU(2)-intertwiner can be constructed as

$\iota(k_{1}k_{2}\cdots k_{n-1})_{m_{1}m_{2}m_{3}\cdots
m_{n}}=\left(\begin{array}[]{ccc}j_{1}&j_{2}&k_{1}\\\
m_{1}&m_{2}&l_{1}\end{array}\right)\epsilon^{k_{1}}_{l_{1}l_{1}^{\prime}}\left(\begin{array}[]{ccc}k_{1}&j_{3}&k_{2}\\\
l_{1}^{\prime}&m_{3}&l_{2}\end{array}\right)\epsilon^{k_{2}}_{l_{2}l_{2}^{\prime}}\cdots\epsilon^{k_{n-2}}_{l_{n-2}l_{n-2}^{\prime}}\left(\begin{array}[]{ccc}k_{n-2}&j_{n-1}&j_{n}\\\
l_{n-2}^{\prime}&m_{n-1}&m_{n}\end{array}\right).$ (146)

Moreover, the Clebsch-Gordan coefficients relates to $3j$-symbol as

$\displaystyle\langle j_{1}m_{1}j_{2}m_{2}|JM\rangle=$
$\displaystyle(-1)^{2j_{2}}\sqrt{2J+1}\,\epsilon^{J}_{MN}\,\left(\begin{array}[]{ccc}J&j_{2}&j_{1}\\\
N&m_{2}&m_{1}\end{array}\right)$ (147) $\displaystyle=$
$\displaystyle(-1)^{j_{1}-j_{2}-J}\sqrt{2J+1}\left(\begin{array}[]{ccc}j_{1}&j_{2}&J\\\
m_{1}&m_{2}&N\end{array}\right)\epsilon^{J}_{NM}$

##  Appendix B the Clebsch-Gordan coefficients with negative parameters

Given $j_{1}$, $m_{1}$, $j_{2}$, $m_{2}$ and $m=m_{1}+m_{2}$, the Clebsch-
Gordan coefficients $\langle
j_{1}m_{1}j_{2}m_{2}|pm\rangle\equiv\left[\begin{array}[]{ccc}j_{1}&j_{2}&p\\\
m_{1}&m_{2}&m\end{array}\right]$ for various $p$ satisfy the difference
equation natura1988nature

$A(p+1)\left[\begin{array}[]{ccc}j_{1}&j_{2}&p+1\\\
m_{1}&m_{2}&m\end{array}\right]+A(p)\left[\begin{array}[]{ccc}j_{1}&j_{2}&p-1\\\
m_{1}&m_{2}&m\end{array}\right]+(A_{0}(p)-m_{1}+m_{2})\left[\begin{array}[]{ccc}j_{1}&j_{2}&p\\\
m_{1}&m_{2}&m\end{array}\right]=0$ (148)

where $\max(|m|,|j_{1}-j_{2}|)\leq p\leq j_{1}+j_{2}$ and

$\displaystyle A(p)=$
$\displaystyle\frac{1}{p}\sqrt{-\frac{(p^{2}-\xi_{1}^{2})(p^{2}-\xi_{2}^{2})(p^{2}-\xi_{3}^{2})}{4p^{2}-1}},$
(149) $\displaystyle A_{0}(p)=$
$\displaystyle\frac{\xi_{1}\xi_{2}\xi_{3}}{p(p+1)}$

with

$\xi_{1}=j_{1}-j_{2},\ \xi_{2}=j_{1}+j_{2}+1,\ \xi_{3}=m.$

With the initial data

$\displaystyle\left[\begin{array}[]{ccc}j_{1}&j_{2}&j_{1}+j_{2}\\\
m_{1}&m_{2}&m\end{array}\right]=\sqrt{\frac{\left(2j_{1}\right)!\left(2j_{2}\right)!\left(j_{1}+j_{2}-m\right)!\left(j_{1}+j_{2}+m\right)!}{\left(2j_{1}+2j_{2}\right)!\left(j_{1}-m_{1}\right)!\left(j_{1}+m_{1}\right)!\left(j_{2}-m_{2}\right)!\left(j_{2}+m_{2}\right)!}}$
(150)

the Clebsch-Gordan coefficients for other values of $p$ are computable with
Eq. (148).

In order to extend the Clebsch-Gordan coefficients to negative parameters, we
define a function

$\mathcal{C}_{0}(x):=\sqrt{\frac{\left(2j_{2}\right)!\,\Gamma(2x+1)\Gamma(x+j_{2}-m+1)\Gamma(x+j_{2}+m+1)}{\left(j_{2}-m_{2}\right)!\left(j_{2}+m_{2}\right)!\Gamma(2x+2j_{2}+1)\Gamma(x-m_{1}+1)\Gamma(x+m_{1}+1)}},$
(151)

with which

$\left[\begin{array}[]{ccc}j_{1}&j_{2}&j_{1}+j_{2}\\\
m_{1}&m_{2}&m\end{array}\right]=\mathcal{C}_{0}(j_{1}).$ (152)

It is remarkable that $\mathcal{C}_{0}(x)$ is well-defined not only for
positive $x$ such that

$x-m_{1}\in\mathbb{Z},\ x\geq|m_{1}|$ (153)

but also for negative $x$ satisfying

$x-m_{1}\in\mathbb{N},\ x\leq\min(-|m_{1}|,-j_{2}-|m|-1)$ (154)

where those gamma functions with negative integers is understood as

$\frac{\Gamma(-m_{1})\cdots\Gamma(-m_{k})}{\Gamma(-n_{1})\cdots\Gamma(-n_{k})}=\lim_{z\to
0}\frac{\Gamma(z-m_{1})\cdots\Gamma(z-m_{k})}{\Gamma(z-n_{1})\cdots\Gamma(z-n_{k})}=(-1)^{n_{1}+\cdots+n_{k}-m_{1}-\cdots-
m_{k}}\frac{n_{1}!\cdots n_{k}!}{m_{1}!\cdots m_{k}!}.$ (155)

By definition, the Clebsch-Gordan coefficients
$\left[\begin{array}[]{ccc}j_{1}&j_{2}&j_{1}+j_{2}-\iota\\\
m_{1}&m_{2}&m\end{array}\right]$ is obtained by applying the recurrence
relation (148) successively for $\iota$ steps with the initial data
$\mathcal{C}_{0}(j_{1},j_{1}+j_{2})$. Then, we defined
$\left[\begin{array}[]{ccc}-j_{1}&j_{2}&-j_{1}+j_{2}-\iota\\\
m_{1}&m_{2}&m\end{array}\right]$, the Clebsch-Gordan coefficients with
negative parameters, as the result by applying the recurrence relation

$\displaystyle\left[\begin{array}[]{ccc}-j_{1}&j_{2}&-q-1\\\
m_{1}&m_{2}&m\end{array}\right]=-\frac{1}{\tilde{A}(-q)}\left(\tilde{A}(-q+1)\left[\begin{array}[]{ccc}-j_{1}&j_{2}&-q+1\\\
m_{1}&m_{2}&m\end{array}\right]+(\tilde{A}_{0}(-q)-m_{1}+m_{2})\left[\begin{array}[]{ccc}-j_{1}&j_{2}&-q\\\
m_{1}&m_{2}&m\end{array}\right]\right)$ (156)

with the initial data $\mathcal{C}_{0}(-j_{1})$, where

$\tilde{A}(q)=A(q)\big{|}_{j_{1}\to-j_{1}},\
\tilde{A}_{0}(q)=A_{0}(q)\big{|}_{j_{1}\to-j_{1}}.$ (157)

This definition extended the Clebsch-Gordan coefficients to negative
parameters. It guarantees that,

$\displaystyle\left[\begin{array}[]{ccc}-j_{1}&j_{2}&-j_{1}+j_{2}-\iota\\\
m_{1}&m_{2}&m\end{array}\right]=\left.\left[\begin{array}[]{ccc}j_{1}&j_{2}&j_{1}+j_{2}-\iota\\\
m_{1}&m_{2}&m\end{array}\right]\right|_{j_{1}\to-j_{1}}.$ (158)

By definition, it has

$\displaystyle\mathcal{C}_{0}(-j_{1})$
$\displaystyle=\sqrt{\frac{\left(2j_{1}-2j_{2}-1\right)!\left(2j_{2}\right)!\left(j_{1}-m_{1}-1\right)!\left(j_{1}+m_{1}-1\right)!}{\left(2j_{1}-1\right)!\left(j_{1}-j_{2}-m-1\right)!\left(j_{1}-j_{2}+m-1\right)!\left(j_{2}-m_{2}\right)!\left(j_{2}+m_{2}\right)!}}$
(159) $\displaystyle=$
$\displaystyle(-1)^{j_{2}+{m_{2}}}\left[\begin{array}[]{ccc}j_{1}-1&j_{2}&j_{1}-j_{2}-1\\\
m_{1}&m_{2}&m\end{array}\right].$

Moreover, $\left[\begin{array}[]{ccc}j_{1}-1&j_{2}&j_{1}-j_{2}+\iota-1\\\
m_{1}&m_{2}&m\end{array}\right]$ can also be obtained by applying successively
the recurrence relation

$\displaystyle\left[\begin{array}[]{ccc}j_{1}&j_{2}&q\\\
m_{1}&m_{2}&m\end{array}\right]=-\frac{1}{B(q)}\left(B(q-1)\left[\begin{array}[]{ccc}j_{1}&j_{2}&q-2\\\
m_{1}&m_{2}&m\end{array}\right]+(B_{0}(q-1)-m_{1}+m_{2})\left[\begin{array}[]{ccc}j_{1}&j_{2}&q-1\\\
m_{1}&m_{2}&m\end{array}\right]\right)$ (160)

with the initial data $\left[\begin{array}[]{ccc}j_{1}&j_{2}&j_{1}-j_{2}-1\\\
m_{1}&m_{2}&m\end{array}\right]$, where

$B(q)=A(q)\big{|}_{j_{1}\to j_{1}-1},\ B_{0}(q):=A_{0}(q)\Big{|}_{j_{1}\to
j_{1}-1}.$ (161)

Furthermore, it can be verified that

$B(q)=-\tilde{A}(-q),\ B_{0}(q-1)=\tilde{A}_{0}(-q).$ (162)

Therefore, according to Eqs. (159), (156) and (160), we finally have

$\left[\begin{array}[]{ccc}-j_{1}&j_{2}&-j_{1}+j_{2}-\iota\\\
m_{1}&m_{2}&m\end{array}\right]=(-1)^{j_{2}+m_{2}-\iota}\left[\begin{array}[]{ccc}j_{1}-1&j_{2}&j_{1}-j_{2}+\iota-1\\\
m_{1}&m_{2}&m\end{array}\right],$ (163)

namely

$\left[\begin{array}[]{ccc}-j_{1}&j_{2}&-j_{1}+\Delta\\\
m_{1}&m_{2}&m\end{array}\right]=(-1)^{\Delta+m_{2}}\left[\begin{array}[]{ccc}j_{1}-1&j_{2}&j_{1}-1-\Delta\\\
m_{1}&m_{2}&m\end{array}\right].$ (164)

##  Appendix C Proof of (42)

One can refer to yang2017graphical ; zhang2018towards ;
makinen2019introduction  for more details on this method. The 2-$j$ symbol is
graphically represented as

$\epsilon^{j}_{mn}=(-1)^{j+n}\delta(m,-n)=\mathord{\vbox{\hbox{\includegraphics[width=43.36464pt]{2j}}}}.$
(165)

The $3j$-symbols is graphically represented as

$\left(\begin{array}[]{ccc}j_{1}&j_{2}&j_{3}\\\
m_{1}&m_{2}&m_{3}\end{array}\right)=\mathord{\vbox{\hbox{\includegraphics[width=65.04034pt]{3j}}}}.$
(166)

The Wigner-D matrix $D^{j}_{mn}(h)$, as a tensor
$h\in\mathcal{H}_{j}\otimes\mathcal{H}_{j}^{*}$, is

$D^{j}_{mn}(h)=\langle
jm|h|jn\rangle=\mathord{\vbox{\hbox{\includegraphics[width=86.72267pt]{wignerd}}}}.$
(167)

For the multiplication operator $D^{\iota}_{ab}(h_{e})$, its action on
$D^{j}_{mn}(h_{e})$ reads

$\displaystyle
D^{\iota}_{ab}(h_{e})D^{j}_{mn}(h_{e})=\sum_{J={j-\iota}}^{j+\iota}d_{J}(-1)^{M-N}\left(\begin{array}[]{ccc}\iota&j&J\\\
a&m&-M\end{array}\right)\left(\begin{array}[]{ccc}\iota&j&J\\\
b&n&-N\end{array}\right)D^{J}_{MN}(h),$ (168)

where $\left(\begin{array}[]{ccc}\iota&j&J\\\ a&m&-M\end{array}\right)$
denotes the Wigner $3j$-symbol. This graphically corresponds to

$\mathord{\vbox{\hbox{\includegraphics[width=86.72267pt]{dhdh}}}}=\sum_{J=|j_{1}-j_{2}|}^{j_{1}+j_{2}}d_{J}\mathord{\vbox{\hbox{\includegraphics[width=86.72267pt]{dhdh2}}}}.$
(169)

For the operators $\hat{p}_{s}^{i}(e)$ and $\hat{p}_{t}^{j}(e)$, by (12) we
have

$\displaystyle\hat{p}_{s}^{\alpha}(e)D^{j}_{mn}(h_{e})$
$\displaystyle=-tw_{j}\epsilon^{j}_{\tilde{m}m}\left(\begin{array}[]{ccc}j&j&1\\\
\tilde{m}&m^{\prime}&\alpha\end{array}\right)D^{j}_{m^{\prime}n}(h_{e}),$
(170) $\displaystyle\hat{p}_{t}^{\alpha}(e)D^{j}_{mn}(h_{e})$
$\displaystyle=tw_{j}\epsilon^{j}_{\tilde{n}n^{\prime}}\left(\begin{array}[]{ccc}j&j&1\\\
\tilde{n}&n&\alpha\end{array}\right)D^{j}_{mn^{\prime}}(h_{e}),$

where we used (145). Thus , one has

$p^{\alpha}_{s}(e)\mathord{\vbox{\hbox{\includegraphics[width=86.72267pt]{dhe}}}}=-tw_{j}\mathord{\vbox{\hbox{\includegraphics[width=86.72267pt]{psdhe}}}}$
(171)
$p^{\alpha}_{t}(e)\mathord{\vbox{\hbox{\includegraphics[width=86.72267pt]{dhe}}}}=tw_{j}\mathord{\vbox{\hbox{\includegraphics[width=86.72267pt]{ptdhe}}}}$
(172)

Given a monomial of holonomies and fluxes $\hat{F}_{e}$. According to Eqs.
(169), (171) and (172), we draw $\hat{F}_{e}$ graphically as

$\int{\rm
d}\mu_{h}\overline{D^{j^{\prime}}_{m^{\prime}n^{\prime}}(h)}\hat{F}_{e}D^{j}_{mn}(h)=\frac{1}{d_{j^{\prime}}}\mathord{\vbox{\hbox{\includegraphics[width=65.04034pt]{FLR}}}},$
(173)

where the dashed lines represents the indices possessed by $\hat{F}_{e}$. Some
dashed lines may carry arrows depending on the form of $\hat{F}_{e}$. With
this formula, the expectation value of $\langle\hat{F}_{e}\rangle_{z_{e}}$ is

$\displaystyle\langle\hat{F}_{e}\rangle_{z_{e}}=\sum_{jj^{\prime}}d_{j}e^{-\frac{t}{2}(j(j+1)+j^{\prime}(j^{\prime}+1))}\mathord{\vbox{\hbox{\includegraphics[width=108.405pt]{expectF4}}}}.$
(174)

To prove (42), we claim that

$\displaystyle\mathord{\vbox{\hbox{\raisebox{0.1pt}{\includegraphics[width=130.08731pt]{pspth1}}}}}=e^{(k-a+\alpha_{1}+\cdots+\alpha_{m})\overline{z_{e}}}\mathord{\vbox{\hbox{\raisebox{0.0pt}{\includegraphics[width=130.08731pt]{pspth0}}}}},$
(175)

which can be verified easily with the fact that
$D^{j}_{mn}(e^{iz\tau_{3}})=e^{zm}\delta_{mn},$ and the non-vanishing
condition of the 3$j$-symbol that
$\left(\begin{array}[]{ccc}j_{1}&j_{2}&j_{3}\\\
m_{1}&m_{2}&m_{3}\end{array}\right)\propto\delta_{m_{1}+m_{2}+m_{3},0}$.

Consequently, it has

$\displaystyle\mathord{\vbox{\hbox{\raisebox{0.0pt}{\includegraphics[width=130.08731pt]{pspth2}}}}}=e^{(\beta_{1}+\cdots\beta_{k})\overline{z_{e}}}\mathord{\vbox{\hbox{\raisebox{0.0pt}{\includegraphics[width=130.08731pt]{pspth1}}}}},$
(176)

which leads to

$\displaystyle\mathord{\vbox{\hbox{\raisebox{0.0pt}{\includegraphics[width=130.08731pt]{pspthexp1}}}}}=e^{-(\beta_{1}+\cdots\beta_{n})\overline{z_{e}}}\mathord{\vbox{\hbox{\raisebox{0.0pt}{\includegraphics[width=130.08731pt]{pspthexp2}}}}}.$
(177)

This graphical equation can be decoded as (42), where the factor $(-1)^{n}$ in
the right hand side is because of the minus sign in the definition of
$\hat{p}_{t}^{\alpha}(e)$.

##  Appendix D Expectation value of
$\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{\iota ab}$

###  D.1  for the case with $\iota=0$

Let us consider

$\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{000}\equiv\hat{F}^{\alpha_{1}\cdots\alpha_{m}}=\hat{p}_{s}^{\alpha_{1}}(e)\cdots\hat{p}_{s}^{\alpha_{m}}(e).$
(178)

Eq. (170) gives us that

$\displaystyle\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}=t^{m}\sum_{j\geq
1/2}d_{j}e^{-tj(j+1)}F_{0}(j,\frac{\partial_{\eta}}{2})\frac{\sinh((2j+1)\eta)}{\sinh(\eta)}$
(179)

where the function $F_{0}$ is given by

$\displaystyle F_{0}(j,k)=$
$\displaystyle(-w_{j})^{m}\epsilon^{j}_{kk_{1}}\left(\begin{array}[]{ccc}j&1&j\\\
k_{1}&\alpha_{1}&n_{1}\end{array}\right)\epsilon^{j}_{n_{1}k_{2}}\left(\begin{array}[]{ccc}j&1&j\\\
k_{2}&\alpha_{2}&n_{2}\end{array}\right)\cdots$ (180)
$\displaystyle\epsilon_{n_{m-2}k_{m-1}}^{j}\left(\begin{array}[]{ccc}j&1&j\\\
k_{m-1}&\alpha_{m-1}&n_{m-1}\end{array}\right)\epsilon_{n_{m-1}k_{m}}^{j}\left(\begin{array}[]{ccc}j&1&j\\\
k_{m}&\alpha_{m}&k\end{array}\right)$

By Eqs. (147) and (164), we achieve

$\displaystyle\left.w_{j}\epsilon^{j}_{\tilde{n}k}\left(\begin{array}[]{ccc}j&1&j\\\
k&\alpha&n\end{array}\right)\right|_{j\to-(j+1)}=(-1)^{\alpha}w_{j}\epsilon^{j}_{\tilde{n}k}\left(\begin{array}[]{ccc}j&1&j\\\
k&\alpha&n\end{array}\right),$ (181)

which results in

$F_{0}(j,k)=F_{0}(-j-1,k).$ (182)

Therefore

$\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}=t^{m}\sum_{n\in\mathbb{Z}}e^{-\frac{t}{4}n^{2}+\frac{t}{4}}nF_{0}(\frac{n-1}{2},\frac{\partial_{\eta}}{2})\frac{e^{n\eta}}{2\sinh(\eta)}$
(183)

where $n\equiv 2j+1$. Substituting the values of the 3-$j$ and 2-$j$ symbols,
we have that

$\displaystyle
F_{0}(j,\frac{\partial_{\eta}}{2})=\delta(\sum_{i=1}^{m}\alpha_{i},0)\left(\prod_{i=1}^{m}\frac{1}{(1+|\alpha_{i}|)^{1/2}}\right)\prod_{n=1}^{m}(\alpha_{n}j-\frac{\partial_{\eta}}{2}+\sum_{k=1}^{n}\alpha_{k})$
(184)

Applying the Poisson summation formula to Eq. (183), we get

$\displaystyle\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}$
(185) $\displaystyle=$
$\displaystyle\delta\left(\sum_{i}\alpha_{i},0\right)t^{m}\prod_{i=1}^{m}\frac{1}{(1+|\alpha_{i}|)^{1/2}}e^{t/4}\int_{-\infty}^{\infty}{\rm
d}x\,e^{-\frac{t}{4}x^{2}}x\prod_{n=1}^{m}(\alpha_{n}\frac{x-1}{2}-\frac{\partial_{\eta}}{2}+\sum_{k=1}^{n}\alpha_{k})\frac{e^{x\eta}}{2\sinh(\eta)}+O(t^{\infty}).$

By the trick

$\partial_{\eta}^{n}\frac{e^{x\eta}}{\sinh(\eta)}=\left.\left((x+\partial_{y})^{n}\frac{e^{x\eta}}{\sinh(y)}\right)\right|_{y\to\eta},$
(186)

Eq. (185) can be finally simplified as

$\displaystyle\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}$
(187) $\displaystyle=$
$\displaystyle\delta\left(\sum_{i}\alpha_{i},0\right)t^{m}\prod_{i=1}^{m}\frac{1}{(1+|\alpha_{i}|)^{1/2}}e^{t/4}\int_{-\infty}^{\infty}{\rm
d}xx\prod_{k=1}^{m}\left(\frac{\alpha_{k}-1}{2}x-\frac{\partial_{y}}{2}+\sum_{i=1}^{k}\alpha_{i}-\frac{\alpha_{k}}{2}\right)\frac{e^{-\frac{t}{4}x^{2}+x\eta}}{2\sinh(y)}\Big{|}_{y\to\eta}+O(t^{\infty})$

The integrand is a polynomial of $x$ multiplied by a Gaussian function
$e^{-tx^{2}/4+x\eta}$. Thus let us integrate $x^{n}e^{-tx^{2}/2+x\eta}$

$\displaystyle\int_{-\infty}^{\infty}{\rm d}xx^{n}e^{-\frac{t}{4}x^{2}+\eta
x}=\frac{t}{e^{t/4}}\frac{\sinh(\eta)}{\eta}\langle
1\rangle\sum_{k=0}^{\lfloor
n/2\rfloor}\binom{n}{2k}\eta^{n-2k}\left(\frac{t}{2}\right)^{k-n}(2k-1)!!$
(188)

where $\lfloor y\rfloor$ denotes the greatest integer less than or equal to
$y$ and Eq. (29) is substituted because we are concerned on the expectation
value with respect to the normalized coherent states, i.e.
$\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}/\langle 1\rangle$.
The leading-order term of the RHS of Eq. (188) is $O(t^{-n+1})$ with $n$ the
power of $x$ in the integrand. Thus it is the highest-power term of $x$ in the
integrand of (187) that leads to the leading-order term of
$\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}$, and the highest-
power and second-highest-power terms of $x$ lead to the next-to-leading order
term of $\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}$, and so
on. Therefore, we finally get

$\displaystyle\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}$
(189) $\displaystyle=$
$\displaystyle\delta\left(\sum_{i}\alpha_{i},0\right)\langle
1\rangle\prod_{i=1}^{m}\frac{1}{(1+|\alpha_{i}|)^{1/2}}\Bigg{\\{}\left(\prod_{k=1}^{m}(\alpha_{k}-1)\right)\eta^{m}+\left(\prod_{k=1}^{m}(\alpha_{k}-1)\right)\frac{(m+1)m}{2}\eta^{m-2}\frac{t}{2}+$
$\displaystyle+\sum_{k=1}^{m}\left(\prod_{i\neq
k}(\alpha_{i}-1)\right)\left(\frac{\coth(\eta)+\alpha_{k}}{2}+\alpha_{1}+\cdots+\alpha_{k-1}\right)\eta^{m-1}t+\cdots\Bigg{\\}}.$

###  D.2 for the case with $\iota\neq 0$

Now let us consider the operator

$\hat{F}_{\iota
ab}^{\alpha_{1}\cdots\alpha_{n}}=\hat{p}_{s}^{\alpha_{1}}(e)\cdots\hat{p}_{s}^{\alpha_{m}}(e)D^{\iota}_{ab}(h_{e}),\
\iota\neq 0.$ (190)

By using (170) and (168), it can be obtained that

$\displaystyle\langle\hat{F}_{\iota
ab}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}=\sum_{jj^{\prime}}d_{j}e^{-\frac{t}{2}(j(j+1)+j^{\prime}(j^{\prime}+1))}d_{j^{\prime}}(-tw_{j^{\prime}})^{m}\sum_{kk^{\prime}}e^{kz_{e}+k^{\prime}\overline{z_{e}}}\epsilon^{j^{\prime}}_{k^{\prime}k_{1}}\left(\begin{array}[]{ccc}j^{\prime}&1&j^{\prime}\\\
k_{1}&\alpha_{1}&n_{1}\end{array}\right)\epsilon^{j^{\prime}}_{n_{1}k_{2}}\left(\begin{array}[]{ccc}j^{\prime}&1&j^{\prime}\\\
k_{2}&\alpha_{2}&n_{2}\end{array}\right)\cdots$ (191)
$\displaystyle\left(\begin{array}[]{ccc}j^{\prime}&1&j^{\prime}\\\
k_{m-1}&\alpha_{m-1}&n_{m-1}\end{array}\right)\epsilon^{j^{\prime}}_{n_{m-1}k_{m}}\left(\begin{array}[]{ccc}j^{\prime}&1&j^{\prime}\\\
k_{m}&\alpha_{m}&n_{m}\end{array}\right)\epsilon^{j^{\prime}}_{n_{m}k_{m+1}}\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k_{m+1}&a^{\prime}&k\end{array}\right)\epsilon^{\iota}_{aa^{\prime}}\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k^{\prime}&b&k^{\prime\prime}\end{array}\right)\epsilon^{j}_{k^{\prime\prime}k}$

Because the products of these $3j$-symbols vanish unless $k=k^{\prime}+b$, the
last equation can be simplified as

$\displaystyle\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{\iota
ab}\rangle_{z_{e}}=t^{m}e^{bz_{e}}\sum_{j,j^{\prime}}e^{-\frac{t}{2}(j(j+1)+j^{\prime}(j^{\prime}+1))}F_{\iota}(j,j^{\prime},\frac{\partial_{\eta}}{2})\frac{\sinh((2j^{\prime}+1)\eta)}{\sinh(\eta)}$
(192)

where $F_{\iota}$ is given by

$\displaystyle F_{\iota}(j,j^{\prime},k^{\prime}):=$ $\displaystyle
d_{j}d_{j^{\prime}}(-w_{j^{\prime}})^{m}\epsilon^{j^{\prime}}_{k^{\prime}k_{1}}\left(\begin{array}[]{ccc}j^{\prime}&1&j^{\prime}\\\
k_{1}&\alpha_{1}&n_{1}\end{array}\right)\epsilon^{j^{\prime}}_{n_{1}k_{2}}\left(\begin{array}[]{ccc}j^{\prime}&1&j^{\prime}\\\
k_{2}&\alpha_{2}&n_{2}\end{array}\right)\cdots\left(\begin{array}[]{ccc}j^{\prime}&1&j^{\prime}\\\
k_{m-1}&\alpha_{m-1}&n_{m-1}\end{array}\right)\epsilon^{j^{\prime}}_{n_{m-1}k_{m}}$
(193) $\displaystyle\left(\begin{array}[]{ccc}j^{\prime}&1&j^{\prime}\\\
k_{m}&\alpha_{m}&n_{m}\end{array}\right)\epsilon^{j^{\prime}}_{n_{m}k_{m+1}}\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k_{m+1}&a^{\prime}&k^{\prime}+b\end{array}\right)\epsilon^{\iota}_{aa^{\prime}}\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k^{\prime}&b&k^{\prime\prime}\end{array}\right)\epsilon^{j}_{k^{\prime\prime}\,k^{\prime}+b}.$

The possible values of $j$ are $|j^{\prime}-\iota|,\ |j^{\prime}-\iota|+1,\
\cdots,\ j^{\prime}+\iota$. By applying Eqs. (147) and (164) again, we have

$\displaystyle\epsilon^{j^{\prime}}_{n_{m}k_{m+1}}\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k_{m+1}&a^{\prime}&k^{\prime}+b\end{array}\right)\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k^{\prime}&b&k^{\prime\prime}\end{array}\right)\epsilon^{j}_{k^{\prime\prime}\,k^{\prime}+b}\Bigg{|}_{j=j^{\prime}-\Delta}$
(194) $\displaystyle=$
$\displaystyle(-1)^{b-a+1}\left[\epsilon^{j^{\prime}}_{n_{m}k_{m+1}}\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k_{m+1}&a^{\prime}&k^{\prime}+b\end{array}\right)\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k^{\prime}&b&k^{\prime\prime}\end{array}\right)\epsilon^{j}_{k^{\prime\prime}\,k^{\prime}+b}\Bigg{|}_{j=j^{\prime}+\Delta}\right]_{j^{\prime}\to-j^{\prime}-1}.$

Together with Eq. (181), it leads to

$F_{\iota}(j^{\prime}+\Delta,j^{\prime},k^{\prime})=-F_{\iota}(j^{\prime}-\Delta,j^{\prime},k^{\prime})\Big{|}_{j^{\prime}\to-j^{\prime}-1}$
(195)

where we use
$F_{\iota}(j,j^{\prime},k)\propto\delta(\sum_{i=1}^{m}\alpha_{i}-a+b,0)$. Thus
for $j^{\prime}=j\pm d$ with $d>0$, one has

$\displaystyle\sum_{2j^{\prime}+d\geq\iota}f(d_{j^{\prime}}^{2})F_{\iota}(j^{\prime}+d,j^{\prime},\frac{\partial\eta}{2})\frac{\sinh(d_{j^{\prime}}\eta)}{\sinh(\eta)}+\sum_{2j^{\prime}-d\geq\iota}f(d_{j^{\prime}}^{2})F_{\iota}(j^{\prime}-d,j^{\prime},\frac{\partial\eta}{2})\frac{\sinh(d_{j^{\prime}}\eta)}{\sinh(\eta)}$
$\displaystyle=$
$\displaystyle\sum_{-2j^{\prime}-2+d\geq\iota}f(d_{j^{\prime}}^{2})F_{\iota}(-j^{\prime}-1+d,-j^{\prime}-1,\frac{\partial\eta}{2})\frac{\sinh(-d_{j^{\prime}}\eta)}{\sinh(\eta)}+\sum_{2j^{\prime}-d\geq\iota}f(d_{j^{\prime}}^{2})F_{\iota}(j^{\prime}-d,j^{\prime},\frac{\partial\eta}{2})\frac{\sinh(d_{j^{\prime}}\eta)}{\sinh(\eta)}$
$\displaystyle=$ $\displaystyle\sum_{\begin{subarray}{c}n\in\mathbb{Z}\\\
n\notin[d-\iota,d+\iota]\end{subarray}}f(d_{j^{\prime}}^{2})F_{\iota}(\frac{n-1}{2}-d,\frac{n-1}{2},\frac{\partial_{\eta}}{2})\frac{\sinh(n\eta)}{\sinh(\eta)}.$

Here the conditions $2j^{\prime}\pm d\geq\iota$ are from the the triangle
condition i.e. $j^{\prime}+j\geq\iota$, of the $3j$-symbol. The second
equality is because of the replacement $j^{\prime}\to-j^{\prime}-1$ in the
first summation, $d_{j^{\prime}}=2j^{\prime}+1\equiv n$ and $f$ is an
arbitrary function. Therefore,

$\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{\iota
ab}\rangle_{z_{e}}=t^{m}e^{bz_{e}}\sum_{\begin{subarray}{c}0\leq d\leq\iota\\\
d+\iota\in\mathbb{Z}\end{subarray}}\frac{2-\delta(d,0)}{2}e^{-\frac{t}{4}\left(2d^{2}-1\right)}\sum_{\begin{subarray}{c}n\in\mathbb{Z}\\\
n\notin[d-\iota,d+\iota]\end{subarray}}e^{-\frac{1}{4}t\left(n^{2}-2dn\right)}F_{\iota}(\frac{n-1}{2}-d,\frac{n-1}{2},\frac{\partial_{\eta}}{2})\frac{\sinh(n\eta)}{\sinh(\eta)}$
(196)

To calculate $F_{\iota}(j,j^{\prime},k^{\prime})$, we notice that

$\displaystyle-
w_{j^{\prime}}\epsilon^{j^{\prime}}_{n_{i-1}k_{i}}\left(\begin{array}[]{ccc}j^{\prime}&1&j^{\prime}\\\
k_{i}&\alpha_{i}&n_{i}\end{array}\right)=\delta(k_{i}+\alpha_{i}+n_{i},0)\left(\alpha_{i}j^{\prime}-n_{i}\right)\sqrt{\frac{\left(j^{\prime}-n_{i-1}\right)!\left(n_{i-1}+j^{\prime}\right)!}{\left(1-\alpha_{i}\right)!\left(\alpha_{i}+1\right)!\left(j^{\prime}-n_{i}\right)!\left(n_{i}+j^{\prime}\right)!}}.$
(197)

Moreover, by applying Eq. (147) and the expression of the Clebsch-Gordan
coefficients (Eq. (5), Section 8.2.1 in varshalovich1988quantum ) , we have

$\displaystyle\epsilon^{j^{\prime}}_{n_{m}k_{m+1}}\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k_{m+1}&a^{\prime}&k\end{array}\right)\left(\begin{array}[]{ccc}j^{\prime}&\iota&j\\\
k^{\prime}&b&k^{\prime\prime}\end{array}\right)\epsilon^{j}_{k^{\prime\prime}\,k}=(-1)^{j^{\prime}-\iota-j+a^{\prime}+b}\delta_{a^{\prime}+k,n_{m}}\delta_{k^{\prime}+b,k}\times$
(198)
$\displaystyle\frac{(j+\iota-j^{\prime})!(j-\iota+j^{\prime})!(-j+\iota+j^{\prime})!}{(j+\iota+j^{\prime}+1)!}\left(\frac{(j^{\prime}+n_{m})!(j^{\prime}-n_{m})!}{(j^{\prime}+k^{\prime})!(j^{\prime}-k^{\prime})!(\iota+a^{\prime})!(\iota-a^{\prime})!(\iota+b)!(\iota-b)!}\right)^{1/2}$
$\displaystyle\sum_{z}\frac{(-1)^{z}(j^{\prime}+\iota+k-z)!(j-k+z)!}{z!(j^{\prime}-j+\iota-z)!(j^{\prime}+k+a^{\prime}-z)!(j-\iota-
k-a^{\prime}+z)!}\sum_{z}\frac{(-1)^{z}(j+\iota+k^{\prime}-z)!(j^{\prime}-k^{\prime}+z)!}{z!(j-j^{\prime}+\iota-z)!(j+k-z)!(j^{\prime}-\iota-k+z)!}.$

Eqs. (197) and (198) lead to

$\displaystyle F_{\iota}(j^{\prime}-d,j^{\prime},k^{\prime})=$
$\displaystyle\delta(\sum_{i=1}^{m}\alpha_{i}-a+b,0)\prod_{i=1}^{m}\frac{1}{(1+|\alpha_{i}|)^{1/2}}\left(\frac{1}{(\iota+a)!(\iota-a)!(\iota+b)!(\iota-b)!}\right)^{1/2}$
(199)
$\displaystyle\frac{(-1)^{d-2a+b}(2j^{\prime}+1)(2j^{\prime}-2d+1)}{(2j^{\prime}-d+\iota+1)_{2\iota+1}}\prod_{n=1}^{m}(\alpha_{i}j^{\prime}-k^{\prime}+\sum_{i=1}^{n}\alpha_{i})$
$\displaystyle\sum_{z}\frac{(-1)^{z}(d+\iota)_{z}(j^{\prime}+k^{\prime}+b-z+\iota)_{\iota+a}(j^{\prime}-d-k^{\prime}-b+z)_{\iota-a}}{z!}$
$\displaystyle\sum_{z}\frac{(-1)^{z}(\iota-d)_{z}(j^{\prime}-d+k^{\prime}-z+\iota)_{\iota-b}(j^{\prime}-k^{\prime}+z)_{\iota+b}}{z!}$

where $(x)_{n}:=x(x-1)\cdots(x-n+1)$ is the falling factorial.

There is subtle issue on the summation over $k$ and $k^{\prime}$ in Eq. (191).
In order to obtain the factor $\sinh((2j^{\prime}+1)\eta)/\sinh(\eta)$ in
(192), one has to consider a summation of $e^{2k^{\prime}\eta}$ over all
values of $k^{\prime}$ in $[-j^{\prime},j^{\prime}]$. However, except the
constraint of $|k^{\prime}|\leq j^{\prime}$ for $k^{\prime}$, there exists
another condition of $|k^{\prime}+b|=|k|\leq j$ implied by $k=k^{\prime}+b$,
which seemingly narrows the range of $k^{\prime}$. In order to preserve the
range of $k^{\prime}$, one can always “artificially engineer”
$F_{\iota}(j,j^{\prime},k^{\prime})$ such that it vanishes for
$|k^{\prime}+b|>j$, which has actually been done by the triangle inequality of
$3j$-symbols in the RHS Eq. (193). Thus the subtlety is now encoded in the
condition that $F_{\iota}(j,j^{\prime},k^{\prime})$ vanishes for
$|k^{\prime}+b|>j$. We need to verify this for the algebraic expression of
$F_{\iota}(j,j^{\prime},k^{\prime})$ in the RHS of Eq. (199). To check it,
assume $k^{\prime}+b=j+\delta=j^{\prime}-d+\delta$ with $0<\delta\leq d+b$.
Then $(j^{\prime}-k^{\prime}+z)_{\iota+b}$, in the second summand over $z$,
vanishes because

$j^{\prime}-k^{\prime}+z\geq d+b-\delta\geq 0,\
j^{\prime}-k^{\prime}+z-(\iota+b)+1\leq\delta<0$ (200)

where $0\leq z\leq\iota-d$ is applied. Therefore, we conclude that the RHS of
Eq.(199) vanishes for $j-b<k^{\prime}\leq j^{\prime}$. Similarly, it can be
checked that $(j^{\prime}-d+k^{\prime}-z+\iota)_{\iota-b}$ vanishes,
indicating the vanishing of the RHS of Eq.(199) for $-j^{\prime}\leq
k^{\prime}+b<-j-b$. We thus claim that the RHS of Eq.(199) vanishes for
$|k^{\prime}+b|>j$, which indicates that the expression of the RHS of Eq.(199)
can be analytically extended to give $F_{\iota}(j,j^{\prime},k^{\prime})$ for
the full range of $k^{\prime}$, i.e. $|k^{\prime}|\leq j$.

Replacing $j^{\prime}$ and $k^{\prime}$ in Eq. (199) by $\frac{n-1}{2}$ and
$\frac{\partial_{\eta}}{2}$ respectively, we have

$\displaystyle
F_{\iota}(\frac{n-1}{2}-d,\frac{n-1}{2},\frac{\partial_{\eta}}{2})=$
$\displaystyle\delta(\sum_{i=1}^{m}\alpha_{i}-a+b,0)\prod_{i=1}^{m}\frac{1}{(1+|\alpha_{i}|)^{1/2}}\left(\frac{1}{(\iota+a)!(\iota-a)!(\iota+b)!(\iota-b)!}\right)^{1/2}$
(201)
$\displaystyle\frac{(-1)^{d-2a+b}n(n-2d)}{(n-d+\iota)_{2\iota+1}}\prod_{k=1}^{m}(\alpha_{i}\frac{n-1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i})$
$\displaystyle\sum_{z=0}^{\iota+d}\frac{(-1)^{z}(\iota+d)_{z}(\frac{n-1}{2}+\frac{\partial_{\eta}}{2}+b-z+\iota)_{\iota+a}(\frac{n-1}{2}-\frac{\partial_{\eta}}{2}-d-b+z)_{\iota-a}}{z!}$
$\displaystyle\sum_{z=0}^{\iota-d}\frac{(-1)^{z}(\iota-d)_{z}(\frac{n-1}{2}+\frac{\partial_{\eta}}{2}-d-z+\iota)_{\iota-b}(\frac{n-1}{2}-\frac{\partial_{\eta}}{2}+z)_{\iota+b}}{z!}.$

In order to apply the Poisson summation formula to (196) to calculate the
summation over $n$, we need to extend the summation to entire $\mathbb{Z}$,
while Eq. (196) sums $n$ over $\mathbb{Z}\setminus[d-\iota,d+\iota]$. However,
because of the term $(n-d+\iota)_{2\iota+1}$ in the denominator of
$F_{\iota}$, this extension is not trivial. We need to prove that
$F_{\iota}(\frac{n-1}{2}-d,\frac{n-1}{2},\frac{\partial_{\eta}}{2})\frac{\sinh(n\eta)}{\sinh(\eta)}$
is well-defined at the points where $F_{\iota}$ itself does not. Once it is
proven, we have

$\displaystyle\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{\iota
ab}\rangle_{z_{e}}=t^{m}e^{bz_{e}}\sum_{\begin{subarray}{c}0\leq d\leq\iota\\\
d+\iota\in\mathbb{Z}\end{subarray}}\frac{2-\delta(d,0)}{2}e^{-\frac{t}{4}\left(2d^{2}-1\right)}\int{\rm
d}x\,e^{-\frac{1}{4}t\left(x^{2}-2dx\right)}F_{\iota}(\frac{x-1}{2}-d,\frac{x-1}{2},\frac{\partial_{\eta}}{2})\frac{\sinh(x\eta)}{\sinh(\eta)}$
(202)
$\displaystyle-\sum_{n\in[d-\iota,d+\iota]\cap\mathbb{Z}}e^{-\frac{1}{4}t\left(n^{2}-2dn\right)}F_{\iota}(\frac{n-1}{2}-d,\frac{n-1}{2},\frac{\partial_{\eta}}{2})\frac{\sinh(n\eta)}{\sinh(\eta)}+O(t^{\infty}).$

Moreover, the integrals in the last equation usually produce a factor
$e^{\eta^{2}/t}$. Thus the second term given by the summation over
$n\in[d-\iota,d+\iota]\cap\mathbb{Z}$ decays exponentially as $t\to 0$ after
normalization. Therefore, we finally obtain

$\displaystyle\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{\iota
ab}\rangle_{z_{e}}=t^{m}e^{bz_{e}}\sum_{\begin{subarray}{c}0\leq d\leq\iota\\\
d+\iota\in\mathbb{Z}\end{subarray}}\frac{2-\delta(d,0)}{2}e^{-\frac{t}{4}\left(2d^{2}-1\right)}\int{\rm
d}x\,e^{-\frac{1}{4}t\left(x^{2}-2dx\right)}F_{\iota}(\frac{x-1}{2}-d,\frac{x-1}{2},\frac{\partial_{\eta}}{2})\frac{\sinh(x\eta)}{\sinh(\eta)}+O(t^{\infty}).$
(203)

In order to show that
$F_{\iota}(\frac{n-1}{2}-d,\frac{n-1}{2},\frac{\partial_{\eta}}{2})\frac{\sinh(n\eta)}{\sinh(\eta)}$
is well-defined at the poles of $F_{\iota}$, we have checked that
$F_{\iota}(\frac{n-1}{2}-d,\frac{n-1}{2},\frac{\partial_{\eta}}{2})$ for at
least $\iota\leq 20$ can be simplified as summation of terms taking the
form101010Here, we conjecture that this statement is true for all $\iota$.

$\frac{1}{x-n}f(x,\partial_{\eta})\frac{\sinh(x\eta)}{\sinh(\eta)}$

where $f(x,\partial_{\eta})$, a polynomial of $x$ and $\partial_{\eta}$,
satisfies that $f(n,\cdot)$ takes $\partial_{\eta}=2k$ for all
$k\in\mathbb{Z}\cap[-\frac{n-1}{2},\frac{n-1}{2}]$ as its zeros. Then, the
following theorem will be helpful.

######  Theorem D.1.

Let $f(x,k)$ be a polynomials of $x$ and $k$ and $n\in\mathbb{Z}$ be some
integer. Let

$h(x):=\frac{1}{x-n}f(x,\partial_{\eta})\frac{\sinh(x\eta)}{\sinh(\eta)}.$
(204)

Then $h(n):=\lim_{x\to n}h(x)$ is well-defined, i.e., $n$ is a removable
singularity of $h(x)$, provided that for $x=n$ $f(x,\partial_{\eta})$
satisfies

$f(n,\partial_{\eta})=g(n,\partial_{\eta})\left(\prod_{k=-\frac{n-1}{2}}^{\frac{n-1}{2}}\left(\partial_{\eta}+2k\right)\right),$
(205)

with some polynomial $g(n,\partial_{\eta})$ of $n$ and $\partial_{\eta}$.
Moreover, if (205) holds, $h(x)$ will take the form

$h(x)=g(x,\partial_{\eta})\left(\tilde{g}_{0}(x,\eta)\frac{\sinh((x-n)\eta)}{x-n}+\sum_{l\geq
1}g_{l}(x,\eta)(x-n)^{l}\right).$ (206)

###### Proof.

We will prove that

$\left(\prod_{k=-\frac{n-1}{2}}^{\frac{n-1}{2}}\left(\partial_{\eta}+2k\right)\right)\frac{\sinh(x\eta)}{\sinh(\eta)}=\sum_{l}g_{l}(x,\eta)(x-n)^{l}$
(207)

with $g_{l}$ some function of $x$ and $\eta$, and $g_{0}$ taking the form

$g_{0}(x,\eta)=\tilde{g}_{0}(x,\eta)\sinh((x-n)\eta).$ (208)

Let $y=x-n$. Then we have

$\frac{\sinh(x\eta)}{\sinh(\eta)}=\frac{1}{\sinh(\eta)}\left(\cosh(y\eta)\sinh(n\eta)+\sinh(y\eta)\cosh(n\eta)\right).$
(209)

Substitute the above expression into the LHS of Eq. (207) and expand the
result. One then expresses the LHS of Eq. (207) as a linear combination of
terms taking the forms of $\partial_{\eta}^{l}(h_{1}(\eta)\sinh(y\eta))$ and
$\partial_{\eta}^{l}(h_{2}(\eta)\cosh(y\eta))$, where $h_{1}$ and $h_{2}$ are
some arbitrary functions. Expanding the action of the differential operator by
Leibniz’s rule, we obtain linear combinations of
$q_{1}(\eta)\cosh(y\eta)y^{m}$ and $q_{2}(\eta)\sinh(y\eta)y^{m}$, where
$q_{1}$ and $q_{2}$ are some arbitrary functions. Thus, (207) is obtained.

To get $g_{0}$, one just act $\partial_{\eta}$ on neither $\cosh(y\eta)$ nor
$\sinh(y\eta)$. Thus

$g_{0}=\left(\prod_{k=-\frac{n-1}{2}}^{\frac{n-1}{2}}\left(\partial_{\eta}+2k\right)\frac{\sinh(n\eta)}{\sinh(\eta)}\right)\cosh(y\eta)+\left(\prod_{k=-\frac{n-1}{2}}^{\frac{n-1}{2}}\left(\partial_{\eta}+2k\right)\frac{\cosh(n\eta)}{\sinh(\eta)}\right)\sinh(y\eta)$
(210)

Since $\frac{\sinh(n\eta)}{\sinh(\eta)}=\sum_{k=-j}^{j}e^{-2k\eta}$ with
$2j+1=n$ and $(\partial_{\eta}+k)e^{-k\eta}=0$, we have

$\prod_{k=-\frac{n-1}{2}}^{\frac{n-1}{2}}\left(\partial_{\eta}+2k\right)\frac{\sinh(n\eta)}{\sinh(\eta)}=0.$
(211)

Therefore

$g_{0}=\left(\prod_{k=-\frac{n-1}{2}}^{\frac{n-1}{2}}\left(\partial_{\eta}+2k\right)\frac{\cosh(n\eta)}{\sinh(\eta)}\right)\sinh(y\eta)$
(212)

∎

##  Appendix E $\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{\iota
ab}\rangle_{z_{e}}$ for the case of $\iota=1$

By setting $\iota=1$ in (196), we can obtain

$\displaystyle\langle\hat{F}_{1ab}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}$
(213) $\displaystyle=$
$\displaystyle-\delta\left(\sum_{i=1}^{m}\alpha_{i}-a+b\right)\frac{(-1)^{a}}{2}t^{m}e^{bz_{e}}\left(\prod_{i=1}^{m}\frac{1}{\sqrt{1+|\alpha_{i}|}}\right)\frac{1}{\sqrt{(1+|a|)(1+|b|)}}\sum_{\begin{subarray}{c}n\in\mathbb{Z}\\\
n\notin[-1,1]\end{subarray}}e^{-\frac{t}{4}(n^{2}-1)}$
$\displaystyle\prod_{k=1}^{m}\left(\alpha_{k}\frac{n-1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\frac{n\left(a\frac{n-1}{2}+b+\frac{\partial_{\eta}}{2}\right)\left(b\frac{n-1}{2}-\frac{\partial_{\eta}}{2}\right)}{\frac{n^{2}-1}{4}}\frac{\sinh(n\eta)}{\sinh(\eta)}$
$\displaystyle+$
$\displaystyle\delta\left(\sum_{i=1}^{m}\alpha_{i}-a+b\right)(-1)^{b-a}t^{m}e^{bz_{e}}\left(\prod_{i=1}^{m}\frac{1}{\sqrt{1+|\alpha_{i}|}}\right)\frac{1}{\sqrt{(1+|a|)(1+|b|)}}\sum_{\begin{subarray}{c}n\in\mathbb{Z}\\\
n\notin[-2,0]\end{subarray}}e^{-\frac{t}{4}(n+1)^{2}}$
$\displaystyle\prod_{k=1}^{m}\left(-\alpha_{k}\frac{n+1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\frac{\left(\frac{n+1}{2}-(1-|b|-b)\frac{\partial_{\eta}}{2}\right)\left(\frac{n+1}{2}-(|b|-1-b)\frac{\partial_{\eta}}{2}+|b|\right)}{\frac{n+1}{2}}\frac{\sinh(n\eta)}{\sinh(\eta)}.$

where in the second term we replaced $n$ by $-n$ for the further convenience.
It is easy to check that Theorem D.1 can be applied to extend the summation in
the last equation to all $n\in\mathbb{Z}$. Therefore, we obtain that

$\displaystyle\langle\hat{F}_{1ab}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}$
(214) $\displaystyle=$
$\displaystyle-\frac{(-1)^{a}}{2}\delta\left(\sum_{i=1}^{m}\alpha_{i}-a+b\right)t^{m}e^{bz_{e}}\left(\prod_{i=1}^{m}\frac{1}{\sqrt{1+|\alpha_{i}|}}\frac{1}{\sqrt{1+|a|}\sqrt{1+|b|}}\right)\int_{-\infty}^{\infty}{\rm
d}x\,e^{-\frac{t}{4}(x^{2}-1)}$
$\displaystyle\prod_{k=1}^{m}\left(\alpha_{k}\frac{x-1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(a\frac{x-1}{2}+\frac{\partial_{\eta}}{2}+b\right)\frac{2x\left(b(x-1)-\partial_{\eta}\right)}{x^{2}-1}\frac{\sinh(x\eta)}{\sinh(\eta)}$
$\displaystyle+$
$\displaystyle(-1)^{b-a}\delta\left(\sum_{i=1}^{m}\alpha_{i}-a+b\right)t^{m}e^{bz_{e}}\left(\prod_{i=1}^{m}\frac{1}{\sqrt{1+|\alpha_{i}|}}\frac{1}{\sqrt{1+|a|}\sqrt{1+|b|}}\right)\int_{-\infty}^{\infty}{\rm
d}xe^{-\frac{t}{4}(x+1)^{2}}$
$\displaystyle\prod_{k=1}^{m}\left(-\alpha_{k}\frac{n+1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(\frac{n+1}{2}-(|b|-1-b)\frac{\partial_{\eta}}{2}+|b|\right)\frac{\left(x+1-(1-|b|-b)\partial_{\eta}\right)}{x+1}\frac{\sinh(n\eta)}{\sinh(\eta)}$
$\displaystyle+O(t^{\infty}).$

To show the algorithms to compute Eq. (214), $m\neq 0$ will be assumed without
the loss of generality. For the case of $m=0$, the algorithm can be applied
very similarly.

To begin with, we will rewrite Eq. (214) as

$\displaystyle\langle\hat{F}_{1ab}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}$
(215) $\displaystyle=$
$\displaystyle-\frac{(-1)^{a}}{2}\delta\left(\sum_{i=1}^{m}\alpha_{i}-a+b\right)t^{m}e^{bz_{e}}\left(\prod_{i=1}^{m}\frac{1}{\sqrt{1+|\alpha_{i}|}}\frac{1}{\sqrt{1+|a|}\sqrt{1+|b|}}\right)\int_{-\infty}^{\infty}{\rm
d}x\,e^{-\frac{t}{4}(x^{2}-1)}$
$\displaystyle\prod_{k=2}^{m}\left(\alpha_{k}\frac{x-1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(a\frac{x-1}{2}+\frac{\partial_{\eta}}{2}+b\right)x\left(b-\frac{\partial_{\eta}}{x-1}\right)\left(\alpha_{1}-\frac{\partial_{\eta}}{x+1}\right)\frac{\sinh(n\eta)}{\sinh(\eta)}$
$\displaystyle+$
$\displaystyle(-1)^{b-a}\delta\left(\sum_{i=1}^{m}\alpha_{i}-a+b\right)t^{m}e^{bz_{e}}\left(\prod_{i=1}^{m}\frac{1}{\sqrt{1+|\alpha_{i}|}}\frac{1}{\sqrt{1+|a|}\sqrt{1+|b|}}\right)\int_{-\infty}^{\infty}{\rm
d}x\,e^{-\frac{t}{4}(x+1)^{2}}$
$\displaystyle\prod_{k=1}^{m}\left(-\alpha_{k}\frac{x+1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{k}\right)\left(\frac{x+1}{2}-(|b|-1-b)\frac{\partial_{\eta}}{2}+|b|\right)\left(1-(1-|b|-b)\frac{\partial_{\eta}}{x+1}\right)\frac{\sinh(x\eta)}{\sinh(\eta)}.$

Then for the first integral, it is

$\displaystyle I_{1}=$ $\displaystyle\int_{-\infty}^{\infty}{\rm
d}x\,e^{-\frac{t}{4}(x^{2}-1)}\prod_{k=2}^{m}\left(\alpha_{k}\frac{x-1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(a\frac{x-1}{2}+\frac{\partial_{\eta}}{2}+b\right)x\left(b-\frac{\partial_{\eta}}{x-1}\right)\left(\alpha_{1}-\frac{\partial_{\eta}}{x+1}\right)\frac{\sinh(x\eta)}{\sinh(\eta)}.$
(216)

Because

$\displaystyle\left(b-\frac{\partial_{\eta}}{x-1}\right)\left(\alpha_{1}-\frac{\partial_{\eta}}{x+1}\right)\frac{\sinh(x\eta)}{\sinh(\eta)}$
(217) $\displaystyle=$
$\displaystyle\Bigg{(}b\alpha_{1}\frac{\sinh(x\eta)}{\sinh(\eta)}-\alpha_{1}\left(\frac{\cosh(\eta
x)}{\sinh(\eta)}-\frac{\sinh(\eta(x-1))}{\sinh^{2}(\eta)(x-1)}\right)-b\left(\frac{\cosh(\eta
x)}{\sinh(\eta)}-\frac{\sinh(\eta(x+1))}{\sinh^{2}(\eta)(x+1)}\right)+$
$\displaystyle\left(\frac{\sinh(\eta
x)}{\sinh(\eta)}+\frac{\cosh(\eta)}{\sinh(\eta)^{3}}\frac{\sinh(\eta(x-1))}{x-1}-\frac{\cosh(\eta)}{\sinh(\eta)^{3}}\frac{\sinh((x+1)\eta)}{x+1}\right)\Bigg{)},$

$I_{1}$ can be deduced further as the following

$\displaystyle I_{1}=$ $\displaystyle\int_{-\infty}^{\infty}{\rm
d}x\,e^{-\frac{t}{4}(x^{2}-1)}\prod_{k=2}^{m}\left(\alpha_{k}\frac{x-1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(a\frac{x-1}{2}+\frac{\partial_{\eta}}{2}+b\right)x\frac{(b\alpha_{1}+1)\sinh(x\eta)-(\alpha_{1}+b)\cosh(\eta
x)}{\sinh(\eta)}$ (218) $\displaystyle+$
$\displaystyle\int_{-\infty}^{\infty}{\rm
d}y\,e^{-\frac{t}{4}(y^{2}+2y)}\prod_{k=2}^{m}\left(\alpha_{k}\frac{y}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(a\frac{y}{2}+\frac{\partial_{\eta}}{2}+b\right)(y+1)\frac{\alpha_{1}\sinh(\eta)+\cosh(\eta)}{\sinh(\eta)^{3}}\frac{\sinh(\eta
y)}{y}$ $\displaystyle-$ $\displaystyle\int_{-\infty}^{\infty}{\rm
d}y\,e^{-\frac{t}{4}(y^{2}-2y)}\prod_{k=2}^{m}\left(\alpha_{k}\frac{y}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k-1}\alpha_{i}\right)\left(a\frac{y}{2}+\frac{\partial_{\eta}}{2}+b-a\right)(y-1)\frac{\cosh(\eta)-b\sinh(\eta)}{\sinh(\eta)^{3}}\frac{\sinh(y\eta)}{y}.$

Then by using the trick

$\partial_{\eta}^{n}e^{\pm x\eta}f(\eta)=\left.(\pm x+\partial_{z})^{n}e^{\pm
x\eta}f(z)\right|_{z\to\eta},$ (219)

we have

$\displaystyle I_{1}=$ $\displaystyle\sum_{s=\pm 1}\int_{-\infty}^{\infty}{\rm
d}x\,e^{-\frac{t}{4}(x^{2}-1)}\prod_{k=2}^{m}\left(\frac{\alpha_{k}-s}{2}x-\frac{\alpha_{k}+\partial_{z}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(\frac{a+s}{2}x+\frac{-a+\partial_{z}}{2}+b\right)$
(220) $\displaystyle
x\frac{(-(\alpha_{1}+b)+s(b\alpha_{1}+1))e^{sx}}{2\sinh(z)}\Big{|}_{z\to\eta}$
$\displaystyle+$ $\displaystyle\sum_{s=\pm 1}\int_{-\infty}^{\infty}{\rm
d}y\,e^{-\frac{t}{4}(y^{2}+2y)}\prod_{k=2}^{m}\left(\frac{\alpha_{k}-s}{2}y-\frac{\partial_{z}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(\frac{a+s}{2}y+\frac{\partial_{z}}{2}+b\right)$
$\displaystyle(y+1)\frac{se^{s\eta
y}}{2y}\frac{\alpha_{1}\sinh(z)+\cosh(z)}{\sinh(z)^{3}}\Big{|}_{z\to\eta}$
$\displaystyle-$ $\displaystyle\sum_{s=\pm 1}\int_{-\infty}^{\infty}{\rm
d}y\,e^{-\frac{t}{4}(y^{2}-2y)}\prod_{k=2}^{m}\left(\frac{\alpha_{k}-s}{2}y-\frac{\partial_{z}}{2}+\sum_{i=1}^{k-1}\alpha_{i}\right)\left(\frac{a+s}{2}y+\frac{\partial_{z}}{2}+b-a\right)$
$\displaystyle(y-1)\frac{se^{s\eta
y}}{y}\frac{\cosh(z)-b\sinh(z)}{\sinh(z)^{3}}\Big{|}_{z\to\eta}$

By defining

$\displaystyle C_{1}=$
$\displaystyle\prod_{k=2}^{m}\left(-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(\frac{\partial_{z}}{2}+b\right)\frac{\alpha_{1}\sinh(\eta)+\cosh(\eta)}{\sinh(\eta)^{3}}$
(221) $\displaystyle C_{2}=$
$\displaystyle-\prod_{k=2}^{m}\left(-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k-1}\alpha_{k}\right)\left(\frac{\partial_{\eta}}{2}+b-a\right)\frac{\cosh(\eta)-b\sinh(\eta)}{\sinh(\eta)^{3}}$

$I_{1}$ can be finally simplified as

$\displaystyle I_{1}=$ $\displaystyle\sum_{s=\pm 1}\int_{-\infty}^{\infty}{\rm
d}x\,e^{-\frac{t}{4}(x^{2}-1)}\prod_{k=2}^{m}\left(\frac{\alpha_{k}-s}{2}x-\frac{\alpha_{k}+\partial_{z}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(\frac{a+s}{2}x+\frac{-a+\partial_{z}}{2}+b\right)$
(222) $\displaystyle
x\frac{(-(\alpha_{1}+b)+s(b\alpha_{1}+1))e^{sx}}{2\sinh(z)}\Big{|}_{z\to\eta}$
$\displaystyle+$ $\displaystyle\sum_{s=\pm 1}\int_{-\infty}^{\infty}{\rm
d}y\,e^{-\frac{t}{4}(y^{2}+2y)}\Bigg{(}\prod_{k=2}^{m}\left(\frac{\alpha_{k}-s}{2}y-\frac{\partial_{z}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(\frac{a+s}{2}y+\frac{\partial_{z}}{2}+b\right)$
$\displaystyle(y+1)\frac{\alpha_{1}\sinh(z)+\cosh(z)}{\sinh(z)^{3}}\Big{|}_{z\to\eta}-C_{1}\Bigg{)}\frac{se^{s\eta
y}}{2y}$ $\displaystyle-$ $\displaystyle\sum_{s=\pm
1}\int_{-\infty}^{\infty}{\rm
d}y\,e^{-\frac{t}{4}(y^{2}-2y)}\Bigg{(}\prod_{k=2}^{m}\left(\frac{\alpha_{k}-s}{2}y-\frac{\partial_{z}}{2}+\sum_{i=1}^{k-1}\alpha_{i}\right)\left(\frac{a+s}{2}y+\frac{\partial_{z}}{2}+b-a\right)$
$\displaystyle(y-1)\frac{\cosh(z)-b\sinh(z)}{\sinh(z)^{3}}\Big{|}_{z\to\eta}-C_{2}\Bigg{)}\frac{se^{s\eta
y}}{y}$
$\displaystyle+C_{1}\int_{-\infty}^{\infty}e^{-\frac{t}{4}(y^{2}+2y)}\frac{\sinh(\eta
y)}{y}{\rm
d}y-C_{2}\int_{-\infty}^{\infty}e^{-\frac{t}{4}(y^{2}-2y)}\frac{\sinh(\eta
y)}{y}{\rm d}y.$

In this expression, the integrands of first three terms involve only
polynomials of the integral variables, which can be computed by using the same
procedure as that to compute
$\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{\frac{1}{2}ab}\rangle$. For the
last two terms, by using the formula

$\displaystyle\int_{-\infty}^{\infty}{\rm
d}ye^{-\frac{t}{4}(y^{2}+by)}\frac{\sinh(y\eta)}{y}=\frac{1}{2}\pi\left(\text{erfi}\left(\frac{1+\frac{tb}{4\eta}}{2\sqrt{\frac{t}{4\eta^{2}}}}\right)+\text{erfi}\left(\frac{1-\frac{tb}{4\eta}}{2\sqrt{\frac{t}{4\eta^{2}}}}\right)\right),$
(223)

and

$\displaystyle\text{erfi}\left(\frac{1}{x}\right)=e^{\frac{1}{x^{2}}}\frac{x}{\pi}\sum_{n=0}^{\infty}\Gamma(\frac{1}{2}+n)x^{2n},$
(224)

a straightforward calculation can be performed as the following

$\displaystyle\int_{-\infty}^{\infty}{\rm d}y\,e^{-\frac{t}{4}(y^{2}\pm
2y)}\frac{\sinh(y\eta)}{y}=\langle
1\rangle\frac{1}{\sqrt{\pi}}\frac{\sinh(\eta)}{\eta^{2}}\sum_{m=0}^{\infty}\sum_{n=0}^{m}\Gamma(\frac{1}{2}+n)\frac{2^{n}}{\eta^{n+m}}\binom{n+m}{m-n}\Bigg{\\{}e^{-\eta}+e^{\eta}\left(-1\right)^{m-n}\Bigg{\\}}\left(\frac{t}{2}\right)^{m+2},$
(225)

where Eq. (29) is substitute because we are concerned about the expectation
value with respect to the normalized coherent state, i.e.
$\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{\iota ab}/\langle 1\rangle$.
This complete our computation for $I_{1}$.

Secondly, $I_{2}$ can be simplified as the following

$\displaystyle I_{2}=$ $\displaystyle\int_{-\infty}^{\infty}{\rm
d}xe^{-\frac{t}{4}(x+1)^{2}}\prod_{k=1}^{m}\left(-\alpha_{k}\frac{x+1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(\frac{x+1}{2}-(|b|-1-b)\frac{\partial_{\eta}}{2}+|b|\right)\frac{\sinh(x\eta)}{\sinh(\eta)}$
(226) $\displaystyle-$ $\displaystyle\int_{-\infty}^{\infty}{\rm
d}xe^{-\frac{t}{4}(x+1)^{2}}\prod_{k=1}^{m}\left(-\alpha_{k}\frac{x+1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right)\left(\frac{x+1}{2}+(|b|-1-b)\frac{\partial_{\eta}}{2}+|b|\right)\frac{(1-|b|-b)\partial_{\eta}}{x+1}\frac{\sinh(x\eta)}{\sinh(\eta)}.$

Because

$\frac{\partial_{\eta}}{x+1}\frac{\sinh(x\eta)}{\sinh(\eta)}=\frac{1}{\sinh(\eta)^{2}}\left(\sinh(\eta)\cosh(x\eta)-\frac{\sinh((x+1)\eta)}{x+1}\right)$
(227)

we have

$\displaystyle I_{2}=$ $\displaystyle\int_{-\infty}^{\infty}{\rm
d}xe^{-\frac{t}{4}(x+1)^{2}}\prod_{k=1}^{m}\left[-\alpha_{k}\frac{x+1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right]\left[\frac{x+1}{2}-(|b|-1-b)\frac{\partial_{\eta}}{2}+|b|\right]\times$
(228) $\displaystyle\frac{\sinh(x\eta)-(1-|b|-b)\cosh(x\eta)}{\sinh(\eta)}$
$\displaystyle+$ $\displaystyle\int_{-\infty}^{\infty}{\rm
d}ye^{-\frac{t}{4}y^{2}}\prod_{k=1}^{m}\left[-\alpha_{k}\frac{y}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right]\left[\frac{y}{2}-(|b|-1-b)\frac{\partial_{\eta}}{2}+|b|\right]\frac{(1-|b|-b)}{\sinh(\eta)^{2}}\frac{\sinh(y\eta)}{y}.$

By defining

$C=\left(\alpha_{1}-\frac{\partial_{\eta}}{2}\right)\left(\alpha_{1}+\alpha_{2}-\frac{\partial_{\eta}}{2}\right)\cdots\left((\alpha_{1}+\cdots+\alpha_{m}-\frac{\partial_{\eta}}{2}\right)\left(|b|-(|b|-1-b)\frac{\partial_{\eta}}{2}\right)\frac{1-|b|-b}{\sinh(\eta)^{2}},$
(229)

we obtain

$\displaystyle I_{2}=$ $\displaystyle\int_{-\infty}^{\infty}{\rm
d}xe^{-\frac{t}{4}(x+1)^{2}}\prod_{k=1}^{m}\left[-\alpha_{k}\frac{x+1}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{i}\right]\times$
(230)
$\displaystyle\left[\frac{x+1}{2}-(|b|-1-b)\frac{\partial_{\eta}}{2}+|b|\right]\frac{\sinh(x\eta)-(1-|b|-b)\cosh(x\eta)}{\sinh(\eta)}$
$\displaystyle+$ $\displaystyle\int_{-\infty}^{\infty}{\rm
d}ye^{-\frac{t}{4}y^{2}}\Bigg{\\{}\prod_{k=1}^{m}\left[-\alpha_{k}\frac{y}{2}-\frac{\partial_{\eta}}{2}+\sum_{i=1}^{k}\alpha_{k}\right]\left[\frac{y}{2}-(|b|-1-b)\frac{\partial_{\eta}}{2}+|b|\right]-C\Bigg{\\}}\frac{(1-|b|-b)}{\sinh(\eta)^{2}}\frac{\sinh(y\eta)}{y}$
$\displaystyle+$ $\displaystyle C\int_{-\infty}^{\infty}{\rm
d}ye^{-\frac{t}{4}y^{2}}\frac{\sinh(y\eta)}{y}.$

Similar as $I_{1}$, the first two integrals are computable because the
integrands therein are just polynomials multiplied by the Gaussian functions.
For the last term, we have

$\displaystyle C\int_{-\infty}^{\infty}{\rm
d}ye^{-\frac{t}{4}y^{2}}\frac{\sinh(y\eta)}{y}=C\mathrm{sgn}(\eta)\pi\text{erfi}\left(\frac{1}{2\sqrt{\frac{t}{4\eta^{2}}}}\right)$
(231) $\displaystyle=$ $\displaystyle\langle 1\rangle
C\frac{t^{2}}{2\sqrt{\pi}e^{t/4}}\frac{\sinh(\eta)}{\eta^{2}}\sum_{n=0}^{\infty}\Gamma(\frac{1}{2}+n)\left(\frac{t}{\eta^{2}}\right)^{n}.$

where Eq. (29) is substitute because we are concerned about the expectation
value with respect to the normalized coherent state, i.e.
$\langle\hat{F}^{\alpha_{1}\cdots\alpha_{m}}_{\iota ab}\rangle/\langle
1\rangle$. Therefore, $I_{2}$ is computable.

In summary, we obtain

$\displaystyle\langle\hat{F}_{1ab}^{\alpha_{1}\cdots\alpha_{m}}\rangle_{z_{e}}=-\frac{(-1)^{a}}{2}\delta\left(\sum_{i=1}^{m}\alpha_{i}-a+b\right)t^{m}e^{bz_{e}}\left(\prod_{i=1}^{m}\frac{1}{\sqrt{1+|\alpha_{i}|}}\frac{1}{\sqrt{1+|a|}\sqrt{1+|b|}}\right)I_{1}$
(232) $\displaystyle+$
$\displaystyle(-1)^{b-a}\delta\left(\sum_{i=1}^{m}\alpha_{i}-a+b\right)t^{m}e^{bz_{e}}\left(\prod_{i=1}^{m}\frac{1}{\sqrt{1+|\alpha_{i}|}}\frac{1}{\sqrt{1+|a|}\sqrt{1+|b|}}\right)I_{2}+O(t^{\infty})$

where $I_{1}$ and $I_{2}$ can be computed by using the algorithm introduced
above. Taking the operator $\hat{F}_{1ab}=D^{1}_{ab}(h_{e})$ as an example, we
can have

$\langle\hat{F}_{1ab}\rangle_{z_{e}}/\langle
1\rangle=\left\\{\begin{array}[]{cc}e^{i\xi}-\left(\frac{e^{i\xi}\tanh\left(\frac{\eta}{2}\right)}{2\eta}+\frac{1}{4}e^{i\xi}\right)t+O(t^{2}),&\
a=b=1,\\\ 1-\frac{\tanh\left(\frac{\eta}{2}\right)}{\eta}t+O(t^{2}),&\
a=b=0,\\\
e^{-i\xi}-\left(\frac{e^{-i\xi}\tanh\left(\frac{\eta}{2}\right)}{2\eta}+\frac{1}{4}e^{-i\xi}\right)t+O(t^{2}),&\
a=b=-1,\end{array}\right.$ (233)

which is compatible with the corresponding result given in Dapor:2017gdk .

##  Appendix F Mathematical supports for Sec. IV

To begin with, we need to study the matrix element of the holonomy and flux.
According to the results in thiemann2001gaugeIII , we have

$\displaystyle\langle\psi_{g_{e}}|\hat{p}_{s}^{\alpha}|\psi_{g_{e}^{\prime}}\rangle=i\langle\psi_{g_{e}}|\psi_{{g_{e}^{\prime}}}\rangle\frac{\mathrm{tr}(\tau^{\alpha}g_{e}^{\prime}g_{e}^{\dagger})}{\sinh(\zeta_{e}^{(0)})}\left(-\zeta_{e}^{(0)}+(\coth(\zeta_{e}^{(0)})-\frac{1}{\zeta_{e}^{(0)}})\frac{t}{2}\right)+O(t^{\infty}).$
(234)

where $\zeta_{e}^{(0)}$ is given by
$2\cosh(\zeta_{e}^{(0)})=\mathrm{tr}(g_{e}^{\prime}g_{e}^{\dagger})$. For
$\hat{p}_{t}^{\alpha}(e)$, one has

$\displaystyle\langle\psi_{g_{e}}|\hat{p}_{t}^{\alpha}|\psi_{g_{e}^{\prime}}\rangle=i\langle\psi_{g_{e}}|\psi_{{g_{e}^{\prime}}}\rangle\frac{\mathrm{tr}(g_{e}^{\dagger}g_{e}^{\prime}\tau^{\alpha})}{\sinh(\zeta_{e}^{(0)})}\left(\zeta_{e}^{(0)}-(\coth(\zeta_{e}^{(0)})-\frac{1}{\zeta_{e}^{(0)}})\frac{t}{2}\right)+O(t^{\infty}).$
(235)

For the holonomy operator $D^{\iota}_{ab}(h_{e})$, we have
thiemann2001gaugeIII

$\displaystyle\langle\psi_{g_{e}}|D^{\frac{1}{2}}_{ab}(h_{e})|\psi_{g_{e}^{\prime}}\rangle$
(236) $\displaystyle=$
$\displaystyle\langle\psi_{g_{e}}|\psi_{g_{e}^{\prime}}\rangle\left(\frac{2\mathrm{tr}(\tau^{k}g_{e}^{\prime}g_{e}^{\dagger})D^{\frac{1}{2}}_{ab}(\tau^{k}g_{e}^{\prime})}{\left(e^{\zeta_{e}^{(0)}}+1\right)^{2}\zeta_{e}^{(0)}}e^{\frac{\zeta_{e}^{(0)}}{2}}+\frac{D^{\frac{1}{2}}_{ab}(g_{e}^{\prime})}{2\zeta_{e}^{(0)}}e^{-\frac{\zeta_{e}^{(0)}}{2}}\right)e^{\frac{-t}{16}}\left(\left(e^{\zeta_{e}^{(0)}}+1\right)\zeta_{e}^{(0)}+(-e^{\zeta_{e}^{(0)}}+1)\frac{t}{4}\right)+O(t^{\infty}).$

According to Eqs. (234), (235) and (236), the matrix elements of the fluxes
and holomomies are of a form described below

$\langle\psi_{g_{e}}|\hat{O}_{i}|\psi_{g_{e}^{\prime}}\rangle=\langle\psi_{g_{e}}|\psi_{g_{e}^{\prime}}\rangle\left(E_{0}(g_{e},g_{e}^{\prime})+tE_{1}(g_{e},g_{e}^{\prime})+O(t^{\infty})\right).$
(237)

Recalling Eq. (56), we are going to consider integrals containing
$\frac{\langle\psi_{g_{e}}|\psi_{g_{e}^{\prime}}\rangle}{\|\psi_{g_{e}}\|\|\psi_{g_{e}^{\prime}}\|}$.
These integrals can be analyzed with the generalized stationary phase
approximation guided by Hörmander’s theorem 7.7.5 in hormander2015analysis .

######  Theorem F.1.

Let $K$ be a compact subset in $\mathbb{R}^{n}$, X be an open neighborhood of
K, and k be a positive integer. If (1) the complex functions $u\in
C^{2k}_{0}(K)$, $f\in C^{3k+1}(X)$ and $\Im(f)\geq 0$ in X, with $\Im(f)$
being the imaginary part of $f$; (2) there is a unique point $x_{0}\in K$
satisfying $\Im(f(x_{0}))=0$, $f^{\prime}(x_{0})=0$, and
$\det(f^{\prime\prime}(x_{0}))\neq 0$ ($f^{\prime\prime}$ denotes the Hessian
matrix), $f^{\prime}\neq 0$ in $K\backslash\left\\{x_{0}\right\\}$ then we
have the following estimation:

$\left|\int_{K}u(x)e^{i\lambda f(x)}dx-e^{i\lambda
f(x_{0})}\left[\mathrm{det}\left(\dfrac{\lambda f^{\prime\prime}(x_{0})}{2\pi
i}\right)\right]^{-\frac{1}{2}}\sum_{s=0}^{k-1}\left(\dfrac{1}{\lambda}\right)^{s}L_{s}u(x_{0})\right|\leq
C\left(\dfrac{1}{\lambda}\right)^{k}\sum_{|\alpha|\leq
2k}\mathrm{sup}\left|D^{\alpha}u\right|.$ (238)

Here the constant C is bounded when f stays in a bounded set in $C^{3k+1}(X)$.
We have used the standard multi-index notation
$\alpha=\langle\alpha_{1},...,\alpha_{n}\rangle$ and

$D^{\alpha}=(-i)^{\alpha}\frac{\partial^{|\alpha|}}{\partial
x_{1}^{\alpha_{1}}...\partial
x_{n}^{\alpha_{n}}},\quad\text{where}\quad|\alpha|=\sum_{i=1}^{n}\alpha_{i}$
(239)

$L_{s}u(x_{0})$ denotes the following operation on u:

$L_{s}u(x_{0})=i^{-s}\sum_{l-m=s}\sum_{2l\geq
3m}\frac{(-1)^{l}2^{-l}}{l!m!}\left[\sum_{a,b=1}^{n}H_{ab}^{-1}(x_{0})\dfrac{\partial^{2}}{\partial
x_{a}\partial x_{b}}\right]^{l}\left(g_{x_{0}}^{m}u\right)\left(x_{0}\right),$
(240)

where $H(x)=f^{\prime\prime}(x)$ denotes the Hessian matrix and the function
$g_{x_{0}}(x)$ is given by

$g_{x_{0}}(x)=f(x)-f(x_{0})-\frac{1}{2}H^{ab}\left(x_{0}\right)\left(x-x_{0}\right)_{a}\left(x-x_{0}\right)_{b}$
(241)

such that
$g_{x_{0}}\left(x_{0}\right)=g^{\prime}_{x_{0}}\left(x_{0}\right)=g^{\prime\prime}_{x_{0}}\left(x_{0}\right)=0$.
For each s, $L_{s}$ is a differential operator of order 2s acting on
$u\left(x\right)$.

###  F.1 Analysis of integrals containing
$\frac{\langle\psi_{g_{e}}|\psi_{g_{e}^{\prime}}\rangle}{\|\psi_{g_{e}}\|\|\psi_{g_{e}^{\prime}}\|}$

Define

$G(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)}):=\frac{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}{\|\psi_{g^{(1)}}\|\|\psi_{g^{(2)}}\|}$
(242)

with $g^{(k)}$ parameterized by

$g^{(k)}=e^{i\vec{p}^{(k)}\cdot\vec{\tau}}e^{\vec{\theta}^{(k)}\cdot\vec{\tau}},\
k=1,2.$ (243)

According to Eq. (28) and Eq. (29),
$G(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)})$ reads

$G(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)})=\frac{\zeta\sqrt{\sinh(p^{(1)})\sinh(p^{(2)})}}{\sqrt{p^{(1)}p^{(2)}}\sinh(\zeta)}e^{-\frac{-2\zeta^{2}+(p^{(1)})^{2}+(p^{(2)})^{2}}{2t}}$
(244)

where $p^{(i)}=\sqrt{\vec{p}^{(i)}\cdot\vec{p}^{(i)}}$ and $\zeta$ is given by

$2\cosh(\zeta)=\mathrm{tr}(g_{1}^{\dagger}g_{2}).$ (245)

Denote

$S(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)}):=-2\zeta^{2}+(p^{(1)})^{2}+(p^{(2)})^{2}.$

We first claim that

######  Lemma F.1.

$\Re(S(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)}))$,
the real part of
$S(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)})$, is
non-negative and vanishes iff $\vec{p}^{(2)}=\vec{p}^{(1)}$ and
$\vec{\theta}^{(2)}=\vec{\theta}^{(1)}$.

To prove this lemma, let us introduce the following proposition which is given
in thiemann2001gaugeII .

######  Proposition 3.

Let $g=e^{i\vec{p}\cdot\vec{\tau}}e^{\vec{\theta}\cdot\vec{\tau}}$.
Considering $\zeta=s+i\phi\in\mathbb{C}$ with $s\in\mathbb{R}$ and
$\phi\in[0,\pi]$ determined by $\cosh(\zeta)=\mathrm{tr}(g)$, we have, with
denoting $p:=\sqrt{\vec{p}\cdot\vec{p}}$,

$\delta=\frac{p^{2}}{4}-s^{2}+\phi^{2}\geq 0$ (246)

where the equality occurs iff
$\theta:=\sqrt{\vec{\theta}\cdot\vec{\theta}}=0$.

Thanks to this proposition, we prove Lemma F.1 as follows.

###### Proof of Lemma F.1.

$g_{1}^{\dagger}g_{2}$ can be decomposed as

$g_{1}^{\dagger}g_{2}=e^{i\vec{x}\cdot\vec{\tau}}e^{\vec{y}\cdot\vec{\tau}}.$

Denote $x=\sqrt{\vec{x}\cdot\vec{x}}$ and $y=\sqrt{\vec{y}\cdot\vec{y}}$. Then

$\Re(S(\vec{p},\vec{\theta},\vec{p},\vec{\theta}))=2\delta-\frac{x^{2}}{2}+(p^{(1)})^{2}+(p^{(2)})^{2}.$
(247)

where $\delta:=x^{2}/4-\Re(\zeta)^{2}$ is non-negative according to the
proposition 3. Thus we only need to prove that
$-x^{2}/2+(p^{(1)})^{2}+(p^{(2)})^{2}\geq 0$.

By definition, we have
$2\cosh(x)=\mathrm{tr}(g_{1}^{\dagger}g_{2}g_{2}^{\dagger}g_{1})$, which leads
to

$2\cosh(x)=\mathrm{tr}(e^{2i\vec{p}^{(1)}\cdot\vec{\tau}}e^{2i\vec{p}^{(2)}\cdot\vec{\tau}}).$
(248)

Since
$e^{i\vec{\mu}\cdot\vec{\tau}}=\cosh(\frac{\mu}{2})\mathbb{I}+2i\frac{\vec{\mu}\cdot\vec{\tau}}{\mu}\sinh(\frac{\mu}{2})$
with $\mu=\sqrt{\vec{\mu}\cdot\vec{\mu}}$, we have

$\displaystyle\cosh(x)=\frac{1-\beta}{2}\cosh(p^{(1)}-p^{(2)})+\frac{1+\beta}{2}\cosh(p^{(1)}+p^{(2)})\leq\cosh(p^{(1)}+p^{(2)})$
(249)

where $\beta=\frac{\vec{p}^{(1)}\cdot\vec{p}^{(2)}}{p^{(1)}p^{(2)}}\in[-1,1]$
and $\cosh(p^{(1)}+p^{(2)})\geq\cosh(p^{(1)}-p^{(2)})$ is used. Moreover,
because of $\sqrt{2(p^{(1)})^{2}+2(p^{(2)})^{2}}\geq p^{(1)}+p^{(2)}\geq 0$,
it has

$\cosh(\sqrt{2(p^{(1)})^{2}+2(p^{(2)})^{2}}\,)\geq\cosh(p^{(1)}+p^{(2)}).$
(250)

Combining the results of (249) and (250), one finally have

$-\frac{x^{2}}{2}+(p^{(1)})^{2}+(p^{(2)})^{2}\geq 0$ (251)

where the equality occurs only if $\vec{p}^{(1)}=\vec{p}^{(2)}$.

In summary we have

$\Re(S(\vec{p},\vec{\theta},\vec{p},\vec{\theta}))\geq 0$ (252)

and $\Re(S(\vec{p},\vec{\theta},\vec{p},\vec{\theta}))=0$ only if
$\vec{p}^{(1)}=\vec{p}^{(2)}$ and $\delta=0$ which means
$\vec{\theta}^{(1)}=\vec{\theta}^{(2)}$. ∎

It turns out below that the integrand of Eq. (56) consists of a Gaussian-like
function

$e^{-\frac{1}{2t}\left(S(\vec{p},\vec{\theta},\vec{p}^{(1)},\vec{\theta}^{(1)})+S(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)})+\cdots+S(\vec{p}^{(k)},\vec{\theta}^{(k)},\vec{p},\vec{\theta})\right)}.$
(253)

Lemma F.1 suggests us to do the stationary phase approximation analysis at
$\vec{p}^{(i)}=\vec{p}$ and $\vec{\theta}^{(i)}=\vec{\theta}$. Notice that
$\vec{p}$ and $\vec{\theta}$ are given to parameterize $g$ as
$g:=e^{i\vec{p}\cdot\vec{\tau}}e^{\vec{\theta}\cdot\vec{\tau}}$, and $g$
labels the coherent state $|\psi_{g}\rangle$ with respect to which the
expectation value of $\hat{O}$ is computed. Thus, rather than considering all
values of $\vec{p}$ and $\vec{\theta}$, it is sufficient to set

$\displaystyle\vec{p}_{o}=(0,0,p),\ \vec{\theta}_{o}=(0,0,\theta)$ (254)

according to (33). Denote

$f_{k;\vec{p},\vec{\theta}}(\vec{p}^{(1)},\vec{\theta}^{(1)},\cdots,\vec{p}^{(k)},\vec{\theta}^{(k)})\equiv
S(\vec{p},\vec{\theta},\vec{p}^{(1)},\vec{\theta}^{(1)})+S(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)})+\cdots+S(\vec{p}^{(k)},\vec{\theta}^{(k)},\vec{p},\vec{\theta}),$
(255)

we have the following result:

######  Theorem F.2.

  * (i)

$\Re(f_{k;\vec{p},\vec{\theta}})\geq 0$ and the equality occurs only when all
$g^{(i)}$ coincide, namely $\vec{p}^{\,(i)}=\vec{p}$ and
$\vec{\theta}^{(i)}=\vec{\theta}$.

  * (ii)

At $\vec{p}^{(i)}=\vec{p}_{o}=(0,0,p)$ and
$\vec{\theta}^{(i)}=\vec{\theta}_{o}=(0,0,\theta)$, it has

$\nabla_{\vec{p}^{(i)}}f_{k;\vec{p}_{o},\vec{\theta}_{o}}=0=\nabla_{\vec{\theta}^{(i)}}f_{k;\vec{p}_{o},\vec{\theta}_{o}},\
\forall i=1,\cdots,k.$ (256)

  * (iii)

The Hessian matrix $f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}$ of
$f_{k;\vec{p}_{o},\vec{\theta}_{o}}$ at $\vec{p}^{(i)}=\vec{p}_{o}=(0,0,p)$
and $\vec{\theta}^{(i)}=\vec{\theta}_{o}=(0,0,\theta)$ is non-degenerate with
the determinant

$\det\left(f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}\right)\Big{|}_{\vec{p}^{(i)}=\vec{p}_{o},\vec{\theta}^{(i)}=\vec{\theta}_{o}}=\left(\frac{1024\sin^{4}\left(\frac{\theta}{2}\right)}{\theta^{4}}\right)^{k}$
(257)

###### Proof.

The first statement is true by using Lemma F.1. For the second statement, let
us consider

$S(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)})=-2\zeta^{2}+(\vec{p}^{(1)})^{2}+(\vec{p}^{(2)})^{2}$
(258)

where $\zeta$ is given by

$\cosh(\zeta)=\frac{1}{2}\mathrm{tr}(e^{\vec{\theta}^{(2)}\cdot\vec{\tau}}e^{-\vec{\theta}^{(1)}\cdot\vec{\tau}}e^{i\vec{p}^{(1)}\cdot\vec{\tau}}e^{i\vec{p}^{(2)}\cdot\vec{\tau}}).$
(259)

Then it has that

$\displaystyle\frac{\partial\cosh(\zeta)}{\partial
p_{j}^{(1)}}\Bigg{|}_{\vec{p}_{o},\vec{\theta}_{o}}=\frac{\partial\cosh(\zeta)}{\partial
p_{j}^{(2)}}\Bigg{|}_{\vec{p}_{o},\vec{\theta}_{o}}$
$\displaystyle=\frac{\delta_{j,3}}{2}\sinh(p)$ (260)
$\displaystyle\frac{\partial\cosh(\zeta)}{\partial\theta_{j}^{(1)}}\Bigg{|}_{\vec{p}_{o},\vec{\theta}_{o}}=-\frac{\partial\cosh(\zeta)}{\partial\theta_{j}^{(2)}}\Bigg{|}_{\vec{p}_{o},\vec{\theta}_{o}}$
$\displaystyle=\frac{i\delta_{j,3}}{2}\sinh(p)$

where the subscript $\vec{p}_{o},\vec{\theta}_{o}$ indicates to take values at
$\vec{p}^{(1)}=\vec{p}_{o}=\vec{p}^{(2)}$ and
$\vec{\theta}^{(1)}=\vec{\theta}_{o}=\vec{\theta}^{(2)}$. According to Eq.
(260), we have

$\displaystyle\nabla_{\vec{p}^{(1)}}S(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)})\Big{|}_{\vec{p}_{o},\vec{\theta}_{o}}=0=\nabla_{\vec{p}^{(2)}}S(\vec{p}^{(1)},\vec{\theta}^{(1)},\vec{p}^{(2)},\vec{\theta}^{(2)})\Big{|}_{\vec{p}_{o},\vec{\theta}_{o}}$
(261)

and, thus,
$\nabla_{\vec{p}^{(i)}}f_{k;\vec{p}_{o},\vec{\theta}_{o}}\big{|}_{\vec{p}_{o},\vec{\theta}_{o}}=0$
for all $i=1,2,\cdots,k$. For
$\nabla_{\vec{\theta}^{(i)}}f_{k;\vec{p}_{o},\vec{\theta}_{o}}\big{|}_{\vec{p}_{o},\vec{\theta}_{o}}$,
the second equation in Eq. (260) gives us

$\nabla_{\vec{\theta}^{(i)}}f_{k;\vec{p}_{o},\vec{\theta}_{o}}\big{|}_{\vec{p}_{o},\vec{\theta}_{o}}=\nabla_{\vec{\theta}^{(i)}}S(\vec{p}^{(i-1)},\vec{\theta}^{(i-1)},\vec{p}^{(i)},\vec{\theta}^{(i)})\Big{|}_{\vec{p}_{o},\vec{\theta}_{o}}+\nabla_{\vec{\theta}^{(i)}}S(\vec{p}^{(i)},\vec{\theta}^{(i)},\vec{p}^{(i+1)},\vec{\theta}^{(i+1)})\Big{|}_{\vec{p}_{o},\vec{\theta}_{o}}=0,$
(262)

which completes the proof of the second statement.

For the last statement, using the conclusion of the second statement, we can
immediately get that

$\frac{\partial f_{k;\vec{p}_{o},\vec{\theta}_{o}}}{\partial
x_{m}^{(i)}\partial y_{n}^{(j)}}\Bigg{|}_{\vec{p}_{o},\vec{\theta}_{o}}=0,\
|i-j|>1$ (263)

where $x_{m}^{(i)}$ and $y_{n}^{(j)}$ represent $\vec{\theta}_{m}^{(i)}$ or
$\vec{p}_{m}^{(i)}$. Therefore, if we order the arguments
$\vec{p}^{(i)},\vec{\theta}^{(i)}$ as

$p^{(1)}_{1},p^{(1)}_{2},p^{(1)}_{3},\theta^{(1)}_{1},\theta^{(1)}_{2},\theta^{(1)}_{3},p^{(2)}_{1},p^{(2)}_{2},p^{(2)}_{3},\theta^{(2)}_{1},\theta^{(2)}_{2},\theta^{(2)}_{3},\cdots.$
(264)

to arrange the matrix elements of
$f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}(\vec{p}_{o},\vec{\theta}_{o})\equiv
f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}\Big{|}_{\vec{p}_{o},\vec{\theta}_{o}}$,
the resulting matrix is block-tridiagonal matrix. Moreover, since all of the
$p$-arguments, as well as the $\theta$-arguments, in
$f_{k;\vec{p}_{o},\vec{\theta}_{o}}$ are symmetric, we conclude that

$\frac{\partial f_{k;\vec{p}_{o},\vec{\theta}_{o}}}{\partial
x_{m}^{(i)}\partial
y_{l}^{(j)}}\Bigg{|}_{\vec{p}_{o},\vec{\theta}_{o}}=\frac{\partial
f_{k;\vec{p}_{o},\vec{\theta}_{o}}}{\partial x_{m}^{(i^{\prime})}\partial
y_{l}^{(j^{\prime})}}\Bigg{|}_{\vec{p}_{o},\vec{\theta}_{o}}.$ (265)

Consequently
$f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}(\vec{p}_{o},\vec{\theta}_{o})$
takes the form

$f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}(\vec{p}_{o},\vec{\theta}_{o})=\left(\begin{array}[]{cccccc}A&B&0&0&\cdots&0\\\
B^{T}&A&B&0&\cdots&0\\\ 0&B^{T}&A&B&\cdots&0\\\
\vdots&\vdots&\vdots&\ddots&\vdots&\vdots\\\ 0&0&\cdots&B^{T}&A&B\\\
0&0&\cdots&0&B^{T}&A\end{array}\right)$ (266)

where $A$ and $B$ are $6\times 6$ matrix. The matrix $A$ and $B$ can be
calculated by considering the case with $k=2$, which gives us

$A=\left(\begin{array}[]{cccccc}\frac{4\tanh\left(\frac{p}{2}\right)}{p}&0&0&-\frac{4\sin^{2}\left(\frac{\theta}{2}\right)\tanh\left(\frac{p}{2}\right)}{\theta}&-\frac{2\sin(\theta)\tanh\left(\frac{p}{2}\right)}{\theta}&0\\\
0&\frac{4\tanh\left(\frac{p}{2}\right)}{p}&0&\frac{2\sin(\theta)\tanh\left(\frac{p}{2}\right)}{\theta}&-\frac{4\sin^{2}\left(\frac{\theta}{2}\right)\tanh\left(\frac{p}{2}\right)}{\theta}&0\\\
0&0&2&0&0&0\\\
-\frac{4\sin^{2}\left(\frac{\theta}{2}\right)\tanh\left(\frac{p}{2}\right)}{\theta}&\frac{2\sin(\theta)\tanh\left(\frac{p}{2}\right)}{\theta}&0&-\frac{4p(\cos(\theta)-1)\coth(p)}{\theta^{2}}&0&0\\\
-\frac{2\sin(\theta)\tanh\left(\frac{p}{2}\right)}{\theta}&-\frac{4\sin^{2}\left(\frac{\theta}{2}\right)\tanh\left(\frac{p}{2}\right)}{\theta}&0&0&-\frac{4p(\cos(\theta)-1)\coth(p)}{\theta^{2}}&0\\\
0&0&0&0&0&2\\\ \end{array}\right)$

and

$B=\left(\begin{array}[]{cccccc}-\frac{2\tanh\left(\frac{p}{2}\right)}{p}&0&0&\frac{2\sin^{2}\left(\frac{\theta}{2}\right)\tanh\left(\frac{p}{2}\right)+i\sin(\theta)}{\theta}&\frac{\sin(\theta)\tanh\left(\frac{p}{2}\right)+i(\cos(\theta)-1)}{\theta}&0\\\
0&-\frac{2\tanh\left(\frac{p}{2}\right)}{p}&0&\frac{-i\cos(\theta)-\sin(\theta)\tanh\left(\frac{p}{2}\right)+i}{\theta}&\frac{2\sin^{2}\left(\frac{\theta}{2}\right)\tanh\left(\frac{p}{2}\right)+i\sin(\theta)}{\theta}&0\\\
0&0&-1&0&0&i\\\
\frac{2\sin^{2}\left(\frac{\theta}{2}\right)\tanh\left(\frac{p}{2}\right)-i\sin(\theta)}{\theta}&\frac{i\left(\cos(\theta)+i\sin(\theta)\tanh\left(\frac{p}{2}\right)-1\right)}{\theta}&0&\frac{2p(\cos(\theta)-1)\coth(p)}{\theta^{2}}&-\frac{2ip(\cos(\theta)-1)}{\theta^{2}}&0\\\
\frac{-i\cos(\theta)+\sin(\theta)\tanh\left(\frac{p}{2}\right)+i}{\theta}&\frac{2\sin^{2}\left(\frac{\theta}{2}\right)\tanh\left(\frac{p}{2}\right)-i\sin(\theta)}{\theta}&0&\frac{2ip(\cos(\theta)-1)}{\theta^{2}}&\frac{2p(\cos(\theta)-1)\coth(p)}{\theta^{2}}&0\\\
0&0&-i&0&0&-1\\\ \end{array}\right).$

With the expression of $A$ and $B$, it can be verified that

$BA^{-1}B^{T}=B^{T}A^{-1}B=0.$ (267)

To calculate
$\det(f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}(\vec{p}_{o},\vec{\theta}_{o}))$,
we define matrices of $\tilde{B}$, $\tilde{C}$ and D of dimensions $6\times
6(k-1)$, $6(k-1)\times 6$ and $6(k-1)\times 6(k-1)$ respectively such that

$f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}(\vec{p}_{o},\vec{\theta}_{o})=\left(\begin{array}[]{cc}A&\tilde{B}\\\
\tilde{C}&D\end{array}\right).$ (268)

Then by using the property of the Schur complement, it has

$\det(f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}(\vec{p}_{o},\vec{\theta}_{o}))=\det(A)\det(D-\tilde{C}A^{-1}\tilde{B})=\det(A)\det(D)$
(269)

where we used $\tilde{C}A^{-1}\tilde{B}=0$ because of Eq. (267) and that
$f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}$ is block-tridiagonal
matrix. Because $D$ is the Hessian matrix
$f^{\prime\prime}_{k-1;\vec{p}_{o},\vec{\theta}_{o}}(\vec{p}_{o},\vec{\theta}_{o})$,
we finally have

$\det(f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}(\vec{p}_{o},\vec{\theta}_{o}))=\det(A)^{k}=\left(\frac{1024\sin^{4}(\frac{\theta}{2})}{\theta^{4}}\right)^{k}.$
(270)

∎

By these results, the stationary phase approximation method introduced in
Theorem F.1 can be applied to calculate the integral (56).

Taking advantage of the above results, we now can come to the proof of Theorem
IV.1.

###  F.2 proof of Theorem IV.1

As in Eq. (56), it has

$\frac{\langle\Psi_{\vec{w}}|\hat{O}|\Psi_{\vec{w}}\rangle}{\langle\Psi_{\vec{w}}|\Psi_{\vec{w}}\rangle}=\int\prod_{j=1}^{k-1}\prod_{e\in
E(\gamma)}\frac{2}{\pi t^{3}}{\rm d}\mu_{H}(u^{(j)}_{e}){\rm
d}^{3}\vec{p}_{e}^{(j)}\prod_{i=1}^{k}\frac{\langle\Psi_{\vec{g}^{(i-1)}}|\hat{O}_{i}|\Psi_{\vec{g}^{(i)}}\rangle}{\|\Psi_{\vec{g}^{(i-1)}}\|\|\Psi_{\vec{g}^{(i)}}\|}$
(271)

where we denoted
$|\Psi_{\vec{g}^{(0)}}\rangle=|\Psi_{\vec{g}^{(k)}}\rangle:=|\Psi_{\vec{w}}\rangle$,
applied Eq. (31) and used the decomposition

$g_{e}^{(i)}=e^{i\vec{p}^{(i)}_{e}\cdot\vec{\tau}}e^{\vec{\theta}^{(i)}_{e}\cdot\vec{\tau}}=e^{i\vec{p}^{(i)}_{e}\cdot\vec{\tau}}u_{e}^{(i)}.$
(272)

By the assumption, we have

$\frac{\langle\Psi_{\vec{g}^{(i-1)}}|\hat{O}_{i}|\Psi_{\vec{g}^{(i)}}\rangle}{\|\Psi_{\vec{g}^{(i-1)}}\|\|\Psi_{\vec{g}^{(i)}}\|}=\frac{\langle\Psi_{\vec{g}^{(i-1)}}|\Psi_{\vec{g}^{(i)}}\rangle}{\|\Psi_{\vec{g}^{(i-1)}}\|\|\Psi_{\vec{g}^{(i)}}\|}\left(E_{0}^{(i)}(\vec{g}^{(i-1)},\vec{g}^{(i)})+tE_{1}^{(i)}(\vec{g}^{(i-1)},\vec{g}^{(i)})+O(t^{\infty})\right).$
(273)

Thus Eq. (271) takes the form

$\displaystyle\frac{\langle\Psi_{\vec{w}}|\hat{O}|\Psi_{\vec{w}}\rangle}{\langle\Psi_{\vec{w}}|\Psi_{\vec{w}}\rangle}=\int\prod_{j=1}^{k-1}\prod_{e\in
E(\gamma)}\frac{2}{\pi t^{3}}{\rm d}\mu_{H}(u^{(j)}_{e}){\rm
d}^{3}\vec{p}_{e}^{(j)}\prod_{i=1}^{k}\frac{\langle\Psi_{\vec{g}^{(i-1)}}|\Psi_{\vec{g}^{(i)}}\rangle}{\|\Psi_{\vec{g}^{(i-1)}}\|\|\Psi_{\vec{g}^{(i)}}\|}\prod_{l=1}^{k}E^{(l)}(\vec{g}^{(l-1)},\vec{g}^{(l)})$
(274)

with $P(\\{\vec{g}^{(i)}\\}_{i=1}^{k})$ denoting the function

$E^{(l)}(\vec{g}^{(l-1)},\vec{g}^{(l)}):=\left(E_{0}^{(l)}(\vec{g}^{(l-1)},\vec{g}^{(l)})+tE_{1}^{(l)}(\vec{g}^{(l-1)},\vec{g}^{(l)})+O(t^{\infty})\right).$
(275)

Eq. (274) can be analyzed with the stationary phase approximation according to
Theorem F.2. It should be noticed that the Haar measure ${\rm d}\mu_{H}$ can
be expressed as

${\rm
d}\mu_{H}(u)=\frac{\sin^{2}\left(\frac{1}{2}\sqrt{\vec{\theta}\cdot\vec{\theta}}\right)}{4\pi^{2}\,(\vec{\theta}\cdot\vec{\theta})}{\rm
d}^{3}\vec{\theta}$ (276)

where $u\in\mathrm{SU}(2)$ is coordinatized as
$u=e^{\vec{\theta}\cdot\vec{\tau}}$ and ${\rm d}^{3}\vec{\theta}$ is the
Lebesgue measure on $\mathbb{R}^{3}$. Substituting the last equation into Eq.
(274) and applying Eq. (238) as well as Eq. (257), we finally obtain

$\displaystyle\frac{\langle\Psi_{\vec{w}}|\hat{O}|\Psi_{\vec{w}}\rangle}{\langle\Psi_{\vec{w}}|\Psi_{\vec{w}}\rangle}=\prod_{e\in
E(\gamma)}\left(\frac{\theta_{e}^{2}}{\sin^{2}(\theta_{e}/2)}\right)^{k-1}\,\sum_{s=0}^{l-1}\left(2t\right)^{s}\mathfrak{D}^{(s)}\Big{|}_{g_{e}^{(k)}=e^{iw_{e}\tau_{3}},\forall
k,e}+O(t^{l})$ (277)

where $\mathfrak{D}^{(s)}$ takes the form

$\mathfrak{D}^{(s)}=(-1)^{s}\sum_{l-j=s}\sum_{2l\geq
3j}\frac{(-1)^{l}2^{-l}}{l!j!}\left[\sum_{a,b=1}^{n}H_{ab}^{-1}(x_{0})\dfrac{\partial^{2}}{\partial
x_{a}\partial
x_{b}}\right]^{l}\left(G^{(j)}\prod_{i=1}^{k}E^{(i)}\right)\Bigg{|}_{g_{e}^{(l)}=e^{iw_{e}\tau_{3}},\forall
l,e}$ (278)

with $G^{(j)}$ being defined as the following,

$G^{(j)}\left(\\{\vec{g}^{(i)}\\}_{i=k}^{k-1}\right)=\left(\prod_{e\in
E(\gamma)}\prod_{n=1}^{k-1}\frac{\sin^{2}(|\vec{\theta}_{e}^{(n)}|/2)}{|\vec{\theta}_{e}^{(n)}|^{2}}\right)\left(\prod_{e\in
E(\gamma)}\prod_{n=0}^{k-1}\frac{\zeta_{e}^{(n,n+1)}\sqrt{\sinh(p_{e}^{(n)})\sinh(p_{e}^{(n+1)})}}{\sqrt{p_{e}^{(n)}p_{e}^{(n+1)}}\sinh(\zeta_{e}^{(n,n+1)})}\right)g_{x_{0}}(\\{\vec{g}^{(i)}\\}_{i=1}^{k-1})^{j}.$
(279)

Here $g_{x_{0}}$ is some function defined by applying (241) to the current
case.

If the leading order term of $\mathfrak{D}^{(s)}$ is claimed to be
$O(t^{d_{s}})$, then each term in the summation over $s$ of Eq. (277) is
$O(t^{s+d_{s}})$. Moreover, for each $l$ in Eq. (278), the derivative acting
on $G^{(j)}\prod E^{(i)}$ is of order $2l$ in total. Because of the properties
given by Eq.(241), the non-vanishing result appears when there are at least
$3j$ derivatives acting on $G^{(j)}$, which indicates that the order of
derivative that acts on the term $\prod E^{(i)}$ in Eq. (278) is $2l-3j$.
Because $2l-3j=s$ and $l=s+j$, there are at most $2s-j$ derivatives acting on
$\prod E^{(i)}$. Further, since $j\geq 0$, the maximum order of derivative
acting on the term $\prod E^{(i)}$ is $2s$, which only occurs for $l=s$. Then,
let us count the leading order of $\mathfrak{D}^{(s)}$ for a given $s$.
According to the expression of $\mathfrak{D}^{(s)}$, once it is evaluated at
the critical point given by $g_{e}^{(k)}=e^{i\omega_{e}\tau_{3}}$, those
$E^{(m)}$ contributed by operators $\hat{O}_{m}$ satisfying (61) will
inevitably increase the power of its leading order term if they are not acted
by any derivative operators. For a fixed $s$ there are at least $n_{s}$ of
these “non-acted” terms with

$n_{s}=\frac{(N_{0}-2s)+|N_{0}-2s|}{2}.$ (280)

Therefore,
$\dfrac{\langle\Psi_{\vec{w}}|\hat{O}|\Psi_{\vec{w}}\rangle}{\langle\Psi_{\vec{w}}|\Psi_{\vec{w}}\rangle}$
is of order of $t^{n}$ with

$n\geq\min\limits_{s\in\mathbb{Z}^{+}}(s+n_{s})=\lfloor\frac{N_{0}+1}{2}\rfloor.$
(281)

###  F.3 proof of Theorem IV.2

At first, denote $\hat{p}^{0}(e)$,
$D^{\frac{1}{2}}_{\frac{1}{2}\frac{1}{2}}(h_{e})$ and
$D^{\frac{1}{2}}_{-\frac{1}{2}-\frac{1}{2}}(h_{e})$ by $\hat{O}_{a}^{\rm d}$
with $a=1,2,3$ respectively, and $p^{\pm 1}(e)$,
$D^{\frac{1}{2}}_{\frac{1}{2}-\frac{1}{2}}(h_{e})$ or
$D^{\frac{1}{2}}_{-\frac{1}{2}\frac{1}{2}}(h_{e})$ by $\hat{O}_{a}^{\rm nd}$
with $a=1,2,3,4$ respectively. According to (234), (235) and (236), the matrix
elements of $\hat{O}_{a}^{i}$ take the form

$\langle\psi_{g^{(1)}}|\hat{O}_{a}^{i}|\psi_{g^{(2)}}\rangle=\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle
E^{i}_{a;0}(g^{(1)},g^{(2)})+O(t),\ \forall i={\rm d,nd}.$ (282)

Then by formulae listed in Sec. F.4, we obtain:

######  Lemma F.2.

Given
$g^{(i)}=e^{i\vec{p}^{(i)}\cdot\vec{\tau}}e^{i\vec{\theta}^{(i)}\cdot\vec{\tau}}$,
we have

$\displaystyle(\partial_{x^{(j)}_{k}}E_{a;0}^{\rm
d})(e^{iw\tau_{3}},e^{iw\tau_{3}})$ $\displaystyle=0,\ j=1,2,\text{ and
}k=1,2$ (283) $\displaystyle(\partial_{x^{(j)}_{3}}E^{\rm
nd}_{a;0})(e^{iw\tau_{3}},e^{iw\tau_{3}})$ $\displaystyle=0,\ j=1,2$

for all $w=p-i\theta\in\mathbb{C}$, where $x^{(i)}_{j}$ denotes $p^{(i)}_{j}$
or $\theta^{(i)}_{j}$. Moreover, consider the matrix $A$ and $B$ in Eq. (266).
$E^{\rm nd}_{a;0}$ satisfies that

$\displaystyle\nabla_{\vec{x}^{(1)}}^{T}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}})A^{-1}B$
$\displaystyle=0=\nabla_{\vec{x}^{(2)}}^{T}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}})A^{-1}B^{T},$ (284)
$\displaystyle\nabla_{\vec{x}^{(1)}}^{T}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}})A^{-1}B^{T}$
$\displaystyle=-\nabla_{\vec{x}^{(1)}}^{T}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}}),$
$\displaystyle\nabla_{\vec{x}^{(2)}}^{T}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}})A^{-1}B$
$\displaystyle=-\nabla_{\vec{x}^{(2)}}^{T}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}}),$ $\displaystyle
B^{T}A^{-1}\nabla_{\vec{x}^{(1)}}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}})$
$\displaystyle=0=BA^{-1}\nabla_{\vec{x}^{(2)}}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}}),$ $\displaystyle
BA^{-1}\nabla_{\vec{x}^{(1)}}E^{\rm nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}})$
$\displaystyle=-\nabla_{\vec{x}^{(1)}}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}}),$ $\displaystyle
B^{T}A^{-1}\nabla_{\vec{x}^{(2)}}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}})$
$\displaystyle=-\nabla_{\vec{x}^{(2)}}E^{\rm
nd}_{a;0}(e^{iw\tau_{k}},e^{iw\tau_{3}}),$

where
$\nabla_{\vec{x}^{(i)}}=(\partial_{\vec{p}_{1}^{(i)}},\partial_{\vec{p}_{2}^{(i)}},\partial_{\vec{p}_{3}^{(i)}},\partial_{\vec{\theta}_{1}^{(i)}},\partial_{\vec{\theta}_{2}^{(i)}},\partial_{\vec{\theta}_{3}^{(i)}})^{T}$
and $\nabla_{\vec{x}^{(i)}}^{T}$ is its transpose.

The second lemma is on the inverse of the Hessian matrix (266). Recalling Eq.
(268), one has that the inverse of
$f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}$, denoted as
$H_{k,\vec{p}_{o},\vec{\theta}_{o}}^{-1}$, is

$H_{k,\vec{p}_{o},\vec{\theta}_{o}}^{-1}=\left(\begin{array}[]{cc}A^{-1}+A^{-1}\tilde{B}(D-\tilde{C}A^{-1}\tilde{B})^{-1}\tilde{C}A^{-1}&-A^{-1}\tilde{B}(D-\tilde{C}A^{-1}\tilde{B})^{-1}\\\
-(D-\tilde{C}A^{-1}\tilde{B})^{-1}\tilde{C}A^{-1}&(D-\tilde{C}A^{-1}\tilde{B})^{-1}\end{array}\right).$
(285)

Since $\tilde{C}A^{-1}\tilde{B}=0$ and
$D=f^{\prime\prime}_{k-1;\vec{p}_{o},\vec{\theta}_{o}}$, we have

$H_{k,\vec{p}_{o},\vec{\theta}_{o}}^{-1}=\left(\begin{array}[]{cc}A^{-1}+A^{-1}\tilde{B}H_{k-1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\tilde{C}A^{-1}&-A^{-1}\tilde{B}H_{k-1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\\\
-H_{k-1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\tilde{C}A^{-1}&H_{k-1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\end{array}\right).$
(286)

For $k=1$, $H_{1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}=A^{-1}$. Thus one has that
$\tilde{B}H_{1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\tilde{C}=0$ and

$H_{2,\vec{p}_{o},\vec{\theta}_{o}}^{-1}=\left(\begin{array}[]{cc}A^{-1}&-A^{-1}\tilde{B}H_{1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\\\
-H_{1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\tilde{C}A^{-1}&H_{1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\end{array}\right).$

Doing this successively, one has
$\tilde{B}H_{k-1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\tilde{C}=0$, and

$H_{k,\vec{p}_{o},\vec{\theta}_{o}}^{-1}=\left(\begin{array}[]{cc}A^{-1}&-A^{-1}\tilde{B}H_{k-1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\\\
-H_{k-1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\tilde{C}A^{-1}&H_{k-1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}\end{array}\right).$
(287)

Finally, $H_{k;\vec{p}_{o},\vec{\theta}_{o}}^{-1}$ can be obtained with this
recurrence relation and the initial data
$H_{1,\vec{p}_{o},\vec{\theta}_{o}}^{-1}=A^{-1}$. The result is as follows:

######  Lemma F.3.

$H_{k;\vec{p}_{o},\vec{\theta}_{o}}^{-1}$ satisfies that

$(H_{k;\vec{p}_{o},\vec{\theta}_{o}}^{-1})_{mn}=\left\\{\begin{array}[]{cc}(-1)^{|m-n|}A^{-1}(BA^{-1})^{|m-n|},&m<n\\\
A^{-1},&m=n\\\ (-1)^{|m-n|}A^{-1}(B^{T}A^{-1})^{|m-n|},&m>n\\\
\end{array}\right.$ (288)

where $H_{k,\vec{p}_{o},\vec{\theta}_{o}}^{-1}$ is arranged as a block matrix
as $f^{\prime\prime}_{k;\vec{p}_{o},\vec{\theta}_{o}}$ in (266)], with
$(H_{k;\vec{p}_{o},\vec{\theta}_{o}}^{-1})_{mn}$ as a block.

Now the theorem (IV.2) can be proven.

###### Proof of Theorem IV.2.

For convenience, we define $s_{o}=M_{+}+N_{+}$. By Theorem IV.1,
$\langle\mathcal{M}\rangle_{z_{e}}$ is of order $t^{s_{o}}$ or higher.
Adopting the result from equation (277), $O(t^{s_{o}})$ only occurs when
$s=s_{o}$ and $\mathfrak{D}^{(s_{o})}$ is of $O(t^{0})$.

Set $s=s_{o}$ in the definition (278) of $\mathfrak{D}^{(s)}$. For given $l$
and $j$, if there is one $E^{(m)}$ taking the form of $E^{\rm nd}_{a;0}+O(t)$
not being acted by derivatives, then the eventual evaluation at the critical
point will vanish. Moreover, in $\mathfrak{D}^{(s_{o})}$, on one hand it
contains $2s_{o}$ of $E^{\rm nd}_{a;0}+O(t)$, and on the other hand the
maximum order of derivative that can act on $\prod E^{(i)}$ is $2s_{o}$, which
only occurs when $l=s_{o}$. Therefore, only when $l=s_{o}$ and $m=0$, all
$E^{\rm nd}_{a;0}+O(t)$ are acted by derivatives. Finally,

$\displaystyle\mathfrak{D}^{(s_{o})}=\left(\frac{\sin^{2}(\frac{\theta}{2})}{\theta^{2}}\right)^{|\mathcal{M}|-1}\frac{2^{-s_{o}}}{s_{o}!}\left[\sum_{a,b=1}^{6(|\mathcal{M}|-1)}H^{-1}_{ab}\frac{\partial^{2}}{\partial
x_{a}\partial
x_{b}}\right]^{s_{o}}\prod_{i=1}^{|\mathcal{M}|-1}E^{(i)}\Bigg{|}_{g^{(l)}=e^{iw\tau_{3}},\forall
l},$ (289)

where $|\mathcal{M}|$ denotes the number of factors in the monomial
$\mathcal{M}$.

To calculate (289), we employ the notion introduced in Lemma F.2 and F.3 to
treat $H^{-1}$ as a block matrix. Then
$\sum_{a,b=1}^{6(|\mathcal{M}|-1)}H^{-1}_{ab}\frac{\partial^{2}}{\partial
x_{a}\partial x_{b}}$ is rewritten as

$\sum_{a,b=1}^{6(|\mathcal{M}|-1)}H^{-1}_{ab}\frac{\partial^{2}}{\partial
x_{a}\partial
x_{b}}=\sum_{m,n=1}^{|\mathcal{M}|-1}\nabla_{\vec{x}^{(m)}}^{T}\left(\hat{H}^{-1}_{|\mathcal{M}|-1;\vec{p}_{o},\vec{\theta}_{o}}\right)_{mn}\nabla_{\vec{x}^{(n)}}.$
(290)

We then expand
$\left[\sum_{m,n=1}^{|\mathcal{M}|-1}\nabla_{\vec{x}^{(m)}}^{T}\left(\hat{H}^{-1}_{|\mathcal{M}|-1;\vec{p}_{o},\vec{\theta}_{o}}\right)_{mn}\nabla_{\vec{x}^{(n)}}\right]^{s_{o}}$
and let each individual term of the expansion act on
$\prod_{i=1}^{|\mathcal{M}|-1}E^{(i)}$. In each individual term of the
expansion, it contains derivative with respect to certain $\vec{x}^{(q)}$.
Because each $E^{\rm nd}_{a;0}+O(t)$ only depends on certain $\vec{x}^{(q)}$,
we only consider the case when all derivatives are paired with all $E^{\rm
nd}_{a;0}+O(t)$ with the same argument. For the other terms of the expansion,
they give vanishing results because of the evaluation at the critical point.

The procedure mentioned above is equivalent to the follows. We first partition
these $E^{\rm nd}_{a;0}+O(t)$ into ordering pairs. Denote all possibilities of
the partition as $\mathcal{P}$. Given a pair $(E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)})+O(t),E^{\rm
nd}_{a_{n};0}(\vec{x}^{(n-1)},\vec{x}^{(n)})+O(t))$ in a partition
$p\in\mathcal{P}$. It can be acted by

$\displaystyle\nabla_{\vec{x}^{(m)}}^{T}(H_{\mathcal{M}-1;\vec{p}_{o},\vec{\theta}_{o}}^{-1})_{mn}\nabla_{\vec{x}^{(n)}},$
(291)
$\displaystyle\nabla_{\vec{x}^{(m-1)}}^{T}(H_{\mathcal{M}-1;\vec{p}_{o},\vec{\theta}_{o}}^{-1})_{m-1,n}\nabla_{\vec{x}^{(n)}},$
(292)
$\displaystyle\nabla_{\vec{x}^{(m)}}^{T}(H_{\mathcal{M}-1;\vec{p}_{o},\vec{\theta}_{o}}^{-1})_{m,n-1}\nabla_{\vec{x}^{(n-1)}},$
(293)
$\displaystyle\nabla_{\vec{x}^{(m-1)}}^{T}(H_{\mathcal{M}-1;\vec{p}_{o},\vec{\theta}_{o}}^{-1})_{m-1,n-1}\nabla_{\vec{x}^{(n-1)}}$
(294)

According to Lemmas F.3 and F.2,

  * (1)

If $m<n$, only the operator (293) gives non-vanishing results which reads

$\displaystyle\nabla_{\vec{x}^{(m)}}^{T}E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)})(H_{\mathcal{M}-1;\vec{p}_{o},\vec{\theta}_{o}}^{-1})_{m,n-1}\nabla_{\vec{x}^{(n-1)}}E^{\rm
nd}_{a_{n};0}(\vec{x}^{(n-1)},\vec{x}^{(n)})$ (295) $\displaystyle=$
$\displaystyle\nabla_{\vec{x}^{(m)}}^{T}E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)})A^{-1}\nabla_{\vec{x}^{(n-1)}}E^{\rm
nd}_{a_{n};0}(\vec{x}^{(n-1)},\vec{x}^{(n)})$

  * (2)

If $m>n$, only the operator (292) gives non-vanishing results which reads

$\displaystyle\nabla_{\vec{x}^{(m-1)}}^{T}E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)})(H_{\mathcal{M}-1;\vec{p}_{o},\vec{\theta}_{o}}^{-1})_{m-1,n}\nabla_{\vec{x}^{(n)}}E^{\rm
nd}_{a_{n};0}(\vec{x}^{(n-1)},\vec{x}^{(n)})$ (296) $\displaystyle=$
$\displaystyle\nabla_{\vec{x}^{(m-1)}}^{T}E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)})A^{-1}\nabla_{\vec{x}^{(n)}}E^{\rm
nd}_{a_{n}}(\vec{x}^{(n-1)},\vec{x}^{(n)})$ $\displaystyle=$
$\displaystyle\nabla_{\vec{x}^{(n)}}^{T}E^{\rm
nd}_{a_{n};0}(\vec{x}^{(n-1)},\vec{x}^{(n)})A^{-1}\nabla_{\vec{x}^{(m-1)}}E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)}),$

where in the last step we used $A=A^{T}$.

It should be reminded that an evaluation at the critical point has been done
in Eqs.(295) and (296).

According to the results in Eqs. (295) and (296), rather than partitioning the
$E^{\rm nd}_{a;0}+O(t)$ into ordering pairs, we can identify the partitions
$p_{1}$ and $p_{2}$ if $p_{1}$ can be the same as $p_{2}$ by reordering each
of its pairs. The set with this identification will be denoted by
$\tilde{\mathcal{P}}$. Then we finally have

$\displaystyle\left[\sum_{m,n=1}^{|\mathcal{M}|-1}\nabla_{\vec{x}^{(m)}}^{T}\left(\hat{H}^{-1}_{|\mathcal{M}|-1;\vec{p}_{o},\vec{\theta}_{o}}\right)_{mn}\nabla_{\vec{x}^{(n)}}\right]^{s_{o}}\prod_{i=1}^{|\mathcal{M}|-1}E^{(i)}\Bigg{|}_{g^{(l)}=e^{iw\tau_{3}},\forall
l}$ (297) $\displaystyle=$ $\displaystyle
2^{s_{o}}\sum_{p\in\tilde{\mathcal{P}}}\prod_{\begin{subarray}{c}(m,n)\in p\\\
\text{with }m<n\end{subarray}}\nabla_{\vec{x}^{(m)}}^{T}E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)})A^{-1}\nabla_{\vec{x}^{(n-1)}}E^{\rm
nd}_{a_{n};0}(\vec{x}^{(n-1)},\vec{x}^{(n)})\prod_{l\in\mathcal{P}^{\mathbb{C}}}E^{\rm
d}_{a_{l};0}(e^{iw\tau_{3}},e^{iw\tau_{3}})+O(t)$

where $\mathcal{P}^{\mathbb{C}}$ denotes the set of those $E^{(i)}$ of the
form $E^{\rm d}_{a;0}$. Substituting the results, we finally have

$\displaystyle\frac{\langle\psi_{z_{e}}|\mathcal{M}|\psi_{z_{e}}\rangle}{\langle\psi_{z_{e}}|\psi_{z_{e}}\rangle}=$
$\displaystyle\left(\left(2t\right)^{s_{o}}\frac{1}{s_{o}!}\sum_{p\in\tilde{\mathcal{P}}}\prod_{\begin{subarray}{c}(m,n)\in
p\\\ \text{with }m<n\end{subarray}}\nabla_{\vec{x}^{(m)}}^{T}E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)})A^{-1}\nabla_{\vec{x}^{(n-1)}}E^{\rm
nd}_{a_{n};0}(\vec{x}^{(n-1)},\vec{x}^{(n)})\right)\times$ (298)
$\displaystyle\prod_{l\in\mathcal{P}^{\mathbb{C}}}E^{\rm
d}_{a_{l};0}(e^{iw\tau_{3}},e^{iw\tau_{3}})+O(t^{s_{o}+1})$

Note that in

$\prod_{\begin{subarray}{c}(m,n)\in p\\\ \text{with
}m<n\end{subarray}}\nabla_{\vec{x}^{(m)}}^{T}E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)})A^{-1}\nabla_{\vec{x}^{(n-1)}}E^{\rm
nd}_{a_{n};0}(\vec{x}^{(n-1)},\vec{x}^{(n)}),$

the first argument in the most left factor and the second argument in the most
right factor cannot be acted by the derivative operator. Therefore, by the
similar discussion as above, one verifies that

$\displaystyle\frac{\langle\psi_{z_{e}}|\mathcal{M}^{\prime}|\psi_{z_{e}}\rangle}{\langle\psi_{z_{e}}|\psi_{z_{e}}\rangle}=$
$\displaystyle\left(\left(2t\right)^{s_{o}}\frac{1}{s_{o}!}\sum_{p\in\tilde{\mathcal{P}}}\prod_{\begin{subarray}{c}(m,n)\in
p\\\ \text{with }m<n\end{subarray}}\nabla_{\vec{x}^{(m)}}^{T}E^{\rm
nd}_{a_{m};0}(\vec{x}^{(m-1)},\vec{x}^{(m)})A^{-1}\nabla_{\vec{x}^{(n-1)}}E^{\rm
nd}_{a_{n};0}(\vec{x}^{(n-1)},\vec{x}^{(n)})\right)+O(t^{s_{o}+1}).$ (299)

Moreover, $E^{\rm d}_{a_{l};(0)}(e^{i\omega\tau_{3}},e^{i\omega\tau_{3}})$
represents the leading order of the expectation value of its corresponding
operator. Since these $E^{\rm
d}_{a_{l};(0)}(e^{i\omega\tau_{3}},e^{i\omega\tau_{3}})$ do not vanish, we
conclude that $\langle\mathcal{M}\rangle_{z_{e}}$ is of order $t^{s_{o}}$ if
and only if $\langle\mathcal{M}^{\prime}\rangle_{z_{e}}$ is. Thus Eq. (65) is
true for the case when both sides are of order $t^{s_{o}}$. Moreover, both
sides of Eq. (65) are of $O(t^{M_{+}+N_{+}})$ or higher, and if the leading
order terms of both sides are not $O(t^{M_{+}+N_{+}})$, then they both vanish
at $O(t^{M_{+}+N_{+}})$.

∎

###  F.4 derivative of the matrix element of $p^{i}(e)$ and
$D^{\frac{1}{2}}_{ab}(h_{e})$

Note that the results show below ignore the terms of order $t$ and higher. We
denote
$\nabla_{\vec{x}^{(i)}}:=\left(\nabla_{\vec{p}^{(i)}},\nabla_{\vec{\theta}^{(i)}}\right)$
and $x^{*}$ is the complex conjugate of $x$.

$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{1}_{s}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(-\frac{1}{2},\frac{1}{2}i\tanh\left(\frac{p}{2}\right),0,-\frac{ip\sin(\theta)\text{csch}(p)}{2\theta},\frac{ip\sin^{2}\left(\frac{\theta}{2}\right)\text{csch}(p)}{\theta},0\right)^{T}$
(300) $\displaystyle=$
$\displaystyle\left(\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{1}_{s}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}\right)^{*},$
$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{2}_{s}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(-\frac{1}{2}i\tanh\left(\frac{p}{2}\right),-\frac{1}{2},0,-\frac{ip\sin^{2}\left(\frac{\theta}{2}\right)\text{csch}(p)}{\theta},-\frac{ip\sin(\theta)\text{csch}(p)}{2\theta},0\right)^{T}$
$\displaystyle=$
$\displaystyle\left(\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{2}_{s}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}\right)^{*},$
$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{3}_{s}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(0,0,-\frac{1}{2},0,0,-\frac{i}{2}\right)^{T}=\left(\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{3}_{s}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}\right)^{*}.$

$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{1}_{t}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\Bigg{(}\frac{e^{i\theta}+e^{p-i\theta}}{2e^{p}+2},\frac{ie^{-i\theta}\left(e^{p}-e^{2i\theta}\right)}{2\left(e^{p}+1\right)},0,\frac{e^{-i\theta}\left(-1+e^{i\theta}\right)p\left(e^{2p}+e^{i\theta}\right)}{2\theta\left(e^{2p}-1\right)},$
(301)
$\displaystyle\frac{ie^{-i\theta}\left(-1+e^{i\theta}\right)p\left(e^{2p}-e^{i\theta}\right)}{2\theta\left(e^{2p}-1\right)},0\Bigg{)}^{T}=\left(\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{1}_{t}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}\right)^{*},$
$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{2}_{t}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\Bigg{(}\frac{ie^{-i\theta}\left(-e^{p}+e^{2i\theta}\right)}{2\left(e^{p}+1\right)},\frac{e^{i\theta}+e^{p-i\theta}}{2e^{p}+2},0,\frac{ie^{-i\theta}\left(-1+e^{i\theta}\right)p\left(-e^{2p}+e^{i\theta}\right)}{2\theta\left(e^{2p}-1\right)},$
$\displaystyle\frac{e^{-i\theta}\left(-1+e^{i\theta}\right)p\left(e^{2p}+e^{i\theta}\right)}{2\theta\left(e^{2p}-1\right)},0\Bigg{)}^{T}=\left(\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{2}_{t}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}\right)^{*},$
$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{3}_{t}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(0,0,\frac{1}{2},0,0,\frac{i}{2}\right)^{T}=\left(\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|\hat{p}^{3}_{t}(e)|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}\right)^{*}.$
$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|D^{\frac{1}{2}}_{\frac{1}{2}\frac{1}{2}}(h_{e})|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(0,0,-\frac{1}{4}e^{-\frac{i\theta}{2}},0,0,-\frac{1}{4}ie^{-\frac{i\theta}{2}}\right)^{T},$
(302)
$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|D^{\frac{1}{2}}_{\frac{1}{2}-\frac{1}{2}}(h_{e})|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(-\frac{e^{\frac{i\theta}{2}}\tanh\left(\frac{p}{2}\right)}{2p},\frac{ie^{\frac{i\theta}{2}}\tanh\left(\frac{p}{2}\right)}{2p},0,-\frac{i\sin\left(\frac{\theta}{2}\right)}{\theta+\theta
e^{p}},-\frac{\sin\left(\frac{\theta}{2}\right)}{\theta+\theta
e^{p}},0\right)^{T},$
$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|D^{\frac{1}{2}}_{-\frac{1}{2}\frac{1}{2}}(h_{e})|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(-\frac{e^{-\frac{i\theta}{2}}\tanh\left(\frac{p}{2}\right)}{2p},-\frac{ie^{-\frac{i\theta}{2}}\tanh\left(\frac{p}{2}\right)}{2p},0,-\frac{ie^{p}\sin\left(\frac{\theta}{2}\right)}{\theta+\theta
e^{p}},\frac{e^{p}\sin\left(\frac{\theta}{2}\right)}{\theta+\theta
e^{p}},0\right)^{T},$
$\displaystyle\nabla_{\vec{x}^{(1)}}\frac{\langle\psi_{g^{(1)}}|D^{\frac{1}{2}}_{\frac{1}{2}\frac{1}{2}}(h_{e})|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(0,0,\frac{1}{4}e^{\frac{i\theta}{2}},0,0,\frac{1}{4}ie^{\frac{i\theta}{2}}\right)^{T},$
$\displaystyle\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|D^{\frac{1}{2}}_{\frac{1}{2}\frac{1}{2}}(h_{e})|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(0,0,\frac{1}{4}e^{-\frac{i\theta}{2}},0,0,-\frac{1}{4}ie^{-\frac{i\theta}{2}}\right)^{T},$
$\displaystyle\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|D^{\frac{1}{2}}_{\frac{1}{2}-\frac{1}{2}}(h_{e})|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(\frac{e^{\frac{i\theta}{2}}\tanh\left(\frac{p}{2}\right)}{2p},-\frac{ie^{\frac{i\theta}{2}}\tanh\left(\frac{p}{2}\right)}{2p},0,-\frac{ie^{p}\sin\left(\frac{\theta}{2}\right)}{\theta+\theta
e^{p}},-\frac{e^{p}\sin\left(\frac{\theta}{2}\right)}{\theta+\theta
e^{p}},0\right)^{T},$
$\displaystyle\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|D^{\frac{1}{2}}_{-\frac{1}{2}\frac{1}{2}}(h_{e})|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(\frac{e^{-\frac{i\theta}{2}}\tanh\left(\frac{p}{2}\right)}{2p},\frac{ie^{-\frac{i\theta}{2}}\tanh\left(\frac{p}{2}\right)}{2p},0,-\frac{i\sin\left(\frac{\theta}{2}\right)}{\theta+\theta
e^{p}},\frac{\sin\left(\frac{\theta}{2}\right)}{\theta+\theta
e^{p}},0\right)^{T},$
$\displaystyle\nabla_{\vec{x}^{(2)}}\frac{\langle\psi_{g^{(1)}}|D^{\frac{1}{2}}_{-\frac{1}{2}-\frac{1}{2}}(h_{e})|\psi_{g^{(2)}}\rangle}{\langle\psi_{g^{(1)}}|\psi_{g^{(2)}}\rangle}=\left(0,0,-\frac{1}{4}e^{\frac{i\theta}{2}},0,0,\frac{1}{4}ie^{\frac{i\theta}{2}}\right)^{T}.$

## References

  * (1) T. Thiemann, _Modern canonical quantum general relativity_ , Cambridge University Press (2007).
  * (2) C. Rovelli and F. Vidotto, _Covariant Loop Quantum Gravity: An Elementary Introduction to Quantum Gravity and Spinfoam Theory_ , Cambridge Monographs on Mathematical Physics, Cambridge University Press (2014).
  * (3) M. Han, Y. Ma and W. Huang, _Fundamental structure of loop quantum gravity_ , _International Journal of Modern Physics D_ 16 (2007) 1397.
  * (4) A. Ashtekar and J. Lewandowski, _Background independent quantum gravity: a status report_ , _Classical and Quantum Gravity_ 21 (2004) R53.
  * (5) I. Mäkinen, _Dynamics in canonical models of loop quantum gravity_ , 2020\.
  * (6) K. Giesel and T. Thiemann, _Algebraic Quantum Gravity (AQG). I. Conceptual Setup_ , _Class. Quant. Grav._ 24 (2007) 2465 [gr-qc/0607099].
  * (7) K. Giesel and T. Thiemann, _Algebraic Quantum Gravity (AQG). II. Semiclassical Analysis_ , _Class. Quant. Grav._ 24 (2007) 2499 [gr-qc/0607100].
  * (8) K. Giesel and T. Thiemann, _Algebraic quantum gravity (AQG). III. Semiclassical perturbation theory_ , _Class. Quant. Grav._ 24 (2007) 2565 [gr-qc/0607101].
  * (9) M. Domagala, K. Giesel, W. Kaminski and J. Lewandowski, _Gravity quantized: Loop Quantum Gravity with a Scalar Field_ , _Phys. Rev. D_ 82 (2010) 104038 [1009.2445].
  * (10) J. Lewandowski, M. Domagala and M. Dziendzikowski, _The dynamics of the massless scalar field coupled to lqg in the polymer quantization_ , _Proc. Sci. QGQGS_ 25 (2011) .
  * (11) E. Alesci, M. Assanioussi, J. Lewandowski and I. Mäkinen, _Hamiltonian operator for loop quantum gravity coupled to a scalar field_ , _Phys. Rev. D_ 91 (2015) 124067.
  * (12) M. Assanioussi, J. Lewandowski and I. Mäkinen, _New scalar constraint operator for loop quantum gravity_ , _Phys. Rev. D_ 92 (2015) 044042.
  * (13) J. Yang and Y. Ma, _New hamiltonian constraint operator for loop quantum gravity_ , _Physics Letters B_ 751 (2015) 343.
  * (14) E. Alesci, G. Botta and G.V. Stagno, _Quantum reduced loop gravity effective hamiltonians from a statistical regularization scheme_ , _Phys. Rev. D_ 97 (2018) 046011.
  * (15) A. Dapor and K. Liegener, _Cosmological Effective Hamiltonian from full Loop Quantum Gravity Dynamics_ , _Phys. Lett. B_ 785 (2018) 506 [1706.09833].
  * (16) M. Han and H. Liu, _Effective dynamics from coherent state path integral of full loop quantum gravity_ , _Physical Review D_ 101 (2020) 046003.
  * (17) T. Thiemann, _Gauge field theory coherent states (gcs): I. general properties_ , _Classical and Quantum Gravity_ 18 (2001) 2025\.
  * (18) T. Thiemann and O. Winkler, _Gauge field theory coherent states (gcs): Ii. peakedness properties_ , _Classical and Quantum Gravity_ 18 (2001) 2561.
  * (19) T. Thiemann and O. Winkler, _Gauge field theory coherent states (gcs): Iii. ehrenfest theorems_ , _Classical and Quantum Gravity_ 18 (2001) 4629.
  * (20) A. Dapor and K. Liegener, _Cosmological coherent state expectation values in loop quantum gravity I. Isotropic kinematics_ , _Class. Quant. Grav._ 35 (2018) 135011 [1710.04015].
  * (21) K. Liegener and E.-A. Zwicknagel, _Expectation values of coherent states for $SU(2)$ Lattice Gauge Theories_, _JHEP_ 02 (2020) 024 [2001.00032].
  * (22) M. Han and H. Liu, _Semiclassical limit of new path integral formulation from reduced phase space loop quantum gravity_ , _Phys. Rev. D_ 102 (2020) 024083 [2005.00988].
  * (23) F. Conrady and L. Freidel, _On the semiclassical limit of 4d spin foam models_ , _Phys. Rev. D_ 78 (2008) 104023 [0809.2280].
  * (24) J.W. Barrett, R. Dowdall, W.J. Fairbairn, F. Hellmann and R. Pereira, _Lorentzian spin foam amplitudes: Graphical calculus and asymptotics_ , _Class.Quant.Grav._ 27 (2010) 165009 [0907.2440].
  * (25) M. Han and M. Zhang, _Asymptotics of spinfoam amplitude on simplicial manifold: Lorentzian theory_ , _Class.Quant.Grav._ 30 (2013) 165012 [1109.0499].
  * (26) E. Bianchi, L. Modesto, C. Rovelli and S. Speziale, _Graviton propagator in loop quantum gravity_ , _Class. Quant. Grav._ 23 (2006) 6989 [gr-qc/0604044].
  * (27) M. Han, Z. Huang and A. Zipfel, _Emergent four-dimensional linearized gravity from a spin foam model_ , _Phys. Rev. D_ 100 (2019) 024060 [1812.02110].
  * (28) K. Liegener and L. Rudnicki, _Algorithmic approach to Cosmological Coherent State Expectation Values in LQG_ , 12, 2020 [2012.07813].
  * (29) M. Assanioussi, J. Lewandowski and I. Mäkinen, _Time evolution in deparametrized models of loop quantum gravity_ , _Phys. Rev. D_ 96 (2017) 024043 [1702.01688].
  * (30) E. Alesci, K. Liegener and A. Zipfel, _Matrix elements of Lorentzian Hamiltonian constraint in loop quantum gravity_ , _Phys. Rev. D_ 88 (2013) 084043 [1306.0861].
  * (31) R. Borissov, R. De Pietri and C. Rovelli, _Matrix elements of Thiemann’s Hamiltonian constraint in loop quantum gravity_ , _Class. Quant. Grav._ 14 (1997) 2793 [gr-qc/9703090].
  * (32) M. Han, Z. Huang, H. Liu and D. Qu, _Numerical computations of next-to-leading order corrections in spinfoam large- $j$ asymptotics_, _Phys. Rev. D_ 102 (2020) 124010 [2007.01998].
  * (33) P. Donà, M. Fanizza, G. Sarno and S. Speziale, _Numerical study of the Lorentzian Engle-Pereira-Rovelli-Livine spin foam amplitude_ , _Phys. Rev. D_ 100 (2019) 106003 [1903.12624].
  * (34) H. Nicolai, K. Peeters and M. Zamaklar, _Loop quantum gravity: An Outside view_ , _Class. Quant. Grav._ 22 (2005) R193 [hep-th/0501114].
  * (35) J. Yang, Y. Ding and Y. Ma, _Alternative quantization of the hamiltonian in loop quantum cosmology_ , _Physics Letters B_ 682 (2009) 1\.
  * (36) W. Kamiński and K. Liegener, _Symmetry restriction and its application to gravity_ , 2009.06311.
  * (37) K. Giesel and T. Thiemann, _Algebraic quantum gravity (aqg): Iv. reduced phase space quantization of loop quantum gravity_ , _Classical and Quantum Gravity_ 27 (2010) 175009.
  * (38) K. Giesel and T. Thiemann, _Scalar Material Reference Systems and Loop Quantum Gravity_ , _Class. Quant. Grav._ 32 (2015) 135015 [1206.3807].
  * (39) K.V. Kuchar and C.G. Torre, _Gaussian reference fluid and interpretation of quantum geometrodynamics_ , _Phys. Rev._ D43 (1991) 419.
  * (40) T. Thiemann, _Canonical Quantum Gravity, Constructive QFT and Renormalisation_ , 2003.13622.
  * (41) T. Thiemann, _Quantum spin dynamics (qsd)_ , _Classical and Quantum Gravity_ 15 (1998) 839.
  * (42) J. Lewandowski and H. Sahlmann, _Symmetric scalar constraint for loop quantum gravity_ , _Phys. Rev. D_ 91 (2015) 044022 [1410.5276].
  * (43) C. Zhang, J. Lewandowski, H. Li and Y. Ma, _Bouncing evolution in a model of loop quantum gravity_ , _Physical Review D_ 99 (2019) 124012\.
  * (44) M. Han and H. Liu, _Improved $\overline{\mu}$-scheme effective dynamics of full loop quantum gravity_, _Phys. Rev. D_ 102 (2020) 064061 [1912.08668].
  * (45) A. Dapor and K. Liegener, _Modifications to Gravitational Wave Equation from Canonical Quantum Gravity_ , _Eur. Phys. J. C_ 80 (2020) 741 [2002.00834].
  * (46) M. Han, H. Li and H. Liu, _Manifestly gauge-invariant cosmological perturbation theory from full loop quantum gravity_ , _Phys. Rev. D_ 102 (2020) 124002 [2005.00883].
  * (47) M. Peskin and D. Schroeder, _An Introduction To Quantum Field Theory_ , Frontiers in Physics, Avalon Publishing (1995).
  * (48) C. Zhang, S. Song and M. Han, _First-order quantum correction in coherent state expectation value of loop-quantum-gravity hamiltonian: I. overview and results_ , _arXiv preprint arXiv:2012.14242_ (2020) .
  * (49) https://github.com/czhangUW/LQG_expectation_values.
  * (50) J.F. Barbero G., _Real ashtekar variables for lorentzian signature space-times_ , _Phys. Rev. D_ 51 (1995) 5507.
  * (51) T. Thiemann, _Quantum spin dynamics (QSD): 7. Symplectic structures and continuum lattice formulations of gauge field theories_ , _Class. Quant. Grav._ 18 (2001) 3293 [hep-th/0005232].
  * (52) C. Rovelli and L. Smolin, _Discreteness of area and volume in quantum gravity_ , _Nuclear Physics B_ 442 (1995) 593.
  * (53) A. Ashtekar and J. Lewandowski, _Quantum theory of geometry: I. area operators_ , _Classical and Quantum Gravity_ 14 (1997) A55.
  * (54) A. Ashtekar and J. Lewandowski, _Quantum theory of geometry ii: Volume operators_ , _Advances in Theoretical and Mathematical Physics_ 1 (1997) 388.
  * (55) _Private communication with Thomas Thiemann_ , .
  * (56) K. Liegener and E.-A. Zwicknagel, _Expectation values of coherent states for su (2) lattice gauge theories_ , _Journal of High Energy Physics_ 2020 (2020) 1.
  * (57) E. Alesci and F. Cianfrani, _Quantum reduced loop gravity: Semiclassical limit_ , _Physical Review D_ 90 (2014) 024006.
  * (58) M. Han and H. Liu, _Loop Quantum Gravity on Dynamical Lattice and Improved Cosmological Effective Dynamics with Inflaton_ , 2101.07659.
  * (59) T. Thiemann, _QSD 5: Quantum gravity as the natural regulator of matter quantum field theories_ , _Class. Quant. Grav._ 15 (1998) 1281 [gr-qc/9705019].
  * (60) H. Sahlmann and T. Thiemann, _Towards the QFT on curved space-time limit of QGR. 1. A General scheme_ , _Class. Quant. Grav._ 23 (2006) 867 [gr-qc/0207030].
  * (61) M. Han, C.-H. Hsiao, S. Song and C. Zhang, _Matter Couplings in Path Integral Formulation of Reduced Phase Space Loop Quantum Gravity_ , in preparation.
  * (62) Y.A.G.A.S. Zhedanov, _Nature of the symmetry group of the 6 $j$-symbol_, _Zh. Eksper. Teoret. Fiz_ 94 (1988) 49.
  * (63) J. Yang and Y. Ma, _Graphical calculus of volume, inverse volume and hamiltonian operators in loop quantum gravity_ , _The European Physical Journal C_ 77 (2017) 235.
  * (64) C. Zhang, J. Lewandowski and Y. Ma, _Towards the self-adjointness of a hamiltonian operator in loop quantum gravity_ , _Physical Review D_ 98 (2018) 086014.
  * (65) I. Mäkinen, _Introduction to su(2) recoupling theory and graphical methods for loop quantum gravity_ , 2019.
  * (66) D.A. Varshalovich, A.N. Moskalev and V.K. Khersonskii, _Quantum theory of angular momentum_ , World Scientific (1988).
  * (67) L. Hörmander, _The analysis of linear partial differential operators I: Distribution theory and Fourier analysis_ , Springer (2015).