TITLE: Definitions of absolute and essential Borel space?
QUESTION [0 upvotes]: All spaces are separable metrizable.
My first question will be about the Borel class $F_\sigma$, but I really want to know if they generalize for all Borel types ($F_\sigma$, $G_\delta$, $F_{\sigma\delta}$, $G_{\delta\sigma}$, etc.)
1. Does "$X$ is an absolute $F_{\sigma}$-space" usually mean the same thing as "$X$ is an $F_{\sigma}$-space"?  And does this mean $X$ is an $F_\sigma$-subset of every metric space in which it is embedded? If so, then $F_\sigma$-spaces would be precisely the $\sigma$-compact spaces.
2. What is the definition of "$X$ is an essential $F_{\sigma\delta}$-space"? Does this mean $X$ is an $F_{\sigma\delta}$-space but does not belong to any lower class in the hierarchy, e.g. $X$ is neither an $F_\sigma$-space nor a $G_\delta$-space?

REPLY [1 votes]: A space $X$ is an absolute Borel set of class (type) $K$ (e.g. $F_\sigma$) iff whenever we embed $X$ into a metric space $M$ via a homeomorphic embedding $i :X \to M$, $i[X]$ is of that class $K$ in $M$. This notion is due to Kuratowksi, see this more recent paper for a reference.

an absolutely closed set is compact, as is easily seen, and an absolute $F_\sigma$ is a $\sigma$-compact space.
an absolute $G_\delta$ is completely metrisable (also a classic fact).

Absolute $F_{\sigma\delta}$ and $G_{\delta\sigma}$ spaces are more complex, e.g. see the aforementioned paper e.g. , also see this paper by Stone.
As to the "essential" part, your hypothesis on its meaning (of that absolute class but no simpler one) makes sense; I could not find an authorative reference for it, though.


TITLE: How to find indicial equation
QUESTION [3 upvotes]: How can I find the indicial equation of $x(x-1)y''+3y'-2y=0$? I tried the method of Frobenius but I keep getting lost in the algebra. Is there any other way to get the indicial equation?

REPLY [13 votes]: $y''+\dfrac{3}{x(x-1)}y'-\dfrac{2}{x(x-1)}y=0$ then $p(x)=\dfrac{3}{x(x-1)}$ and $q(x)=-\dfrac{2}{x(x-1)}$. The equation has two regular singular points $x=0$ and $x=1$. For $x=0$ we see
$$p_0=\lim_{x\to0}xp(x)=\lim_{x\to0}\dfrac{3}{x-1}=-3$$
and
$$q_0=\lim_{x\to0}x^2q(x)=\lim_{x\to0}\dfrac{-2x}{x-1}=0$$
then the indicial equation is $r(r-1)+p_0r+q_0=0$ or $r^2-4r=0$ shows $r=0$ and $r=4$. In this case $4-0\in\mathbb{Z}$ so $r=4$ gives a solution and let $y=x^4(a_0+a_1x+a_2x^2+\cdots)$ to find the answer. Do like this with point  $x=1$.


text:
\section{Analysis of \fedavg under \cref{asm2}}
\label{sec:fedavg}
In this section we study the convergence of \fedavg under \cref{asm2}.
We provide a complete, non-asymptotic version of \cref{fedavg:a2} and provide the proof.
We formally define \fedavg in \cref{alg:fedavg} for reference.

Formally we use $\mathcal{F}_t$ to denote the $\sigma$-algebra generated by $\{w_\tau^m\}_{\tau \leq t, m \in [M]}$. Since \fedavg is Markovian, conditioning on $\mathcal{F}_t$ is equivalent to conditioning on $\{w_t^m\}_{m \in [M]}$.

\begin{algorithm}
  \caption{Federated Averaging (a.k.a. Local SGD, Parallel SGD)}
  \begin{algorithmic}[1]
      \label{alg:fedavg}
      \Procedure{\fedavg }{$\eta$}
      \State Initialize $= w_0^{m} = w_0$ for all $m \in [M]$
      \For{$t = 0, \ldots, T-1$}
      \For{every worker $m\in [M]$ {\bf in parallel}}
      \State $g_t^m \gets \nabla f(w_t^{m} ; \xi_t^m)$
      \Comment{Query gradient at $w_t^{m}$}
      \State $v_{t+1}^{m} \gets w_t^{m} - \eta \cdot g_t^m$
      \Comment{Compute next iterate candidate $v_{t+1}^{m}$}
      \If{sync}
      \State $w_{t+1}^m \gets \frac{1}{M} \sum_{m=1}^M v_{t+1}^m$
      \Comment{Average and broadcast}
      \Else
      \State $w_{t+1}^m \gets  v_{t+1}^m$
      \Comment{Candidates assigned to be the next iterates }
      \EndIf
      \EndFor
      \EndFor
      \EndProcedure
  \end{algorithmic}
\end{algorithm}

\subsection{Main theorem and lemma: Complete version of \cref{fedavg:a2}}
\begin{theorem}
  \label{fedavg:a2:full}
  Let $F$ be $\mu > 0$-strongly convex, and assume \cref{asm2}, then for
  \begin{equation}
    \eta := \min \left\{ \frac{1}{4L}, \frac{2}{\mu T} \log \left( \euler + \min \left\{ \frac{\mu^2 M T^2 D_0^2}{\sigma^2}, \frac{\mu^6 T^5 D_0^2}{Q^2 K^2 \sigma^4}\right\} \right)\right\},
  \end{equation}
  \fedavg yields
  \begin{align}
    &  \expt \left[ F \left( \sum_{t=0}^{T-1} \frac{\rho_t}{S_T} \overline{w_t} \right)  \right] - F^*
    + \frac{\mu}{2} \expt [\|\overline{w_T} - w^*\|^2]
    \\
    \leq &
    \exp \left( -\frac{\mu T}{8 L} \right) 4L D_0^2
    +
    \frac{3 \sigma^2}{ \mu MT } \log \left( \euler + \frac{\mu^2 M T^2 D_0^2}{\sigma^2}  \right)
    +
    \frac{3073 Q^2 K^2 \sigma^4}{\mu^5 T^4} \log^4 \left( \euler + \frac{\mu^6 T^5 D_0^2}{Q^2 K^2 \sigma^4}  \right).
  \end{align}
  where $\rho_t := (1 - \frac{1}{2} \eta \mu)^{T-t-1}$, $S_T := \sum_{t=0}^{T-1}\rho_t$, and $D_0 = \|\overline{w_0}-w^*\|$.
\end{theorem}
The proof of \cref{fedavg:a2:full} is based on the following two lemmas regarding the convergence and \nth{4}-order stability of \fedavg. The averaging technique applied here is similar to \citep{Stich-arXiv19}.
\begin{lemma}[Perturbed iterate analysis for \fedavg under \cref{asm2}]
  \label{fedavg:a2:conv:main}
  Let $F$ be $\mu > 0$-strongly convex, and assume \cref{asm2}, then for $\eta \in (0, \frac{1}{4L}]$, \fedavg satisfies
  \begin{align}
         & \expt \left[ F \left( \sum_{t=0}^{T-1} \frac{\rho_t}{S_T} \overline{w_t} \right)  \right] - F^*
    + \frac{\mu}{2} \expt [\|\overline{w_T} - w^*\|^2]
    \\
    \leq & \frac{1}{\eta} \exp \left( - \frac{1}{2} \eta \mu T \right) D_0^2
    +
    \frac{1}{M} \eta \sigma^2 + \frac{Q^2}{\mu} \left( \max_{0 \leq t < T}  \frac{1}{M} \sum_{m=1}^M \expt \left[\|\overline{w_t} - w_t^m\|^4 \right] \right).
  \end{align}
  where $\rho_t, S_T$ are defined in the statement of \cref{fedavg:a2:full}.
\end{lemma}
The proof of \cref{fedavg:a2:conv:main} is deferred to \cref{sec:fedavg:a2:conv}.
\begin{lemma}[\nth{4}-order discrepancy overhead bound for \fedavg]
  \label{fedavg:a2:stab:main}
  In the same settings of \cref{fedavg:a2:conv:main}, \fedavg satisfies (for any $t$)
  \begin{align}
    \expt \left[ \frac{1}{M} \sum_{m=1}^M \| \overline{w}_t - w_t^m \|^4 \right]
    \leq
    192 \eta^4 K^2 \sigma^4.
  \end{align}
\end{lemma}
The proof of \cref{fedavg:a2:stab:main} is deferred to \cref{sec:fedavg:a2:stab}.

Combining \cref{fedavg:a2:conv:main,fedavg:a2:stab:main} gives
\begin{lemma}[Convergence of \fedavg under \cref{asm2} for general $\eta$]
  \label{fedavg:a2:general:eta}
  In the same settings of \cref{fedavg:a2:conv:main},
  \fedavg yields
  \begin{equation}
      \expt \left[ F \left( \sum_{t=0}^{T-1} \frac{\rho_t}{S_T} \overline{w_t} \right)  \right] - F^*
 + \frac{\mu}{2} \expt [\|\overline{w_T} - w^*\|^2]
    \leq
 \frac{1}{\eta} \exp \left( - \frac{1}{2} \eta \mu T \right) D_0^2
 +
 \frac{1}{M} \eta \sigma^2 + \frac{192 \eta^4 Q^2 K^2 \sigma^4 }{\mu}.
 \label{eq:fedavg:a2:general:eta}
  \end{equation}
\end{lemma}
\begin{proof}[Proof of \cref{fedavg:a2:general:eta}]
  Immediate from \cref{fedavg:a2:conv:main,fedavg:a2:stab:main}.
