text:
\begin{document}
\bibliographystyle{IEEEtran}

\title{Charge-then-Forward: Wireless Powered Communication for Multiuser Relay Networks}
\IEEEoverridecommandlockouts
\author{Mengyu Liu and Yuan Liu,~\IEEEmembership{Senior Member,~IEEE}

\thanks{

The authors are  with the School of Electronic and Information Engineering,
South China University of Technology, Guangzhou 510641, China (e-mail: liu.mengyu@mail.scut.edu.cn,  eeyliu@scut.edu.cn).
}
}

\maketitle

\vspace{-1.5cm}
\begin{abstract}

This paper studies a  relay-assisted wireless powered  communication network (R-WPCN) consisting of multiple source-destination pairs and a hybrid relay node (HRN).

We consider  a ``charge-then-forward" protocol at the HRN, in which the HRN with constant energy supply first acts as an energy transmitter to charge the sources, and then forwards the information from the sources to their destinations through time division multiple access (TDMA) or frequency division multiple access (FDMA).  Processing costs at the wireless-powered sources are taken into account.

Our goal is to maximize the sum-rate of all transmission pairs by jointly optimizing the time, frequency  and power resources.

The formulated optimization problems for both TDMA and FDMA are non-convex.

For the TDMA scheme, by  appropriate transformation, the problem is reformulated as a convex problem and be optimally solved.

For the FDMA case, we find the asymptotically optimal solution in the dual domain.

Furthermore, suboptimal algorithms are proposed for both schemes to tradeoff the complexity and performance.

Finally, the simulation results validate the effectiveness of the proposed schemes.

\end{abstract}

\textbf{\textit{Index Terms---} powered  communication network (WPCN), resource allocation, cooperative relay.}

\section{Introduction}

The rapid growth of high-speed data and multimedia services  increases energy consumption for better quality-of-services. However, conventional battery-powered communications have to replace or recharge batteries manually to extend their lifetime, which is inconvenient, unsafe, and costly. Recently, radio-frequency (RF) signal enabled   wireless power transfer (WPT) has drawn great attention as it essentially provides more cost-effective and  green energy supplies for wireless devices, where RF signals are used as the carriers to convey wireless energy to low-power wireless devices.

There are two main directions of WPT among  the current related researches. One line of WPT focuses on so-called simultaneous wireless information and power transfer (SWIPT), where the  same RF signal carries both energy and information at the same time \cite{MIMO_zhangrui,SWIPT_morden}.

Due to the practical limitation of receivers that the received signals cannot be used to perform energy harvesting  and information decoding  simultaneously, two practical receiver  architectures, namely time switching (TS) and power splitting (PS), were proposed in \cite{MIMO_zhangrui}. For TS, the received signal is either used for energy harvesting or information decoding, whereas for PS, the received signal is split into two separate streams with one stream for energy harvesting and the other for information decoding at the same time.

SWIPT has been investigated extensively in different systems, e.g., the fading channels \cite{fading_channel}, relay channels \cite{relay_channel_2,relay_channel_3,12-relay,15-ts-relay} and orthogonal frequency division multiple access (OFDMA) channels \cite{OFDMA_channels_1,OFDMA_channels_2,OFDMA_channels_4,me,R1_x1}.

The newly emerging wireless powered communication network (WPCN) is another line of WPT where ambient RF signals are used to power wireless devices \cite{opportunitiesWPCN}.
There are two basic applications about WPCN.

One is that energy transmitters and information access points (APs) are located separately where energy transmitters transmit energy to wireless devices and then wireless devices transmit their information to APs using their harvested energy from energy transmitters \cite{7-WPCN}.

Another application is that a hybrid AP (HAP) performs the roles of energy transmitter and AP integrally.

For instance, a ``harvest-then-transmit'' protocol was proposed in \cite{1-wpcn-he}, where HAP first broadcasts wireless energy to all wireless devices in the downlink and then wireless devices utilize the harvested energy to transmit their independent information to HAP in the uplink based on TDMA.

Different from the HAP in SWIPT  that coordinates wireless energy and information, the HAP in WPCN only broadcasts wireless energy.

An important application for WPCN lies in relay-assisted  WPCN (R-WPCN),  where relays are used to assist information transmission in R-WPCN.

There are two categories among the current related works about R-WPCN: one is source powering relay \cite{R1_x2,4-user-co-relay,3-relay-full} and the second is relay powering source \cite{R3_x1,kwan_relaycapacity }.

For the first category, i.e., source powering relay,

the authors in \cite{R1_x2} proposed a ``harvest-then-cooperate" protocol, where both source and relay can harvest energy from the RF signals from a base-station.

A  two-user R-WPCN was studied in \cite{4-user-co-relay}  where a nearer user to HAP harvests energy sent by HAP and relays information of the farther user in half-duplex.

In \cite{3-relay-full}, the full-duplex relay not only is powered by the source but also harvests energy from itself by energy recycling.

As for the category of relay powering source,  the throughput maximization problem was investigated in \cite{R3_x1}, where the source can  harvest energy from the access point and/or relay before information transmission.

The authors in \cite{kwan_relaycapacity} studied the channel capacity subject to an additional energy transmission cost at the energy harvesting sources.

Note that both \cite{R3_x1} and \cite{kwan_relaycapacity} considered a single source-destination pair.

\begin{figure}[t]
\begin{centering}
\centering
\includegraphics[scale = 0.3]{Figure1.eps}
\vspace{-0.1cm}
\caption{System model of the considered multiuser R-WPCN. }\label{fig:ff1}
\end{centering}
\vspace{-0.1cm}
\end{figure}

In this paper, we consider a new R-WPCN consisting of multiple source-destination pairs assisted by a single hybrid relay node (HRN), as shown in Fig. \ref{fig:ff1}.  We assume that the HRN in this paper has constant energy supply, while the sources nodes have no embedded power supply so that they have to be powered by the HRN before information transmission, i.e., ``harvest-then-transmit" protocol is applied at the sources. The HRN thus acts double roles, one for an energy transmitter and the other for an information helper. That is, the HRN first charges the sources and then forwards their information, i.e., ``charge-then-forward" protocol is considered at the HRN.

As the considered HRN has double roles, i.e., energy transmitter and information helper, energy charging and information forwarding of the HRN are mutually influenced and restricted since the HRN's total energy is fixed. That is, encouraging the energy charging will increase the transmit power of sources at the first hop but decrease the information forwarding at the second hop. How to find the optimal tradeoff that maximizes the system sum-rate is non-trivial. In addition, based on TDMA and FDMA for multiuser information transmission, the network resources, like time, power, and frequency, are highly coupled and the formulated optimization problems are non-convex and difficult to solve. Moreover, we consider the processing cost at the wireless-powered  sources, which further complicates the problems.

The main contributions of this paper are summarized as follows:
\begin{itemize}
  \item  We consider a new multiuser R-WPCN based on a ``charge-then-forward" relaying protocol, where the HRN first powers the energy-free sources and then forwards the information from the sources to their destinations by TDMA and  FDMA.  Processing cost is considered at the wireless-powered sources.

  \item  Depending on whether TDMA or FDMA  is adopted, we formulate two optimization problems respectively for sum-rate maximization, which are both non-convex. We propose efficient algorithms to find the  optimal solutions. In addition, suboptimal algorithms are proposed for both schemes to tradeoff the complexity and performance.

  \item  We provide some useful insights into the R-WPCN. For example,  the time of WPT should be as small as possible so that the time for wireless information transmission (WIT) can be maximized for sum-rate maximization. In addition,  due to the doubly distance-dependent signal attenuation for both WPT and the first hop of WIT, it is shown that the sum-rate decreases when the HRN moves from the sources to the destinations.

\end{itemize}

The remainder of this paper is organized as follows. In Section II, we introduce the system model  of multiuser R-WPCN and problem formulations based on TDMA and FDMA, respectively.

Section III presents the optimal and suboptimal resource allocation algorithms for the TDMA based problem. In  the next, the  asymptotically optimal and suboptimal algorithms for the FDMA based problem are presented in Section IV.

In Section V, we evaluate the performance of proposed algorithms by simulations. Finally, Section VI concludes the paper.

\section{System Model and Problem Formulation}

As shown in Fig. \ref{fig:ff1}, we consider a general two-hop R-WPCN with multiple source-destination pairs as well as a HRN which not only transfers energy to the sources but also forwards information from the sources to the destinations.  All nodes are equipped with a single antenna. We assume that the HRN is half-duplex due to the practical consideration, and there is no direct link between each source-destination pair due to the shielding effect caused by obstacles. As a result, each pair needs the assistance  of the HRN to forward information.

In this paper, we consider that the HRN has a constant energy supply while the source nodes have no embedded energy and thus have to harvest energy for information transmission. In addition, we assume that each source has the energy harvesting function to store the energy.

In particular, we consider the ``charge-then-forward" relaying protocol to coordinate power and information transfer, in which the HRN first acts as a wireless power beacon to charge the sources then as a helper for forwarding their information.

Specifically, the whole transmission is divided into two continuous phases. The first phase is used for WPT conducted by the HRN. The second phase is WIT, i.e.,  the sources use the harvested energy to transmit their independent information to their destinations via the assistance of the HRN in the second phase based on TDMA or FDMA. The sources do not store the harvested energy for future, i.e., all the energy harvested during WPT phase is used for WIT.

The global channel state information (CSI) of the network is assumed to be known at the HRN where the central processing task is embedded.

RF power transfer crucially depends on the available CSI of the nodes, which needs additional resources to acquire and the straightforward way is channel estimation via pilot signals, similar to conventional wireless communication systems. In our R-WPCN, whether the sources transmit pilots and the HRN estimates CSI, or the HRN transmits pilots and the sources estimate CSI, the sources are required to have initial energy to transmit/decode the pilot signals at the beginning of the training phase (before WPT in transmission phase). Thus it is reasonable to assume that the wireless-powered sources reserve some circuit power at the beginning for channel estimation, since the energy used to channel estimation is much smaller than that of information transmission in practice. CSI acquisition in WPT systems is very important but seems to be beyond of the scope of this paper.

In this paper, we consider a block fading wireless environment so that the channel impulse response can be treated as time invariant in the block duration. As a result, the channel gains within the block duration remain unchanged  (but can vary in different block durations).

For convenience, we assume that the transmission time of each block is normalized to be unit.

