Guo and Sun

RESEARCH

On the Performance of Cooperative Non-orthogonal Multiple Access With Energy Harvesting Over Nakagami-m Fading Channels

Weidong Guo^{1,2}^{*} and Lin Sun^{2}

*Correspondence: gud2001@qfnu.edu.cn 1School of Information Science and Engineering, Shandong University, 72 Binhai Road, 266237 Qingdao, China Full list of author information is available at the end of the article

Abstract In this paper, we analyze the performance of a full-duplex (FD) cooperative non-orthogonal multiple access (NOMA) model in multi-relay multiuser networks. In contrast to conventional cooperative networks, relays in the considered network have no embedded energy supply, but they can harvest energy from the radio-frequency signals broadcasted by the source to facilitate the cooperative NOMA transmission. Based on this structure, we propose a novel relay selection scheme in Nakagami-m fading channels by considering both the channel state information (CSI) and the battery statuses of relays. Assuming each relay has finite energy storage for the energy accumulation, we adpot a finite Markov chain to capture the evolution of relay batteries and simplify the performance analysis by making some reasonable assumptions. General closed-form expressions of the outage probability and the ergodic sum-rate are derived. All the theoretical results are validated by Monte-Carlo simulations. The impacts of various system parameters, such as the number of relays and the self-interference (SI) at the involved relay, on the performance are extensively investigated. It is shown that the FD relaying achieves better performance in the moderate and high SNR regions and the usage of NOMA with FD relaying outperforms the conventional orthogonal multiple access (OMA) network in terms of terminal fairness.

Keywords: Non-orthogonal multiple access; Full-duplex; Energy harvesting; Markov chain

1 Introduction

Guo and Sun Page 2 of 18

With the aim of increasing the capacity of wireless systems, non-orthogonal multiple access (NOMA) has been recognized as a promising technique for the fifth generation (5G) networks due to its superior spectral efficiency [1]. The key idea of NOMA is serving multiple users in the same time solt over the same frequency band, but with different power levels to achieve a tradeoff between system throughput and user fairness [2]. Lan et al. [3], explored the impact of error propagation of successive interference cancellation (SIC) and user velocity on the NOMA performance. Their results showed that even in the worst error propagation scenario, NOMA can outperform the conventional orthogonal multiple access (OMA) for various user mobilities. Chen et al. [4], studied NOMA over the downlink of a wireless system. Traditional minimum-mean-squared-error (MMSE) precoding matrices have been used. Considering that wireless relaying is an effective means to combat the effects of channel fading and shadowing on the transmitted signal, it is reasonable to integrate NOMA into relaying networks. In [5], a new cooperative NOMA scheme was proposed and analyzed in terms of outage probability and diversity gain. The authors in [6] studied the impact of relay selection on the performance of cooperative NOMA. In particular, a two-stage relay selection strategy is proposed, which can achieve the minimal outage probability among all the possible selection schemes, and achieve the full diversity gain. Yu et al. [7] propose a cooperative truncated automatic repeat-request (T-ARQ)-based NOMA scheme with relay selection. The numerical results show that the proposed scheme can achieve much higher throughput than the T-ARQ-based time-division multiple access schemes.

Although cooperative NOMA is capable of improving the performance for farusers by using a half-duplex relay, it results in additional bandwidth cost for the system. To avoid this issue, full duplex (FD) is a promising technique for the next generation wireless systems, because it can double the spectrum efficiency by realizing transmission and reception on the same carrier frequency, simultaneously [8]. Two main types of FD relay techniques, namely FD amplify-and-forward (AF) relaying and FD decode-and-forward (DF) relaying, have been discussed in [9]-[11]. The outage probability of FD-AF relaying were analyzed in [9], which considers the processing delay of relaying in practical scenarios. Osorio et al. [10], analyzed the outage probability of FD-AF relaying with direct link available. The authors in [11] characterized the outage performance of FD-DF relaying. It demonstrated that the optimal duplex mode can be selected according to the outage probability.

In addition to improving the spectral efficiency which is the key motivation of NOMA, another critical goal of 5G networks is to maximize the energy efficiency. Energy harvesting (EH), has recently received considerable attention as a sustainable solution to alleviate the energy shortage of wireless networks. Radio frequency (RF) signals radiated by ambient transmitters can be treated as a viable new source for EH. A fundamental limitation of EH-based wireless communications lies in the restricted transmission range. Recently, relays with EH capabilities have attracted lots of research interests, where the relay can use the energy harvested from the source signals to perform the data forwarding [12]. The relay can be wirelessly powered and the communication coverage can be expanded. In [13], a wireless cooperative network was considered, in which multiple source-destination pairs communicate with each other via an EH relay. In [14], a best cooperative mechanism was proposed to facilitate the wireless EH and spectrum sharing, which can greatly improve both the energy and spectral efficiencies. Nevertheless, the amount of energy harvested from RF radiation is often restricted and it is desirable for relays to accumulate the harvested energy in the energy storage such as super-capacitors or rechargeable batteries [15]. In [16], a threshold-based energy storage scheme was employed at relays, the stored energy in the battery of each relay actually forms a Markov chain over time. By investigating the properties of this Markov chain, the asymptotic average throughput was derived. In [17], to support an efficient utilization of harvested energy to improve the system throughput, a harvest-use-store relaying strategy with distributed beamforming were proposed. Zhai et al. [18], proposed a high-efficient EH-based cooperative transmission protocol using either power-splitting or time-switching technique, where both the source and the relay are capable of harvesting energy from a power beacon and then using the harvested energy for the cooperative DF or AF transmission.

In order to improve the spectral efficiency and energy efficiency, we consider the FD-NOMA relaying for the downlink transmission over Nakagami-m fading channels and provide a comprehensive performance analysis. Our analysis provides more general results, as the Nakagami-m fading can incorporate the most commonly used Rayleigh fading and Rice fading as special cases, which can be obtained by adjusting the fading parameter. The main contributions of this paper are summarized as follows:

• We consider a NOMA-based downlink FD-AF relaying network over Nakagamim fading channels. The relays have no other energy supplies, but they are equipped with chargeable battery and can harvest and store the wireless energy broadcasted by the source. Based on this structure, a new relay selection scheme in Nakagami-m fading channels is proposed considering both the channel state information (CSI) and the battery statuses of relays.

• This model combines NOMA, cooperative communication, FD and EH together reasonably and can improve spectral efficiency and energy efficiency simultaneously.

• We model the amount of energy stored in the battery of each relay as finite states and use a finite Markov chain to model the variation of energy at each relay. Both the transition probability and the steady-state probability are derived for the next performance analysis.

• The asymptotic expressions of outage probability and the ergodic rate for both near and far users are derived. Finally, simulation results are presented to validate our theoretical analysis and verify the superiority of FD-NOMA relaying over the traditional half-duplex relaying based OMA.