\end{proof}
\cref{fedavg:a2:full} then follows by plugging an appropriate $\eta$ to \cref{fedavg:a2:general:eta}.
\begin{proof}[Proof of \cref{fedavg:a2:full}]
  To simplify the notation, denote the terms on the RHS of \cref{eq:fedavg:a2:general:eta} as
  \begin{equation}
    \varphi_{\downarrow}(\eta) := \frac{1}{\eta} \exp \left( - \frac{1}{2} \eta \mu T \right) D_0^2,
    \qquad
    \varphi_{\uparrow}(\eta) :=
    \frac{1}{M} \eta \sigma^2 + \frac{192 \eta^4 Q^2 K^2 \sigma^4 }{\mu}.
  \end{equation}
  Let
  \begin{equation}
    \eta_0 := \frac{2}{\mu T} \log \left( \euler + \min \left\{ \frac{\mu^2 M T^2 D_0^2}{\sigma^2}, \frac{\mu^6 T^5 D_0^2}{Q^2 K^2 \sigma^4}\right\} \right),
    \qquad
    \text{then }
    \eta = \min \left\{ \frac{1}{4L}, \eta_0\right\}.
  \end{equation}
  Therefore $\varphi_{\downarrow}(\eta) \leq \varphi_{\downarrow}(\frac{1}{4L}) + \varphi_{\downarrow}(\eta_0)$, where
  \begin{equation}
    \varphi_{\downarrow} \left( \frac{1}{4L} \right) = \exp \left( -\frac{\mu T}{8 L} \right) 4L D_0^2,
    \label{eq:fedavg:a2:proof:1}
  \end{equation}
  and
  \begin{equation}
    \varphi_{\downarrow}(\eta_0) \leq \frac{\mu T}{2} D_0^2 \cdot \left( \min \left\{ \frac{\mu^2 M T^2 D_0^2}{\sigma^2}, \frac{\mu^6 T^5 D_0^2}{Q^2 K^2 \sigma^4}\right\}  \right)^{-1}
    \leq
    \frac{\sigma^2}{2 \mu MT} + \frac{Q^2 K^2 \sigma^4}{2 \mu^5 T^4}.
    \label{eq:fedavg:a2:proof:2}
  \end{equation}
  On the other hand
  \begin{align}
    \varphi_{\uparrow}(\eta) \leq \varphi_{\uparrow}(\eta_0)
    \leq
    \frac{2 \sigma^2}{\mu MT} \log \left( \euler + \frac{\mu^2 M T^2 D_0^2}{\sigma^2}  \right)
    +
    \frac{3072 Q^2 K^2 \sigma^4}{\mu^5 T^4} \log^4 \left( \euler + \frac{\mu^6 T^5 D_0^2}{Q^2 K^2 \sigma^4}  \right).
    \label{eq:fedavg:a2:proof:3}
  \end{align}
  Combining \cref{fedavg:a2:general:eta,eq:fedavg:a2:proof:1,eq:fedavg:a2:proof:2,eq:fedavg:a2:proof:3} gives
  \begin{align}
    &  \expt \left[ F \left( \sum_{t=0}^{T-1} \frac{\rho_t}{S_T} \overline{w_t} \right)  \right] - F^*
    + \frac{\mu}{2} \expt [\|\overline{w_T} - w^*\|^2]
    \\
    \leq &
    \exp \left( -\frac{\mu T}{8 L} \right) 4L D_0^2
    +
    \frac{3 \sigma^2}{ \mu MT } \log \left( \euler + \frac{\mu^2 M T^2 D_0^2}{\sigma^2}  \right)
    +
    \frac{3073 Q^2 K^2 \sigma^4}{\mu^5 T^4} \log^4 \left( \euler + \frac{\mu^6 T^5 D_0^2}{Q^2 K^2 \sigma^4}  \right).
  \end{align}
\end{proof}

\subsection{Perturbed iterative analysis for \fedavg: Proof of \cref{fedavg:a2:conv:main}}
\label{sec:fedavg:a2:conv}
We first state and proof the following proposition on one-step analysis.
\begin{proposition}
  \label{fedavg:a2:conv:onestep}
  Under the same assumption of \cref{fedavg:a2:conv:main}, for all $t$, the following inequality holds
  \begin{equation}
    \expt \left[ \|\overline{w_{t+1}} - w^*\|^2  |\mathcal{F}_t \right]
    \leq
    \left( 1 - \frac{1}{2}\eta \mu  \right)  \|\overline{w_t} - w^* \|^2 -  \eta (F(\overline{w_t} ) - F^*)
    +  \frac{\eta Q^2}{\mu M}  \sum_{m=1}^M \|\overline{w_t} - w_t^m\|^4  + \frac{\eta^2\sigma^2}{M}.
  \end{equation}
\end{proposition}

\begin{proof}[Proof of \cref{fedavg:a2:conv:onestep}]
  By definition of the $\fedavg$ procedure (see \cref{alg:fedavg}), for all $m \in [M]$, $v_{t+1}^m = w_t^m - \eta \nabla f(w_t^m; \xi_t^m)$. Taking average over $m = 1,\ldots, M$ gives
  \begin{equation}
    \overline{w_{t+1}} - w^* = w_t - \eta \cdot \frac{1}{M} \sum_{m=1}^M \nabla f(w_t^m; \xi_t^m) - w^*.
  \end{equation}
  Taking conditional expectation, by bounded variance \cref{asm1}(c),
  \begin{equation}
    \expt \left[ \|\overline{w_{t+1}} - w^*\|^2  |\mathcal{F}_t \right]
    =
    \left\| w_t - \eta \cdot \frac{1}{M} \sum_{m=1}^M \nabla F(w_t^m) - w^*  \right\|^2 + \frac{1}{M} \eta^2 \sigma^2.
    \label{eq:fedavg:a2:conv:onestep:0}
  \end{equation}
  Now we analyze the $\left\| w_t - \eta \cdot \frac{1}{M} \sum_{m=1}^M \nabla F(w_t^m) - w^*  \right\|^2$ term as follows
  \begin{align}
    & \left\| w_t - \eta \cdot \frac{1}{M} \sum_{m=1}^M \nabla F(w_t^m) - w^*  \right\|^2
    \\
    =    &
    \left\| \overline{w_t} - \eta \cdot \nabla F (\overline{w_t}) - w^* + \eta \left( \nabla F (\overline{w_t}) -   \frac{1}{M} \sum_{m=1}^M \nabla F(w_t^m)  \right) \right\|^2
    \\
    \leq & \left( 1 + \frac{1}{2} \eta \mu  \right) \left\| \overline{w_t} - \eta \nabla F (\overline{w_t}) - w^* \right\|^2
    + \eta^2 \left( 1 + \frac{2}{\eta \mu} \right) \left\|  \nabla F (\overline{w_t}) -   \frac{1}{M} \sum_{m=1}^M \nabla F(w_t^m)   \right\|^2
    \tag{apply \cref{helper:unbalanced:ineq} with $\zeta = \frac{1}{2}\eta \mu$}
    \\
    \leq & \left( 1 + \frac{1}{2} \eta \mu  \right) \left\| \overline{w_t} - \eta \nabla F (\overline{w_t}) - w^* \right\|^2
    + \eta^2 \left( 1 + \frac{2}{\eta \mu} \right) \frac{Q^2}{4M} \sum_{m=1}^M \|\overline{w_t} - w_t^m\|^4
    \tag{by \cref{helper:3rd:Lip}}
    \\
    \leq & \left( 1 + \frac{1}{2} \eta \mu  \right) \left\| \overline{w_t} - \eta \nabla F (\overline{w_t}) - w^* \right\|^2
    + \frac{\eta Q^2}{\mu M} \sum_{m=1}^M \|\overline{w_t} - w_t^m\|^4.
    \label{eq:fedavg:a2:conv:onestep:1}
  \end{align}
  where the last inequality is due to $1 + \frac{2}{\eta \mu} \leq \frac{4}{\eta\mu}$ since $\eta \mu \leq \eta L \leq \frac{1}{4}$.

  The first term of the RHS of \cref{eq:fedavg:a2:conv:onestep:1} is bounded as
  \begin{align}
         & \left\| \overline{w_t} - \eta \nabla F (\overline{w_t}) - w^* \right\|^2
    \\
    =    & \left\|  \overline{w_t} - w^* \right\|^2 - 2 \eta \left\langle \nabla F (\overline{w_t}),   \overline{w_t} - w^* \right\rangle + \eta^2 \|\nabla F(\overline{w_t})\|^2
    \tag{expansion of squared norm}
    \\
    \leq & \left\|  \overline{w_t} - w^* \right\|^2 - \eta \left( \mu \|\overline{w_t} - w^* \|^2 - 2 (F(\overline{w_t} ) - F^*) \right) + \eta^2 \cdot (2 L (F(\overline{w_t} ) - F^*))
    \tag{$\mu$-strongly convexity and $L$-smoothness by \cref{asm1}}
    \\
    =    &
    (1 - \eta \mu )  \|\overline{w_t} - w^* \|^2 - 2\eta (1 - \eta L) (F(\overline{w_t} ) - F^*)
    \\
    \leq & (1 - \eta \mu )  \|\overline{w_t} - w^* \|^2 - \eta (F(\overline{w_t} ) - F^*).
    \tag{since $\eta \leq \frac{1}{2L}$}
  \end{align}
  Multiplying $(1 + \frac{1}{2} \eta \mu)$ on both sides gives (note that $(1 + \frac{1}{2} \eta \mu)(1 - \eta \mu) \leq (1 - \frac{1}{2} \eta \mu)$)
  \begin{align}
    & \left( 1 + \frac{1}{2} \eta \mu  \right) \left\| \overline{w_t} - \eta_t \nabla F (\overline{w_t}) - w^* \right\|^2
    \\
    \leq &
    \left( 1 + \frac{1}{2} \eta \mu  \right)\left( 1 - \eta \mu \right)  \|\overline{w_t} - w^* \|^2 -  \eta \left( 1 + \frac{1}{2} \eta \mu  \right) \left( F (\overline{w_t}) - F^* \right)
    \\
    \leq & \left( 1 - \frac{1}{2} \eta \mu \right) \|\overline{w_t} - w^* \|^2 -  \eta  \left( F (\overline{w_t}) - F^* \right).
    \label{eq:fedavg:a2:conv:onestep:2}
  \end{align}
  Combining \cref{eq:fedavg:a2:conv:onestep:0,eq:fedavg:a2:conv:onestep:1,eq:fedavg:a2:conv:onestep:2}
  completes the proof of \cref{fedavg:a2:conv:onestep}.