\subsection{TDMA Case}

 \begin{figure}
  \centering
  \subfigure[TDMA]{
    \label{fa}
    \includegraphics[width=0.45\textwidth]{Figure2a.eps}}
  \hspace{1in}
  \subfigure[FDMA]{
    \label{fb}
    \includegraphics[width=0.45\textwidth]{Figure2b.eps}}
  \caption{The ``charge-then-forward" relaying protocol based on (a) TDMA and (b) FDMA information transmission.}
  \label{fig:subfig}
\end{figure}

We first consider the case of TDMA-based information transmission as shown in Fig. 2(a). The total transmission time is divided into $K+1$ time slots where the first time slot, say slot 0 with time duration $\alpha_0$, is allocated for WPT and all source nodes harvest energy from the HRN,  while the rest $K$ slots are assigned to WIT of the $K$ pairs where each pair $k$ is allocated $\alpha_k$ time duration. Moreover, for each pair's information transmission, $\alpha_{k}$ is further divided into two equal sub-slots with $\alpha_k/2$ for the first hop and the rest $\alpha_k/2$ for the second hop.

By normalizing the whole time to be unit, we have
\begin{equation}
\sum\limits_{k=0}^K\alpha_{k}\leq1,\quad 0\leq\alpha_{k}\leq1, \quad k=0,1,\cdots,K.  \label{alpha_k_TDMA}
\end{equation}
In addition, the power of the HRN used at slot $k$ is denoted as $p_{k}$. Besides, we consider that there is a peak power constraint on $p_k$, i.e., $0 \leq p_k \leq P_{\rm{peak}}.$ Denote the maximum transmit power of the HRN as $P$,

then the energy constraint at the HRN is given by

\begin{equation}
\alpha_{0}p_{0}+\sum\limits_{k=1}^K\frac{\alpha_{k}}{2}p_k\leq P,  \label{power_C_TDMA}
\end{equation}

where  $\alpha_{k}/2$ is the transmission time of the HRN in the second hop for each pair $k$.  Note that the terms of power and energy are interchangeably used here since the duration of each block is normalized to be $T=1$ unit.

We consider energy accumulation for TDMA case, i.e., source $k$ harvests and accumulates energy from the previous slots, i.e., slot 0 to slot $k-1$.

The channel power gain from HRN to source $k$ for WPT and the channel power gain from source $i$ to source $k$ are denoted by $g_{{\rm r},k}$ and $g_{i,k}$, respectively.

Denote source $i$'s transmit power at slot $i$ as $q_i$.

Then the harvested energy of source $k$ can be expressed as

\begin{align}  \label{E_k_TDMA}
E_k = \left\{
        \begin{array}{ll}
          \eta\alpha_0p_0g_{{\rm r},k}, &\hbox{$k=1$} \\
          \eta \bigg( \alpha_0 p_0 g_{{\rm r},k} + \sum\limits_{i=1}^{k-1} \frac{\alpha_i }{2}  p_i g_{{\rm r},k} \\ \quad\quad\quad\quad + \sum\limits_{i=1}^{k-1} \frac{\alpha_i }{2} q_i g_{i,k}  \bigg), &\hbox{$k=2,\cdots,K,$}
        \end{array}
      \right.
\end{align}

which comprises three parts: the first term is the energy harvested in WPT phase, and the last two terms correspond to the energy harvested from the HRN and sources in the previous WIT phase. Here $0<\eta<1$ is the energy conversion efficiency at the sources.

As a result, the energy causality constraint at source $k$ is given by
\begin{align}
\frac{\alpha_k}{2} q_k + E_k^c \leq  E_k, \quad  k=1,\cdots,K, \label{C_energysource_TDMA}
\end{align}
where $\alpha_k/2$ is the time of the first hop for the information transmission of pair $k$ and $E_k^c$ is the non-zero energy processing cost at source $k$.

Moreover, the additional White Gaussian noise (AWGN) at  each  node is modeled as circularly symmetric complex Gaussian (CSCG) random variable with zero mean and variance $\sigma^{2}$. Denote the channel power gains for the first and second hops of WIT for pair $k$ as $h_{1,k}$ and $h_{2,k}$, respectively. Using decode-and-forward (DF) relaying strategy, the achievable rate for each pair $k$ is given by

\begin{align}  \label{R_k_TDMA}
R_{k}=
\frac{\alpha_{k}}{2}\min\left\{\log_{2}\left(1+\frac{q_k h_{1,k}}{\sigma^{2}}\right),\log_{2}\left(1+\frac{p_{k} h_{2,k}}{\sigma^{2}}\right)\right\}.
\end{align}

Our objective is to maximize the sum-rate of all pairs by jointly optimizing the time allocation, the transmit power of sources and HRN. Let  $\boldsymbol{p} = \{p_k\} $, $\boldsymbol{q} = \{q_k\} $,
$\boldsymbol{\alpha} = \{\alpha_k\} $ and $\boldsymbol{R} = \{R_k\} $,

the problem can be mathematically formulated as

\begin{subequations}
\begin{align}
{\rm (P1):}
\max_{\{\boldsymbol{p},\boldsymbol{q},\boldsymbol{\alpha},\boldsymbol{R}\}}\quad\quad&\sum\limits_{k=1}^K R_{k} \label{TDMA_c1}  \\
{\rm s.t.}\quad\quad&\eqref{alpha_k_TDMA}, \eqref{power_C_TDMA}, \eqref{C_energysource_TDMA},\eqref{R_k_TDMA}\nonumber \\
\quad\quad&0\leq p_{k}\leq P_{\rm{peak}},  k=0,1,\cdots,K.  \label{tdma_peak_c}
\end{align}
\end{subequations}

Problem (P1) is non-convex since the rate expression \eqref{R_k_TDMA} is not jointly concave in the variables. We will optimally solve this problem in Section III.

\subsection{FDMA Case}

We also consider the case of FDMA-based information transmission for the multiple pairs as shown in Fig. 2(b), where the total time is divided into two time slots, i.e., slot 0 and slot 1 utilized for WPT (using energy signals) and WIT, respectively.

The time duration of slot 0 and slot 1 are denoted by $\alpha_{0}$ and $\alpha_{1}$ with

\begin{equation}
\alpha_{0}+\alpha_{1}\leq1,\quad 0\leq \alpha_{0},\alpha_1\leq 1.  \label{x_k,n_OFDMA}
\end{equation}

We assume that the HRN broadcasts energy signals over the entire bandwidth in the phase of WPT, while information signals are conveyed by using FDMA over $N$ subcarriers (SCs) in the next phase of WIT.

For information transmission, we define a binary SC allocation variable $x_{k,n}$ with $x_{k,n}=1$ representing that SC $n$ is allocated to pair $k$ for WIT and $x_{k,n}=0$ otherwise. Each SC is allocated  to at most one pair at slot 1 for WIT to avoid interference.   The  constraint can be expressed as

\begin{eqnarray}
\sum\limits_{k=1}^Kx_{k,n}\leq 1, \forall n,\quad x_{k,n} \in \{0,1\},\quad \forall n, k=1,\cdots,K. \label{x_k,n}
\end{eqnarray}

The channel power gains for WPT of source $k$, the first and second hops of pair $k$ over SC $n$ for WIT are denoted as $g_{{\rm r},k}$, $h_{1,k,n}$ and $h_{2,k,n}$, respectively.

The transmit power of HRN for WPT at slot 0 is denoted as $p_{0}$, and the power of the HRN for forwarding pair $k$'s information on SC $n$ at slot 1 is $p_{k,n}$.

Note that WPT is conducted over the entire bandwidth, thus there is no index $n$ for $g_{{\rm r},k}$ and $p_0$.

The total transmit energy constraint of the HRN is thus given by

\begin{equation}
\alpha_{0}p_{0}+\sum\limits_{n=1}^N\sum\limits_{k=1}^K\frac{\alpha_{1}}{2} p_{k,n} \leq P,  \label{FDMA_P}
\end{equation}

where $\alpha_{1}/2$ represents the transmission time of  the HRN in the second hop.

Moreover, we define  source $k$'s transmit power on SC $n$ at slot 1 as $q_{k,n}$. Different from TDMA case, since all sources transmit their information at the same time in slot 1, the harvested energy of sources are only from HRN during WPT phase. Therefore, the energy  constraint at source $k$ is given by

\begin{equation}
\sum\limits_{n=1}^N\frac{\alpha_{1}}{2} q_{k,n} + E_k^c \leq \eta \alpha_0 p_0 g_{{\rm r},k},  \quad  k=1,\cdots,K,  \label{energy_C_OFDMA}
\end{equation}

where the $\alpha_{1}/2$ represents the time of the first hop during information transmission.

The achievable rates of the first and second hops for pair $k$ over SC $n$ can be respectively written as:
\begin{align}
R_{1,k,n} &= \frac{\alpha_{1}}{2N}\log_{2}\left(1+\frac{q_{k,n}h_{1,k,n}}{\sigma^{2}}\right),\forall n,k,\\
R_{2,k,n} &= \frac{\alpha_{1}}{2N} \log_{2}\left(1+\frac{p_{k,n}h_{2,k,n}}{\sigma^{2}}\right), \forall n,k.
\end{align}
The achievable rate of pair $k$ by using DF relaying strategy is the minimum of the rates achieved in the two hops, which can be expressed as

\begin{eqnarray} \label{FDMA_minrate}
R_{k,n} = \min\left\{ R_{1,k,n} ,R_{2,k,n}\right\},\quad \forall n,k=1,\cdots,K.
\end{eqnarray}

Our goal is  maximizing the sum-rate of all transmission pairs by varying the transmit power of the sources and HRN,  SC assignment and time allocation.

Let $\boldsymbol{p} = \{p_0,p_{k,n}\}$, $\boldsymbol{q}= \{q_{k,n}\}$, $\boldsymbol{x}=\{x_{k,n}\}$, $\boldsymbol{\alpha}=\{\alpha_0,\alpha_1\}$ and $\boldsymbol{R}=\{R_{k,n}\}$,  the optimization problem can be mathematically formulated as

\begin{subequations}
\begin{align}
&{\rm (P2):}

\max_{\{\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{R}\}}\quad\sum\limits_{k=1}^K\sum\limits_{n=1}^Nx_{k,n}R_{k,n} \label{P2}\\
&{\rm s.t.}\quad\eqref{x_k,n_OFDMA}, \eqref{x_k,n}, \eqref{FDMA_P}, \eqref{energy_C_OFDMA}\nonumber\\
&\quad\quad R_{k,n}\leq R_{1,k,n},R_{k,n}\leq R_{2,k,n}, \forall n,k=1\cdots,K,\label{r_k,n_10}\\
& \quad\quad 0\leq p_{0}, p_{k,n}\leq P_{\rm{peak}},  \forall n, k=1,\cdots,K. \label{peak_FDMA}
\end{align}
\end{subequations}

Problem (P2) is also non-convex  since both binary and continuous variables are involved, which is a mixed-integer programming problem. The asymptotically optimal solution for Problem (P2) will be obtained in Section IV.

\section{Resource Allocation in TDMA Case}

In this section, we study the TDMA case by solving Problem (P1). Problem (P1) is not convex and cannot be solved in its original form. Therefore, to make the problem tractable, we introduce a set of new variables $\boldsymbol{m} = \{\alpha_k q_k /2\}$ and $\boldsymbol{s}=\{\alpha_0p_0,\alpha_k p_k /2\}$. Clearly, $\boldsymbol{m}$ and $\boldsymbol{s}$ can be viewed as the actual transmit energy of the sources and HRN, respectively.

Problem (P1) is equivalent to the following problem:
\begin{subequations}
\begin{align}
&{\rm (P1'):} \label{TDMA_problem}
\max_{\{\boldsymbol{\alpha},\boldsymbol{m},\boldsymbol{s},\boldsymbol{R}\}}\quad\sum_{k=1}^K R_k  \\
&{\rm s.t.}\sum\limits_{k=0}^K\alpha_{k}\leq1,\quad 0\leq\alpha_{k}\leq1, \quad k=0,1,\cdots,K,  \label{alpha_s_TDMA}\\
&\quad \sum_{k=0}^{K} s_k \leq P, \label{energy_Q}\\
&\quad R_{k} \leq R_{1,k},  R_{k} \leq R_{2,k}, \quad  k=1,\cdots,K,\label{C_rate2_TDMA}\\
& \quad m_k + E_k^c \leq  \eta s_0 g_{{\rm r},k}, \quad k=1, \label{energy_source_TDMAk1}\\
& \quad m_k + E_k^c \leq  \eta \left( \sum_{i=0}^{k-1} s_i g_{{\rm r},k} + \sum_{i=1}^{k-1} m_i g_{i,k}  \right), k=2,\cdots,K, \label{energy_source_TDMA} \\
& \quad 0 \leq s_0 \leq \alpha_0 P_{\rm peak},  0 \leq s_k \leq \frac{\alpha_k}{2} P_{\rm peak},  k=1,\cdots,K. \label{alpha_k_peak}
\end{align}
\end{subequations}
where
\begin{align}
R_{1,k} &= \frac{\alpha_k}{2} \log_2 \left(1 + \frac{2 m_k h_{1,k}}{\alpha_k \sigma^2} \right),\quad  k=1,\cdots,K,\\
R_{2,k} &= \frac{\alpha_k}{2} \log_2 \left(1 + \frac{2 s_k h_{2,k}}{\alpha_k \sigma^2} \right),\quad  k=1,\cdots,K.
\end{align}

Since constraint \eqref{C_rate2_TDMA} is convex and the other constraints of Problem (P1$'$) are affine, Problem (P1$'$) is convex in its current form.

In the literature \cite{first_order_1,first_order_2,first_order_3,first_order_4}, the first-order method can be used to solve these non-convex problems by approximating the non-convex objective functions and constraints into convex ones. However, in this paper, by appropriate variable transformation, Problem (P1) is reformulated to be convex, which can thus be optimally solved by applying the Lagrange duality method, as will be shown next.

We first introduce non-negative Lagrangian multipliers $\boldsymbol{\lambda} = \{\lambda_k\}\succeq0$ and $\boldsymbol{\beta}=\{\beta_k\}\succeq 0$ associated with the rate constraint \eqref{C_rate2_TDMA}, $\boldsymbol{\nu}=\{\nu_k\}\succeq0$ associated with the energy  causality constraints \eqref{energy_source_TDMAk1} and \eqref{energy_source_TDMA}. In addition, non-negative Lagrangian multipliers $\mu\geq 0$ and $\xi\geq 0$ are associated with the total time constraint \eqref{alpha_s_TDMA} and total energy constraint at HRN \eqref{energy_Q}.

Then, the Lagrangian of Problem (P1$'$) is given by
\begin{align}
&\mathcal{L}\left(\boldsymbol{s},\boldsymbol{m},\boldsymbol{\alpha},\boldsymbol{R},\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi\right) \nonumber \\
= & \sum_{k=1}^{K} \left[ R_k +   \lambda_k \left( R_{1,k} - R_k \right) +  \beta_k \left( R_{2,k} - R_k \right)\right]\nonumber \\
& +\mu \left( 1-\sum_{k=0}^{K} \alpha_k \right) + \xi \left( P - \sum_{k=0}^{K} s_k \right) \nonumber \\
 & + \nu_1(\eta s_0 g_{{\rm r},1}-E_1^c-m_1) \nonumber \\
 & + \sum_{i=2}^{K}\nu_i \left(\sum_{k=0}^{i-1} \eta s_k g_{{\rm r},i} + \sum_{k=1}^{i-1} \eta m_k g_{k,i} - E_i^c - m_i \right). \label{L_TDMA_R}
\end{align}

Denote $\mathcal{D}$ as the set of $\{\boldsymbol{s},\boldsymbol{m},\boldsymbol{\alpha},\boldsymbol{R}\}$ satisfying the primary constraints, then the dual function of Problem (P1$'$) is given by
\begin{align}
g(\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi) = \max_{\{\boldsymbol{s},\boldsymbol{m},\boldsymbol{\alpha},\boldsymbol{R}\} \in \mathcal{D}} \mathcal{L}\left(\boldsymbol{s},\boldsymbol{m},\boldsymbol{\alpha},\boldsymbol{R},\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi\right). \label{dual_function_TDMA}
\end{align}

To compute the dual function $g(\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi)$, we need to find the optimal $\{\boldsymbol{s}^*,\boldsymbol{m}^*,\boldsymbol{\alpha}^*,\boldsymbol{R}^*\}$ to maximize the Lagrangian under the
given dual variables $\{\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi\}$.  In the following we present the derivations in detail.

\subsection{Optimizing $\{\boldsymbol{s},\boldsymbol{m},\boldsymbol{\alpha},\boldsymbol{R}\}$ for Given  $\{\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi\}$}

\subsubsection{Maximizing Lagrangian over $\{R_{k}\}$} The part of the dual function with respect to the rate variable $\{R_{k}\}$ is given by
\begin{align}
g_R(\boldsymbol{\lambda},\boldsymbol{\beta}) = \max_{\boldsymbol{R}\succeq 0} \sum_{k=1}^{K} (1-\lambda_k-\beta_k)R_k.
\end{align}
To make sure that the dual function is bounded, we have $1-\lambda_k-\mu_k \equiv 0$. In such case, $g_R(\boldsymbol{\lambda},\boldsymbol{\beta})\equiv0$ \cite{Mojianhua} and we obtain that $\beta_k = 1-\lambda_k$. Note that $0\leq \lambda_k \leq 1$ such that $\beta_k$ is non-negative. By substituting these results above into \eqref{L_TDMA_R}, the Lagrangian can be rewritten as:
\begin{align}
&\mathcal{L}\left(\boldsymbol{s},\boldsymbol{m},\boldsymbol{\alpha},\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\right) \nonumber \\
= & \sum_{k=1}^{K}  [\lambda_kR_{1,k} + (1-\lambda_{k})R_{2,k}]+ \mu \left( 1-\sum_{k=0}^{K} \alpha_k \right) \nonumber \\
&+ \xi \left( P - \sum_{k=0}^{K} s_k \right)  + \nu_1(\eta s_0 g_{{\rm r},1}-E_1^c-m_1) \nonumber \\
&+ \sum_{i=2}^{K}\nu_i \left(\sum_{k=0}^{i-1} \eta s_k g_{{\rm r},i} + \sum_{k=1}^{i-1} \eta m_k g_{k,i} - E_i^c - m_i \right). \label{L_TDMA_NoR}
\end{align}

\subsubsection{Maximizing Lagrangian over $\{m_k\}_{k=1}^{K},$ $\{s_k\}_{k=1}^{K}$ and $\{\alpha_k\}_{k=1}^{K}$} Observing the Lagrangian in \eqref{L_TDMA_NoR}, we find that the dual function in \eqref{dual_function_TDMA} can be decomposed into $K+1$ independent functions:
\begin{flalign}
g(\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi) = \sum_{k=0}^{K} g_k (\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi) + \mu + \xi P - \sum_{i=1}^{K} \nu_i E_i^c, \label{dual_zuhe}
\end{flalign}
where
\begin{align} \label{g_k_TDMA}
g_k (\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi) \triangleq \max_{\{\boldsymbol{s},\boldsymbol{m},\boldsymbol{\alpha}\} \in \mathcal{D}}
\mathcal{L}_k \left(\boldsymbol{s},\boldsymbol{m},\boldsymbol{\alpha},\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\right)
\end{align}
with
\begin{align} \label{L_k_TDMA}
& \mathcal{L}_k \left(\boldsymbol{s},\boldsymbol{m},\boldsymbol{\alpha},\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\right) = \nonumber\\
& \left\{
  \begin{array}{ll}
    -\mu_k \alpha_k - \xi s_k + \sum\limits_{i=1}^{K} \eta \nu_i g_{{\rm r},i} s_k,  \quad\quad\quad\quad\quad \quad \quad \quad \hbox{$k=0$,} \\
\lambda_k R_{1,k} + (1-\lambda_k) R_{2,k} \\ - \mu \alpha_k - \xi s_k - \nu_k m_k \\ + \sum\limits_{i=k+1}^{K} \eta \nu_i ( g_{{\rm r},i} s_k +  g_{k,i} m_k) , \quad\quad\quad \hbox{$k=1,\cdots,K-1$,} \\
    \lambda_k R_{1,k} + (1-\lambda_k) R_{2,k} - \mu \alpha_k - \xi s_k - \nu_k m_k, \quad \hbox{$k=K$.}
  \end{array}
\right.
\end{align}

For given dual point $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$, maximizing \eqref{L_TDMA_NoR} over $\{m_k\}_{k=1}^{K},\{s_k\}_{k=1}^{K},\{\alpha_k\}_{k=1}^{K}$ is equivalent to solving  \eqref{g_k_TDMA} for $k=1,\cdots,K.$  From \eqref{L_k_TDMA}, the partial derivatives of $\mathcal{L}_k$ with respect to $s_k$ and $m_k$ can be given by \eqref{piandaos} and \eqref{piandaom} on the top of the next page.
\begin{figure*}[ht]
\begin{align}
\frac{\partial \mathcal{L}_k}{\partial s_k} &=\left\{
                                               \begin{array}{ll}
                                                 \frac{(1-\lambda_k)\alpha_k h_{2,k}}{(\alpha_k \sigma^2 + 2s_k h_{2,k})\ln2} + \sum\limits_{i=k+1}^{K}\eta \nu_i g_{{\rm r},i} - \xi, & \hbox{$k=1,\cdots,K-1$} \\
                                                 \frac{(1-\lambda_k)\alpha_k h_{2,k}}{(\alpha_k \sigma^2 + 2s_k h_{2,k})\ln2}  - \xi, & \hbox{$k=K$}
                                               \end{array}
                                             \right. \label{piandaos}\\
\frac{\partial \mathcal{L}_k}{\partial m_k} &= \left\{
                                                \begin{array}{ll}
                                                  \frac{\lambda_k\alpha_k h_{1,k}}{(\alpha_k \sigma^2 + 2m_k h_{1,k})\ln2} + \sum\limits_{i=k+1}^{K}\eta \nu_i g_{k,i} - \nu_k, & \hbox{$k=1,\cdots,K-1$} \\
                                                  \frac{\lambda_k\alpha_k h_{1,k}}{(\alpha_k \sigma^2 + 2m_k h_{1,k})\ln2} - \nu_k, & \hbox{$k=K$}
                                                \end{array}
                                              \right. \label{piandaom}
\end{align}
\hrulefill
\end{figure*}
Given $\alpha_k,k=1,\cdots,K$, the optimal energy variables $s_k$ and $m_k$ that maximize $\mathcal{L}_k$ can be obtained by setting $\frac{\partial \mathcal{L}_k}{\partial s_k}=0$ and $\frac{\partial \mathcal{L}_k}{\partial m_k} =0$ and are given by \eqref{s_k_TDMA} and \eqref{m_k_TDMA} on the top of the next page.
\begin{figure*}[ht]
\begin{align}
s_k &= \left\{
         \begin{array}{ll}
            \frac{\alpha_k}{2} \min \left\{ \left[ \frac{1-\lambda_k}{\left(\xi - \sum_{i=k+1}^{K} \eta \nu_i g_{{\rm r},i}\right)\ln2} - \frac{\sigma^2}{h_{2,k}} \right]^+,P_{\rm peak} \right\}, & \hbox{$k=1,\cdots,K-1$} \label{s_k_TDMA} \\
           \frac{\alpha_k}{2} \min \left\{ \left( \frac{1-\lambda_k}{\xi \ln2} - \frac{\sigma^2}{h_{2,k}} \right)^+,P_{\rm peak} \right\}, & \hbox{$k=K$}
         \end{array}
       \right.\\
m_k &= \left\{
         \begin{array}{ll}
           \frac{\alpha_k}{2} \left[ \frac{\lambda_k}{\left( \nu_k - \sum_{i=k+1}^{K} \eta\nu_ig_{k,i} \right)\ln2} - \frac{\sigma^2}{h_{1,k}}, \right]^+, & \hbox{$k=1,\cdots,K-1$}  \label{m_k_TDMA}\\
           \frac{\alpha_k}{2} \left( \frac{\lambda_k}{ \nu_k \ln2} - \frac{\sigma^2}{h_{1,k}}, \right)^+, & \hbox{$k=K$}
         \end{array}
       \right.
\end{align}
\hrulefill
\end{figure*}

With given $s_k$ and $m_k$,

we can easily prove that $\frac{\partial \mathcal{L}_k}{\partial \alpha_k}$ is a decreasing function of $\alpha_k$. As a result, the optimal $\alpha_k$ with given $s_k$ and $m_k$ can be found by a simple bisection search over $0\leq \alpha_k \leq 1.$

To summarize, for $k=1,\cdots,K$, Problem \eqref{g_k_TDMA} can be solved by iteratively optimizing between $\{s_k,m_k\}$ and $\alpha_k$ with one of them fixed at one time, which is known as block-coordinate descent (BCD) method.

\subsubsection{Maximizing Lagrangian over $s_0$ and $\alpha_0$} Next, we study the solution of Problem \eqref{g_k_TDMA} for $k=0$, which is a linear programming problem (LP). From \eqref{L_k_TDMA}, to maximize $\mathcal{L}_0$ we have
\begin{align}
s_0 &= \left\{
        \begin{array}{ll}
          \alpha_0 P_{\rm peak}, & \hbox{if $-\xi + \sum\limits_{i=0}^{K} \eta \nu_i g_{{\rm r},i} > 0$,} \label{s0_TDMA}\\
          0, & \hbox{otherwise.}
        \end{array}
      \right. \\

\alpha_0 &= \left\{
        \begin{array}{ll}
          1, & \hbox{if $-\mu -\xi P_{\rm peak} + \sum\limits_{i=0}^{K} \eta \nu_i g_{{\rm r},i} P_{\rm peak}> 0$,}  \label{alpha0_TDMA}\\
          0, & \hbox{otherwise.}
        \end{array}
      \right.
\end{align}

\subsection{Optimizing Dual Variables $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$}

As a dual function is always convex \cite{convexbook}, we adopt the ellipsoid method to simultaneously iterate the dual variables $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$ to the optimal ones by using the defined subgradients as follows:

\begin{equation}     \label{citidu}
\boldsymbol{\triangle}=\left[
  \begin{array}{lcl}
    \Delta \lambda_{k} = R_{1,k}-R_{2,k}, k=1,\cdots,K\\
    \Delta \mu  = 1-\sum\limits_{k=0}^{K}\alpha_k\\
\Delta \xi =P-\sum\limits_{k=0}^{K}s_k\\
\Delta \nu_1 = \eta s_0 g_{{\rm r},1} - E_1^c - m_1\\
\Delta \nu_k = \sum\limits_{i=0}^{k-1}\eta s_i g_{{\rm r},k} - \sum\limits_{i=1}^{k-1}\eta m_i g_{i,k}\\ \quad\quad\quad- E_k^c - m_k, k=2,\cdots,K
  \end{array}
\right]
\end{equation}

\subsection{Discussion on Optimality and Complexity}

The optimal $s_k^*$, $m_k^*$ and $\alpha_k^*$ for $k=1,\cdots,K$ are obtained at optimal $\{\boldsymbol{\lambda}^*,\boldsymbol{\nu}^*,\mu^*,\xi^*\}$, then the optimal $\alpha_0^*$ is given by $\alpha_0^* = 1-\sum_{k=1}^{K}\alpha_k^*$.  With $\{\alpha_k^*\}_{k=0}^{K}$, $\{s_k^*\}_{k=1}^{K}$ and $\{m_k^*\}_{k=1}^{K}$, Problem (P1$'$) becomes a LP with variable $s_0$. The optimal value of $s_0^*$ is obtained by solving this LP.

To summarize, the algorithm to solve Problem (P1$'$) is given in Algorithm 1. The time complexity of steps 3-7 is of order $K^2$. The complexity of step 9 is $\mathcal{O}(K^2)$. Therefore, the complexity of steps 3-9 is given by $\mathcal{O}(K^2)$. Note that step 10 iterates $\mathcal{O}(q^2)$ to converge, where $q$ is the number of dual variables and $q=2K+2$ in our case.
Thus the complexity of steps 1-10 is $\mathcal{O}(q^2K^2)$. The time complexity of the LP is $\mathcal{O}(K)$. Therefore, the complexity of Algorithm 1 is $\mathcal{O}(q^2K^2+K)$.

\begin{algorithm}[tb]
\caption{Optimal Algorithm for Problem (P1$'$)}
\begin{algorithmic}[1]
\STATE Initialize $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$.
\REPEAT
\STATE Initialize $\alpha_k = 1/K, k = 1,\cdots,K$.
\REPEAT
\STATE Compute $\{s_k\}_{k=1}^{K}$ and $\{m_k\}_{k=1}^{K}$ by \eqref{s_k_TDMA} and \eqref{m_k_TDMA}, respectively.
\STATE Obtain $\{\alpha_k\}_{k=1}^{K}$  with given $\{s_k\}$ and $\{m_k\}$ by bisection search.
\UNTIL improvement of $\mathcal{L}_k,k=1,\cdots,K$ converges to a prescribed accuracy.
\STATE Compute $s_0$ and $\alpha_0$ by \eqref{s0_TDMA} and \eqref{alpha0_TDMA}, respectively.
\STATE Update $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$ according to the ellipsoid method via \eqref{citidu}.
\UNTIL $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$ converge to a prescribed accuracy.
\STATE Set $s_k^*=s_k$, $m_k^*=m_k$,  $\alpha_k^*=\alpha_k$ for $k=1,\cdots,K$, and $\alpha_0^* = 1-\sum\limits_{k=1}^{K}\alpha_k^*$.
\STATE Obtain $s_0^*$ by solving Problem (P1$'$) with $\{s_k^*\}_{k=1}^{K}$, $\{m_k^*\}_{k=1}^{K}$ and$\{\alpha_k^*\}_{k=0}^{K}$.
\end{algorithmic}
\end{algorithm}

\begin{proposition} \label{TDMA_peak_infty}
For the TDMA case with $K$ source-destination pairs and $P_{\rm peak}\to\infty$, the maximum sum-rate by solving Problem (P1$'$) is  achieved by $\alpha_0^*\to0$.
\end{proposition}

\textit{Proof:} Clearly, we have $\alpha_0>0$ and $s_0>0$; otherwise, no energy will be harvested at the sources. Since the objective function of Problem (P1$'$) is an increasing function of $\alpha_k$ for $k=1,\cdots,K$ from constraint \eqref{C_rate2_TDMA},  when it comes to the extreme case with $P_{\rm peak}\to\infty$, for any given $s_k$ and $m_k$ satisfying constraints \eqref{energy_Q}, \eqref{energy_source_TDMAk1} and \eqref{energy_source_TDMA}, the optimal solution must be achieved by $\sum_{k=1}^{K}\alpha_k\to 1$  according to constraint \eqref{alpha_s_TDMA}. In this case, $\alpha_0^*\to0$ and $p_0^*\to \infty$ are required to guarantee positive harvested energy at the sources. The proof is thus completed.

\begin{proposition} \label{TDMA_peak_fty}
For the TDMA case with $K$ source-destination pairs and finite $P_{\rm{peak}}$,  the maximum sum-rate for Problem (P1$'$) is achieved by $p_{0}^*=P_{\rm{peak}}$.
\end{proposition}

\textit{Proof:} Please refer to Appendix A.

By Proposition \ref{TDMA_peak_infty} and Proposition \ref{TDMA_peak_fty}, it can be inferred that Problem (P1$'$) is actually a problem of energy and time allocation at the HRN, i.e., allocating energy and time for WPT and each WIT. Therefore,  for any given energy  allocated for WPT (i.e., $s_0=\alpha_0 p_0$), the HRN should charge the sources at its maximum available power (i.e., $p_0=P_{\rm peak}$), so that the time used for WPT $\alpha_0=s_0/p_0$ can be as small as possible and more time $1-\alpha_0$ can be allocated to WIT due to the  sum-rate maximization goal. In particular, when $P_{\rm peak}\to \infty$, the portion of transmission time $\alpha_0$ for WPT should asymptotically go to zero, which means that the sources can harvest sufficient energy in a sufficiently small time and  almost whole time is allocated to  WIT.

\subsection{Suboptimal Algorithm}

The complexity of the optimal algorithm becomes high as the number of pairs increases, mainly due to the dual updates. By simplifying the system model and eliminating the dual updates, in this section, we present an efficient suboptimal algorithm which significantly reduces the complexity.

At first, in WIT phase, the received power at each source in other periods is from the relay and other sources, which are both small. Specifically, the received energy from other sources is negligible due to the double energy
decay, i.e., the energy decay of relay-to-source and then source-to-source. As DF relaying protocol is adopted, the transmission power of the relay could relatively match the source's transmit power, and thus the relay's transmit power for forwarding is also small. As a result, in this section, we consider that the harvested energy at the sources is only from the WPT phase. With give $\alpha_0$, the transmit power of source $k$ can be given by
\begin{align}
q_k = \frac{2(\eta\alpha_0p_0g_{{\rm r},k}-E_k^c)^+}{\alpha_k},\quad  k=1,\cdots,K. \label{q_k+}
\end{align}

Second, due to Proposition \ref{TDMA_peak_fty}, we let $p_0=P_{\rm peak}$. Moreover, we assume that the equal power allocation (EPA) at the HRN in the WIT phase, the transmit power at the HRN for pair $k$ is thus given by
\begin{align} \label{average_P_TDMA}
p_k = \min\left\{\frac{2(P-\alpha_0P_{\rm peak})}{1-\alpha_0},P_{\rm peak}\right\}, k =1,\cdots,K.
\end{align}

Third, due to the energy decay in the WPT phase, the transmit power of sources may be small, thus the performance of this considered dual-hop relaying system may depend on the rate of first hop under most cases. As a result, in this section, we only focus on maximizing the sum rate of the first hop. Therefore, we have the following problem:

\begin{align}

\max_{\boldsymbol{\alpha}\in \mathcal{D}}\quad\quad&\sum\limits_{k=1}^K \frac{\alpha_k}{2}\log_2\left(1+\frac{2A_k }{\alpha_k\sigma^2}\right) \label{sumrate_first_TDMA}

\end{align}
where $A_k \triangleq (\eta\alpha_0P_{\rm peak}g_{{\rm r},k}-E_k^c)^+h_{1,k}$.

\begin{proposition} \label{sub_TDMA}
The optimal solution of Problem \eqref{sumrate_first_TDMA} with given $\alpha_0$ is given by
\begin{align}
q_k &= \frac{2}{(1-\alpha_0)h_{1,k}}\sum_{k=1}^{K} A_k,\quad  k=1\cdots,K,\label{q_sub_TDMA}\\
\alpha_k &= \frac{ A_k}{\sum_{k=1}^{K} A_k } (1-\alpha_0),\quad  k=1\cdots,K.\label{alpha_sub_TDMA}
\end{align}
\end{proposition}

\textit{Proof: }Please refer to Appendix B.

With given $\alpha_0$, we can obtain a set of $\{\boldsymbol{\alpha},\boldsymbol{p},\boldsymbol{q}\}$ by \eqref{average_P_TDMA}, \eqref{q_sub_TDMA} and \eqref{alpha_sub_TDMA}. Then, the optimal $\alpha_0$ maximizing the sum-rate can be found by the one-dimensional search.

\begin{algorithm}[tb]
\caption{Suboptimal Algorithm for Problem (P1)}
\begin{algorithmic}[1]
\STATE Divide $\alpha_{0}$ in $[0,1]$ with fixed step $\epsilon$.
\FOR{ each $\alpha_0 P_{\rm peak} \leq  P$ }

\STATE Compute the time allocation for WIT $\{\alpha_k\}_{k=1}^{K}$ according to \eqref{alpha_sub_TDMA}.
\STATE Compute the power allocation for WIT $\{p_k\}_{k=1}^{K}$ and $\{q_k\}_{k=1}^{K}$ according to  \eqref{average_P_TDMA} and \eqref{q_sub_TDMA}, respectively.
\STATE Compute the sum-rate according to \eqref{R_k_TDMA} with given $\alpha_0$.

\ENDFOR
\STATE Choose the optimal $\alpha_0^*$ that has the maximum sum-rate.

\end{algorithmic}
\end{algorithm}

To summarize, the above suboptimal algorithm is given in Algorithm 2.  The complexity of steps 3-5 is $\mathcal{O}(K)$. The complexity for searching $\alpha_0$ is $\mathcal{O}(1/\epsilon)$. Therefore, the whole complexity of Algorithm 2 is $\mathcal{O}(K/\epsilon)$, which is  linear in $K$ and much lower than that of the optimal algorithm in above subsection.

\section{Resource Allocation in FDMA Case}

Problem (P2) is a mixed integer programming and thus is NP-hard and non-convex.  However, it has been shown that the duality gap of the resource allocation problems in FDMA systems becomes zero when the number of SCs goes to large \cite{dual_gap1,dual_gap2}. This means that the optimal solution obtained in dual domain is equivalent to the optimal solution of the original non-convex problem due to the zero duality gap. Thus we solve Problem (P2) in dual domain.

At first, we introduce non-negative Lagrangian multipliers $\boldsymbol{\lambda}=\{\lambda_{k,n}\}\succeq 0$ and $\boldsymbol{\beta}=\{\beta_{k,n}\}\succeq 0$ corresponding to the two rates of the first and second hops in \eqref{r_k,n_10}, and $\boldsymbol{\nu}=\{\nu_k\}\succeq 0$ associated with the energy causality constraint \eqref{energy_C_OFDMA}. Moreover, $\mu\geq 0$, $\xi\geq 0$ are introduced to associate with the total time constraint \eqref{x_k,n_OFDMA} and total energy constraint \eqref{FDMA_P}, respectively. Then the dual function of Problem (P2) can be defined as
\begin{align}
g(\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi) \triangleq \max_{\{\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{R}\}\in\mathcal{D}} \mathcal{L}(\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{R},\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi),
\end{align}
where $\mathcal{D}$ is the set of all primal variables $\{\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{R}\}$ satisfying  the constraints, and the Lagrangian of Problem (P2) is
\begin{align} \label{L_FDMA}
&\mathcal{L}(\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{R},\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi) \nonumber\\
= &\sum_{k=1}^{K}\sum_{n=1}^{N} x_{k,n} [ R_{k,n}+ \lambda_{k,n}(R_{1,k,n}-R_{k,n}) \nonumber\\
&+ \beta_{k,n}(R_{2,k,n}-R_{k,n}) ] + \mu (1-\alpha_0-\alpha_1) \nonumber\\
&+ \xi\left( P -\alpha_0 p_0 - \sum_{n=1}^{N} \sum_{k=1}^{K} \frac{\alpha_1}{2}p_{k,n}\right) \nonumber\\
&+ \sum_{k=1}^{K} \nu_k \left(\eta \alpha_0 p_0 g_{{\rm r},k}  - \sum_{n=1}^{N} \frac{\alpha_1}{2}q_{k,n}- E_k^c\right).
\end{align}

Computing the dual function $g(\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi)$ requires to determine the optimal $\{\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{R}\}$ for given dual variables $\{\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi\}$. In the following we present the derivations in detail.

\subsection{Optimizing $\{\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{R}\}$ for Given  $\{\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi\}$}

\subsubsection{Maximizing Lagrangian over $\{R_{k,n}\}$} Similar to TDMA case, the part of dual function with respect to $\{R_{k,n}\}$ is given by
\begin{align}
g_R(\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi) = \max_{\boldsymbol{R}\succeq 0} \sum_{k=1}^{K} \sum_{n=1}^{N} (1-\lambda_{k,n}-\beta_{k,n})x_{k,n}R_{k,n}.
\end{align}

To make sure that the dual function is bounded, we have $1-\lambda_{k,n}-\beta_{k,n} \equiv 0$. In such case, $g_R(\boldsymbol{\lambda},\boldsymbol{\beta},\boldsymbol{\nu},\mu,\xi) \equiv 0$ and we obtain that $\beta_{k,n} = 1 - \lambda_{k,n}.$
Note that $0\leq\lambda_{k,n}\leq1$ to make sure that $\beta_{k,n}$ is non-negative. By substituting the result above into \eqref{L_FDMA}, we have
\begin{align} \label{L_FDMA_xiao}
&\mathcal{L}(\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi) \nonumber\\
=&  \sum_{k=1}^{K} \sum_{n=1}^{N} x_{k,n} \left[ \lambda_{k,n}R_{1,k,n} + (1-\lambda_{k,n})R_{2,k,n} \right] \nonumber\\
&+ \mu (1-\alpha_0-\alpha_1) \nonumber\\
&+ \xi\left( P -\alpha_0 p_0 - \sum_{n=1}^{N} \sum_{k=1}^{K} \frac{\alpha_1}{2}p_{k,n}\right) \nonumber\\
& + \sum_{k=1}^{K} \nu_k \left(\eta \alpha_0 p_0 g_{{\rm r},k}  - \sum_{n=1}^{N} \frac{\alpha_1}{2}q_{k,n}- E_k^c\right).
\end{align}

\subsubsection{Maximizing Lagrangian over $\{p_{k,n}\}$, $\{q_{k,n}\}$, $\{x_{k,n}\}$ and $\alpha_1$ } Observing the Lagrangian in \eqref{L_FDMA_xiao}, we can rewrite \eqref{L_FDMA_xiao} as follows:
\begin{align}
&\mathcal{L}(\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi) \nonumber\\
=& \sum_{n=1}^{N}\mathcal{L}_n (\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{\lambda},\boldsymbol{\nu},\xi) + \mu(1-\alpha_0-\alpha_1) \nonumber\\
 & + \xi(P -\alpha_0p_0) + \sum_{k=1}^{K}\nu_k\left(\eta\alpha_0p_0g_{{\rm r},k}-E_k^c\right),
\end{align}
where
\begin{align} \label{Ln_FDMA}
&\mathcal{L}_n (\boldsymbol{p},\boldsymbol{q},\boldsymbol{x},\boldsymbol{\alpha},\boldsymbol{\lambda},\boldsymbol{\nu},\xi)\nonumber\\
=&\sum_{k=1}^{K} x_{k,n} \left[ \lambda_{k,n}R_{1,k,n} + (1-\lambda_{k,n})R_{2,k,n} \right]\nonumber\\
&- \sum_{k=1}^{K} \frac{\alpha_1}{2}\left( \xi p_{k,n} + \nu_k q_{k,n} \right).
\end{align}
Maximizing $\mathcal{L}$ over $\{p_{k,n}\}$ and $\{q_{k,n}\}$ is equivalent to maximizing each $\mathcal{L}_n$ over $p_{k,n}$ and $q_{k,n}$, which is shown in the following. Each SC $n$ should be allocated to at most one pair and we can apply exhaustive method for SC to obtain the optimal $k^*$. Particulary, assume that SC $n$ is selected for pair $k$, then \eqref{Ln_FDMA} is equivalent to
\begin{align}  \label{Ln_fix}
\mathcal{L}_n (p_{k,n},q_{k,n},\alpha_1,\lambda_{k,n},\nu_k,\xi)
=& \lambda_{k,n} R_{1,k,n} +  (1-\lambda_{k,n}) R_{2,k,n}\nonumber\\
-&  \frac{\alpha_1}{2}\left( \xi p_{k,n} + \nu_k q_{k,n} \right).
\end{align}
By differentiating \eqref{Ln_fix} with respect to $p_{k,n}$ and $q_{k,n}$, and letting them to zero, the $p_{k,n}$ and $q_{k,n}$ maximizing $\mathcal{L}_n$ are given by
\begin{align}
p_{k,n} &= \min\left\{\left( \frac{1-\lambda_{k,n}}{\xi N \ln2} - \frac{\sigma^2}{h_{2,k,n}}\right)^+,P_{\rm peak}\right\}, \label{p}\\
q_{k,n} &= \left( \frac{\lambda_{k,n}}{\nu_k N \ln2} - \frac{\sigma^2}{h_{1,k,n}}\right)^+.\label{q}
\end{align}

After computing the power allocations $p_{k,n}$ and $q_{k,n}$ for WIT, we can obtain the SC allocation  maximizing each $\mathcal{L}_n$ as

\begin{eqnarray}  \label{eq:x(k,n)}
x_{k,n}=
\begin{cases}
1,       &\mbox{if} \quad k=\mathop{\argmax}_{k}\mathcal{L}_{n},\\
0,   &\mbox{otherwise}.
\end{cases}
\end{eqnarray}

With given $p_{k,n}$, $q_{k,n}$ and $x_{k,n}$, the Lagrangian \eqref{L_FDMA_xiao} becomes a linear function of $\alpha_1$. From \eqref{L_FDMA_xiao}, to maximize $\mathcal{L}$, we have
\begin{align} \label{alpha1_FDMA}
\alpha_1 = \left\{
             \begin{array}{ll}
               1, & \hbox{if $\Psi  > 0$,} \\
               0, & \hbox{otherwise.}
             \end{array}
           \right.
\end{align}
where $\Psi$ is given by

\begin{align}
\Psi =& \sum_{k=1}^{K}\sum_{n=1}^{N} x_{k,n} \left[ \lambda_{k,n}R_{1,k,n} + (1-\lambda_{k,n})R_{2,k,n} \right] \nonumber\\
&-  \frac{1}{2}\sum_{k=1}^{K}\sum_{n=1}^{N} x_{k,n}(\xi p_{k,n} - \nu_kq_{k,n}).
\end{align}
Here, $p_{k,n}$, $q_{k,n}$ and $x_{k,n}$ have been obtained by \eqref{p}, \eqref{q} and \eqref{eq:x(k,n)}, respectively.

\subsubsection{Maximizing Lagrangian over $\alpha_0$ and $p_0$} From \eqref{L_FDMA_xiao}, to maximize $\mathcal{L}$, we have
\begin{align}
\alpha_0 &= \left\{
             \begin{array}{ll}
               1, & \hbox{if $\sum\limits_{k=1}^{K} \eta \nu_k g_{{\rm r},k} P_{\rm peak} -\xi P_{\rm peak} - \mu > 0$,} \\
               0, & \hbox{otherwise.}
             \end{array}
           \right.\label{alpha_0FDMA}\\
p_0 &= \left\{
        \begin{array}{ll}
          P_{\rm peak}, & \hbox{if $\sum\limits_{k=1}^{K} \eta \nu_k g_{{\rm r},k} - \xi > 0$,} \\
          0, & \hbox{otherwise.}
        \end{array}
      \right. \label{p0FDMA}
\end{align}

\subsection{Optimizing Dual Variables $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$}

Similar to TDMA case, the ellipsoid method can be employed to update $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$ toward optimal $\{\boldsymbol{\lambda}^*,\boldsymbol{\nu}^*,\mu^*,\xi^*\}$ with global convergence \cite{convexbook}, the subgradients required for which are

\begin{equation}     \label{subtidu}
\boldsymbol{\triangle}=\left[
  \begin{array}{lcl}
    \Delta \lambda_{k,n} = x_{k,n}(R_{1,k,n}-R_{2,k,n}) , \forall k,n\\
    \Delta \mu  = 1-\alpha_0-\alpha_1\\
\Delta \xi = P -\alpha_0p_0 - \sum\limits_{k=1}^{K}\sum\limits_{n=1}^{N}\frac{\alpha_1}{2}p_{k,n}\\
\Delta \nu_k = \eta\alpha_0p_0g_{{\rm r},k} - \sum\limits_{n=1}^{N}\frac{\alpha_1}{2}q_{k,n} - E_k^c, \forall k
  \end{array}
\right]
\end{equation}

\subsection{Discussions on Optimality and Complexity}

It is worth noting that the optimal $\boldsymbol{p}^*$, $\boldsymbol{q}^*$ and $\boldsymbol{x}^*$ are obtained at optimal $\{\boldsymbol{\lambda}^*,\boldsymbol{\nu}^*,\mu^*,\xi^*\}$. With given optimal $\boldsymbol{p}^*$, $\boldsymbol{q}^*$ and $\boldsymbol{x}^*$, Problem (P2) becomes a LP with $\boldsymbol{\alpha}$. The optimal $\boldsymbol{\alpha}^*$ can be obtained by solving this LP.

\begin{algorithm}[tb]
\caption{Optimal Algorithm for Problem (P2)}
\begin{algorithmic}[1]
\STATE \textbf{initialize} $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$.
\REPEAT
\FOR{each SC $n$}
\STATE Compute  $\{p_{k,n}\}$  and $\{q_{k,n}\}$ according to  \eqref{p} and \eqref{q}, respectively.
\STATE Obtain the  subcarrier allocation $\{x_{k,n}\}$ according to \eqref{eq:x(k,n)}.
\ENDFOR
\STATE Compute $\alpha_1,\alpha_0,p_0$ according to \eqref{alpha1_FDMA}, \eqref{alpha_0FDMA} and \eqref{p0FDMA}, respectively.
\STATE Update $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$  by the ellipsoid method  using the subgradients defined in \eqref{subtidu}.
\UNTIL $\{\boldsymbol{\lambda},\boldsymbol{\nu},\mu,\xi\}$ converge to a prescribed accuracy.
\STATE Set $\boldsymbol{p}^* = \boldsymbol{p}$, $\boldsymbol{q}^* = \boldsymbol{q}$ and $\boldsymbol{x}^*=\boldsymbol{x}$.
\STATE Obtain $\boldsymbol{\alpha}^*$ by solving Problem (P2) with $\boldsymbol{p} = \boldsymbol{p}^*$, $\boldsymbol{q} = \boldsymbol{q}^*$ and $\boldsymbol{x}=\boldsymbol{x}^*$.
\end{algorithmic}
\end{algorithm}

To summarize, the algorithm for Problem (P2) is given by Algorithm 3. The complexity of steps 4-5 is  $\mathcal{O}(K N)$. The complexity of steps 3-6 is $\mathcal{O}(K N^2)$.
The ellipsoid method needs complexity of $\mathcal{O}(q)$, where $q$ is the number of dual variables and $q=KN+K+2$ in our case. The time complexity of the LP is $\mathcal{O}(K)$. Therefore, the whole complexity of Algorithm 3 is $\mathcal{O}(q^2KN^2 + K)$.

\begin{proposition}\label{FDMA_0}
In the case of the FDMA case with $K$ source-destination pairs and $P_{\rm peak} \to \infty$, the  maximum sum-rate obtained by solving Problem (P2) is given by $\alpha_0^*\to 0$ and $\alpha_1^* \to 1$.
\end{proposition}

\textit{Proof:} Clearly, we can easily obtain that $\alpha_0>0$; otherwise, no energy is harvested at the sources  in the considered R-WPCN.  Thus,  $\alpha_1<1$.  Denote $s_{0}=\alpha_0 p_{0}$, $s_{k,n}=\alpha_1 p_{k,n}/2$ and $m_{k,n}=\alpha_1 q_{k,n}/2$ as energy variables. For ang given $s_0,s_{k,n}$, $m_{k,n}$ and $x_{k,n}$ satisfying the primary constraints \eqref{x_k,n_OFDMA}, \eqref{x_k,n}, \eqref{FDMA_P}, \eqref{energy_C_OFDMA}, the objective function of Problem (P2) is an increasing function of $\alpha_1$ according to \eqref{r_k,n_10}. Thus, the sum-rate is maximized when $\alpha_1^* \to 1$, which follows $\alpha_0^*\to 0$.  The proof is thus completed.

By Proposition \ref{FDMA_0},  the positive harvested energy at the sources is achieved under the assumption that the HRN is able to transmit an infinite power due to $\alpha_0 \to 0$. For a finite $P_{\rm peak}$, a nonzero time ratio should be scheduled to the WPT phase to harvest sufficient energy for WIT.

Similar to the TDMA case, when it comes to the more general case with $P_{\rm peak}<\infty$, we have the following proposition,

\begin{proposition} \label{FDMA_dingli_peak}
In the case of the FDMA case with $K$ source-destination pairs and $P_{\rm peak} < \infty$, the  maximum sum-rate is achieved by $p_0^*=P_{\rm peak}$.
\end{proposition}

\textit{Proof:} The proof is similar as the proof of Proposition \ref{TDMA_peak_fty}, and thus is omitted here.

\subsection{Suboptimal Algorithm}

The main complexity of the above optimal algorithm is resulted from the ellipsoid method. In this section, we propose a suboptimal algorithm by assuming equal power allocation (EPA) over SCs, which can eliminate the dual updates, while the optimal $\alpha_0$ is  obtained by the one-dimensional search.

First, due to the energy decay in the WPT phase, the transmit power of sources may be small, thus the performance of this considered dual-hop relaying system may depend on the rate of first hop under most cases. As a result, we heuristically  choose the pair having the maximum $h_{1,k,n}$ to occupy SC $n$,  which is given by

\begin{eqnarray}  \label{sub_x}
x_{k,n}=
\begin{cases}
1,       &\mbox{if} \quad k=k^{\ast}=\mathop{\argmax}_{k} h_{1,k,n},\\
0,   &\mbox{otherwise}.
\end{cases}
\end{eqnarray}

With given $\alpha_0$, we let $p_0=P_{\rm peak}$ according to Proposition \ref{FDMA_dingli_peak}. Then, from  \eqref{FDMA_P} and \eqref{energy_C_OFDMA}, the power allocation by assuming EPA can be given by
\begin{align}
p_{k,n} &=\left\{
  \begin{array}{ll}
     \min \left\{\left[ \frac{2(P-\alpha_0 P_{\rm peak})}{\alpha_1 N}\right ]^+,P_{\rm peak}\right\}, & \hbox{if $x_{k,n}=1$,} \label{p_sub}\\
    0, & \hbox{otherwise.}
  \end{array}
\right.\\
q_{k,n}&=\left\{
  \begin{array}{ll}
     \frac{2 (\eta\alpha_0 P_{\rm peak}g_{{\rm r},k} - E_k^c)^+ }{\alpha_1 M_k}, & \hbox{if $x_{k,n}=1$,} \label{q_sub} \\
    0, & \hbox{otherwise.}
  \end{array}
\right.
\end{align}
where $M_k$ is the number of SCs assigned to source $k$, which can be determined via \eqref{sub_x}.

\begin{algorithm}[tb]
\caption{Subptimal Algorithm for Problem (P2)}
\begin{algorithmic}[1]

\STATE Divide $\alpha_{0}$ in $[0,1]$ with fixed step $\epsilon$.
\FOR{ each $\alpha_{0}P_{\rm{peak}}\leq P $}
\STATE Obtain the SC allocation $\{x_{k,n}\}$ by \eqref{sub_x}.
\STATE Compute the power allocations $\{p_{k,n}\}$ and $\{q_{k,n}\}$ according to \eqref{p_sub} and \eqref{q_sub}, respectively.
\STATE Compute the sum-rate according to \eqref{FDMA_minrate} for given $\alpha_0$.
\ENDFOR
\STATE Choose the optimal $\alpha_0^*$ that has the maximum sum-rate.
\end{algorithmic}
\end{algorithm}

The above suboptimal algorithm is summarized in Algorithm 4. The complexity of  steps 3-5 is $\mathcal{O}(K N)$. The complexity for searching $\alpha_0$ is $\mathcal{O}(1/\epsilon)$,
therefore the whole complexity of the algorithm is $\mathcal{O}(K N/\epsilon)$, which is  much lower than that of Algorithm 3.

\section{Numerical Results}

\begin{figure}[t]
\begin{centering}
\centering
\includegraphics[scale=1.2]{Figure3.eps}
\vspace{-0.1cm}
\caption{The location of R-WPCN in the simulations.}\label{location}
\end{centering}
\vspace{-0.1cm}
\end{figure}

In this section, we provide extensive numerical results to evaluate the performance of the proposed algorithms.

In the simulation setup, we assume that the total bandwidth is 10 MHz for both TDMA and FDMA cases and  the noise spectral density is assumed to be $-174$ dBm/Hz. We set  the energy conversion
efficiency as $\eta=0.8$  and the processing cost as $E_k^c = 1 \times 10^{-7}$ Joule   at all sources.

For all simulations, we set $K=4$ pairs of sources-destinations  unless otherwise noted.

Moreover,  we consider a two-dimensional plane of node location as shown in Fig. \ref{location}, where the source nodes and destination nodes are randomly but uniformly distributed in the corresponding square regions, and the HRN can move along with the x-axis from $-$5 m to 5 m. The HRN is assumed to locate at (0,0) unless otherwise noted.

In this paper, the pass-loss exponent is 3 and we adopt Richan fading channel model for the small-scale fading, where the Richan factor is set to be 3.

\begin{figure}[t]
\begin{centering}
\includegraphics[scale=0.65]{Figure4.eps}
\vspace{-0.1cm}
\caption{Sum-rate versus  total transmit power $P$ in  TDMA-based  R-WPCN, where $P_{\rm peak}=2P$.  }\label{fig:TDMA_performance_compare}
\end{centering}
\vspace{-0.1cm}
\end{figure}

Fig. \ref{fig:TDMA_performance_compare} demonstrates the sum-rate versus the total transmit power $P$ in the TDMA case with $P_{\rm peak}=2P$.

We introduce  the following two benchmark schemes for the purpose of performance comparison. First, the equal energy allocation (EEA) scheme is considered, where the energy allocated for WPT at the HRN is fixed as $\alpha_0 P_{\rm peak}=P/2$, while the optimal time allocation is still obtained as Algorithm 1.

In addition, we consider the equal resource allocation (ERA) for WIT with given $\alpha_0$, where the equal time and power allocations for WIT are assumed and the optimal $\alpha_0^*$ is obtained as Algorithm 2.

For all  schemes, the sum-rate is observed to increase with the total transmit power $P$.

Compared with EEA and ERA, we can see that the proposed optimal and suboptimal algorithms achieve  better performance.

And the suboptimal algorithm is observed to perform very closely to the optimal algorithm, which demonstrates the effectiveness of the proposed suboptimal resource allocations.

\begin{figure}
\begin{centering}
\includegraphics[scale=0.63]{Figure5.eps}
\vspace{-0.1cm}
\caption{Duality gap versus number of SCs, where $P=30$ dBm, $P_{\rm peak}=2P$.} \label{fig:dual_gap}
\end{centering}
\vspace{-0.1cm}
\end{figure}

Fig. \ref{fig:dual_gap} illustrates the duality gaps versus the different numbers of SCs $N$. The duality gap is shown to decrease with the SCs number $N$. It can be observed that the duality gap is indeed approximately zero with 64 SCs, which  verifies the effectiveness of the proposed dual-based algorithm of FDMA case.

\begin{figure}[t]
\begin{centering}
\includegraphics[scale=0.65]{Figure6.eps}
\vspace{-0.1cm}
\caption{Sum-rate versus total transmit power $P$ in FDMA-based  R-WPCN, where $N=64$, $P_{\rm{peak}}=2P$.} \label{FDMA_perform}
\end{centering}
\vspace{-0.1cm}
\end{figure}

Fig. \ref{FDMA_perform} demonstrates the sum-rate versus the total transmit power $P$ in the FDMA case. The  peak power constraint $P_{\rm{peak}}$ is set to be $2P$.

For performance comparison, we consider the following three benchmarking schemes.

The first one is subcarrier pairing scheme where the subcarrier allocation in the two hops can be different \cite{pairing1}.

The second  one is equal energy allocation (EEA) where the  WPT energy  is fixed as $\alpha_0  P_{\rm peak} = P/2$, while the optimal power and SC allocations are still obtained as Algorithm 3.

The last benchmark is that the subcarrier assignment is fixed (FSA) while the optimal time and power allocations are also jointly optimized as Algorithm 3.

For all schemes, we can observe that the sum-rate is increasing with the total transmit power $P$ and the proposed optimal scheme achieves considerable gain compared with the other benchmark  schemes.

Compared with the optimal algorithm, it can be observed that the performance gain achieved by subcarrier pairing is limited but requires additional $\mathcal{O}(N^3)$ complexity by Hungarian algorithm.

Besides, the suboptimal algorithm with low complexity also has good performance.

Moreover, we can observe that the EEA scheme is only efficient under some particular system setup ($24 < P < 32$ dBm), while the suboptimal algorithm has a good performance over a wide range of transmit power compared with the EEA scheme, which demonstrates the superiority  of the proposed suboptimal algorithm.

The poor performance of FSA compared to the proposed schemes  indicates that dynamic SC allocation  provides significant improvement in terms of the sum-rate.

 Á½ÖÖÄ£ÐÍÏÂ¶Ô±È  ·åÖµ¹¦ÂÊÎªºá×ø±ê
\begin{figure}[t]
\begin{centering}
\includegraphics[scale=0.65]{Figure7.eps}

\vspace{-0.1cm}
\caption{System performance  versus peak power constraint $P_{\rm{peak}}$, where $N=64$ and $P=30$ dBm.}
\label{optimal_a0_value}
\end{centering}
\vspace{-0.1cm}
\end{figure}

Fig. \ref{optimal_a0_value} illustrates  system performance  versus $P_{\rm{peak}}$  with  the fixed total transmit power $P=30$  dBm.

First, it can be observed that the sum-rate increases with $P_{\rm{peak}}$ and reaches a plateau when $P_{\rm peak}>40$ dBm. This is because that the total available power $P$ is fixed and thus the sum-rate must be bounded even if $P_{\rm peak}$ becomes sufficiently large.

Second, for both schemes, we can see that the optimal $\alpha_0^*$ decreases with $P_{\rm peak}$. This is because that less time for WPT is needed to obtain the harvested energy requirement with a larger $P_{\rm peak}$.

In addition, we can also observe that the energy for WPT increases with $P_{\rm peak}$ and then remains fixed for both schemes, where almost all the available energy is used for WPT. This is because that the optimal WPT energy is under the peak power constraint, which becomes infeasible when $P_{\rm peak}$ becomes sufficiently large.

\begin{figure}[t]
\begin{centering}
\includegraphics[scale=0.65]{Figure8.eps}

\vspace{-0.1cm}
\caption{Optimal energy for WPT and sum-rate versus the distance $d_1$, where $N=64$, $P=30$ dBm and $P_{\rm peak}=2P.$}
\label{distance_A_compare_TDMA_FDMA}
\end{centering}
\vspace{-0.1cm}
\end{figure}

Fig. \ref{distance_A_compare_TDMA_FDMA}  examines the effect of the relay position on the energy for WPT and the sum-rate, where the total transmit power is set to be 30 dBm and the HRN moves along with the x-axis from $d=-5$ m to $d=5$ m.

First, we can observe that the energy allocated for WPT is increasing with  $d$ for both cases. This may be because that the channel gains for WPT become worse with a longer distance from sources to HRN, which requires more energy at the HRN allocated to WPT.

Besides, the sum-rate for both schemes is observed to decrease with $d$,  which is different from the traditional relaying systems.

On the one hand, due to the doubly distance-dependent signal attenuation for both WPT and the first hop of WIT, the rate of the first hop cannot be improved though the WPT energy becomes larger.

On the other hand, since the energy for WIT is decreasing with $d$ due to the larger WPT energy, the rate of the second hop is also bottlenecked.  As a result, the HRN should be located in proximity to the sources.

\begin{figure}[t]
\begin{centering}
\includegraphics[scale=0.65]{Figure9.eps}

\vspace{-0.1cm}
\caption{Sum-rate and optimal energy for WPT versus the number of pairs, where $N=64$, $P=30$ dBm and $P_{\rm peak}=2P.$}
\label{Kcompare_TDMA_FDMA}
\end{centering}
\vspace{-0.1cm}
\end{figure}

Fig. \ref{Kcompare_TDMA_FDMA} depicts the sum-rate and the optimal energy for WPT versus the number of pairs $K$ with $P = 30$ dBm.

First, it can be observed that the FDMA scheme achieves higher sum-rate compared to the TDMA scheme.

Besides, we can see that the sum-rate for both schemes increases with the number of pairs and the rate tends to be saturated due to the limited bandwidth and transmission power of the system.

In addition, it can be observed that the energy for WPT decreases with the increasing number of pairs, this may be because that  more pairs can harvest energy from the same energy signals broadcasted by the HRN, which results in  a higher energy efficiency.

\section{Conclusions}

This paper studied a new  R-WPCN based on TDMA and FDMA for multiple source-destination pairs with the help of a HRN which acts as both roles of an energy transmitter and an information forwarder.

By considering  the ``charge-then-forward" protocol at the HRN, we studied the sum-rate maximization problem for both TDMA and FDMA cases. For the TDMA case, we proposed a global optimal solution, while for the FDMA case, we designed an asymptotically optimal solution.

To tradeoff the performance and complexity, the suboptimal algorithms for both cases were also proposed.

Extensive simulations showed that our proposed algorithms significantly outperform the conventional schemes.

The following directions  can be considered for possible future works.  First, the HRN is equipped with large-scale antenna array (or massive MIMO). Specifically, massive MIMO can generate concentrated energy beams to power wireless nodes and thus deal with the challenge of long-distance WPT. Second, multiple HRNs with full-duplex could be considered to improve spectral efficiency, where distributed coordinated beamforming at the HRNs for both WPT and WIT will be interesting. Third, it may also consider that the direct link exists between the sources and destinations.  In this case the resource allocation schemes will be largely different. Last, it will be interesting to extend our work to the multiple frames, where sources can accumulate energy and then transmit information over different frames, and the HRN can dynamically allocate its resources over different frames.

\section*{Acknowledgement}
The authors would like to thank the anonymous reviewers. Their valuable suggestions greatly improved the quality of this paper.

\appendices

\section{Proof of Proposition \ref{TDMA_peak_fty} }

It can be proved by contradiction as follows:

We denote the optimal solution of Problem $\rm (P1')$ as $\{\alpha_k^*\}_{k=0}^{K}$, $\{s_k^*\}_{k=0}^{K}$ and $\{m_k^*\}_{k=1}^{K}$.  Suppose that $s_0^* < \alpha_0^* P_{\rm peak}$, i.e., $p_0^*<P_{\rm peak}$. Then,  we  consider the following solution $s_0^*$ and $\{\tilde{\alpha}_k\}_{k=0}^{K}$, where $s_0^* = \tilde{\alpha}_0 P_{\rm peak}$, i.e., $\tilde{p}_0 = P_{\rm peak}$.

Since $p_0^* < \tilde{p}_0$, we have $\alpha_0^* > \tilde{\alpha}_0$ with the same optimal $s_0^*$. From \eqref{alpha_s_TDMA} and \eqref{alpha_k_peak}, we can obtain that $\tilde{\alpha}_k \geq \alpha_k^*$ for $ k=1,\cdots,K$.

Moreover,  according to constraint  \eqref{C_rate2_TDMA}, as the objective function of Problem $\rm (P1')$ is an increasing function of $\{\alpha_k\}_{k=1}^{K}$ with given $s_k^*$ and $m_k^*$,  the case $\{\tilde{\alpha}_k\}_{k=0}^{K}$ achieves higher sum-rate than the case $\{\alpha_k^*\}_{k=0}^{K}$.

Thus, this contradicts with the assumption that the $\{\alpha_k^*\}_{k=0}^{K}$ is the optimal solution.

Therefore, the optimal solution must be given by $\{\tilde{\alpha}_k\}_{k=0}^{K}$ and $\tilde{p}_0 = s_0^*/\tilde{\alpha}_0 = P_{\rm peak}$, which completes the proof.

\section{Proof of Proposition \ref{sub_TDMA}}

We can easily prove that Problem \eqref{sumrate_first_TDMA}  is convex with given $\alpha_0$, which can be optimally solved by the Lagrangian dual method.

The Lagrangian of Problem \eqref{sumrate_first_TDMA}  is given by
\begin{align}
\mathcal{L}(\boldsymbol{\alpha},\lambda) = \sum\limits_{k=1}^K \frac{\alpha_k}{2}\log_2\left(1+\frac{2A_k }{\alpha_k\sigma^2}\right) + \lambda(1-\sum_{k=0}^{K}\alpha_k),\nonumber
\end{align}
where $\lambda\geq 0$ denotes the Lagrangian multiplier associated with the constraint \eqref{alpha_k_TDMA}. Since the problem is convex, we can find its optimal solution by using Karush-Kuhn-Tucker (KKT) conditions. With given $\alpha_0$, let us denote the primary and dual optimality values of Problem \eqref{sumrate_first_TDMA}  as $\{\alpha_k^*\}_{k=1}^{K}$ and $\lambda^*$. By differentiating $\mathcal{L}(\boldsymbol{\alpha},\lambda)$ with respect to $\alpha_k$ and using KKT stationary conditions, we obtain
\begin{align} \label{daoshu_x}
\frac{1}{2}\log_2(1+2x_k) - \frac{1}{2\ln2} + \frac{1}{2\ln2(1+2x_k)}=\lambda^*,
\end{align}
where $x_k = \frac{A_k }{\alpha_k^*\sigma^2},k=1,\cdots,K.$

To make \eqref{daoshu_x} more clearly, we denote $y_k = \frac{1}{1+2x_k}$, then \eqref{daoshu_x} is equivalent to
\begin{align} \label{daoshu_y}
\log_2y_k = \frac{\ln y_k}{\ln2}= \frac{y_k}{\ln2}-\frac{1}{\ln2}-2\lambda^*.
\end{align}
Thus, we have $-y_k - \delta = -\ln y_k = \ln\frac{1}{y_k}$, where $\delta = -1-2\ln2 \lambda^*$. Then, it is easy for us to get that $e^{-y_k-\delta}=\frac{1}{y_k}$, which is equivalent to $-y_k e^{-y_k}=-e^\delta$. As a result, the solution of \eqref{daoshu_y} is given by $y_k = -W(-e^\delta)$.
Finally, due to $y_k = \frac{1}{1+2x_k}$ and $x_k =\frac{ A_k }{\alpha_k^* \sigma^2}$, we can obtain that the solution of $\frac{\partial\mathcal{L}_k}{\partial\alpha_k}=0$  is given by
\begin{align} \label{TDMA_timeW}
\alpha_k^* = \frac{-2 A_k   W(-e^\delta)}{\sigma^2(1+W(-e^\delta))}, \quad \forall k =1,\cdots,K.
\end{align}

From \eqref{TDMA_timeW}, we can find that $\alpha_k^*,k=1,\cdots,K$ is proportional to $A_k$. Moreover, it can be easily verified that $\sum_{k=1}^{k=K}\alpha_k^*=1-\alpha_0$ must hold for Problem \eqref{sumrate_first_TDMA} with given $\alpha_0$. Thus the optimal $\{\alpha_k^*\}_{k=1}^{K}$ with given $\alpha_0$ is thus given by
\begin{align} \label{alpha_sub}
\alpha_k^* = \frac{ A_k}{\sum_{k=1}^{K} A_k } (1-\alpha_0), \quad \forall k=1\cdots,K.
\end{align}
By plugging \eqref{alpha_sub} into \eqref{q_k+}, then the optimal $\{q_k^*\}_{k=1}^{K}$ with given $\alpha_0$ is given by
\begin{align}
q_k^* = \frac{2}{(1-\alpha_0)h_{1,k}}\sum_{k=1}^{K} A_k, \quad \forall k=1\cdots,K.
\end{align}
The proof is thus completed.

\bibliography{Reference}

\end{document}


TITLE: Why is the conjugate of a function always convex?
QUESTION [9 upvotes]: The conjugate of a function $f$ is given (for some $y \in \operatorname{dom}(f)$) as:
$$f^*(y) = \sup_{x \in \operatorname{dom}(f)}\left(y^Tx - f(x)\right)$$
It is known that $f^*$ is convex even if $f$ is not. I would like to know how to prove this.

REPLY [15 votes]: For any given $x \in \operatorname{dom} f$, the function:
$$y \mapsto y^\top x - f(x)$$
is an affine function, which is convex and lower semicontinuous (and, indeed, continuous). The epigraph of these functions is therefore closed and convex.
Then, $f^*$ is the pointwise supremum of the above functions, and its epigraph is the intersection of the above affine functions' epigraphs. Each are closed and convex, proving the epigraph of $f^*$ is closed and convex, thus $f^*$ is convex and lower semicontinuous.


TITLE: Trapezoidal integration in matlab - code explanation
QUESTION [0 upvotes]: I have a simple matlab script for calculating an Integral using trapezoidal rule.
The idea is calculate Fourier series, and it used for the Fourier coefficients.
function [area] = recArea(func, xl, xr, N, n)
    if nargin < 5
        n = 1;
    end
    w = ((xl - xr) / N);
    y = eval(func);
    area = sum((y(1:end-1) + y(2:end)) / 2*w);
end

Call example:

$a_n$ = recArea(func, xl, xr, N)

in the function:

function [y] = fourier(m, N, func, xl, xr)

Where:
$func$ is the function to calculate it's fourier sereis
$xl, xr$ the integral bounds: $\int_{xl}^{ar}$.
Can you please explain me how this code works?
Thanks in advance.

REPLY [2 votes]: The domain $[x_l, x_r]$ is uniformly partitioned with $N$ sub-intervals, meaning that the width of any sub-interval is $\frac{x_r-x_l}{N}$.
For trapezoid rule integration, you compute the area of the trapezoid formed by the line segment connecting $(x_i,f(x_i))$ to $(x_{i+1},f(x_{i+1}))$ (draw it). The area under this line is basically the area of the rectangle formed by drawing a line segment parallel to the x-axis from $(x_i,f(x_i))$ to $(x_{i+1},f(x_i))$, which can be expressed by:
$$A_1 = \textrm{width}\times\textrm{height} = \frac{x_r-x_l}{N}\times f(x_i).$$
Then, you add the area of the triangle whose base is the line segment just described, and whose hypotenuse is the line from $(x_i,f(x_i))$ to $(x_{i+1},f(x_{i+1}))$. This has area:
$$A_2 = \frac12\textrm{base}\times\textrm{height} = \frac12 \frac{x_r-x_l}{N} \times \left(f(x_{i+1})-f(x_i)\right).$$
Add $A_1$ and $A_2$, and you get the area of a sub-interval as:
$$A_i = A_1+A_2 = \frac{x_r-x_l}{2N}\left(f(x_i)+f(x_{i+1})\right).$$
Sum these over all intervals, and that exactly describes what's happening by the line of code that says area = sum(....


TITLE: Complex Analysis - Argument Principle vs. Rouche's Theorem
QUESTION [2 upvotes]: The Argument Principle
Suppose a function $f$ is meromorphic on an open set that contains a circle $C$ and its interior. Further assume that $f$ has no zeroes on $C$ (but may have zeroes in the interior of $C$). Then,
$$\frac{1}{2\pi i} \oint_C \frac{f'(z)}{f(z)} dz \ =$$
the number of zeroes of $f$ in the interior of $C$ (counted with multiplicity) minus the number of poles of $f$ in the interior of $C$ (counted with multiplicity).
Rouche's Theorem
Suppose that $f$ and $g$ are holomorphic functions on a domain and let $C$ be a circle whose interior also lies in the domain. If $|f(z)-g(z)| \lt |f(z)|$ for all $z$ on $C$, then $f$ and $g$ have the same number of zeroes in the interior of $C$, counting multiplicities.
My Question
Why does the Argument Principle require $f$ to have no zeroes on $C$, whereas this requirement is relaxed for Rouche's Theorem?
Thank you.

REPLY [4 votes]: In the Argument Principle: If $f$ has a zero  $z_0 \in C$, then $z_0$ is a pole for the quotient $f'/f$. Therefore you cannot define $\int_C f'(z)/f(z)\,dz$.
In Rouche's Theorem: The inequality $|f-g|<|f|$ on $C$ implies that $f$ and $g$ have no zeroes on $C$. In fact, if $g$ (resp. $f$) had a zero $z_0 \in C$, then $|f(z_0)|<|f(z_0)|$
(resp. $|g(z_0)|<0$), which is a contradiction.


TITLE: Defining a recurrence relation for number of words of length $n$ formed from the alphabet $\{x, y, z\}$ that do not contain the string $xxx$
QUESTION [0 upvotes]: $a_{n}$ describes the number of words that can be composed of this particular set $\{x,y,z\}$. The sequence $xxx$ must not appear in the word.
Example: $a_{1}=3$,  $a_{2}=9$,  $a_{3}=26$
The answer is: $a_{n}=2a_{n-1}+2a_{n-2}+2a_{n-3}$
How did they reach this solution? Any sort of guidance will be appreciated.

REPLY [0 votes]: You have correctly found the initial conditions.
Notice that any admissible word must begin with $y$, $z$, $xy$, $xz$, $xxy$, or $xxz$.
An admissible word of length $n$ that begins with $y$ or $z$ can be formed by appending an admissible word of length $n - 1$, of which there are $a_{n - 1}$, to a word that begins with $y$ or $z$.  Thus, there are $2a_{n - 1}$ words that begin with $y$ or $z$.
An admissible word of length $n$ that begins with $xy$ or $xz$ can be formed by appending an admissible word of length $n - 2$, of which there are $a_{n - 2}$, to a word that begins with $y$ or $z$.  Thus, there are $2a_{n - 2}$ words that begin with $y$ or $z$.
I will let you write the argument for words that begin with $xxy$ or $xxz$.


TITLE: Analogue of Tate or Hodge conjecture for varieties over $\mathbb Q_p$
QUESTION [17 upvotes]: I've been learning about p-adic Hodge theory recently (I'm a beginner), and I've been wondering about the following question the past couple of weeks.  Sorry for the long setup, it's mainly background; the question is down towards the bottom.
Background: For a field $k$, and a second field $K$ of characteristic 0, we can consider the $K$-linear category $Mot(k)_{K}$ of motives defined by smooth projective varieties over $k$ (pure motives) up to homological or numerical equivalence (let's assume homological equivalence equals numerical equivalence for each Weil cohomology we consider).  Then this is a semi-simple abelian category (due to Jannsen), and there are realization functors to various categories of vector spaces extending the functor sending a variety $X$ to its cohomology groups (e.g. its etale cohomology groups $H^i(X,\mathbb{Q}_\ell)$ for $\ell \not= p$ in which case $K = \mathbb{Q}_\ell$).
In the case when $k = \mathbb{C}$ and we take $K = \mathbb{Q}$, we can consider the singular cohomology realization $Mot(\mathbb{C})_\mathbb{Q} \rightarrow {Vect}_{\mathbb{Q}}$ sending $X$ to $H_{sing}^*(X,\mathbb{Q})$.  One way of stating the Hodge conjecture is that when we consider the extra structure of the Hodge decomposition on $H^*_{sing}(X,\mathbb{Q}) \otimes \mathbb{C}$ and therefore consider the enriched realization $Mot(\mathbb{C})_{\mathbb{Q}} \rightarrow {Hdge}_{\mathbb{Q}}$ to pure $\mathbb{Q}$-Hodge structures, this functor is fully faithful.
Similarly, if $k$ is a finite field or a number field, we can take $K = \mathbb{Q}_\ell$, we can consider the functor $Mot(k)_{\mathbb{Q}_{\ell}} \rightarrow GalRep(\mathbb{Q}_\ell)$ which sends $X$ to its $\ell$-adic cohomology $H^i(X_{\overline{k}},\mathbb{Q}_\ell)$ as an $\ell$-adic Galois representation.  Then the Tate conjecture is that this functor is fully faithful.  In the number field case, there are some other conjecturally fully faithful enriched realization functors discussed in Ch. 7 of Andre's book "Une Introduction Aux Motifs".
Motivation: Now in $p$-adic Hodge theory, we deal with varieties over a complete discrete valuation field (for example $k = \mathbb{Q}_p$), and study the extra structures which one obtains on the de Rham cohomology of such a variety from Hodge theory, crystalline cohomology, and etale cohomology.  To me it seems that one is trying to endow the de Rham cohomology of a variety $X$ with as much "natural" structure as possible, at least enough to determine the Galois representation $H^i(X_{\overline{k}},\mathbb{Q}_p)$, but possibly one can add more.  So I wonder the following:

Let $k$ be a finite extension of $\mathbb{Q}_p$.  Is it reasonable to believe that there exist a characteristic-0 field $K$ and a fully faithful functor $Mot(k)_K \rightarrow \mathcal{C}$, where $\mathcal{C}$ is the category of $K$-vector spaces with some "natural" extra structure, for example that considered in $p$-adic Hodge theory?  If not, why not?  If so, is there a candidate for $\mathcal{C}$?

EDIT: I want to add an additional condition: that when one takes the forgetful functor $\mathcal{C} \rightarrow Vect(K)$, the composition $Mot(k)_K \rightarrow Vect(K)$ is the functor defined by some Weil cohomology (e.g., de Rham cohomology or $p$-adic etale cohomology).
My impression from the example given by user ulrich in the comment to
Should the etale cohomology of a smooth projective variety (over rationals) be semi-simple; why?
is that we CANNOT take $K = \mathbb{Q}_p$ and take the functor sending $X$ to $H^*(X_{\overline{k}},\mathbb{Q}_p)$ , because the the $p$-adic cohomology of $X$ might not be semi-simple, for example when $X$ is an elliptic curve with non-split multiplicative reduction.  But perhaps there's another possible $\mathcal{C}$?
EDIT: As wccanard explained in a comment (now deleted), the example above doesn't actually show that the $p$-adic etale cohomology functor isn't fully faithful, since the subobject of $H^1(X,\mathbb{Q}_p)$ which is not split is not "motivic".  However, there are other reasons why $p$-adic etale cohomology can't work as explained in wccarnard's answer below.

REPLY [7 votes]: This is not an answer, it's just a note that some of the more standard examples of functors don't seem to work because they're not fully faithful. This was originally in the comments but it clogged them up so I deleted them and moved them here.
The observation I want to make is simply that if $k=\mathbb{Q}_p$ then I think that none of the following work in general: $K=\mathbb{Q}_\ell$ and the functor is $\ell$-adic etale cohomology; $K=\mathbb{Q}_p$ and the functor is $p$-adic etale cohomology; and finally $K=\mathbb{Q}_p$, the target category is $K$-vector spaces equipped with a filtration, a Frobenius and a monodromy operator in the usual $p$-adic Hodge theory waym and the functor is $D_{pst}$.
The reason none of these work is the following. Fix an elliptic curve over $\mathbb{F}_p$ with $a_p=0$ (assume such a curve exists). Then there are uncountably many lifts of this curve to $\mathbb{Q}_p$, so there are two non-isogenous lifts. The motives attached to the $h^1$ of these curves will, I imagine, have no non-zero maps between them. However the functors described above are all controlled by the special fibre so produce isomorphic linear algebra objects. In the $D_{pst}$ case what is happening is that the Frobenius has trace zero and the $Fil^1$ can't be either of the eigenspaces so it must be another line; however all other lines give isomorphic filtered $\phi$-modules.
The reason this doesn't contradict the Hodge conjecture is that you can take coefficient field equal to $\mathbb{Q}$ in that setting, which is much smaller, and then the filtration (which is defined once you tensor up to $\mathbb{C}$ gives you far more information.


TITLE: A curve $(t,t^2,t^3,\ldots,t^n) \in\mathbb R^n$ is not contained in any hyperplane
QUESTION [3 upvotes]: Let $C$ the curve given by $\gamma:\mathbb R\to\mathbb R^n$, $\gamma(t)=(t,t^2,t^3,\ldots,t^n)$.

Show that $C$ is contained in no hyperplane.

2.Find the tangent line and osculating plane at $(0,0,\ldots,0)$.
For (1) my idea is to find the curvature of $\gamma$. We have,
$$\gamma'(t)=(1,2t,3t^2,...,nt^{n-1})\quad \gamma''(t)=(0,2,6t,...,n(n-1)t^{n-2}).$$
Then $k(t)=||\gamma''(t)||=\sqrt[]{4+36t^2+\cdots+(n(n-1)t^{n-2})^2}>0$ so it is not contained in any hyperplane. Is it ok? If so, is there another way to solve this without using curvature?
To get the tangent line and osculating plane, I just take the binormal vector $B(t)=\gamma'(t)\times \gamma''(t)$ and them i get the plane.
Is it ok?
Thanks for your help.

REPLY [2 votes]: Forget the curvature.
Let's take $n+1$ different points along the curve, the first of them being at $t=0$, just to make things look nice. Draw $n$ vectors from that point to all other points. Now, the non-(hyper)planarity of the curve is equivalent to rank of that system being $n$. In other words, these vectors must be a basis of all $\mathbb R^n$. The easiest way to verify that is to check the determinant of a matrix composed of these vectors - it must not be 0. But wait, that's the Vandermonde matrix, and as long as all points are different (which I said they are), $\det\ne0$.


TITLE: When is an $L^p(\mathbb{R})$ function bounded everywhere?
QUESTION [0 upvotes]: I have a fairly elementary question regarding Lebesgue spaces that I can't seem to figure out.
As the title suggests, I wish to know what additional condition(s) do we need to impose on a function $f\in L^p(\mathbb{R})$ so that it is bounded everywhere (as oppose to almost everywhere)? Would continuity suffice?

REPLY [1 votes]: Smooth counterexample: Let me do an example of a smooth function that is in $L^p(\mathbb{R})$, but not bounded. Pick your favourite smooth function $\varphi$ supported in $[0;1/2]$ and such that $\Vert \varphi \Vert_p=1$. Then we define
$$ f(x) = \sum_{n\geq 1} n^{1/2}\varphi(n^{2p}(x-n)). $$
Then we have
$$ \Vert f \Vert_p \leq \sum_{n\geq 1} n^{1/2}\Vert \varphi(n^{2p}(\cdot - n) \Vert_p
= \sum_{n\geq 1} \frac{1}{n^{3/2}}. $$
Here we used a couple of substitutions
$$ \Vert \varphi(n^{2p}(\cdot - n) \Vert_p =\left( \int_{\mathbb{R}} \vert \varphi(n^{2p}(x-n)) \vert^p dx \right)^{1/p}
= \left( \int_\mathbb{R} \vert \varphi(n^{2p}y) dy \right)^{1/p}
= \left(\int_{\mathbb{R}} n^{-2p} \vert \varphi(z) dz \right)^{1/p}
= n^{-2} \Vert \varphi \Vert_p = n^{-2}. $$
Thus, $f$ is in $L^p(\mathbb{R})$. As the $\varphi(n^{2p}(x-n))$ have disjoint support and $\varphi$ is smooth, so is $f$. Furthermore, again as the $\varphi(n^{2p}(x-n))$ have disjoint support, we get that $\Vert f \Vert_{L^\infty([0;n+1))}= n^{1/2}\Vert \varphi \Vert_\infty$ and hence $\Vert f \Vert_\infty = \infty$.
Of course, this is just a complicated version of the example of copper.hat from the comments. I just wanted to show that continuity or differentiability is not quite enough.
Lipschitz continuity is enough: As I wrote in the comments, lipschitz continuity would be sufficient. As $f$ is lipschitz continous, is is enough to show that $\lim_{\vert x \vert \rightarrow \infty} f(x) = 0$.
Assume by contradiction that $\lim_{x\rightarrow \infty} f(x)=0$ does not hold. Then there exists $C >0$ and a sequence $(x_n)_{n\geq 1}$ such that $\lim_{n\rightarrow\infty} x_n =\infty$, such that $f(x_n)\geq 2C$ (I assumed wlog that those values are positive and that $\vert x_n - x_m \vert \geq 1$ for $n\neq m$). But our function is lipschitz continuous with, say, lipschitz constant $L>0$. Then we have for $y\in [x_n; x_n +a]$ with $a:=\min \{1/2, C/L \}$
$$ f(y)\geq  f(x_n)- \vert f(y)-f(x_n)\vert \geq 2C - L \vert x_n - y\vert \geq C.$$
Thus, we get
$$\Vert f \Vert_p^p \geq \sum_{n\geq 1} \int_{x_n}^{x_n+a} \vert f(y) \vert^p dy \geq \infty,$$
which contradicts our assumption that $f\in L^p(\mathbb{R})$.
Alternatively, you could assume that $f$ is $C^1(\mathbb{R})$ and impose some integrability for $f'$.
As PhoemueX pointed out in the comments below, uniform continuity works as well (would be essentially the same proof).


TITLE: Area of intersection between Cycloid and square
QUESTION [1 upvotes]: I understand how to compute the area under the whole cycloid. But how do you compute a partial area.
Let say, we want to know the area of the intersection between a cycloid generated by a circle of radius $r$ that stretches from $-r\pi$ to $+r\pi$ and a square of side 1 that stretches from $-1/2$ to $+1/2$.
With: $1/(2\pi) < r < 1/2$
I tried using the parametric equation the same way I did to compute the whole cycloid area. But I am stuck when I switch from an integral over x to an integral over t as t run from and to an angle than depends on $r$.
How do you proceed to compute the area ?

REPLY [1 votes]: Let the square have sidelength $2a$ in order to make things dimensionally correct. Parametrize the cycloid by
$$\gamma:\quad x(t)=r(t+\sin t), \quad y(t)=r(1+\cos t)\qquad(-\pi\leq t\leq\pi)\ .$$
The vertical sides of the square intersect $\gamma$ in two points $P$, $Q$. The area $S$ in question is then bounded by four arcs, among them a part $\gamma'$ of $\gamma$ which belongs to a certain subinterval $-t_0\leq t\leq t_0$. The value $t_0\in\>]0,\pi[\>$ has to be determined numerically from the equation $$x(t)=r(t+\sin t)=a\ .$$
In order to compute the area of $S$ we use Green's theorem
$${\rm area}(S)=-\int_{\partial S} y\>dx=\int_{\gamma'}y\>dx=\int_{-t_0}^{t_0}y(t)\,\dot x(t)\>dt\ .$$
Here we have made use of the fact that $y\,dx=0$ along three of the four  arcs constituting $\partial S$, and that $\gamma'$ has the wrong sense of  direction in $\partial S$.


TITLE: Taylor expansion in C++ libm uses different constants
QUESTION [4 upvotes]: The Taylor expansion of the arcsine is usually given by
$$
\arcsin(x) = x+ \frac{1}{6} x^3 + \frac{3}{40} x^5 + \frac{5}{112} x^7 + \ldots
$$
This (https://github.com/simonbyrne/apple-libm/blob/master/Source/Intel/asin.c) implementation of the arcsine in libm uses constants that are very close to the above but slightly different.  I understand that they can use better constants given the fact they know what taylor expansion they will stop at and to reduce the error over the fixed range they are interested in.  My question is how do they determine those constants?  How could I duplicate it if I wanted to use a different taylor expansion or a taylor expansion of different order?

REPLY [4 votes]: See Computer Approximations by Hart, Cheney, et al. Likely what you're seeing are the result of Chebyshev (uniform) approximation over specific intervals. Have a look at the reference cited in the code for CR-LIBM.


TITLE: Bounded real-valued function on $[0,1]$, non integrable?
QUESTION [2 upvotes]: Does such a function exist? If yes, it must be a very pathological case. I am talking here about Lebesgue integrability.
For instance, if $f(x)=1$ if $x$ is rational and zero otherwise, then $\int_0^1 f(x)dx = 0$. So you need to find an example more pathological than that. A possible example is the following.
Let $f(x)$ be the realization of a Gaussian random variable $Z_x$ with mean equal to $0$ and variance equal to $1$. Let us assume that the $Z_x$'s are identically and independently distributed. Such a function $f(x)$ is nowhere continuous, and can be seen as the realization of a white noise. However, you could argue that its integral on $[0,t]$ is the value $B(t)$ of a realization of a Brownian motion starting with $B(0)=0$, and measured at time $t$. Thus $\int_0^1 f(x) dx = B(1)$. Note that Brownian motions are nowhere differentiable, so maybe there is a contradiction in what I am saying here.
Anyway, I've never found counter-examples: a function bounded on $[0, 1]$ but not integrable in that interval. Can you show an example?

REPLY [3 votes]: Let $f$ be a bounded function on $[0,1]$.
Either $f$ is measurable, and then
$$
\int_0^1 |f| ≤ \sup |f|\ \int_0^1 1\,\mathrm d x = \sup |f| < \infty
$$
so $f$ is integrable.
Either $f$ is not measurable. This exists iff you assume assume the axiom of choice. You can then take any non measurable set $\Omega$ and take $f = \chi_\Omega$ the characteristic function of this set, as suggested by Nate Eldredge. Then by definition, this function is not integrable.


TITLE: An application of Hurwitz's Theorem
QUESTION [0 upvotes]: Suppose $f(z)=\sum_{n=0}a_nz^n$ such that the radius of convergence is 1 and $a_0\neq0$.Let $s_n(z)=a_0+a_1z+\cdots+a_nz^n$ be the partial sums.Let $1>\delta>0$ such that $f\neq 0$ on $|z|=1-\delta$.Prove that there is $K>0$ and integer N such that for all $n\ge N$ the modulus of the product of zeros of $s_n(z)$ inside $|z|=1-\delta$ is at least $K$.
We know that $f(0)\neq0$,what bothers me is that since for sufficiently large $n$,the zeros of $s_n(z)$ should be near $f$,what if $f$ has no zeros inside $|z|=1-\delta$.The reason to choose a smaller circle than the unit one seems to avoid the zeros of $s_n$ goes to the boundary.I sense that I shall invoke the radius of convergence to derive some constraints on the coefficients so as to ensure polynomial $s_n$ will have zeros inside the circle,which is near zero.

REPLY [1 votes]: A Hurwitz-based proof may be intuitively clear, but writing it out in every detail takes a bit of work. I'd say: $f$ has distinct zeros $a_1,\dots,a_m$ with multiplicities $k_1,\dots,k_m$. Let $r>0$ be small enough so that the $r$-neighborhoods of the zeros are disjoint and stay inside of $|z|=1-\delta$. Also make sure that $r<\frac12 \min|a_j|$. By Hurwitz, for large enough $n$ the function $s_n$ has $k_j$ zeros inside of the $r$-neighborhood of $a_j$; and it has no other zeroes inside  of $|z|=1-\delta$. It follows that the product of its zeroes is at least $$\prod_{j=1}^m (|a_j|-r)^{k_j} \ge \prod_{j=1}^m (|a_j|/2)^{k_j} $$
which is the desired lower bound.
An alternative (in my opinion, better) approach is to use Jensen's formula.

Remark. The statement could be made cleaner and stronger: for every $\delta\in (0,1)$ there exist $K>0$ and $N$ such that  for all $n\ge N$ the modulus of the product of zeros of $s_n$ inside $|z|=1-\delta$ is at least $K$.
Indeed, if $\delta$ happens to be such that $f$ has a zero on the circle $|z|=\delta$, you can use some larger value of $\delta$ instead. There are at most countably many values of $\delta$ to avoid, since $f$ has at most countably many zeros in the unit disk.


TITLE: Smallest $\sigma$-algebra
QUESTION [0 upvotes]: I have a measure space $(\mathcal{X},\mathbb{E},\mu)$. I am considering the paving of all $\mu$-nullsets given by:
$\mathbb{N}_{\mu}= \left\{ N \subseteq \mathcal{X} : \exists E \in \mathbb{E} \hspace{0,1cm} \text{such that } N \subseteq E \hspace{0,1cm} \text{and }\mu(E)=0   \right\}$.
I am also considering the paving $E_{\mu}$ which is $\mathbb{E}$ extended by the nullsets:
$\mathbb{E}_{\mu}= \left\{ E \cup N : E \in \mathbb{E} \hspace{0,1cm} \text{and }N \in \mathbb{N}_{\mu}   \right\}$
A paving is understood to be an arbitrary collection of subsets.
I wish to show that $\mathbb{E}_{\mu}$ is the smallest $\sigma$-algebra containing $\mathbb{E}$ and $\mathbb{N}_{\mu}$.
My thoughts
I know that the $\sigma$-algebra generated by a set $\mathcal{G}$ is the smallest $\sigma$-algebra containing $\mathcal{G}$.
Using this I can say that the smallest $\sigma$-algebra containing $\mathbb{E}$ and $\mathbb{N}_{\mu}$ is $\sigma(\mathbb{E} \cup \mathbb{N}_{\mu})$. So I have to show that $\sigma(\mathbb{E} \cup \mathbb{N}_{\mu})= \mathbb{E}_{\mu}$. For this purpose I have to show 1) $\mathbb{E}_{\mu} \subseteq \sigma(\mathbb{E} \cup \mathbb{N}_{\mu})$ and 2) $\sigma(\mathbb{E} \cup \mathbb{N}_{\mu})\subseteq \mathbb{E}_{\mu}$.
1) Assume that $K \in \mathbb{E}_{\mu}$. Then $K$ can be written as $K=E \cup N$. Since $\sigma(\mathbb{E} \cup \mathbb{N}_{\mu})$ consists of countable unions, intersection and complements of the elements $E \cup N$ then $K \in \sigma(\mathbb{E} \cup \mathbb{N}_{\mu})$ and therefore we have that $\mathbb{E}_{\mu} \subseteq \sigma(\mathbb{E} \cup \mathbb{N}_{\mu})$.
2) I have shown that $\mathbb{E} \subseteq \mathbb{E}_{\mu}$ and $\mathbb{N}_{\mu} \subseteq \mathbb{E}_{\mu}$. Therefore we have that $\mathbb{E} \cup \mathbb{N}_{\mu} \in \mathbb{E}_{\mu}$. By the properties of $\sigma$-algebra we have that $\sigma(\mathbb{E} \cup \mathbb{N}_{\mu}) \subseteq \sigma(\mathbb{E}_{\mu})= \mathbb{E}_{\mu}$. The last equality holds since I have shown $\mathbb{E}_{\mu}$ to be a $\sigma$-algebra.
Is my my proof correct?

REPLY [1 votes]: Let $\mathcal{A}$ be the smallest $\sigma$-algebra containing both $\mathbb E$ and $\mathbb N_\mu$.
Let $A\in \Bbb E_\mu$, so $A = E \cup N$ with $E \in \Bbb E$ and $B \in \Bbb N_\mu$. So $E,N \in \mathcal{A}$ by definition, and so (as $\mathcal{A}$ is closed under unions), $E \cup N = A \in \mathcal{A}$. So $\Bbb E_\mu \subseteq \mathcal{A}$.
Secondly, $\Bbb E_\mu$ is a $\sigma$-algebra (check the axioms !, as is done here) and contains $\Bbb E$ (we can take $N=\emptyset \in \Bbb N_\mu$ to see that $E=E \cup \emptyset \in \Bbb E_\mu$ for any $E \in \Bbb E$, and similarly $\Bbb N_\mu \subseteq \Bbb E_\mu$ by the same idea (take $\emptyset \in \Bbb E$ as well). As $\mathcal{A}$ is the minimal such $\sigma$-algebra, $\mathcal{A} \subseteq \Bbb E_\mu$, showing the other inclusion.