The rest of this paper is organized as follows. In section 2, we introduce the system model. The relay selection scheme is illustrated in section 3. The closedform expressions of outage probability and ergodic rate are derived in section 4. Numerical and simulations results are presented in section 5 and conclusions are drawn in section 6.

2 Method

Next, we describe the setting and design of the system and the relay selection rule we provide.

2.1 System model In our system, a source S intends to communicate to a number of potential users D_{j} (j = 1, 2, ..., N) over channels with flat fading. Multiple potential FD relays R_{i} (i = 1, 2, ...M) are willing to amplify and forward the signal from S to the users. We assume that there exist direct links between S and users. S and each D_{j} are configured with one antenna and operate in half-duplex (HD) mode, while relays are equipped with one transmit antenna and one receive antenna that enables a FD operation. The channel between any transmitter and any receiver is assumed to follow Nakagami-m fading. The channels pertaining to the direct link, first-hop and second-hop undergo independent identically (i. i. d.) fading and the channel coefficients are denoted as h_{SD_{j}}, h_{SR_{i}} and h_{R_{i}D_{j}} with parameters m_{0}, m_{1} and m_{2}, respectively. The channel power gain G_{SD_{j}} = |h_{SD_{j}} |^{2}, G_{SR_{i}} = |h_{SR_{i}}|^{2} and G_{R_{i}D_{j}} = |h_{R_{i}D_{j}} |^{2} thus follow the Gamma distribution, with mean Φ_{SD}, Φ_{SR} and Φ_{RD}, respectively.

In a general case, due to the imperfect isolation or cancellation process, FD operation may suffer from residual self-interference (SI). In order to simplify the analysis, we assume each relay with spatially separated transmit and receive antennas is fixed, and the interference channel is non-fading, so the instantaneous interference-to-noise-ratio (INR) equals the corresponding average INR. As shown later in the simulations, this assumption does not degrade the accuracy of our analysis. Without loss of generality, we assume that users’ channel gains are ordered as Ĝ_{SD_{1}} ≤ Ĝ_{SD_{2}} ≤ ⋯ ≤ Ĝ_{SD_{N}} (Ĝ_{SD_{j}} is the sorted random variable of G_{SD_{j}}). In what follows, we focus on the analysis of the group users consisting of D_{p} and D_{q} (G_{SD_{P}} < G_{SD_{q}}). The additive white Gaussian noise (AWGN) of all links each link is assumed to have zero mean and equal variance σ^{2}.

Figure 1 A reference model for multi-relay cooperative FD-NOMA network.

In this paper, we use the conventional AF protocol for the NOMA downlink cooperative transmission and consider a normalized unit block time (i.e., T = 1). Before the transmission, each relay checks their battery at the beginning of a transmission block and judges whether it has enough energy to forward the source information. If the relay does not have enough energy or is not selected to forward the information, it will perform EH in this time block and store the harvested energy into its battery. We assume P_{S} is the transmission power of source and is sufficiently large such that the energy harvested from the noise is negligible. Thus, for a relay node R_{i}, the amount of energy harvested from the source can be expressed as

E_{R_{i}} = ηP_{S}G_{SR_{i}}T,  (1)

where η is the energy conversion efficiency at each relay.

For those relays with sufficient energy, they report their CSI to S for the relay selection. Let E_{i} denote the amount of energy in the battery of relay R_{i}, we define those relays with enough energy as the eligible set

ϕ = {R_{i}| E_{i} ≥ V C, i = 1, 2, ...M} ,  (2)

where C is the EH threshold to activate the EH circuit and V is a positive integer. If we assume the selected relay for forwarding data has the same transmitting power as S, V must satisfy

(V − 1) C < P_{S}T ≤ V C.  (3)

Among this set of relays, Rm̄ is selected and the index of Rm̄ can be expressed as

m̄ = arg max_{i:R_{i}∈ϕ} {G^{SR_{i}}} .  (4)

In (4), the relay with the best channel power gain of the first hop in this set is selected for forwarding information, all other relays perform EH in this block time. We assume the transmitting power of the selected relay (P_{Rm̄}) is the same as S, which means P_{S} = P_{Rm̄} = P. Although the energy consumed of Rm̄ for sending date in one block time is P T, which can be seen from (3), we still consider it as V C, because C is the basic unit of the relay battery. The difference between P T and V C does not severely degrade the accuracy of the performance analysis.

After selection, S broadcasts the signal x_{S}[k] = \sqrt{a_{p}}x_{p}[k]+\sqrt{a_{q}}x_{q}[k] (a_{p} and a_{q} are the power allocation factors and satisfy a_{p}+a_{q} = 1, a_{p} > a_{q}), while the selected relay Rm̄ receives x_{S}[k] and forwards the signal generated by co-channel transmission to D_{p} and D_{q} in FD mode, while the corresponding D_{p} and D_{q} simultaneously receive the signals transmitted from S and forwarded by Rm̄. In the kth time slot, k=1,2,3..., the signals received at Rm̄, D_{p} and D_{q} can be, respectively, expressed as

y_{Rm̄}[k] = h_{SRm̄} \sqrt{P} x_{S}[k] + h_{RRm̄} \sqrt{P} Ay_{Rm̄}[k − T_{d}] + n_{Rm̄}[k],  (5)
y_{D_{p}}[k] = h_{SD_{p}} \sqrt{P} x_{S}[k] + h_{Rm̄D_{p}} \sqrt{P} Ay_{Rm̄}[k − T_{d}] + n_{D_{p}}[k],  (6)
y_{D_{q}}[k] = h_{SD_{q}} \sqrt{P} x_{S}[k] + h_{Rm̄D_{q}} \sqrt{P} Ay_{Rm̄}[k − T_{d}] + n_{D_{q}}[k],  (7)

where A is the amplification coefficient, h_{RRm̄} is the residual SI between transmit antenna and receive antenna of relay, and T_{d} is the processing delay. n_{Rm̄}[k], n_{D_{p}}[k] and n_{D_{q}}[k] are independent circularly-symmetric Gaussian random variables. Without loss of generality and for the overall simplicity, we assume γ = P /σ^{2} is the average signal-to-noise-ratio (SNR). We also set γm̄ = γE (|h_{RRm̄}|^{2}) as an actual measured value, which averages out the fading effect. With unity energy transmit signals and after some simple substitutions, we can obtain [19]

A = (G_{SRm̄}P + |h_{RRm̄}|^{2}P + σ^{2})^{− \frac{1}{2}} .  (8)

In this study, we assume that the two signals from S and Rm̄ are resolvable absolutely (herein the relaying link corresponding to direct link has small time delay for any transmitted signals, the same as in [20]), therefore the instantaneous signalto-interference-plus-noise-ratios (SINRs) at D_{p} from S and Rm̄ are, respectively, expressed by