\end{proof}
With \cref{fedavg:a2:conv:onestep} at hand we are ready to prove \cref{fedavg:a2:conv:main}. The telescoping techniques applied here are similar to \citep{Stich-arXiv19}.
\begin{proof}[Proof of \cref{fedavg:a2:conv:main}]
  Telescoping \cref{fedavg:a2:conv:onestep} yields
  \begin{align}
         & \expt \left[\|\overline{w_T}-w^*\|^2 \right]
    + \eta \sum_{t=0}^{T-1} \left( 1 - \frac{1}{2} \eta \mu \right)^{T-t-1}\left( \expt[F (\overline{w_t})] - F^* \right)
    \\
    \leq &
    \left( 1 - \frac{1}{2} \eta \mu \right)^T \| \overline{w_0} - w^*\|^2
    +\sum_{t=0}^{T-1} \left(1 - \frac{1}{2} \eta \mu \right)^{T-t-1}
    \left(  \frac{1}{M} \eta^2 \sigma^2 + \frac{\eta Q^2}{\mu M} \sum_{m=1}^M \expt \left[ \|\overline{w_t} - w_t^m\|^4 \right]  \right)
    \\
    \leq &
    \left( 1 - \frac{1}{2} \eta \mu \right)^T \| \overline{w_0} - w^*\|^2
    + S_T
    \left(  \frac{1}{M} \eta^2 \sigma^2 + \frac{\eta Q^2}{\mu} \max_{0 \leq t < T} \frac{1}{M} \sum_{m=1}^M \expt \left[\|\overline{w_t} - w_t^m\|^4 \right] \right).
  \end{align}
  Multiplying $\frac{1}{\eta S_T}$ on both sides and rearranging,
  \begin{align}
         & \sum_{t=0}^{T-1} \frac{\rho_t}{S_T} (\expt[F(\overline{w_t})] - F^*)
    + \frac{1}{\eta S_T} \expt [\|\overline{w_T} - w^*\|^2]
    \\
    \leq &
    \frac{(1 - \frac{1}{2}\eta \mu)^T}{\eta S_T} \|\overline{w_0} - w^*\|^2
    +
    \frac{1}{M} \eta \sigma^2 + \frac{Q^2}{\mu} \left( \max_{0 \leq t < T}  \frac{1}{M} \sum_{m=1}^M \expt \left[\|\overline{w_t} - w_t^m\|^4 \right] \right).
    \label{eq:fedavg:a2:conv:main:1}
  \end{align}
  Note that $S_T := \sum_{t=0}^{T-1} \rho_t = \frac{1 - (1 - \frac{1}{2} \eta \mu)^T}{\frac{1}{2} \eta \mu}$, we have
  \begin{equation}
    \frac{1}{\eta S_T} = \frac{\mu}{2\left( 1 - (1 - \frac{1}{2} \eta \mu)^T \right)} \geq \frac{\mu}{2},
    \label{eq:fedavg:a2:conv:main:2}
  \end{equation}
  and
  \begin{equation}
    \frac{(1 - \frac{1}{2}\eta \mu)^T}{\eta S_T} = \frac{\mu (1 - \frac{1}{2} \eta \mu)^T}{2 \left( 1 - (1 - \frac{1}{2} \eta \mu)^{T} \right)}
    \leq \frac{\mu(1 - \frac{1}{2} \eta \mu)^T}{\eta \mu}
    \leq \frac{1}{\eta} \exp \left( - \frac{1}{2} \eta \mu T \right).
    \label{eq:fedavg:a2:conv:main:3}
  \end{equation}
  Also by convexity
  \begin{equation}
    \sum_{t=0}^{T-1} \frac{\rho_t}{S_T} (\expt[F(\overline{w_t})] - F^*) \geq
    \expt \left[ F \left( \sum_{t=0}^{T-1} \frac{\rho_t}{S_T} \overline{w_t} \right) \right] - F^*.
    \label{eq:fedavg:a2:conv:main:4}
  \end{equation}
  Plugging \cref{eq:fedavg:a2:conv:main:2,eq:fedavg:a2:conv:main:3,eq:fedavg:a2:conv:main:4} to \cref{eq:fedavg:a2:conv:main:1} gives
  \begin{align}
         & \expt \left[ F \left( \sum_{t=0}^{T-1} \frac{\rho_t}{S_T} \overline{w_t} \right)  \right] - F^*
    + \frac{\mu}{2} \expt [\|\overline{w_T} - w^*\|^2]
    \\
    \leq & \frac{1}{\eta} \exp \left( - \frac{1}{2} \eta \mu T \right)\|\overline{w_0} - w^*\|^2
    +
    \frac{1}{M} \eta \sigma^2 + \frac{Q^2}{\mu} \left( \max_{0 \leq t < T}  \frac{1}{M} \sum_{m=1}^M \expt \left[\|\overline{w_t} - w_t^m\|^4 \right] \right).
  \end{align}
\end{proof}

\subsection{Discrepancy overhead bound for \fedavg: Proof of \cref{fedavg:a2:stab:main}}
\label{sec:fedavg:a2:stab}
In this subsection we will prove \cref{fedavg:a2:stab:main} regarding the 4th order stability of \fedavg.
We introduce a few more notations to simplify the discussions.
Let $m_1, m_2 \in [M]$ be two arbitrary distinct workers. For any timestep $t$, let $\Delta_t := w_t^{m_1} - w_t^{m_2}$, and $\Delta_t^\varepsilon := \varepsilon_t^{m_1} - \varepsilon_t^{m_2}$ where $\varepsilon_t^{m} = \nabla f(w_t^m; \xi_t^m) - \nabla F(w_t^m)$ be the bias of the gradient oracle of the $m$-th worker evaluated at $w_t$. Let $\Delta_t^\nabla := \nabla F(w_t^{m_1}) - \nabla F(w_t^{m_2})$.

We first state and prove the following proposition on one-step 4th-order stability. The proof is analogous to the \nth{4}-order convergence analysis of \fedavg in \citep{Dieuleveut.Patel-NeurIPS19}.
\begin{proposition}
  \label{fedavg:a2:stab:1}
  In the same setting of \cref{fedavg:a2:stab:main}, for all $t$,
  \begin{equation}
    \sqrt{\expt \|\Delta_{t+1}\|^4}
    \leq
    \sqrt{\expt \|\Delta_t\|^4} + \sqrt{192} \eta^2 \sigma^2.
  \end{equation}
\end{proposition}

\begin{proof}[Proof of \cref{fedavg:a2:stab:1}]
  If $t+1$ is a synchronized step, then the result follows trivially. We assume from now on that $t+1$ is not a synchronized step, then
  \begin{align}
         & \expt [\|\Delta_{t+1}\|^4 |\mathcal{F}_t] = \expt \left[ \|\Delta_t - \eta (\Delta_t^{\nabla} + \Delta_t^\varepsilon) \|^4 |\mathcal{F}_t \right]
    \\
    =    & \expt \left[ \left( \|\Delta_t\|^2 - 2 \eta \langle \Delta_t, \Delta_t^\nabla + \Delta_t^\varepsilon \rangle + \eta^2 \|\Delta_t^\nabla + \Delta_t^\varepsilon\|^2 \right)^2 \middle| \mathcal{F}_t \right]
    \\
    =    & \expt \|\Delta_t\|^4 - 4\eta \|\Delta_t\|^2 \langle \Delta_t, \Delta_t^\nabla \rangle
    + 4 \eta^2 \expt \left[ \langle \Delta_t, \Delta_t^\nabla + \Delta_t^\varepsilon \rangle^2 |\mathcal{F}_t \right]
    + 2 \eta^2 \|\Delta_t\|^2  \expt \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^2 |\mathcal{F}_t \right]
    \\
         & \quad - 4 \eta^3 \expt \left[ \langle \Delta_t, \Delta_t^\nabla + \Delta_t^\varepsilon \rangle  \cdot \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^2 |\mathcal{F}_t  \right]
    + \eta^4 \expt \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^4 |\mathcal{F}_t \right]
    \\
    \leq & \expt \|\Delta_t\|^4 - 4\eta \|\Delta_t\|^2 \langle \Delta_t, \Delta_t^\nabla \rangle + 6 \eta^2 \|\Delta_t\|^2  \expt \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^2 |\mathcal{F}_t \right]
    \\
         & \quad + 4 \eta^3 \|\Delta_t\| \expt  \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^3 |\mathcal{F}_t \right]
    + \eta^4 \expt \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^4 |\mathcal{F}_t \right] \tag{Cauchy-Schwarz inequality}
    \\
    \leq & \expt \|\Delta_t\|^4 - 4\eta \|\Delta_t\|^2 \langle \Delta_t, \Delta_t^\nabla \rangle
    + 8 \eta^2  \|\Delta_t\|^2  \expt \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^2 |\mathcal{F}_t \right]
    + 3 \eta^4 \expt \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^4 |\mathcal{F}_t \right],
    \label{eq:fedavg:a2:stab:prop:1}
  \end{align}
  where the last inequality is due to
  \begin{equation}
    4 \eta^3 \|\Delta_t\| \expt  \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^3 |\mathcal{F}_t \right]
    \leq 2 \eta^2 \|\Delta_t\|^2  \expt  \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^2 |\mathcal{F}_t \right] +
    2 \eta^4 \expt  \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^4 |\mathcal{F}_t \right]
  \end{equation}
  by AM-GM inequality.

  Note that by $L$-smoothness and convexity, we have the following inequality by standard convex analysis (\cf, Theorem 2.1.5 of \citep{Nesterov-18}),
  \begin{align}
    \|\Delta_t^\nabla\|^2
    =
    \| \nabla F(w_t^{m_1}) - \nabla F(w_t^{m_2})\|^2
    \leq
    L \left\langle w_t^{m_1} - w_t^{m_2},  \nabla F(w_t^{m_1}) - \nabla F(w_t^{m_2})  \right\rangle
    =
    L \langle \Delta_t, \Delta_t^\nabla \rangle.
    \label{eq:fedavg:a2:stab:prop:2}
  \end{align}
  Consequently
  \begin{equation}
    \expt \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^2 |\mathcal{F}_t \right] = \| \Delta_t^\nabla\|^2 +  \expt \left[ \|\Delta_t^\varepsilon\|^2 |\mathcal{F}_t \right]
    \leq \| \Delta_t^\nabla\|^2 + 2\sigma^2
    \leq L \langle \Delta_t, \Delta_t^\nabla \rangle + 2\sigma^2.
  \end{equation}
  Similarly
  \begin{align}
         & \expt \left[ \| \Delta_t^\nabla + \Delta_t^\varepsilon \|^4 |\mathcal{F}_t \right]
    \leq
    8 \|\Delta_t^\nabla\|^4 + 8 \expt \left[ \|\Delta_t^\varepsilon\|^4 |\mathcal{F}_t \right]
    \tag{AM-GM inequality}
    \\
    \leq &
    8 \|\Delta_t^\nabla\|^4  + 64 \sigma^4
    \tag{by \cref{helper:diff:4th}}
    \\
    \leq & 8L^2 \|\Delta_t^2\|^2 \|\Delta_t^\nabla\|^2 + 64 \sigma^4
    \tag{by $L$-smoothness}
    \\
    \leq & 8L^3 \|\Delta_t^2\|^2 \langle \Delta_t, \Delta_t^\nabla \rangle + 64 \sigma^4.
    \tag{by \cref{eq:fedavg:a2:stab:prop:2}}
  \end{align}
  Plugging the above two bounds to \cref{eq:fedavg:a2:stab:prop:1} gives
  \begin{align}
    \expt [\|\Delta_{t+1}\|^4 |\mathcal{F}_t]
    \leq
    \|\Delta_t\|^4 - 4\eta (1 - 2 \eta L - 6 \eta^3 L^3) \|\Delta_t\|^2 \langle \Delta_t, \Delta_t^\nabla \rangle + 16 \eta^2 \|\Delta_t\|^2 \sigma^2 + 192 \eta^4 \sigma^4.
    \label{eq:fedavg:a2:stab:prop:3}
  \end{align}
  Since $\eta L \leq \frac{1}{4}$ we have $(1 - 2 \eta L - 6 \eta^3 L^3) > 0$. By convexity $\langle \Delta_t, \Delta_t^\nabla \rangle \geq 0$. Hence the second term on the RHS of \cref{eq:fedavg:a2:stab:prop:3} is non-positive. We conclude that
  \begin{align}
    \expt [\|\Delta_{t+1}\|^4 |\mathcal{F}_t]
    \leq
    \|\Delta_t\|^4 + 16 \eta^2  \sigma^2 \|\Delta_t\|^2 + 192 \eta^4 \sigma^4.
  \end{align}
  Taking expectation gives
  \begin{align}
         & \expt [\|\Delta_{t+1}\|^4] \leq \expt [ \|\Delta_t\|^4 ] + 16 \eta^2 \sigma^2 \expt[\|\Delta_t\|^2] + 192 \eta^4 \sigma^4
    \\
    \leq & \expt [ \|\Delta_t\|^4 ] + 16 \eta^2 \sigma^2 \sqrt{\expt[\|\Delta_t\|^4]} + 192 \eta^4 \sigma^4 = \left( \sqrt{\expt \|\Delta_t\|^4} + \sqrt{192} \eta^2 \sigma^2 \right)^2.
  \end{align}
  Taking square root on both sides completes the proof.
\end{proof}
With \cref{fedavg:a2:stab:1} at hand we are ready to prove \cref{fedavg:a2:stab:main}.
\begin{proof}[Proof of \cref{fedavg:a2:stab:main}]
  Let $t_0$ be the latest synchronized prior to $t$, then telescoping \cref{fedavg:a2:stab:1} yields (note that $\Delta_{t_0} = 0$)
  \begin{equation}
    \sqrt{\expt \|\Delta_{t}\|^4} \leq \sqrt{192} \eta^2 \sigma^2 (t - t_0) \leq \sqrt{192} \eta^2 K \sigma^2,
  \end{equation}
  where the last inequality is because $K$ is the synchronization gap. Thus
  \begin{align}
    \frac{1}{M} \sum_{m=1}^M \expt \left[ \left\| \overline{w_t} - w_t^m \right\|^4  \right] \leq \expt [\|\Delta_t\|^4] \leq 192 \eta^4 K^2 \sigma^4,
  \end{align}
  where the first ``$\leq$'' is due to Jensen's inequality.
\end{proof}

config: arxiv


text:
\begin{document}
\title{
Compute-and-Forward: Optimization Over Multi-Source-Multi-Relay Networks}
\author{Zhi~Chen~\IEEEmembership{Member,~IEEE,}, Pingyi~Fan~\IEEEmembership{Senior Member,~IEEE,}
        and~Khaled~Ben~Letaief~\IEEEmembership{Fellow,~IEEE}
\thanks{Z. Chen and P. Fan are with the Department of Electrical Engineering, Tsinghua University, Beijing, China, 100084. Emails: chenzhi06@mails.tsinghua.edu.cn; fpy@tsinghua.edu.cn.

K. B. Letaief is with Department of Electronic and Computer Engineering,
Hong Kong University of Science and Technology, Hong Kong (e-mail:
eekhaled@ece.ust.hk).}}
\maketitle

\maketitle
\baselineskip 24pt
\begin{abstract}\\
\baselineskip=18pt
In this work, we investigate a multi-source multi-cast network with
the aid of {an arbitrary number of} relays, where it is assumed that no
direct link is available at each S-D pair. {The aim is
to find the fundamental limit on
the maximal common multicast throughput of all source nodes if resource allocations are available.} A transmission protocol employing
the relaying strategy, namely, compute-and-forward (CPF), is proposed.
{We also adjust the
methods in the literature to obtain the integer network-constructed
coefficient matrix (a naive method, a local optimal method as well as a global optimal method)
to fit for the general topology with an arbitrary number of relays.}

Two transmission scenarios are addressed. The first scenario is delay-stringent transmission where each message must be delivered within one slot.
The second scenario is delay-tolerant transmission where
no delay constraint is imposed. The associated optimization problems
to maximize the short-term and long-term common multicast throughput
are formulated and solved, and the optimal allocation of power and
time slots are presented. Performance comparisons show that the CPF
strategy outperforms conventional decode-and-forward (DF) strategy.
{It is also shown that with more relays, the CPF strategy performs even
better due to the increased diversity.}
Finally, by simulation, it is observed that for a large network in relatively
high SNR regime, CPF with the local optimal method for the
network-constructed matrix can perform close to that with the global
optimal method.
\baselineskip=18pt
\end{abstract}
\begin{keywords}
Compute-and-forward, resource allocation, delay-stringent, delay-tolerant, fading
\end{keywords}
\IEEEpeerreviewmaketitle

\section{Introduction}
Network coding, an efficient way to mitigate network interference and improve throughput, was firstly proposed by Yeung \emph{et al} in \cite{ahlswede2000network}\cite{li2003linear} for wireline networks. Employing it, relay nodes capability  is expanding from only simply forwarding messages to forwarding some functions of different messages to multiple destinations nodes. The intended nodes can then detract the message required as long as they have prior knowledge of the rest messages. In this way, wireline network throughput is improved. Furthermore, network coding is shown to be promising in wireless networks in terms of throughput improvement in \cite{fragouli2007wireless}\cite{li2010network}\cite{niu2010cooperative}\cite{zhang2009network}.

However, due to broadcast nature in wireless communications, the performance potential of conventional digital network coding (DNC) is strictly constrained by the interference from transmissions of other irrelevant transmitters in wireless networks. For instance, for a simple three-node, two-way relay network (TWRN), the relay node needs to jointly decode two individual messages in the multi-access phase with DNC \cite{oechtering2008stability}\cite{chen2012energy}\cite{Chen2012two-way5}, whereas the performance is degraded due to the fact that one user's message is regarded as interference to the message of the other user at the relay node in the multi-access uplink. As it is, this interference constrains the achievable rate pair especially in high signal-to-noise ratio (SNR) regime. Other strategies, like amplify-and-forward (AF) and compress-and-forward (CF), has its intrinsic advantages without decoding individual messages, but the noise term will be amplified along with messages delivery throughout the network, resulting in degradation of network performance.

To this end, a smart way of network coding, which was referred to as
physical layer network coding (PNC) \cite{popovski2007physical}
\cite{zhang2006physical}\cite{zhang2009channel}\cite{wang2009channel}
\cite{lu2011asynchronous}, i.e., compute-and-forward (CPF) \cite{nazer2011compute}
for general multi-source multi-relay networks, attracted increasing attention.
Typically for a TWRN, in the uplink phase, the two source nodes simultaneously
transmit their messages to the relay node, and the relay node merely decodes a
linear combination of these two messages with integer coefficients, other than
 employing joint decoding. In the downlink phase, the relay node can transmit to
 the two source nodes the linear combined  message. In this case, the two
 sources can subtract their transmitted messages and then decode the
 intended message. In this way, the relay node mitigates the
 interference coming from joint decoding in the uplink phase by
 employing digital network coding, and also avoids the noise amplification
 when using analog network coding (ANC) \cite{katti2007embracing}
 \cite{zhang2009optimal}. In \cite{narayanan2007joint}, PNC was shown
 to perform close to a capacity upper bound and its performance gain
 over DNC and ANC was demonstrated.

More importantly, PNC is shown to achieve high performance for more general
multi-source multi-relay networks in \cite{nazer2011compute}
\cite{nazer2011reliable}\cite{xiao2009design}\cite{xiao2010multiple}.
In the celebrated work \cite{nazer2011compute}, PNC was referred to
as compute-and-forward (CPF) strategy. With this strategy, each relay
will decode a function message formed by a linear combination of the
messages from all source nodes with a selected integer coefficient vector.
Each destination node hence obtains different function messages from
various relay nodes and decodes all the source messages as long as
the integer coefficients constructed matrix is in full rank. In the literature,
the outage probability performance of CPF is demonstrated to
outperform other relaying strategies, such as decode-and-forward (DF)
and amplify-and-forward (AF). In all these works,
a type of linear codes, lattice code were employed to achieve
the derived CPF capacity region.
On the other hand, Nazer \emph{et al} mainly investigated the outage rate
of the CPF rate over S-R links in \cite{nazer2011compute},
whereas the transmission over R-D links was not discussed.
In \cite{osmane2011compute}, how to obtain the locally optimal integer
coefficient vector at each relay node to maximize its computation rate
was addressed,
however the full rank of the matrix constructed by all the integer
coefficient vectors of all relays was not guaranteed hence the
destination nodes may still not be able to decode all source messages.
In \cite{wei2012compute} it jointly optimized the integer coefficient
vector at each relay to finally obtain the optimal common computation
rate with the guaranteed full rank matrix constructed by these integer
coefficient vectors, at the cost that the achievable computing rates at
some relay nodes may be reduced to satisfy the full rank requirement.
In addition, the outage performance of CPF under some specific network
configurations were addressed in the literature, e.g.,
\cite{osmane2011compute} for a three-transmitter multi-access
network and \cite{wang2012outage} for multi-way relay networks
with the aid of only one relay node.

{Note that all the previous works in the literature for general networks
only considered the outage performance of CPF strategy with
constant transmit power and no optimal resource allocation was
studied. In this work, we therefore aim to investigate the achievable throughput with CPF
with adjustable resource allocation strategies, assuming that
full channel state information (CSI) is available
at all transmitters. By doing so,
the fundamental limit on the maximal common multicast
throughput can be obtained, which is useful in
the system design of multi-source multi-relay networks.}

There are also some works investigating the potential
of lattice codes as a capacity-achieving codes in
\cite{narayanan2007joint}\cite{zamir2002nested}\cite{erez2004achieving}
\cite{gunduz2009multi}, where in \cite{zamir2002nested}
Zamir \emph{et al} discussed nested lattice codes for
structured multi-terminal binning, and in \cite{erez2004achieving}
Zamir \emph{et al} showed the capacity-achievable property of
lattice codes over AWGN channels. In \cite{gunduz2009multi},
lattice codes was employed for a multi-way relay network and
the capacity region to within a half-bit of the cut-set bound
was shown to be achievable with it.

In this work, a general multi-source multi-relay multicast network is considered and the lattice codes are employed to realize CPF transmission.
{The performance with CPF in terms of the fundamental limit on the achievable common
multicast throughput is investigated},
with both the S-R links as well as the R-D links taken into account.
Two cases of interest will be studied. One is a delay-stringent scenario,
where each multicast transmission from all sources to all destinations
must be finished in one slot. The other is a delay-tolerant
transmission scenario, where the multicast transmission from all sources
to all destinations can be finished within arbitrary finite number of slots, i.e., no delay constraints are imposed. The major contributions of this work
are listed as follows.

    \begin{itemize}
    \item We design a CPF based multi-source multicast transmission protocol
    for the topology {with an arbitrary number of relays}.
    \item
    For the delay-stringent scenario, an optimization problem to maximize the common multicast throughput over one block with the specified channel gains employing CPF, is formulated and solved.
    \item
    For the delay-tolerant scenario, an optimization problem to
    maximize the average common multicast throughput employing CPF
    by allocating time and power resources, is formulated and solved analytically.
    {In addition, the convexity of the formulated problem is
    proved analytically.}
    \item
    We find that through simulations,
    \begin{enumerate}
    \item with an arbitrary number of source nodes,
    CPF with global-optimized network-constructed matrix
    outperforms DF strategy, which verifies the superiority of CPF over DF.
    \item CPF performs better with the increasing number of relay nodes.
    \item with a small number of source nodes, CPF with local-optimized network-constructed matrix performs slightly worse than DF for delay-stringent scenario, while slightly better than DF for delay-tolerant scenario, due to higher rank failure probability.
    \item with a relatively large number of source nodes, CPF with local
    optimized network-constructed matrix approximates
    the performance of CPF with global optimized network-constructed
    matrix and outperforms DF, due to the reduced rank failure probability.
    {This finding makes the implementation of CPF more practical and flexible,
    due to the greatly reduced network overhead
    information exchange required by the global optimized method
    for forming matrix.}
    \end{enumerate}

\end{itemize}