γ_{SD_{p}} = \frac{a_{p}Ĝ_{SD_{p}}γ}{a_{q}Ĝ_{SD_{p}}γ + 1},  (9)

γ_{Rm̄D_{p}} = \frac{a_{p}γ^{2}G_{SRm̄}G_{Rm̄D_{p}}}{γ (G_{SRm̄} + G_{Rm̄D_{p}}) + γG_{Rm̄D_{p}} (a_{q}γG_{SRm̄} + γm̄) + γm̄ + 1}. (10)

Using the selection combining (SC), we can obtain the instantaneous SINR of D_{p} as

γ_{D_{p}} = max (γ_{SD_{p}}, γ_{Rm̄D_{p}}) .  (11)

As for D_{q}, the SIC technique is applied to decode information of D_{p}. Particularly, D_{q} first decodes the message of D_{p}, then subtracts this component from the received signal to detect its own information. Therefore, the received SINRs at D_{q} for the decoding of the information of D_{p} are given by

γ_{SD_{q→p}} = \frac{a_{p}Ĝ_{SD_{q}}γ}{a_{q}Ĝ_{SD_{q}}γ + 1},  (12)

γ_{Rm̄D_{q→p}} = \frac{a_{p}γ^{2}G_{SRm̄}G_{Rm̄D_{q}}}{γ (G_{SRm̄} + G_{Rm̄D_{q}}) + γG_{Rm̄D_{q}} (a_{q}γG_{SRm̄} + γm̄) + γm̄ + 1}. (13)

After the SIC decoding, the instantaneous SINRs at D_{q} for the decoding of the information of its own information are expressed by

γ_{SD_{q}} = a_{q}G_{SD_{q}}γ,  (14)

γ_{Rm̄D_{q}} = \frac{a_{q}γ^{2}G_{SRm̄}G_{Rm̄D_{q}}}{γ (G_{SRm̄} + G_{Rm̄D_{q}}) + γG_{Rm̄D_{q}}γm̄ + γm̄ + 1}.  (15)

Similarly to (10), according to the SC, the instantaneous SINR of D_{q} can be expressed as

γ_{D_{q}} = max (γ_{SD_{q}}, γ_{Rm̄D_{q}}) .  (16)

2.2 Relay selection rule In the proposed scheme, each relay checks their battery at the beginning of a transmission block and judges whether it has enough energy to forward the source information. For those relays with sufficient energy, they report their CSI to S for making the relay selection decision. All the other relays apart from Rm̄ perform the EH in this time block. We assume each relay accumulates the harvested energy using a finite energy storage with capacity of BC (B=1, 2, 3...). We also define C as the EH threshold and P T as the energy required for the date forwarding. The relay can harvest uC amount of energy in one time slot and u is given by

uC ≤ E_{R_{i}} < (u + 1) C. (u = 0, 1, 2.....)  (17)

This assumption is closer to the practical scenario, and the evolution of the battery status of each relay can be modeled as a finite-state Markov chain. Using the transition probability matrix of this chain, we can get the steady-state probability vector which can be used for analyzing the performance.

Due to the lack of general form of steady-state probability, the above analysis is computationally intense when M is large. To facilitate the computation, we propose an approximate approach based on two simplified assumptions. Firstly, we denote the energy amount of relay at the selection epoch as a random variable Z. The exact distribution of Z is very complex. To ease the computation, we approximate Z by a uniform random variable over [0, BC]. The adopted approximation is inspired by the fact that the amount of harvested energy in a transmission block follows the geometric distribution with parameter 1/2 [21]. In general, these conditions may not always hold. We will discuss the effectiveness of the assumption through numerical results. Secondly, we found that an arbitrary relay may be either short of enough power to participate in the relay selection or otherwise, so the evolution of energy amount of relay is captured by using two states, either active or inactive. With this simplified two-state Markov chain, a relay is in s_{0} if the relay lacks sufficient energy to transmit, or in s_{1} when the relay has enough energy for transmission. Next, we will explain how to obtain the transition probability matrix of the two-state Markov chain.

The transition from state s_{0} to s_{0} happens when a relay has no enough energy to transmit (i.e., Z < V C) in the current block and the accumulated energy after the EH in a time block remains below the given threshold. The corresponding transition probability is given by

p_{0,0} = Pr (Z + E_{R_{i}} < V C| 0 ≤ Z < V C) = Pr (Ẑ + E_{R_{i}} < V C) ,  (18)

where Ẑ is a truncated random variable defined as