The remainder of this work is organized as follows.
In Section II, we present the system model of a multi-source
multicast network with the aid of multiple relay nodes,
and describes clearly the transmission protocol and the decoding procedure at the destination nodes.
In Section III and Section IV, delay-limited scenario and delay-tolerant
scenario are investigated, respectively.
{The associated problems to find the fundamental limit on the maximal common multicast throughput
of the entire network} are formulated and solved, by jointly allocating
time and energy resources for each transmit phase.
Simulation results are presented in Section V.
Finally, we conclude this work in Section VI.

\begin{figure}[!h]
   \centering
   \includegraphics[width = 6cm]{system_model.eps}
   \caption{System model for a multi-source multicast network with the aid of multiple relay nodes. The direct link between each S-D pair is assumed to be unavailable.} \label{fig:system_model}
   \end{figure}

\section{System Model}
In this work, we mainly focus on a multi-source, multi-relay, multicast network, as shown in Fig. \ref{fig:system_model}, which consists of $M$ sources ($S_1$,\ldots, $S_M$), $K$ relays ($R_1$,\ldots, $R_K$) and $L$ destinations ($D_1$,\ldots, $D_L$). Each source node or relay node is
equipped with one antenna and works in half-duplex mode, i.e., cannot
transmit and receive data simultaneously. No direct link between any
S-D pair is assumed to be available. Hence, transmission must be assisted
by relay nodes.

A block fading channel model is also assumed for each link. The stochastic and instantaneous channel gain of each link are assumed to be known at both the transmitter and the receivers, which can be realized via feedback.
{In this work, lattice coding is adopted, as it can not only achieve close-to capacity rate,
but also preserves the linear property, i.e., the linear combination of lattice
codes is also a lattice point, which is essential for CPF strategies.}

We define $h_{ml}$ as the channel fading coefficient of the link $S_l$-$R_m$ and ${\bf z}_m$ as the i.i.d. additive white Gaussian noise vector, i.e., ${\bf z}_m \sim {\it N}({\bf 0},{\bf I}_n)$. We also denote ${\bf h}_m=[h_{m1} \ldots h_{mM}]^T$ by the channel coefficient vector consisting of the links from all sources to $R_m$. Similarly, we denote ${\bf g}_m=[g_{m1} \ldots g_{mL}]^T$ by the channel coefficient vector consisting of the links from the $m$th relay to all destinations and $g_{m_{\min}}=\min_i |g_{mi}|^2 $ as the minimum channel gain of the links from the $m$th relay to all destinations.
We also assume that all nodes are with the same average power constraint $P_0$
In addition, we assume the integer coefficient vector adopted at the $m$th relay in decoding the function message is ${\bf a}_m \in {\mathbb{Z}}^M$, which is carefully selected to form a
probably decodable integer-combined version. It is interestingly noted that
the $m$th relay can select different ${\bf a}_m$ to form different decodable function messages, albeit at different CPF rates.

As shown in Fig. \ref{fig:transmission}, a CPF based transmission protocol consists of $M+1$ phases. The detailed procedure for the case $K=M$\footnote{It will be discussed in Section \ref{sec:general}
for the transmission procedure for the two cases of $K<M$ and $K>M$.} is illustrated as follows.
\begin{itemize}
\item In Phase $1$, all source nodes transmit their messages simultaneously to all the relay nodes. At the end of this phase, each relay node decodes one or more linear equations of the combination of individual transmitted messages from all sources with selected integer NC coefficient vectors.
\item In Phase $i$ ($i=2, \ldots, M+1$), the $(i-1)$th relay delivers its decoded function message to all destination nodes.
    At the end of Phase $i$, all the destination nodes decode the function message received and store it for source-message decoding at the end of Phase $M$.
\item Source-message decoding at the end of Phase $M+1$: with $M$ decoded function messages
from the relay nodes, each destination nodes tries to recover all original messages.
This would be possible if sufficient amount of equations are received reliably,
i.e., rank $M$ (the number of source nodes) of the matrix constructed by all these integer coefficients is achieved at each destination node.

\end{itemize}

\begin{figure}[!h]
   \centering
   \includegraphics[width = 7cm]{Transmission_protocol.eps}
   \caption{{Transmission protocol for a multi-source multicast network
   with the aid of relays. The direct link between each S-D pair is assumed
    to be unavailable. In this figure, each $S_i$ wants to broadcast a message
    $x_i$ to all destinations nodes and they simultaneously transmit
    in the first slot.
   In the following slots, each relay node forwards a decoded combined
   message
   to all destinations.
   For example, in the $i+1$th slot, $R_i$ broadcasts a message
   $f_i(x_1,x_2,\cdots,x_M)=\sum_{j=1}^M a_{ij}x_j$ to all destination nodes.}} \label{fig:transmission}
   \end{figure}

Note that for the CPF phase (Phase 1), from \cite{nazer2011compute},
with the specified ${\bf a}_m$, the CPF rate at the $m$th relay
for the real-valued AWGN networks, achieved by
lattice coding, is
\begin{align}
R_{\mathrm{CPF}}^m = \max_{{\bf a}_m \in \mathbb{R}}\frac{1}{2}\log^+ \left( ||{\bf a}_m||^2 - \frac{P({\bf h}_m^{\mathrm{T}}{\bf a}_m)^2}{1+P||{\bf h}_m||^2} \right)^{-1}. \label{eqn:optrate}
\end{align}
where $P$ is the transmit power. The achievable common rate of all relays, is hence given by,
\begin{align}
R_{\mathrm{CPF}}=\min_i R^i_{\mathrm{CPF}} \quad i=1,\ldots, K.
\end{align}
Correspondingly, the required transmit power at all source nodes
for the $m$th relay with the common transmit rate $R_{\mathrm{CPF}}$, i.e., $P_{\mathrm{CPF}}^m$, is given by,
\begin{align}
P_{\mathrm{CPF}}^m=\frac{1-2^{2R_{\mathrm{CPF}}}b_m}{2^{2R_{\mathrm{CPF}}}c_m-a_m} \label{eqn:app1_2}
\end{align}
where
$a_m=|{\bf h}_m|^2$, $b_m=|{\bf a}_m|^2$,
$c_m=|{\bf h}_m|^2|{\bf a}_m|^2-|{\bf h}_m^{\mathrm{T}} {\bf a}_m|^2$ and
$d_m=|{\bf h}_m^{\mathrm{T}} {\bf a}_m|^2=a_mb_m-c_m$.
The required transmit power at all source nodes for the CPF phase is given by,
$$P_{\mathrm{CPF}}=\max_m P_{\mathrm{CPF}}^m.$$ At this power level, all relays can enjoy the common computing rate $R_{\mathrm{CPF}}$.

In addition, for the achievable common CPF rate, we would like to show an interesting property of
the coefficients $a_m$, $b_m$ and $c_m$ related to $R_{\mathrm{CPF}}$,
which is summarized in Lemma \ref{appendix_1}.
\begin{lemma} \label{appendix_1}
With positive transmit power at relay nodes, we have
\begin{align}
\max_m (\frac{1}{b_m})<2^{2R_{\mathrm{CPF}}}<\min_m \frac{a_m}{c_m} \quad m=1,\ldots,K
\end{align}
\end{lemma}
The proof is given in Appendix \ref{appen_2} and is omitted here.

{On the other hand, for the relaying phases, the relay nodes forward the function messages
to all destination nodes consecutively. Since all destination nodes need to
successfully receive the function messages, the broadcast rate is determined by
the minimum channel gain of the R-D links for the transmitted relay node, i.e., $g_{m_{\min}}$. Hence, the
broadcast rate of the $i$th relay is determined by,}
 \begin{align}
 R_{r_i}=\frac{1}{2}\log_2(1+P_{r_i}g_{m_{\min}}).
 \end{align}

\subsection{Case Study: A Simple Example} \label{sec:case}
{Consider a three-source, three-relay and three-destination network.
The message transmitted by each source is $x_i$ ($i=1,2,3$) in the first hop.
The first relay decodes a function message of $2x_1+3x_2+5x_3$ ($z_1$) and forwards it to all destinations in the second hop.
The second and the third relay decode function messages of
$x_1+x_2+x_3$ ($z_2$) and $2x_1+x_2+5x_3$ ($z_3$) respectively. Hence they transmit
these function messages in the third hop and the fourth hop consecutively.
Assume that all destination nodes successfully receive the three function messages and attain the coefficient vectors of the three relays.
To obtain the three original source messages, all destination nodes
are then required to solve the following equation below.}
\begin{align}
\begin{bmatrix} 2 & 3 & 5 \\ 1 & 1 & 1\\ 2 & 1 & 5 \end{bmatrix}
\begin{bmatrix} x_1 \\ x_2\\ x_3 \end{bmatrix}
=A\begin{bmatrix} x_1 \\ x_2\\ x_3 \end{bmatrix}=\begin{bmatrix} z_1 \\ z_2\\ z_3 \end{bmatrix}
\end{align}
{To ensure the uniqueness of the solution to this matrix equation,
$A$ must be a full-rank matrix, i.e., $rank(A)=3$ for the given case
and each destination can then decode all the original messages from the sources.

Unfortunately, by letting $z_3=z_1=2x_1+3x_2+5x_3$,
we have $rank(A)=2$ and $x_1$ and $x_3$ can not be decoded at the destinations
and the entire transmission fails.
The importance of selection of the coefficient vectors to guarantee $A$ is
in full rank is therefore verified and we discuss it} in Sec. \ref{sec:vector} for $K=M$
and Sec. \ref{sec:general} for the general topology.

\begin{figure}[!t]
   \centering
   \includegraphics[width = 5.5cm]{decoding_procedure.eps}
   \caption{{Decoding procedure for a multi-source multicast network
   with the aid of relays. The direct link between each S-D pair is
   assumed to be unavailable. In this figure, $S_i$ wants to broadcast
   a message $x_i$ to all destinations nodes. They simultaneously transmit
   in the first slot and the relay node $R_i$ decodes a combined version
   of all source messages, $z_i=\sum_{i=1}^M a_{ij}x_j$.
   In each of the following slots, each destination
   node decodes a combined message from the relays and it
   tries to decode the source messages by solving a matrix
   equation at the end of the transmission procedure as depicted.}} \label{fig:decoding}
   \end{figure}

\subsection{Integer Coefficient Vector for the case $K=M$} \label{sec:vector}

From the discussion above, it is observed that selection of
coefficient vector is crucial for the achievable CPF rates.
Hence for clarity, we shall briefly review
the three methods in the literature to
obtain the integer coefficient vectors. {Note that all these
methods in the literature were only presented for the case $K=M$.}
For a given ${\bf h}_m$ at the $m$th relay, we can find ${\bf a}_m$ by
\begin{itemize}
\item method a) (naive method): obtain the integer coefficient vector ${\bf a}_m$ individually by solving $a_{mj}=\mathrm{round}(h_{mj})$ ($j=1,\ldots,M$) at the $m$th relay, where the function $\mathrm{round}(\cdot)$ returns the closest integer to $\{\cdot\}$.
\item method b) (local optimal method): obtain the locally-optimal integer coefficient vector ${\bf a}_m$ individually as in \cite{osmane2011compute} with the respective ${\bf h}_m$ at the $m$th relay.
\item method c) (global optimal method): obtain all the integer coefficient vectors (${\bf a}_m$) ($m=1,\ldots,M$) at all relays by employing the jointly optimization method in \cite{wei2012compute} to ensure that the full rank requirement of the matrix constructed by the integer coefficient vectors of all relays is satisfied.

\end{itemize}
For methods a) and b), only the local channel information is required. Therefore, they can not guarantee the full rank requirement at the destination nodes, i.e., the corresponding destination node may ultimately fail to decode the source messages and lead to the failure of the entire transmission.

On the other hand, method c) can guarantee the full matrix requirement
at the cost of indispensable additional signaling overhead among the relay nodes,
due to the joint optimal search procedure. It is therefore expected
to perform better than method a) and b) at the sacrifice of higher overhead.

\subsection{General Scenario} \label{sec:general}
{From the case study in} \ref{sec:case}, {for the general M-source multicast network,
we must also have $rank(A)=M$ to decode all source messages at the destinations, as shown in} Fig. \ref{fig:decoding}. {It is readily concluded that at least $M$ transmissions are required by the relay nodes.
However, if more than $M$ relay transmissions are allowed,
it incurs inefficient spectrum resource usage and the throughput is decreased.
In this sense, we only allow $M$ (the number of source nodes) relay phases over the R-D transmissions in this paper.

Therefore, if $K \ge M$, we can always select the best $M$ relays for transmission.
For instance, method b) can be extended in this case by selecting the best $M$ relays
with the highest achievable CPF rates among all relays in an descending order.
For an extended method c), however, we need to select the best
$M$ vectors
with their highest achievable CPF rates among all relays in an descending order,
with the constraint that the corresponding coefficient vectors linearly
independent.

For the case that $K < M$, by letting some relays to decode more linear combined messages in the first phase (selecting more than one ${\bf a}_m$), $M$ relay
transmission phases are also attainable to make the destinations capable of decoding all source messages.
For example, an implementation algorithm employing an extended method c) is given as follows.}
\begin{itemize}
\item At the end of Phase 1, each relay node tries to decode a number of function messages with
the corresponding achievable CPF rates.
\item Each relay node broadcasts its achievable rates and the message index to a central controller and
it lists different CPF rates in an descending order and selects the highest $M$
CPF rates with their coefficient vectors linearly independent.
\item The central controller notifies each relay its selected rates and the associated message index, followed by
the M relaying transmission phases.

\end{itemize}
{Hence, for the case $K<M$, some relay nodes with higher CPF rates
 may need to transmit in multiple slots.}

In the following, {we seek to find the fundamental limit on the optimal common multicast throughput over the delay-stringent and the delay-tolerant scenarios respectively.}
It is emphasized that the considered relaying phases still consists of $M$ phases.
As discussed above, it follows from the full rank requirement of the constructed matrix
if $K<M$ and more efficient usage of spectrum resources
if $K>M$. Hence, in this work, we only allow $M$ relaying phases.

\section{Delay-Stringent Transmission}
In this section, we shall address the achievable common multicast throughput with CPF under stringent delay constraints, i.e., messages must be delivered within one slot from the source nodes to the destination nodes. In this case, the aim is to maximize the achievable common multicast throughput by optimally allocating the time resources to each phase within each slot. This can be applied to some realtime communication applications with the minimum rate requirement within each slot.

Firstly, we denote $f_{\mathrm{CPF}}$ as the time fraction assigned to the first phase for S-R transmission with CPF and $f_i$ as the time fraction allotted to the $(i+1)$th phase or the $i$th relay for transmission.

In addition, $P_{\mathrm{CPF}}$ and $P_{i}$ are denoted
as the transmit powers of each source node at the first phase
(compute-and-forward phase) and the $i$th relay, respectively.
We also denote {$R_{\mathrm{CPF}}$ and $R_{i}$ as the transmit
 rate of the corresponding phase. Hence, the products
 $f_{\mathrm{CPF}}R_{\mathrm{CPF}}$ and $f_iR_{i}$ are the throughput
 of the CPF phase and the $i$th relay phase, respectively. The end to end
 throughput per slot is hence determined by the minimum of the
 throughput of all phases, namely,
 $\min(f_{\mathrm{CPF}}R_{\mathrm{CPF}},f_iR_{i})$ ($i=1,\cdots,M$).}

Based on the above analysis, the optimization problem within one slot for given channel coefficients of all links, termed as {\bf P1}, is formulated as follows,
\begin{align}
\max_{f_i,f_{\mathrm{CPF}}} \min \bigl( f_{\mathrm{CPF}}R_{\mathrm{CPF}}(P_{\mathrm{CPF}}),f_iR_{i}(P_i) \bigr)  \label{eq:opt_delay}
\end{align}
where $i=1,\ldots M$. The associated power constraint and the physical constraint are given as follows.
\begin{align}

P_{\mathrm{CPF}} \leq P_{0} \label{eq:con1_delay}\\
P_i \leq P_{0} \label{eq:con2_delay} \\
f_{\mathrm{CPF}}+\sum_{i=1}^M f_i \leq 1 \label{eq:con3_delay}
\end{align}
where the objective function in (\ref{eq:opt_delay}) is to maximize
the minimum throughput of each phase.
(\ref{eq:con1_delay}) and (\ref{eq:con2_delay}) give the power
 constraint of the CPF phase and the $i$th relaying phase. Hence, by transmitting at
 $P_0$, each phase achieves its optimal rate.

{Generally speaking, this is a max-min optimization problem and seems
difficult in the first glance.
We hence rewrite {\bf P1} in an equivalent form as follows,}
\begin{align}
\max \quad f_{\mathrm{CPF}}R_{\mathrm{CPF}}(P_0)
\end{align}
subject to
\begin{align}
f_{\mathrm{CPF}}R_{\mathrm{CPF}}(P_0)=f_iR_{i}(P_0)
\end{align}
and the physical constraint in (\ref{eq:con3_delay}).
{This follows from the fact that the achievable rate
is an increasing function of the transmit power and hence
to achieve optimality requires transmission with the highest power allowed.}

This is however a simple linear optimization problem and the solution to {\bf P1} is given by,
\begin{align}
f_{\mathrm{CPF}}^* = \frac{L^*}{R_{\mathrm{CPF}}(P_{sm})} \label{eq:time_delay_relay}\\
f_i^* = \frac{L^*}{R_{i}(P_{r_i})} \label{eq:time_delay_cpf}\\
L^* = \frac{1}{\frac{1}{R_{\mathrm{CPF}}(P_{sm})}+\sum_{i=1}^{M}\frac{1}{R_{i}(P_{r_i})}} \label{eq:throughput_delay}
\end{align}
where the asterisk denotes optimality and $L^*$ is the optimal achievable common multicast throughput with given channel coefficients in delay-stringent applications.

Note that here we consider the simple case that the average power
constraint at each node is imposed as the constraint for each node
for the discussed slot, with the aim to optimize the throughput
within each slot.

\newtheorem{remark}{Remark}
\begin{remark}
A typical application of it is the realtime multicast multimedia which requires the minimum data
rate for the discussed slot. Another typical case would be a very slow fading scenario where
channels can be regarded as quasi-static. For other applications without
the minimum rate constraint or in a relatively fast fading environment, a possible extension
of {\bf P1} is to optimize the throughput within multiple slots by allocating different
power resources to each slot due to channel variations while still maintaining the
delay-stringent constraint, such that the throughput can be improved.
\end{remark}

\section{Delay-Tolerant Transmission}
In this section, we shall address the performance achievable by employing compute-and-forward with full channel state information available at the transmitters. We also assume that the derived integer coefficient vectors of all relays are available at all the destination nodes. As it is, all these messages can be obtained via a feedback channel in a time-sharing manner or frequency-sharing manner for each node.

For clarity, we denote $\bar{P}_{\mathrm{CPF}}$ as the average power consumed
at each source node for the CPF phase and $\bar{P}_{i}$ the average power
consumed of the $i$th relay, respectively. Correspondingly,
{$\bar{R}_{\mathrm{CPF}}$ as the average CPF rate achievable
for the CPF phase and $\bar{R}_{i}$ the average rate
at the $i$th relay, respectively.}
Note that all these values are averaged over
the associated channel coefficient distributions. {Hence
$f_i\bar{R}_{\mathrm{CPF}}$ and $f_i\bar{R}_{i}$ are the
average throughput of the CPF phase and the $i$th relaying phase, respectively.
Hence, the throughput of the transmission protocol considered for the delay-tolerant scenario
is given by}
$$\min ( f_{\mathrm{CPF}}\bar{R}_{\mathrm{CPF}}({\bf h}_l), f_i\bar{R}_i({\bf g}_{i})).$$
{To achieve the optimal throughput, we hence need to adjust
power and time resources allocated for each phase.} The associated problem to
maximize the average common multicast throughput, referred to as {\bf P2},
is formulated as follows.

\begin{align}
\max_{R_{\mathrm{CPF}}({\bf h}_l), R_i({\bf g}_i)} \quad \min ( f_{\mathrm{CPF}}\bar{R}_{\mathrm{CPF}}({\bf h}_l), f_i\bar{R}_i({\bf g}_{i})) \label{opt_CSIT}
\end{align}
where $l,i=1,\ldots M$. The associated constraints are given as follows.
\begin{align}
&\bar{P}_{\mathrm{CPF}} \leq P_{0}\label{con1_CSIT}\\
&\bar{P}_i \leq P_{0} \quad i=1,\ldots,M \label{con2_CSIT}\\
&f_{\mathrm{CPF}}+\sum_{i=1}^M f_i \leq 1 \label{con3_CSIT}
\end{align}

The optimization problem above aims to optimally allocate time resources to different phases in order to mitigate the performance degradation caused by bottleneck links. In this sense,
{\bf P2} can also be transformed into an equivalent optimization problem below
\begin{align}
\max_{R_{\mathrm{CPF}}({\bf h}_l), R_i({\bf g}_i)} \quad f_{\mathrm{CPF}}\bar{R}_{\mathrm{CPF}}
\end{align}
subject to the average power constraints in (\ref{con1_CSIT})-(\ref{con2_CSIT}), and the physical constraint in (\ref{con3_CSIT}), and
\begin{align}

&f_{\mathrm{CPF}}\bar{R}_{\mathrm{CPF}}=f_i\bar{R}_i \quad \forall i \label{con4_CSIT}

\end{align}
where (\ref{con4_CSIT}) points out the fact that the product of average rate of each phase and the time resource allotted to that phase should be made equal for each phase for optimality.
Note that it cannot be readily observed that {\bf P2} is a convex optimization problem due to the relationship between $P_{\mathrm{CPF}}$ and $R_{\mathrm{CPF}}$ in (\ref{eqn:app1_2}). Fortunately,
it is verified that $P_{\mathrm{CPF}}$ is indeed a convex function of $R_{\mathrm{CPF}}$ and therefore {\bf P2} is a convex optimization problem. This observation is summarized in Theorem \ref{thm_convex} and the proof is presented in Appendix \ref{appen_1}.
\begin{thm}\label{thm_convex}
{\bf P2} is a convex optimization problem.
\end{thm}
The detailed proof is given in Appendix \ref{appen_1} and omitted here.
According to Theorem \ref{thm_convex},