Ẑ = \{ \begin{matrix}Z  Z < V C \\ 0  Z ≥ V C \end{matrix} .  (19)

Since Z is approximated as uniformly distributed, the probability density function (PDF) of Ẑ can be obtained as

fẐ (z) = \frac{1}{BC} u(V C − z) + (1 − \frac{V}{B} ) δ (z − V C) ,  (20)

where u (·) and δ (·) denote the unit step function and the Dirac delta function, respectively. Considering the channels obey Nakagami-m distribution, a general random variable h_{Λ}, Λ ∈ {SD_{j}, SR_{i}, R_{i}D_{j}}, is subject to a Gamma distribution, so the corresponding PDF is

f_{G_{Λ}}(x) = (\frac{m_{Λ}}{Φ_{Λ}} )^{m_{Λ}} \frac{x^{m_{Λ}−1}}{Γ(m_{Λ})} exp (−\frac{m_{Λ}x}{Φ_{Λ}}) ,  (21)

where Γ (·) is Gamma function [22]. The corresponding cumulative distribution function (CDF) is expressed by

F_{G_{Λ}} (x) = 1 − exp (−\frac{m_{Λ}x}{Φ_{Λ}}) ∑_{n=0}^{m_{Λ}−1} \frac{(m_{Λ}x/Φ_{Λ})^{n}}{n!} .  (22)

If m_{Λ} takes positive integer values, then with the help of [22, (3.351)], (18) can be calculated as

\begin{aligned}p_{0,0} & = ∫_{0}^{V C} Pr (G_{SR_{i}} < \frac{V C − z}{ηP T} )fẐ (z) dz \\ & = 1 − \frac{1}{BC} ∑_{n=0}^{m_{1}−1} (\frac{1}{m_{1}/Ē_{R_{i}}} − exp (−\frac{m_{1}V C}{Ē_{R_{i}}} ) ∑_{k=0}^{n} \frac{n!(V C)^{k}}{k!(m_{1}/Ē_{R_{i}})^{1−k}} ), \end{aligned} (23)

where Ē_{R_{i}} = E (E_{R_{i}}) = P ηΦ_{SR}.

The transition from s_{0} to s_{1} happens when the relay enters the EH mode in the current bock and the accumulated energy exceeds V C. Hence, we have

p_{0,1} = Pr (Z + E_{R_{i}} ≥ V C| 0 ≤ Z < V C) = Pr (Ẑ + E_{R_{i}} ≥ V C) .  (24)

Similar to the derivation of (22), we can obtain

p_{0,1} = \frac{1}{BC} ∑_{n=0}^{m_{1}−1} (\frac{1}{m_{1}/Ē_{R_{i}}} − exp (−\frac{m_{1}V C}{Ē_{R_{i}}} ) ∑_{k=0}^{n} \frac{n!(V C)^{k}}{k!(m_{1}/Ē_{R_{i}})^{1−k}} ). (25)

If the relay which is in state s_{1} is selected as the best relay for the data forwarding and its harvested energy transits from state s_{1} to state s_{0}. In this case, the transition probability p_{1,0} is then given by

p_{1,0} = Pr (Z − V C < V C| V C ≤ Z ≤ BC) Pr (Rm̄ = R_{i})
= Pr (Z′ < 2V C) Pr (Rm̄ = R_{i}) ,

(26)

where Z′ is a truncated random variable defined as

Z′ = \{ \begin{matrix}Z & V C ≤ Z ≤ BC \\ 0 & Z < V C \end{matrix} . (27)

Since Z is uniformly distributed, the PDF of Z′ can be obtained as

f_{Z′} (z) = \frac{1}{BC} (u(z − BC) − u(z − V C)) + \frac{V}{B} δ (z − V C) .  (28)

Using this PDF, the first term on the right-hand side of (26) can be obtained as

Pr (Z′ < 2V C) = \{ \begin{matrix}2V /B & 2V C < BC \\ 1 & 2V C ≥ BC \end{matrix} . (29)

As to the second term, we note that the i. i. d. fading assumption implies each relay in ϕ has an equal chance to be selected as the best relay. To simplify the analysis, we approximate the cardinality of ϕ by its mean such that

Pr (Rm̄ = R_{i}) ≈ \frac{1}{M̄_{e}} = \frac{1}{Mv_{1}} ,  (30)

where M̄_{e} is the average number of relays in ϕ and v_{1} is the steady-state probability of state s_{1}. Combining (29) and (30), we can obtain the closed-form of p_{1,0}.

The transition probability from s_{1} to s_{1} can be derived in the similar manner as the previous case, so we omit the derivation here.

When the two-state Markov chain is formulated, the steady-state probability vector can be directly obtained as

v = (v_{0}, v_{1}) = (\frac{p_{1,0}}{p_{0,1} + p_{1,0}} , \frac{p_{0,1}}{p_{0,1} + p_{1,0}} ) .  (31)

We note that both v_{0} and v_{1} involve p_{1,0}, which is a function of v_{1}. By substituting p_{0,1} and p_{1,0} into (31), v_{1} can be solved explicitly. Take the condition 2V C < BC in (29), for example, v_{1} can be derived in closed form as

v_{1} = 1 − \frac{2V}{BM} (\frac{1}{BC} ∑_{n=0}^{m_{i}−1} (\frac{1}{m_{1}/Ē_{R_{i}}} − exp (−\frac{m_{1}V C}{Ē_{R_{i}}} ) ∑_{k=0}^{n} \frac{n!(V C)^{k}}{k!(m_{1}/Ē_{R_{i}})^{1−k}} ))^{−1}.

(32)

With v_{1} at hand, we can get Pr (|ϕ| = Ω) (Ω = 0, 1, 2...), which follows the binomial distribution with the probability mass function given as

Pr (|ϕ| = Ω) = \binom{M}{Ω}(v_{1})^{Ω}(1 − v_{1})^{M−Ω},  (33)

where \binom{a}{b} = \frac{a!}{(a−b)!b!}. With the close form of Pr (|ϕ| = Ω), we can analyze the performance in the following part.

3 Performance analysis

3.1 Outage behavior According to the total probability law, the outage probability of D_{p} or D_{q} can be expressed as

P_{out}^{Γ} = ∑_{Ω=1}^{M} Pr (|ϕ| = Ω)P_{out}^{Γ,Ω},  (34)

where Γ = {p, q}. As Pr (|ϕ| = Ω) can be got in (33), the next step is to get P_{out}^{Γ,Ω}, which will be given in the following.

The scheme would be in outage if the SINR falls below the pre-set threshold, so the outage probability of D_{p} is given by

P_{out}^{p,Ω} = Pr (γ_{D_{p}} < γ_{thp}) = Pr (max (γ_{SD_{p}}, γ_{Rm̄D_{p}}) < γ_{thp})

= \underbrace{Pr (γ_{SD_{p}} < γ_{thp})}_{J_{1}} \underbrace{Pr (γ_{Rm̄D_{p}} < γ_{thp})}_{J_{2}}, (35)

where γ_{thp} is the threshold of D_{p}. Next, according to order statistics [23], the PDF and CDF of the sorted random variable ĥ_{SD_{p}} are expressed, respectively, by

fĜ_{SDp} (x) = ∑_{j=0}^{N−p} (−1)^{j}Z_{p}^{N} \binom{N − p}{j}f_{G_{SD_{j}}} (x)(F_{G_{SD_{j}}} (x) )^{p+j−1},  (36)

FĜ_{SDp} (x) = ∑_{j=0}^{N−p} \frac{(−1)^{j}Z_{p}^{N}}{p + j} \binom{N − p}{j}_{(F_{G_{SD_{j}}} (x) )}^{p+j},  (37)

where Z_{p}^{N} = N!/((N − p)! (p − 1)!). Based on the above definitions, J_{1} is calculated by

\begin{aligned}J_{1} & = Pr (\frac{a_{p}Ĝ_{SD_{p}}γ}{a_{q}Ĝ_{SD_{p}}γ + 1} < γ_{thp}) = ∫_{0}^{τ} fĜ_{SDp} (x)dx \\ & = ∑_{j=0}^{N−p} Z_{p}^{N} \frac{(−1)^{j}}{p + j} \binom{N − p}{j}(1 − exp (−\frac{m_{0}τ}{Φ_{SD}} ) ∑_{n=0}^{m_{0}−1} \frac{(m_{0}τ/Φ_{SD})^{n}}{n!} )^{p+j}, \end{aligned}

(38)

where τ = γ_{thp}/γ (a_{p} − a_{q}γ_{thp}). It is noteworthy that (38) exists if and only if γ_{thp} < a_{p}/a_{q}, otherwise the probability will be always be one. J_{2} can be simplified as

\begin{aligned}J_{2} = & Pr (G_{Rm̄D_{p}} ≤ τ) \\ & + Pr (G_{Rm̄D_{p}} > τ, h_{SRm̄} ≤ \frac{τ (γG_{Rm̄D_{p}}γm̄ + γG_{Rm̄D_{p}} + γm̄ + 1)}{γ (G_{Rm̄D_{p}} − τ)}) \\ = & ∫_{0}^{τ} f_{G_{Rm̄D_{p}}} (x) dx + ∫_{τ}^{∞} f_{G_{Rm̄D_{p}}} (y) dy ∫0^{\frac{τ(γyγm̄+γy+γm̄+1)}{γ(y−τ)}} f_{G_{SRm̄}} (x)dx. \end{aligned} (39)

To solve (39), we first need to calculate f_{G_{SRm̄}} (x), which can be get according to (36) when p = 1 and N = Ω. Using f_{G_{SRm̄}} (x) and with the help of [22, (3.471.9)], J_{2} can be formulated as

J_{2} =1 − exp (−\frac{m_{2}τ}{Φ_{RD}} ) ∑_{n=0}^{m_{2}−1} \frac{(m_{2}τ/Φ_{RD})^{n}}{n!}

+ \frac{Z_{1}^{Ω}}{Γ(m_{2})}(\frac{m_{2}}{Φ_{RD}} )^{m_{2}} exp (−\frac{m_{2}τ}{Φ_{RD}} )

× ∑_{j=0}^{Ω−1} \frac{(−1)^{j}}{j + 1} \binom{Ω − 1}{j} ∑_{r=0}^{j} (−1)^{r} exp (−\frac{m_{1}rum̄}{Φ_{SR}} )

× ∑_{r_{i}≥0,r_{0}+r_{1}+...+r_{m_{1}−1}=r} ∏_{z=0}^{m_{1}−1} (z!)^{r_{z}}^{\binom{r}{r_{0},r_{1},...,r_{m_{1}−1}}(\frac{m_{1}um̄}{Φ_{SR}} )^{∑{m_{1}−1 \atop z=0} zr_{z}}}

× ∑_{k_{1}=0}^{∑_{z=0}^{m_{1}−1} zr_{z}} ∑_{k_{2}=0}^{m_{2}−1} \binom{∑_{z=0}^{m_{1}−1} zr_{z}}{k_{1}}\binom{m_{2} − 1}{k_{2}}d^{k_{1}}τ^{m_{2}−1−k_{2}}

× 2(\frac{m_{1}rum̄dΦ_{RD}}{Φ_{SR}m_{2}} )^{\frac{k_{2}−k_{1}+1}{2}} K_{k_{2}−k_{1}+1} (\sqrt{\frac{2m_{1}rum̄dm_{2}}{Φ_{SR}Φ_{RD}}} ) ,

(40)

where um̄ = τ (γ_{m̄} + 1), d = τ + 1/γ, Z_{1}^{Ω} = Ω!/(Ω − 1)!, \binom{r}{r_{0},r_{1},...,r_{m_{1}−1}} =

\frac{r!}{r_{0}!,r_{1}!,...,r_{m_{1}−1}!}, ∑_{r_{i}≥0,r_{0}+r_{1}+...+r_{m_{1}−1}=r} means the sum of all the possible values of

which satisfies r_{0} + r_{1} + ... + r_{m_{1}−1} = r and K_{v} (·) expresses the vth-order modified Bessel function of the second kind. By substituting (33), (38) and (40) into (34), the outage probability of D_{p} can be obtained.

Next, the analysis of an outage event at D_{q} is presented. Owing to the mechanism of SIC, the outage probability of D_{q} is given as

P_{out}^{q,Ω} = Pr (γ_{D_{q}} < γ_{thq})
= Pr (γ_{SD_{q}} < γ_{thq}) Pr (γ_{RD_{q}} < γ_{thq})
= \underbrace{(1 − Pr (γ_{SD_{q→p}} ≥ γ_{thp}) Pr (γ_{SD_{q}} ≥ γ_{thq}))}_{J_{3}}

× \underbrace{(1 − Pr (γ_{Rm̄D_{q→p}} ≥ γ_{thp}) Pr (γ_{Rm̄D_{q}} ≥ γ_{thq}))}_{J_{4}},

(41)

where γ_{thq} is the threshold of D_{q}. J_{3} can be expanded as

J_{3} =1 − Pr (\frac{a_{p}γĜ_{SD_{q}}}{a_{q}γĜ_{SD_{q}} + 1} ≥ γ_{thp}) Pr (γa_{q}Ĝ_{SD_{q}} ≥ γ_{thq})
=1 − Pr (Ĝ_{SD_{q}} ≥ max (\frac{γ_{thp}}{γ (a_{p} − a_{q}γ_{thp})}, \frac{γ_{thq}}{γa_{q}} ) =^{Δ} β)
= ∑_{j=0}^{N−q} Z_{p}^{N} \frac{(−1)^{j}}{q + j} \binom{N − q}{j}(1 − exp (−\frac{m_{0}β}{Φ_{SD}} ) ∑_{n=0}^{m_{0}−1} \frac{(m_{0}β/Φ_{SD})^{n}}{n!} )^{q+j}.

(42)

Next, J_{4} can be written as

J_{4} =1 − Pr (\frac{a_{p}γ^{2}G_{SRm̄}G_{Rm̄D_{q}}}{γ (G_{SRm̄} + G_{Rm̄D_{q}}) + γG_{Rm̄D_{q}} (a_{q}γG_{SRm̄} + γ_{Rm̄}) + γ_{m̄} + 1} ≥ γ_{thp})

× Pr (\frac{a_{q}γ^{2}G_{SRm̄}G_{Rm̄D_{q}}}{γ (G_{SRm̄} + G_{Rm̄D_{q}}) + γG_{Rm̄D_{q}}γ_{m̄} + γ_{m̄} + 1} ≥ γ_{thq})

=1 − ∫_{β}^{∞} ∫ ∞_{\frac{β(γyγm̄ +γy+γm̄ +1)}{γ(y−β)}} f_{G_{SRm̄}} (x) f_{G_{Rm̄D_{q}}} (y) dxdy

= ∫_{0}^{β} f_{G_{Rm̄D_{q}}} (y) dy+ ∫_{β}^{∞} f_{G_{Rm̄D_{q}}} (y) dy ∫0^{\frac{β(γyγm̄ +γy+γm̄ +1)}{γ(y−β)}} f_{G_{SRm̄}} (x) dx.

(43)

The derivation of (43) is the same as (39), so we omit it here. So far, we have derived the outage probability of D_{q} by substituting (33), (42) and (43) into (34).

3.2 Ergodic sum-rate In this section, we analyse the ergodic sum-rate for a pair of users, which is given by

R_{sum} = R_{ave}^{p} + R_{ave}^{q} = E (log_{2} (1 + γ_{D_{p}})) + E (log_{2} (1 + γ_{D_{q}})) .  (44)

The rate of D_{p} can be calculated in the following, where approximation is attained under the condition of γ → ∞

R_{ave}^{p} = E (log_{2} (1 + γ_{D_{p}}))

= E (log_{2} (1 + max (γ_{SD_{p}}, γ_{Rm̄D_{p}})))

≈ E (log_{2} (1 + max (\frac{a_{p}}{a_{q}} , \frac{a_{p}γG_{SRm̄}}{a_{q}γG_{SRm̄}} ))) = log_{2} (\frac{a_{p}}{a_{q}} ) .

(45)

The rate of D_{q} can be calculated as

R_{ave}^{q} =E (log_{2} (1 + γ_{D_{q}}))
= Pr (γ_{SD_{q→p}} ≥ γ_{thp}) Pr (γ_{Rm̄D_{q→p}} ≤ γ_{thp}) \underbrace{E (log_{2} (1 + γ_{SD_{q}}))}_{J_{5}}

+ Pr (γ_{SD_{q→p}} ≤ γ_{thp}) Pr (γ_{Rm̄D_{q→p}} ≥ γ_{thp}) \underbrace{E (log_{2} (1 + γ_{Rm̄D_{q}}))}_{J_{6}}

+ Pr (γ_{SD_{q→p}} ≥ γ_{thp}) Pr (γ_{Rm̄D_{q→p}} ≥ γ_{thp})

× \underbrace{E (log_{2} (1 + max (γ_{SD_{q}}, γ_{Rm̄D_{q}})))}_{J_{7}} .

(46)

All the probabilities in the (46) can be obtained according to the above analysis, we only need derive J_{5}, J_{6} and J_{7} in the following. Firstly, J_{5} can be written as

J_{5} = E (\frac{1}{ln 2} ln (1 + a_{q}Ĝ_{SD_{q}}γ))

= \frac{1}{ln 2} ∫_{0}^{∞} ln (1 + a_{q}γx) fĜ_{SDq} (x) dx

= \frac{1}{ln 2} ∫_{0}^{∞} \frac{1 − FĜ_{SDq} (x/a_{q}γ)}{1 + x}dx.

(47)

Using [22, (3.351.2)] and [22, (3.353.5)], the closed form of J_{5} can be get

J_{5} =\frac{Z_{q}^{N}}{ln 2Γ(m_{0})} ∑_{j=0}^{N−q} (−1)^{j}\binom{N − q}{j} ∑_{r=0}^{q+j−1}

× ∑_{r_{i}≥0,r_{0}+r_{1}+...+r_{m_{2}−1}=r} (−1)^{r} \begin{matrix}( & {r \atop r_{0},r_{1},...,r_{m_{2}−1}} & ) \\ ∏ & {m_{2}−1 \atop z=0} & (z!)^{r_{z}} \end{matrix} (∑_{z=0}^{m_{0}−1} zr_{z} + m_{0} − 1)!

× ∑_{t=1}^{∑_{z=0}^{m_{0}−1} zr_{z}+m_{0}−1} \frac{Φ_{SD}^{−t}}{t!(a_{q}γ)^{t}m_{0}^{−t}(r + 1)^{∑{m_{0}−1 \atop z=0} zr_{z}+m_{0}−t}}

× ((−1)^{t−1} exp (\frac{m_{0} (r + 1)}{a_{q}γΦ_{SD}}) Ei (−\frac{m_{0} (r + 1)}{a_{q}γΦ_{SD}})

+ ∑_{o=1}^{t} (o − 1)!(−1)^{t−o}(\frac{m_{0} (r + 1)}{a_{q}γΦ_{SD}})^{−o}) .

(48)

where E_{i}(x) = ∫_{x}^{∞} t^{−1} exp(−t)dt denotes the exponential integral function [22, (5.1.1)]. Secondly, In order to calculate J_{6}, we need to obtain the asymptotic ex-

Guo and Sun Page 13 of 18

pression as

J_{6} = E (log_{2} (1 + \frac{a_{q}γ^{2}G_{SRm̄}G_{Rm̄D_{q}}}{γ (G_{SRm̄} + G_{Rm̄D_{q}}) + γG_{Rm̄D_{q}}γ_{m̄} + γ_{m̄} + 1}))

= E (log_{2} (1 + \frac{a_{q}γ^{2} \frac{G_{SRm̄}}{γ_{m̄} +1} G_{Rm̄D_{q}}}{\frac{γG_{SRm̄}}{γ_{m̄} +1} + γG_{Rm̄D_{q}} + 1}))

< E (log_{2} (1 + a_{q}γ min (G_{Rm̄D_{q}}, \frac{G_{SRm̄}}{γ_{m̄} + 1}))) .

(49)

Now denoting a random variable

L = a_{q}γ min (G_{Rm̄D_{q}}, \frac{G_{SRm̄}}{γ_{m̄} + 1}) .  (50)

Thus an upper bound of J_{6} can be rewritten as

E (log_{2} (1 + γ_{Rm̄D_{q}})) < E (\frac{1}{ln 2} ln (1 + l)) = \frac{1}{ln 2} ∫_{0}^{∞} \frac{1 − F_{L} (l)}{1 + l}dl.  (51)

Submitting the PDFs of G_{SRm̄} and G_{Rm̄D_{q}} into (51) and using [22, (3.353.5)], we can get

J_{6} <\frac{Z_{1}^{Ω}}{ln 2Γ(m_{1})} ∑_{n=0}^{m_{2}−1} \frac{(m_{2}/Φ_{RD}a_{q}γ)^{n}}{n!} ∑_{j=0}^{Ω−1} (−1)^{j}\binom{Ω − 1}{j}

× ∑_{r=0}^{j} ∑_{r_{i}≥0,r_{0}+r_{1}+...+r_{m_{1}−1}=r} \begin{matrix}(−1) & r \\ ∏ & {m_{1}−1 \atop z=0} & (z!)^{r_{z}} \end{matrix}

× ∑_{o=0}^{∑{m_{1}−1 \atop z=0} zrz +m_{1}−1} \frac{(∑_{z=0}^{m_{1}−1} zr_{z} + m_{1} − 1)!((γ_{m̄} + 1))^{o}Φ_{SR}^{−o}}{o!(a_{q}γ)^{o}((r + 1))^{∑_{z=0}^{m_{1}−1} zr_{z}+m_{1}−o}m_{1}^{−o}}

× ((−1)^{n+o−1} exp (\frac{m_{1} (γ_{m̄} + 1) (r + 1)}{Φ_{SR}a_{q}γ} + \frac{m_{2}}{a_{q}γΦ_{RD}} )

× Ei (−\frac{m_{1} (γ_{m̄} + 1) (r + 1)}{Φ_{SR}a_{q}γ} − \frac{m_{2}}{a_{q}γΦ_{RD}} )

+ ∑_{ε=1}^{n+o} (ε − 1)!(−1)^{n+o−ε}(\frac{m_{1} (γ_{m̄} + 1) (r + 1)}{Φ_{SR}a_{q}γ} + \frac{m_{2}}{a_{q}γΦ_{RD}} )^{−ε}) .

(52)

Thirdly, J_{7} can be rewritten as

J_{7} = E (log_{2} (1 + max (a_{q}γG_{SD_{q}}, \frac{a_{q}γ^{2} \frac{G_{SRm̄}}{γ_{Rm̄Rm̄} +1}G_{Rm̄D_{q}}}{γ \frac{G_{SRm̄}}{γ_{Rm̄Rm̄} +1} + γG_{Rm̄D_{q}} + 1})))

< E (log_{2} (1 + max (a_{q}γG_{SD_{q}}, a_{q}γ min (\frac{G_{SRm̄}}{γ_{Rm̄Rm̄} + 1}, G_{Rm̄D_{q}})))) .

(53)

As the derivation process of (49), we denote the random variable

W = max (a_{q}γh_{SD_{q}}, a_{q}γ min (\frac{G_{SRm̄}}{γ_{Rm̄Rm̄} + 1}, G_{Rm̄D_{q}})) .  (54)

Then, (53) can be obtained as

E (log_{2} (1 + max (γ_{SD_{q}}, γ_{Rm̄D_{q}}))) < E (\frac{1}{ln 2} ln (1 + W)) = \frac{1}{ln 2} ∫_{0}^{∞} \frac{1 − F_{W} (w)}{1 + w}dw.

(55)

The derivation of (55) is the same as (49), so we omit it here. By submitting (48), (52), (55) and all the probabilities into (46), we can get the close form expression of R_{ave}^{q}. With the R_{ave}^{p} and R_{ave}^{q}, we can get the expression of R_{sum}.

4 Results and discussion

In this section, Monto-Carlo simulations are performed to validate our theoretical analysis. In all the simulations, the energy conversion efficiency η = 0.5, the fixed transmission rate of the source is 1 bit per channel use (bpcu). We set C as the multiple of the source transmission energy, i.e., C = δP T, where δ > 0 is the scaling factor. Unless stated otherwise, the following sets of parameters are used: N = 3, σ^{2} = 1, m_{0} = m_{1} = m_{2} = 2, Φ_{SD} = Φ_{SR} = Φ_{RD} = 1, a_{p} = 8/9, a_{q} = 1/9, γ_{thp} = 1dB, γ_{thq} = 2dB, p = 1, q = 3.

Figure 2 Outage probability with different scheme (SI=-10dB, δ = 0.5, M=6, B = 10).

4.1 Results 4.1.1 Numerical and simulation results of outage behavior Fig. 2 shows the performance of our proposed scheme versus the transmitting SNR of source. We can see that the approximate analysis of finite Markov chain is accurate in the whole SNRs. Compared with the HD and conventional OMA relaying scheme, the outage behavior of each terminal in the FD relaying system is better, even in the existence of SI at the relay. We can also find that the diversity order for each user is equal to that in the HD cooperative NOMA system. Since the theoretical analyses agree well with the simulations in medium and high SNR ranges, we will only plot the analytical results in the remaining outage probability figures.

In Fig. 3, we investigate the impacts of SI on the outage probability. We can see that that the outage performance of the FD relaying (only take one user, i.e. D_{p}, e.g. and vice versa) system will become worse than that in the HD relaying network for high-SNR values when SI is extremely severe. But both two schemes can have better outage performance than using the conventional OMA scheme. Also, two users can be served simultaneously. Therefore, the user fairness of the FD cooperative NOMA system can be improved.

Figure 3 Outage probability with different SI (δ = 0.5, M=6, B= 10).

In Fig. 4, we illustrate the impacts of δ on the performance of the proposed scheme with different SI in medium SNR conditions (P /σ^{2}= 20dB). For all the curves, the trends are the same for all schemes, namely, the probability first decreases then increases as δ varies from 0 to 1. This means that, when the other parameters are determined, there must be an optimal value of δ. However, the values of the inflection points are not always the same for different users and SI. We can see that for the user with better channel condition (D_{q}), the optimal value is smaller than user with worse channel condition (D_{p}), this is because D_{p} has the worse direct link, so it needs more power to enhance the relay link. We can also see that for the FD cooperative NOMA system, the optimal value shifts to the left with the increase in the level of SI and when δ is too big, the outage behavior is nearly the same under different values of SI. More importantly, compared to the HD cooperative NOMA system, the outage behavior of D_{p} is getting inferior since the value of SI is large. However, the optimal value of δ can easily be obtained by a one-dimensional exhaustive search, with this optimal value of δ, the system can resist fading more effectively.

4.1.2 Numerical and simulation results of ergodic sum-rate Fig. 5 presents the rates of the two paired users varies with system SNR. Firstly, the sum-rate increases with the growth of transmitting SNR, and asymptotic curves of the ergodic sum-rate and simulation results are almost overlapping in the medium-

Figure 4 Outage probability against (M=6, B = 10).

and high-SNR regions, we will only plot the analytical results in the remaining outage rate figures. Secondly, the results clearly show that the FD cooperative NOMA system can achieve a larger sum-rate compared with the HD relaying system and conventional OMA system, but the double rate cannot be achieved owing to the effect of SI. More importantly, with the increase of system SNR, the rate of D_{p} approaches a certain value, whereas D_{q} can obtain improved performance. Both the two conclusions are consistent with the derivations in (45) and (46).

Figure 5 Rate with different user and different scheme (SI=-10dB, δ = 0.5, M=6, B = 10).

In Fig. 6, we investigate the impact of battery size on the proposed finite battery scheme by varying the relay number and SI. we can observe that the performance increases as B increases. However, the gain provided by a larger battery size does not increase when exceeds a certain value and this value is approximate to 6. From this figure, we can also see that the battery can not offer enough energy for forwarding the date when B = 1, so the performances are the same for all curves and the performance gap for different SI increases as the relay number increases.

Figure 6 Sum-rate against B (δ = 0.5).

4.2 Discussion Compared with existing literatures, the channel in this paper we use is more in line with the actual communication scenario and the model we propose combines NOMA, cooperative communication, FD and EH together reasonably and can improve spectral efficiency and energy efficiency simultaneously. However, there are still some shortages in our proposed scheme. The joint optimization of multiple parameters such as the transmission power, power splitting factor and user pairing has not been mentioned in this paper, this joint optimization can make the proposed schene more useful in practical application which we will study in the future.

5 Conclusions

In this paper, we have proposed a relay selection scheme for the FD-NOMA networks with EH over Nakagami-m fading channels. We modeled the amount of harvested energy at each relay using a finite Markov chain and then derived the approximate closed-form expression of the outage probability and ergodic sum-rate. Simulations were carried out to verify the effectiveness of our theoretical analysis. We concluded that through carefully choosing the system parameters of the network, (e.g., EH threshold), FD-NOMA performs much better than the HD-NOMA or the conventional OMA.

Abbreviations FD: full-duplex; NOMA: Non-orthogonal multiple access; CSI: Channel state information; SI: Self-interference; OMA: Orthogonal multiple access; SIC: Successive interference cancellation; MMSE: Minimum-mean-squared-error; T-ARQ: Truncated automatic repeat-request; AF: Amplify-and-forward; DF: Decode-and-forward; EH: Energy harvesting; RF: Radio frequency; HD: Half-duplex; INR: Interference-to-noise ratio; AWGN: Additive white Gaussian noise; SNR: Signal-to-noise-ratio; SINR: Signal-to-interference-plus-noise-ratio; SC: selection combining; PDF: Probability density function; CDF: Cumulative distribution function Acknowledgements Not applicable.

Availability of data and materials Not available online. Please contact corresponding author for data requests.

Funding This work was supported by the Shandong Provincial Natural Science Foundation, China (Grant No. ZR2017MD019).

Author’s contributions WG proposed the system model, analyzed the performance and drafted the paper. LS participated in the simulations and helped to draft the manuscript. All authors read and approved the final manuscript.

Competing interests The authors declare that they have no competing interests.

Author details 1School of Information Science and Engineering, Shandong University, 72 Binhai Road, 266237 Qingdao, China. 2School of Physics and Physical Engineering, Qufu Normal University, 57 Jingxuanxi Road, 273165 Qufu, China.

References

1. Y. Saito, A. Benjebbour, Y. Kishiyama, T. Nakamura, System level performance evaluation of downlink non-orthogonal multiple access (NOMA). in Proceedings of the IEEE Annual International Symposium on Personal, In door and Mobile Radio Communications (PIMRC), London, United Kingdom, 2013.

2. Z. Ding, Z. Yang, P. Fan, H. V. Poor, On the performance of non-orthogonal multiple access in 5G systems with randomly deployed users. IEEE Signal Process. Lett. 21(12), 1501-1505 (2014).

3. Y. Lan, A. Benjebboiu, X. Chen, A. Li, H. Jiang, Considerations on downlink non-orthogonal multiple access (NOMA) combined with closed-loop SU-MIMO. in Proceedings of the IEEE International Conference on Signal Process and Communication System (ICSPCS), Gold Coast, Australia, 2014.

4. X. Chen, A. Benjebbour, Y. Lan, A. Li, H. Jiang, Impact of rank optimization on downlink non-orthogonal multiple access (NOMA) with SU-MIMO. in Proceedings of the IEEE International Conference on Communication System (ICCS), Macau, China, 2014.

5. Z. Ding, M. Peng, H. V. Poor, Cooperative non-orthogonal multiple access in 5G systems. IEEE Commun. Lett. 19(8),1462-1465 (2015).

6. Z. Ding, H. Dai, H. V. Poor, Relay selection for cooperative NOMA. IEEE Wireless Commun. Lett. 5(4), 416-419 (2016).

7. Z. Yu, C. Zhai, W. Ni, D. Wang, Non-orthogonal multiple access with cooperative truncated ARQ and relay selection. IEEE Access, 7, 56228-56243 (2019).

8. Z. Zhang, X. Chai, K. Long, A. V. Vasilakos, L. Hanzo, Full duplex techniques for 5G networks: Self-interference cancellation, protocol design, and relay selection. IEEE Commun. Mag. 53(5), 128-137 (2015).

9. Q. Wang, Y. Dong, X. Xu, X. Tao, Outage probability of full-duplex AF relaying with processing delay and residual self-interference. IEEE Commun. Lett. 19(5), 783-786 (2015).

10. D. P. M. Osorio, E. E. B. Olivo, H. Alves, J. C. S. S. Filho, M. Latva-Aho, Exploiting the direct link in full-duplex amplify-and-forward relaying networks. IEEE Signal Processing Lett. 22(10), 1766-1770 (2015).

11. T. Kwon, S. Lim, S. Choi, D. Hong, Optimal duplex mode for DF relay in terms of the outage probability. IEEE Trans. Veh. Technol. 59(7), 3628-3634 (2010).

12. A. A. Nasir, X. Zhou, S. Durrani, R. A. Kennedy, Relaying protocols for wireless energy harvesting and information processing. IEEE Trans. Wireless Commun. 12(7), 3622-3636 (2013).

13. Z. Ding, S. M. Perlaza, I. Esnaola, H. V. Poor, Power allocation strategies in energy harvesting wireless cooperative networks. IEEE Trans Wireless Commun. 13(2), 846-860 (2014).

14. H. Gao, W. Ejaz, M. Jo, Cooperative Wireless Energy Harvesting and Spectrum Sharing in 5G Networks. IEEE Access. 4, 3647-3658 (2016).

15. S. Sudevalayam, P. Kulkarni, Energy harvesting sensor nodes: Survey and implications. IEEE Commun. Surveys Tuts. 13(3), 443-461 (2011).

16. C. Huang, J. Zhang, P. Zhang, S. Cui, Threshold-based transmissions for large relay networks powered by renewable energy. in Proceedings of the IEEE Global Communications Conferences (ICC), Budapest, Hungary, 2013.

17. Z. Zhou, M. Peng, Z. Zhao, W. Wang, R. S. Blum, Wireless-powered cooperative communications: Power-splitting relaying with energy accumulation. IEEE J. Select. Areas Commun. 34(4), 969-982 (2016).

18. C. Zhai, Z. Yu, X. Wang, High-efficient cooperative relaying with wireless powered source and relay. Computer Networks, 152, 199-209 (2019).

19. Z. Zhang, Z. Ma, Z. Ding, M. Xiao, G. K. Karagiannidis, Full-duplex two-way and one-way relaying: average rate, outage probability, and tradeoffs. IEEE Trans. Wireless Commun. 15(6), 3920-3933 (2016).

20. X. Yue, Y. Liu, S. Kang, A. Nallanathan, Z. Ding, Exploiting full/Half-duplex user relaying in NOMA systems. IEEE Trans. Commun. 66(2), 560-575 (2018).

21. V. Rego, Characterization of equilibrium queue length distributions in M/GI/I queues. Computers and Operations Research, 15(1), 7-17 (1988).

22. I. Gradshteyn, I. Ryzhik, Table of integrals, series and products. Academic Press: San Diego, CA (1994).

23. H. A. David, H. N. Nagaraja, Order statistics. 3rd edn, John Wiley, New York (2003).