{\bf P2} can be solved by the Lagrangian multiplier method. The associated Lagrangian multiplier function is given by,
\begin{align}
&F(R_{\mathrm{CPF}},\gamma_i,\beta_i)=f_{\mathrm{CPF}}\bar{R}_{\mathrm{CPF}}
-\beta_0 \left(\bar{P}_{\mathrm{CPF}}-P_{0} \right) \nonumber\\
&- \sum_{i=1}^M \beta_i \left(\bar{P}_{i}-P_{0}  \right) - \sum_{i=1}^M \gamma_i \left( f_{\mathrm{CPF}}\bar{R}_{\mathrm{CPF}}-f_i\bar{R}_{i} \right)\nonumber\\
&-\alpha( f_{\mathrm{CPF}}+\sum_{i=1}^M f_i -1)
\end{align}
where $\beta_0$, $\beta_i$ are the Lagrangian
multipliers with respect to the power constraints of the CPF phase
and the relaying phases. $\gamma_i$ is the Lagrangian multiplier
with respect to the rate constraint. $\alpha$ is the Lagrangian multiplier for the physical
constraint.

The KKT conditions are given by,
\begin{align}
&f_{\mathrm{CPF}}-\beta_0 \cdot \frac{dP_{\mathrm{CPF}}}{dR_{\mathrm{CPF}}} -\sum_{i=1}^L \gamma_i f_{\mathrm{CPF}} = 0 \label{kkt_1}\\
&-2\beta_i\ln 2\frac{2^{2R_{i}}}{g_{r_i}}+\gamma_if_i=0 \label{kkt_2} \\
&\bar{R}_{\mathrm{CPF}} - \sum_{i=1}^M \gamma_i\bar{R}_{\mathrm{CPF}}  - \alpha = 0 \label{kkt_3} \\
&\gamma_i\bar{R}_i-\alpha = 0\label{kkt_4}\\
&f_{\mathrm{CPF}}+\sum_{i=1}^M f_i = 1 \label{kkt_5}
\end{align}

From (\ref{kkt_2}), we then obtain the optimal rate as well as the allotted power for the relaying phases.

\begin{align}
R^*_{r_i} &= \frac{1}{2}\log_2[\frac{f_i^*\gamma_i^*}{2\beta_i^*\ln 2}g_{r_i}]^+, \quad i=1, \cdots, M \label{eqn:rate_relay}\\
P^*_{r_i}&=[\frac{f_i^*\gamma_i^*}{2\beta_i^*\ln 2 }-\frac{1}{g_{r_i}}]^+ \quad i = 1,\cdots, M \label{eqn:power_relay}
\end{align}
where (\ref{eqn:power_relay}) can be directly obtained from the Shannon capacity formula $P=(2^{2R}-1)/g$ and (\ref{eqn:rate_relay}).

For the computation rate, by deriving the first derivative of $R_{\mathrm{CPF}}$ to $P_{\mathrm{CPF}}$ in (\ref{eqn:app1_2}) and insert it into (\ref{kkt_1}) and employing some arithmetic operations, we arrive at the follow equality.
\begin{align}
\frac{l}{\beta_0}=\frac{dx}{(cx-a)^2} \label{equality}
\end{align}
where $l=f_{\mathrm{CPF}}(1-\sum_{i=1}^M \gamma_i)/(2\ln2)$ and $x=2^{2R_{\mathrm{CPF}}}$.
Here we omit the subscripts of $a_m$, $b_m$ and $c_m$ for simplicity and we assume that
$R_{\mathrm{CPF}}^m=\min_i R_{\mathrm{CPF}}$, i.e., the $m$th offers the lowest
CPF rate.

By applying some arithmetic operations on (\ref{equality}), we arrive at a mono basic quadratic equation,
\begin{align}
lc^2x^2-(2acl+\beta_0d)x+a^2l=0
\end{align}
Therefore, the possible solutions to this quadratic equation are given by,
\begin{align}
x=\frac{2acl+ \beta_0 d \pm \sqrt{\beta_0^2d^2+4acld\beta_0} }{2lc^2}
\end{align}

Furthermore, it can be shown that only one out of these two possible solutions to this quadratic equation is feasible.

Firstly, we consider a possible solution, $x_1=\frac{2acl+ \beta_0 d + \sqrt{\beta_0^2d^2+4acld\beta_0} }{2lc^2}$. We then have
\begin{align}
x_1=&\frac{2acl+ \beta_0 d + \sqrt{\beta_0^2d^2+4acld\beta_0} }{2lc^2}\\
&>\frac{2acl}{2lc^2}=\frac{a}{c}
\end{align}
Note that from Lemma \ref{appendix_1}, it is observed that $x<a/c$. Hence this solution is infeasible.

Secondly, for the other possible solution, $x_2=\frac{2acl+ \beta_0 d - \sqrt{\beta_0^2d^2+4acld\beta_0} }{2lc^2}$, we have
\begin{align}
x_2&=\frac{2acl+ \beta_0 d - \sqrt{\beta_0^2d^2+4acld\beta_0} }{2lc^2} \nonumber \\
&<\frac{2acl}{2lc^2}=\frac{a}{c}
\end{align}
On the other hand, since
$$2acl+ \beta_0 d - \sqrt{\beta_0^2d^2+4acld\beta_0}>0$$ we have $x_2>0$.
Hence the feasibility of $x_2$ is demonstrated.
We then arrive at the unique feasible solution to $x=2^{2R_{\mathrm{CPF}}}$, which is,
\begin{align}
x=\frac{2acl+ \beta_0 d - \sqrt{\beta_0^2d^2+4acld\beta_0} }{2lc^2} \label{eqn:kkt_solution}
\end{align}

Therefore, by replacing $x$ with $2^{2R_{\mathrm{CPF}}}$,
the optimal CPF rate can be obtained accordingly from (\ref{eqn:kkt_solution}).
we hence arrive at Theorem \ref{thm_solution}.

\begin{thm} \label{thm_solution}
The solution to {\bf P2} is given by,
\begin{align}
R^*_{\mathrm{CPF}}&=\frac{1}{2}\log_2 [\frac{2a_mc_ml+\beta_0^*d_m- \sqrt{(\beta_0^*)^2d_m^2+4a_mc_md_m l \beta_0^*} }{2lc_m^2}]^+ \label{eqn:rate_CDF}\\
R^*_{r_i} &= \frac{1}{2}\log_2[\frac{f_i^*\gamma_i^*}{2\beta_i^*\ln 2}g_{r_i}]^+, \quad i=1, \cdots, M \label{eqn:rate_relay_1}
\end{align}
where we set $l=f_{\mathrm{CPF}}^*(1-\sum_{i=1}^M \gamma_i^*)/(2\ln2)$ and assume that $$R^m_{\mathrm{CPF}}=\min_i R^i_{\mathrm{CPF}}  = R^*_{\mathrm{CPF}}.$$
 Note again that the CPF rate of the $m$th relay is assumed to be the common CPF rate of all relays with the specified channel gains\footnote{The procedure to determine the common CPF rate and the specific index of the associated relay is given at the end of this page.}.

The associated allocated power with respect to the specified channel gains are given as follows,
\begin{align}
P^*_{\mathrm{CPF}}&=[\frac{1-b_m \cdot 2^{2R^*_{\mathrm{CPF}}}}{c_m \cdot 2^{2R^*_{\mathrm{CPF}}}-a_m}]^+\label{eqn:power_CDF}\\
P^*_{r_i}&=[\frac{f_i^*\gamma_i^*}{2\beta_i^*\ln 2 }-\frac{1}{g_{r_i}}]^+ \quad i = 1,\cdots, M \label{eqn:power_relay_1}
\end{align}

\end{thm}

Note that in Theorem \ref{thm_solution} we simply assume the $m$th relay enjoys the minimum CPF rate, i.e., the common CPF rate. However, it may not be so in all cases. A procedure to determine the common CPF rate and the index of the associated relay is therefore given as follows.
\begin{enumerate}
\item Initialization: set $m=1$.
\item Assuming the $m$th relay enjoys the minimum CPF rate, compute $R^m_{\mathrm{CPF}}$ by (\ref{eqn:rate_CDF}) and then compute the CPF power required at all relay nodes, i.e., $P^i_{\mathrm{CPF}}$ ($i=1,\ldots, M$) by (\ref{eqn:power_CDF}) with $R^m_{\mathrm{CPF}}$ as the common CPF rate. If $P^m_{\mathrm{CPF}}=\min_i P^i_{\mathrm{CPF}}$, go to 4) and terminates. Otherwise, go to 3).
\item Set $m=m+1$, go to 2).
\item Output the common CPF rate, the common CPF power required as well as the index of the associated $m$th relay.
\end{enumerate}

According to Theorem \ref{thm_solution} and the algorithm above,
a multi-dimensional bisection search method can be implemented to
obtain the optimal solution to {\bf P2}, with the convergence conditions
in (\ref{kkt_3}) and (\ref{kkt_4}). It is noted that, in each
iteration with the estimated $\beta_i$, $\gamma_i$ ($i=0,\ldots,M$),
one can readily obtain $\bar{R}_i$, $\bar{R}_{\mathrm{CPF}}$ and the
associated average power consumption, we then compute $f_i$ and
$f_{\mathrm{CPF}}$ from (\ref{con4_CSIT}) and (\ref{kkt_5}) and
hence $\alpha$ from (\ref{kkt_4}). After getting all parameters, we
check the sign of the left hand side of (\ref{kkt_3}) and
update $\beta_i$, $\gamma_i$
accordingly. Following this procedure, the optimal solution to {\bf P2}
can be obtained.

{It should be noted that the bisection method adopted can guarantee global convergence
at a very slow convergence rate of $1/2$} \cite{Cheney12}. {However, by carefully selecting the
start intervals of the parameters, the adopted bisection algorithm
is shown to be able to converge to the
global optimal point in tens of iterations in simulation, which is acceptable in
practical implementations.}

\section{Simulation}
We now present some simulation results to compare the achievable rates by employing CPF and DF.
Full channel state information is assumed to be available at the associated
transmitters. The average power constraint at all nodes are assumed to be the same in the simulation setting for simplicity. The channels are assumed to be real valued fading channels and their gains are modeled by zero mean and unit gain Gaussian variables. In addition, the noises at all nodes are additive white Gaussian variables with zero mean and unit variance. For ease of computation, we shall mainly focus on a multicast network with two source nodes, two relay nodes and two destination nodes ($L=M=2$), if not specified.

For comparison, here we briefly give a DF protocol for a multi-source multi-relay multicast network with $2M$ phases, which is,
\begin{itemize}
\item In the first $i$th ($i=1,\ldots,M$) phase, the $i$th source node transmits its data to the $i$th relay nodes at rate $R_i$.
\item In the $(M+i)$th phase ($i=1,\ldots,M$), the $i$th relay node broadcasts the data from the $i$th source node to all destination nodes at rate $R_{M+i}$.
\end{itemize}
For delay-stringent applications, the problem to maximize the common multicast throughput within one slot with DF, is formulated as {\bf P3} below.
\begin{align}
&\max_{f_i} \quad \min_i(f_{i}R_{i}) \label{eqn:DF_DS_opt} \\

s.t.\quad &  \nonumber\\

&P_{i} \leq P_{s_i} \quad i=1,\cdots,M\\
&P_i \leq P_{r_i} \quad i=M+1,\cdots,2M\\

&\sum_{i=1}^{2M} f_i \leq 1
\end{align}
The solution to {\bf P3} is similar to that of {\bf P1} and is hence omitted.
In addition, the average throughput for delay-stringent applications can be readily obtained by averaging over the channel distributions.

Similarly, the problem to optimize the averaged common multicast throughput for delay-tolerant applications is formulated as {\bf P4} below.
\begin{align}
\max_{f_i,R_i({\bf h}_l),R_i({\bf g}_i)} \quad f_{i}\bar{R}_{i} \label{eqn:DF_DT_opt}
\end{align}
subject to the following constraints
\begin{align}
&\bar{P}_{l} \leq P_{s_l}\quad l=1,\cdots,M\\
&\bar{P}_{M+i} \leq P_{r_i} \quad i=1,\cdots,M\\
&f_{l}\bar{R}_{l}=f_i\bar{R}_i \quad \forall l,i \\%\label{con4_CSIT}\\
&\sum_{i=1}^{2M} f_i \leq 1
\end{align}
Note that {\bf P4} is a standard convex optimization problem and can be readily solved by KKT conditions. The solution to it is however omitted for brevity.

For clarity, a table describing the associated transmit strategies linked to different applications is presented below.
\begin{table}[h!]
\caption{list of Different Optimization problems with their adopted strategies as well as the related applications}
\centering
\begin{tabular}{|c|c|m{4.7cm}|}
ine
Problem Index &  Strategy   &  Detailed Description  \\
ine
{\bf P1} & CPF-DS & compute-and-forward strategy employed in the delay-stringent application \\
ine
{\bf P2} & CPF-DT & compute-and-forward strategy employed in the delay-tolerant application \\
ine
{\bf P3} & DF-DS & decode-and-forward strategy employed in the delay-stringent application \\
ine
{\bf P4} & DF-DT & decode-and-forward strategy employed in the delay-tolerant application\\
ine
\end{tabular}
\label{Tab1}
\end{table}

Before presenting the performance of the proposed strategies, we would like to show the probability that the network integer vector constructed matrix is not in full rank by employing the proposed methods. It is worth mentioning that, if this matrix is not in full rank, each destination node can not recover all messages from the source nodes. Hence, full rank requirement of the network integer vectors  plays a crucial role in applying CPF strategy.

In Fig.\ \ref{fig:rank}, the probability of rank failure of each method is shown. It is seen that with global optimization, method c) satisfies full rank requirement and the advantage of applying method c) is verified.
Both method a) and method b) have non-zero failure probabilities,
among which method a) has a constant failure probability since its
determination criterion is independent of transmit power. Method b)
is with a decreasing failure probability with the increase of
transmit SNR. Interestingly, with more users at high SNR regime ($K=L=M=4$ at the SNR level over 20dB), it is observed the rank failure probability of method b) is negligible.
On the other hand, it is worth mentioning that using method c) requires additional overhead cost, which would be a critical obstacle for a large-scale multicast network in implementation, as each node needs to exchange some control information on how to construct a global-optimal full network coding system matrix by using CPF.

In Fig.\ \ref{fig:DS_comp}, the optimal averaged common multicast throughput by using CPF strategy as well as DF strategy with delay-stringent constraints are compared over $10000$ randomly generated channel realizations.
It is shown that CPF-DS with the integer network channel coefficient vectors found by the global optimal method, i.e., method c), outperforms that with the local optimal method (method b)) and the naive method (method a)). It is also observed that CPF-DS employing method c) outperforms DF in terms of achievable throughput. However, it is seen that CPF-DS with method a) or b) performs worse than DF strategy, due to their non-negligible rank failure probabilities as shown in Fig.\ \ref{fig:rank}.

In Fig.\ \ref{fig:CPF_DS_time_comp}, the advantage of optimal time allocation is shown for delay-stringent case.
It is observed that the performance of CPF-DS can be greatly improved with optimal time resource allocation. For instance, in the regime of $30$dB transmit power for CPF-DS with method c), an additional $0.2$ bit/s/Hz throughput improvement is achieved  by using optimal time resource allocation, compared with that using equal time-resource allocation.

In Fig.\ \ref{fig:DT_comp}, the optimal common multicast throughput by using different strategies for delay-tolerant applications are shown. It is shown that CPF-DT with method c) outperforms DF strategy in terms of common throughput. For instance, in the regime of 30dB transmit power, employing CPF-DT with method c), over 10\% throughput improvement is achievable, compared with DF-DT. It is also interesting to note that CPF-DT with method b) is slightly better than DF, whereas CPF-DT with method a) is worse than DF,
which is due to the high probability of rank failure and the far-from-optimal
integer coefficient vectors sorted by adopting method a).

Similar to Fig.\ \ref{fig:CPF_DS_time_comp}, the advantage of optimal time allocation is shown for the delay-tolerant application in Fig.\ \ref{fig:CPF_DT_time_comp}. With transmit power at $30$dB, throughput is increased from roughly $0.9$ bit/s/Hz to over $1$ bit/s/Hz for {\bf P2} with optimal time splitting.

In Fig.\ \ref{fig:DS_DT}, the achievable throughput by using CPF strategy for both the delay-stringent and delay-tolerant cases are compared. It is observed that without delay constraints, higher throughput is expected to be achievable. With transmit power at $30$dB, an additional $0.15$ bit/s/Hz throughput improvement is achieved for CPF-DT compared with CPF-DS, where both of them employ method c).

In Fig.\ \ref{fig:DS_DT_less_relay},{ we are interested in the topology
with
arbitrary number of relay nodes ($K=3$, $K=2$ and $K=1$). For the case $K=1$,
the sole relay node needs to decode two
function messages for successful source-message decoding at the destination nodes.
It is observed that with less relay nodes, the optimal
achievable common rate with CPF is decreased, for both the delay-stringent and
delay-tolerant scenarios. This is intuitively due to the reduced
cooperative diversity coming from the decreased number of relays. For the case that
$K=3$, it is observed that the optimal throughput of CPF
is further improved than that with $K=2$, which comes from the increased cooperative diversity.

On the other hand,
it is also observed that with the single relay node ($K=1$), the
CPF strategy
performs slightly worse than DF strategy for delay-tolerant case and roughly
as good as DF for delay-stringent case, due to the
reduced relaying diversity.
Taken into account
that more overhead information is required for CPF strategy,
it is intuitively concluded that
DF is still a good choice for transmission in a small-scale network with
less potential relay nodes
than source nodes.}

In Fig.\ \ref{fig:four_DT},
the optimal common multicast throughput for delay-tolerant
application using different strategies is shown for a four-source,
four-relay and four-user multicast network ($K=M=L=4$).
The performance gain of employing CPF over DF is hence verified for
a larger network. It is observed that CPF employing method b) performs
only slightly worse than that with method c) in terms of throughput,
since method b) in a larger network has a lower rank failure probability.
Hence, it is intuitively learned that, it may be worthwhile to employ
method b) for CPF in large networks in the medium to high SNR regime.
In this way, we can not only achieve close to optimal performance as given
by employing method c), but as well avoid the overhead cost incurred
by employing method c).

\begin{figure}[t]
   \centering
   \includegraphics[width = 7.8cm]{rank.eps}
   \caption{Rank failure probability of CPF employing different methods to
   obtain the integer channel vectors ($K=L=M$).} \label{fig:rank}
   \end{figure}

\begin{figure}[t]
   \centering
   \includegraphics[width = 7.8cm]{DS_comp.eps}
   \caption{Optimal common throughput for delay-stringent case by
   using different strategies ($K=L=M=2$).} \label{fig:DS_comp}
   \end{figure}

\begin{figure}[t]
   \centering
   \includegraphics[width = 7.8cm]{CPF_DS_time_comp.eps}
   \caption{Optimal common throughput for delay-stringent case with
   optimal time allocation and with equal time splitting ($K=L=M=2$). For equal time splitting case, we set $f_{\mathrm{CDF}}=f_1=f_2=1/3$.} \label{fig:CPF_DS_time_comp}
   \end{figure}

   \begin{figure}[t]
   \centering
   \includegraphics[width = 7.8cm]{DT_comp.eps}
   \caption{Optimal common throughput for delay-tolerant case by using
   different strategies ($K=L=M=2$).} \label{fig:DT_comp}
   \end{figure}

   \begin{figure}[t]
   \centering
   \includegraphics[width = 7.8cm]{CPF_DT_time_comp.eps}
   \caption{Optimal common throughput for delay-tolerant case
   with optimal time allocation and with equal time splitting ($K=L=M=2$). For equal time splitting case, we set $f_{\mathrm{CDF}}=f_1=f_2=1/3$.} \label{fig:CPF_DT_time_comp}
   \end{figure}

   \begin{figure}[t]
   \centering
   \includegraphics[width = 7.8cm]{DS_DT.eps}
   \caption{Optimal common throughput for delay-tolerant case
   and delay-stringent case ($K=L=M=2$).} \label{fig:DS_DT}
   \end{figure}

   \begin{figure}[t]
   \centering
   \includegraphics[width = 7.8cm]{DS_DT_less_relay.eps}
   \caption{{Optimal common throughput for delay-tolerant case and
    delay-stringent case for two-source, two-destination network. The
    topologies with three relays ($K=3$), two relays ($K=2$) and one relay
    ($K=1$) are evaluated and compared.} } \label{fig:DS_DT_less_relay}
   \end{figure}

   \begin{figure}[t]
   \centering
   \includegraphics[width = 7.8cm]{four_DT.eps}
   \caption{Optimal common throughput for delay-tolerant case
   in a network consisted of four source nodes, four relay nodes
   and four destination nodes ($K=L=M=4$).} \label{fig:four_DT}
   \end{figure}

\section{Conclusion}
In this work, we considered a multi-source multicast network
with the aid of {an arbitrary number of} relay nodes. {We tried to find the fundamental limit on the maximal common
multicast throughput of all S-D pairs}. To this end, a transmission
protocol employing compute-and-forward strategy was proposed {for an arbitrary number of relays}.
Delay-stringent transmission as well as delay-tolerant transmission
applications were both investigated. The associated optimization
problems were formulated and solved, through the allocation of time
and energy resources. {Various simulation was done for validation of the performance improvement of CPF over the conventional DF in terms of throughput.}
{It was shown that with the increasing number of relay nodes, the CPF strategy
can perform better due to the increased diversity.}
Finally, it was intuitively shown that, using CPF with method b) was a good choice
for large networks in medium to high SNR regime, as it not only provided
performance quite close to CPF with method c), but also avoided the
additional communication cost incurred by using method c).

\section*{Acknowledgments}
This work was partially supported by NSFC grant No. 61171064 and The
National 973 Project of China grant No. 2012CB316102.