{'fraction_non_alphanumeric': 0.04185530646011746, 'fraction_numerical': 0.023164057528318696, 'mean_word_length': 4.062126092959042, 'pattern_counts': {'":': 0, '<': 17, '<?xml version=': 0, '>': 15, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 17, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Optimal Resource Allocation for Full-Duplex IoT Systems Underlaying Cellular Networks with Mutual SIC NOMA

Device-to-device (D2D) and non-orthogonal multiple access (NOMA) are promising technologies to meet the challenges of the next generations of mobile communications in terms of network density and diversity for internet of things (IoT) services. This paper tackles the problem of maximizing the D2D sum-throughput in an IoT system underlaying a cellular network, through optimal channel and power allocation. NOMA is used to manage the interference between cellular users and full-duplex (FD) IoT devices. To this aim, mutual successive interference cancellation (SIC) conditions are identified to allow simultaneously the removal of the D2D devices interference at the level of the base station and the removal of the cellular users (CU) interference at the level of D2D devices. To optimally solve the joint channel and power allocation (PA) problem, a time-efficient solution of the PA problem in the FD context is elaborated. By means of graphical representation, the complex non-convex PA problem is efficiently solved in constant time complexity. This enables the global optimal resolution by successively solving the separate PA and channel assignment problems. The performance of the proposed strategy is compared against the classical state-of-the-art FD and HD scenarios, where SIC is not applied between CUs and IoT devices. The results show that important gains can be achieved by applying mutual SIC NOMA in the IoT-cellular context, in either HD or FD scenarios.

Abstract-Device-to-device (D2D) and non-orthogonal multiple access (NOMA) are promising technologies to meet the challenges of the next generations of mobile communications in terms of network density and diversity for internet of things (IoT) services. This paper tackles the problem of maximizing the D2D sum-throughput in an IoT system underlaying a cellular network, through optimal channel and power allocation. NOMA is used to manage the interference between cellular users and full-duplex (FD) IoT devices. To this aim, mutual successive interference cancellation (SIC) conditions are identified to allow simultaneously the removal of the D2D devices interference at the level of the base station and the removal of the cellular users (CU) interference at the level of D2D devices. To optimally solve the joint channel and power allocation (PA) problem, a time-efficient solution of the PA problem in the FD context is elaborated. By means of graphical representation, the complex non-convex PA problem is efficiently solved in constant time complexity. This enables the global optimal resolution by successively solving the separate PA and channel assignment problems. The performance of the proposed strategy is compared against the classical state-ofthe-art FD and HD scenarios, where SIC is not applied between CUs and IoT devices. The results show that important gains can be achieved by applying mutual SIC NOMA in the IoT-cellular context, in either HD or FD scenarios. communication, a node can send and receive simultaneously using the same frequency resource. The achieved gain can go up to a virtual two-fold increase in spectral efficiency (SE) compared to half-duplex (HD) send-then-receive systems. In return, a self-interference (SI) is incurred due to the transmitted signal looping back into the receiver, thus limiting its appeal compared to HD to the point where this latter may even outperform FD in some cases. Nonetheless, the improvement in antenna architecture and in SI cancellation circuitry dramatically reduces the residual self-interference (RSI) [3][4][5], advocating for the use of FD for future communication standards.
The increasing demand for higher spectral efficiency and massive connectivity for internet of things (IoT) steered research towards non-orthogonal multiple access (NOMA) techniques. The sharing of multiple devices over the same time and frequency resource enables important SE gains, lower latency communications, and increased number of connected IoT devices [6][7][8][9][10][11][12]. In power-domain NOMA, signals are differentiated in the power dimension [13], where superposition coding of users signals is used at the transmitter, and successive interference cancellation (SIC) is performed at the receiver side. At the level of a receiver, the message with the highest power is decoded first and then subtracted from the total received signal, then the message with the second highest power is extracted and so on until the user decodes its own message. Works such [14][15][16] used NOMA and FD for cooperative relaying as well as secrecy provision, calling on to virtual MIMO as physical layer security enabler for fifth generation centric IoT applications. In this paper, we study the resource allocation problem for D2D communications systems underlaying cellular networks using the NOMA technique coupled with FD transmission scenarios. This being said, the broader problem of resource optimization in IoT could be tackled via multiple other tools such as next generation reconfigurable intelligent surfaces (RIS) [17], quantum computing inspired metaheuristics [18], and much more.
A. Related Works is used in the D2D pairs, and CUs are grouped in NOMA clusters. The influence of the HD-D2D users over the SIC decoding orders of CUs is accounted for in both the block assignment and the PA phase, because the interference they generate may change the decoding order. However, NOMA SIC is not used to decode the interfering signals of the collocated D2D pairs. The same is true for [20], but additional power constraints are introduced on the D2D pairs to maintain the same SIC decoding orders at CUs as for the case of D2Ddisabled systems. The work in [21] introduces the concept of D2D group, where a D2D transmitter communicates with multiple D2D receivers via NOMA. To maximize the network sum-throughput, sub-channel allocation is conducted using many-to-one matching for CU-D2D grouping and optimal PA is approximated iteratively via successive convex approximation. When limiting the number of multiplexed D2Ds to one per CU user, the work in [22] provides a joint D2D-CU grouping and PA strategy for energy efficiency maximization: the Kuhn-Munkres technique is applied successively for channel allocation, while optimal PA is obtained using the Karush-Kuhn-Tucker conditions. In all the preceding studies, NOMA is applied either between the CU users [20], or between users of the same D2D group [21], [22], but the interference cancellation of the D2D signals at the level of CU users (and inversely) is not considered. At most, attention is given towards managing the SIC decoding order at the level of the CUs in [19], [20], or at the level of the D2D receivers in [21], [22].
The work in [23] tackles the problem of HD-D2D throughput maximization in an uplink system where NOMA is used between D2D and CU users. If the D2D causes strong interference on the BS, its signal can be decoded then subtracted before retrieving the CU signal. However, FD-D2D is not studied and SIC occurs only at the level of the BS, i.e. not at the devices levels. Besides, the information-theoretic conditions for SIC feasibility are not considered in the study. In [24], an efficient graph-based scheme is proposed to maximize the D2D sum-rate of an uplink system. To that end, an interlay mode is introduced to HD-D2D communication where a D2D pair can join a NOMA group to remove the interference between it and the cellular NOMA users. However, the conditions for applying SIC -and thus for determining the SIC decoding order -are only conditioned by the ascending order of channel gains between the senders and the receivers. In other words, the interfering signals that can be canceled are the ones that are attributed channel gains better than that of the useful signal, regardless of their power level at reception. This may lead to outage probabilities of one if no PA measures are taken to guarantee SIC stability as shown in [10]. The work in [25] incorporates NOMA into D2D cellular networks to maximize system connectivity. Unlike [24], the D2D NOMAaided modes are defined according to the SIC orders at the level of the D2D and the BS. The SIC decoding orders are governed by the strong interfering signal which is bound to the channel conditions as well as the used PA. The optimal PA and mode selection are solved in the presence of decoding signalto-interference-plus-noise ratios (SINRs) threshold constraints, then the user pairing problem is turned into a min-cost max-flow problem which is solved by the Ford-Fulkerson algorithm. However, the case of FD-D2D NOMA-aided network was not addressed neither in this study, nor in the entire literature combining NOMA and D2D.

B. Contributions
In this paper, we study the combination of NOMA with HD-D2D and FD-D2D systems using mutual SIC. We introduced the concept of mutual SIC in our previous works [26][27][28] where we showed that the signals of two or more users multiplexed in NOMA, and powered by distributed antennas, can be decoded and removed at the level of every user in the NOMA cluster. In a D2D cellular scenario, this translates into removing the interference of the D2D devices at the level of the BS, and the interference of CU users at the level of D2D pairs. The objective of the paper is to maximize the D2D sum-throughput, through joint optimal channel and power allocation, while maintaining the quality of service (QoS) requirement for all CU users. The main contributions of the paper are summarized as follows: • We derive the power multiplexing conditions (PMC) and SIC conditions allowing for the interference cancellation between D2D and CU users. The SIC constraints are the set of conditions that make mutual SIC feasible from the information theory perspective, i.e. the conditions on achievable rates at the respective levels of the users. The PMCs are the set of conditions that make the mutual SIC technique feasible from a practical implementation perspective, i.e. they guarantee that the power level of the message to be decoded, , at a given SIC iteration, is greater than the sum of the remaining messages, so that can be detected and discerned from the noise plus background interference. This guarantees SIC stability, since every signal is ensured to be the dominant signal during its decoding [9][10]. • We show that the PMCs imply the SIC conditions for both HD and FD transmission modes, which greatly reduces the PA problem complexity for the case of FD-SIC. • We solve analytically the PA problem for all the transmission strategies, especially for the case of FD-SIC where an efficient procedure is provided to optimally solve the D2D rate maximization problem with constant time complexity. • We show that the optimal solution of the joint PA and channel allocation problem can be achieved by successively resolving the PA problem and then the channel allocation problem. The remainder of the paper is organized as follows: section II presents the system model and formulates the joint channel and power allocation problem, decomposing the resource allocation into separate PA and channel allocation problems. The PA problems of FD and HD without SIC (FD-NoSIC, HD-NoSIC) are solved in section III, while the PA problem with SIC is reformulated for HD and FD (HD-SIC, FD-SIC) in section IV. Mutual SIC PA is solved for the case of HD transmission in section V. In sections VI, VII, and VIII the conditions of mutual SIC for FD-D2D are derived, the problem constraint reduction is performed, and then the proposed geometrical resolution is exposed, allowing for a cost-effective resolution of the FD-SIC PA problem. The channel allocation is discussed in section IX. Simulation results are presented in section X, and conclusions are drawn in section XI.

Direct Link
Interference Link The interference channel gains between a CU , on the one hand, and ,1 and ,2 on the other hand, are denoted by ℎ respectively. The direct link between the CU and the BS has a squared channel gain denoted by ℎ , . The message , transmitted by with power , reaches the BS with a power level ℎ , , and causes an interference level of ℎ ,1 , and ℎ ,2 , at ,1 and ,2 respectively. Each device , ( ∈ {1, 2}) of the th D2D pair can transmit a message , of power , to the other D2D user and suffers from both the interference of user and its RSI power , , , with , denoting the SI cancellation capability. The D2D inter-user channel gain is denoted by ℎ and the interference channel gains from ,1 and ,2 to the BS are denoted by ℎ , ,1 and ℎ , ,2 respectively. In this study, a frequency-non-selective channel is assumed, so that the channel gains are independent from the sub-band frequency and account only for large-scale fading including path-loss and shadowing. Table I contains the main notations used in the paper.
In this work, it is assumed that, prior to resource allocation and data exchange, a D2D discovery phase [29] takes place in the system, during which the D2D devices inform the BS about their desire to initiate a D2D link, and forward to the BS their estimates of the D2D-CU links (ℎ ,1 , , ℎ ,2 , ), as well as the D2D links (ℎ ). Therefore, the BS is assumed to have perfect knowledge of the long-term evolution of the different channel gains, through signaling exchange between the different entities. The BS then performs resource allocation based on these estimated channel gains to optimally pair the D2Ds to CUs and to instruct D2D-CU pairs of the required Maximum achievable D2D rate for every D2D-CU couple P L 1 , P L 2 P L 3 , P L 4 Planes associated to 1 , 2 , 3 , and 4 resp. transmit powers on their collocated channels, according to the selected transmission scenario.

A. Formulation of the Joint Channel and Power Allocation Problem
Let be the channel allocation matrix, with the element ( , ) at the th row and th column equaling one if D2D pair is collocated with CU and zero otherwise. Also, let R 2 ( , ) be the maximum achievable D2D rate of the pair when collocated with . Channel allocation is performed such that a D2D pair is multiplexed over a single UL channel, on the one hand, and such that a maximum of one D2D pair is multiplexed over a UL channel, on the other. The joint channel and power allocation problem for the maximization of the total D2D throughput can be cast as: with R 2 ( , ) the solution to: with , the minimum target rate of , its achieved rate, and 2 ( , ) the D2D rate, sum of achieved rates by ,1 ( ,1 ) and ,2 ( ,2 ). The in the subscripts of , , ,1, , ,2, refers to the maximum transmit powers of , ,1 and ,2 respectively. From the structure of Problem (1), and since CU users are allocated orthogonal channels, the performance of a given D2D-CU pair is independent from the network activity over the remaining channels in the system. Therefore, one can optimize the throughput of all possible D2D-CU pairs, constructing a × table of achievable rates, and then proceed to the optimal allocation of channels to D2Ds in a second phase, thus pairing D2Ds to CUs based on their achievable rate to maximize the D2D sum-throughput in the system. The aim of the following sections is to obtain the optimal PAs of the four transmission methods FD-NoSIC, HD-NoSIC, HD-SIC and FD-SIC in order to build their corresponding tables of achievable rates R − 2 , R − 2 , R − 2 , and R − 2 respectively. Based on these tables, the optimal channel allocation is conducted in Section IX. The overall procedure spanning from the PA problems formulations in section III, to the optimal channel assignment in section IX, is depicted in Fig. 2.

III. POWER ALLOCATION FOR NO-SIC SCENARIOS
From hereinafter, since the optimal D2D rate of all ( , ) couples is to be computed and because the resolution of the PAs is independent of the elected D2D-CU couple, we drop the indices relative to a specific D2D pair and CU user. Hence, user designates the CU user at hand, and 1 and 2 are the corresponding D2D pair. The involved channels gains are therefore denoted as ℎ , ℎ , 1 , ℎ , 2 , ℎ 1 , , ℎ 2 , and ℎ , , and the transmit powers of 1 , 2 and are 1 , 2 , , with their power limits 1, , 2, , , .

A. FD-NoSIC
In FD, 1 and 2 transmit simultaneously, thus they both suffer from RSI. Since, in this method, SIC is not attempted at the levels of 1 , 2 and the BS, the SINRs at the level of the BS and the D2D users are given by: with 2 the additive Gaussian noise power. The achieved rates are expressed according to the Shannon capacity theorem : with the bandwidth of each UL channel resource. Due to the interference terms in (3), Problem (2) is non-convex. To solve it, a geometrical representation can be used, leading to the analytical global solution in [30]. This method is adopted in our work to derive the results of the FD-NoSIC scenario in the performance assessment section.

B. HD-NoSIC
The time slot is now divided into two equal half-time slots where 1 and 2 alternately transmit and receive information. To maximize the total D2D rate, the optimization is conducted separately in the two half-time slots. In the first half, 1 transmits information ( 2 = 0). In Problem (2), the objective function and CU rate are now: Also, Problem (2) is constrained only by eqs. (2a) to (2c). Note that ,1 is the transmit power of during the first halftime slot.
2 ,1 is strictly increasing with 1 and decreasing with ,1 ; therefore, to maximize 2 ,1 = 2 , 1 should be increased and ,1 decreased as long as ,1 satisfies the minimum rate condition of the CU. Consequently, 1 should be increased as much as possible and then ,1 is obtained as a function of 1 ( ,1 = ( 1 )) by enforcing an equality between ,1 and , . If for 1 = 1, , ( 1, ) ≤ , , the couple ( 1, , ( 1, )) is retained as the ( 1 , ,1 ) solution; otherwise, the couple ( −1 ( , ), , ) delivers the best solution. The same reasoning is applied for the second half-time slot (where 1 = 0) to maximize 2 ,2 = 1 . The total user and D2D rates are given by: IV. PA PROBLEM MODIFICATION FOR HD AND FD WITH MUTUAL SIC (HD-SIC AND FD-SIC) Using a SIC receiver at the level of the BS and the D2D users, interfering messages can be decoded then subtracted from the received message, canceling thereby the interference in both FD and HD scenarios. Let 1 and 2 be the messages transmitted by the devices 1 and 2 . In the case of FD, the BS can decode and subtract successively 1 then 2 , or inversely, before proceeding to the decoding of (the message transmitted by the CU); hence, two decoding orders are possible. Users 1 and 2 can also remove the interference of , leading to the following SINR expressions: The SINRs are replaced in eqs. (4) and (5) to obtain and 2 = 1 + 2 that will be used in Problem (2). For the case of HD, the SINRs in the first half-time slot are: In the second half-time slot, is the same and 1 = 2 ℎ / 2 . Problem (2) is now reformulated in each time slot by expressing the rates using the present SINRs. However, additional constraints relative to the SIC feasibility must be added to the problem. They are derived in the next sections for HD-SIC, then for FD-SIC. Afterwards, Problem (2) is solved for the two SIC scenarios.

V. PA FOR HD-SIC SCENARIO
Consider the first half-time slot, where and 1 are transmitting and and 2 are receiving. Hereafter, we develop the mutual SIC constraints between and 2 (as a receiver). Let be the SINR of the message at the level of user ( is either 1 , 2 or , and is either 1, 2 or ). For to successfully decode the message 1 transmitted by 1 to 2 , the received rate of 1 at the level of must be greater than the rate of 1 at the level of 2 . Thus, we must have: 1 >  messages to two separate receivers, which was studied in [26]. It was shown that the SINR conditions lead to: In addition to condition (6), the PMCs must be verified, in order to ensure that the message to be decoded first at the level of a receiver has a higher power level than that of the remaining message [28]. The PMCs for the decoding of and 1 at the level of 2 and are given by: Note that (6) is satisfied if (7) is satisfied, since (6) is equivalent to < . Therefore, the PMCs encompass the rate conditions while being more restrictive. Problem (2) now only includes the additional constraint (7) for the first time slot. The HD-SIC rate expressions are as follows: ).
Maximizing 2 lies in the increase of 1 . Also, guaranteeing the CU rate to , = (2 , − 1) 2 /ℎ , . However, due to the PMCs, the increase in 1 is very likely to increase ,1 according to the range of allowed values in (7), leading to an excess of CU rate. Since maximization of network throughput (i.e. sum of D2D and CU rates) is not the objective of this study, we select from the range of admissible ,1 values, the one closest to , . With that criterion in mind, the power allocation problem for D2D rate maximization is solved by increasing 1 as much as possible (possibly until 1, ) and adjusting ,1 accordingly. The proposed PA procedure, illustrated in problem constraints. Note that if 1, is too low (< , / ), the problem is not feasible even when (6) is verified. For the second time slot, the same methodology is followed, where the PMCs and the new necessary and sufficient channel conditions are given by: As a conclusion, in the HD-SIC scenario, the system checks for the validity of the channel condition corresponding to the half-time slot before going through the procedure described above. If the channel condition is not favorable or if no solution exists (i.e. , > , or 1, < , / for the first half, and 2, < , / for the second half), the system reverts to the HD-NoSIC solution of section III-B. This leads to four combinations of SIC/NoSIC procedures, two for every half-time slot, and they are all included in the HD-SIC algorithm.

VI. DERIVATION OF THE SIC CONDITIONS FOR FD
MUTUAL SIC In this scenario, we are looking for the conditions that allow 1 to decode , 2 to decode , and to decode 1 and 2 . As already mentioned, two decoding orders are possible at the level of .
A. First decoding order: decodes 2 then 1 We first start by studying the mutual SIC constraints between and 1 (as a receiver). For to successfully decode the message 2 transmitted by 2 to 1 , we must have: Since practical systems are interference-limited [31], [32], the noise power is negligible compared to the interfering terms which yields the SIC condition: In addition to condition (10), the PMCs must be verified. Since decodes 2 first, then we have the following PMC for the decoding of 2 : For 1 to be able to remove the interference of prior to retrieving 2 , we must have 1 > , which leads to: (12) and the corresponding PMC is: Regarding the mutual SIC between the receivers and 2 , the decoding of 1 at the level of requires 1 to be greater than 1 2 : Note that 1 does not include 2 since 2 is decoded and canceled prior to 1 . The corresponding PMC is given by: At the level of 2 , 2 must be greater than to decode and subtract before retrieving 1 . This yields the following condition: Finally, the PMC at the level of 2 is given by: B. Second decoding order: decodes 1 then 2 Following the same reasoning as in section VI-A, for the case where 1 is decoded before 2 at the level of , the PMC and rate constraints for a full SIC between 1 and , and 2 and , are obtained and listed below: In addition to constraints eqs. (2a) to (2d), Problem (2) now includes eight new constraints that express the full SIC feasibility (either equations (10) to (17) or (18) to (25), depending on the decoding order). Solving this optimization problem with inequality constraints by means of the standard Karush-Kuhn-Tucker conditions implies exploring all the possible combinations of active/inactive constraints (an inequality constraint is active if it is verified with equality). This results in a total of 2 12 − 1 combinations to be considered. To reduce this exorbitant complexity, the interplay between SIC rate conditions and PMCs is analyzed in the next section, targeting the removal of redundant constraints.

VII. PA PROBLEM SIMPLIFICATION OF FD-SIC BY CONSTRAINT REDUCTION
Consider the first decoding order at the level of where 2 is decoded before 1 . The PMCs for the decoding of 1 at the level of and of at the level of 2 are given by (15) and (17). By multiplying (15) by ℎ 2 , and adding it to (17) multiplied by ℎ , , one can eliminate to obtain: which is the SIC condition (16) introduced to remove at the level of 2 . Also, eliminating 1 from the two PMCs by means of adding (15) multiplied by ℎ to (17) multiplied by ℎ , 1 yields (14). Consequently, the PMCs for the decoding of 1 at the level of , and at the level of 2 imply their counterpart rate conditions. Moreover, it is noted from (16) that the same necessary condition (6) that is found in HD-SIC between and 2 as receivers, is obtained for the application of FD-SIC between 2 and : Note that if (6) is not true, (16) becomes impossible to satisfy no matter 1 and 2 ; however, when (6) is true, (16) can be satisfied under an adequate power play between 1 and 2 .
We now move to the PMC and SIC conditions for the decoding of 2 and at the level of and 1 respectively, i.e. (11), (13), (10) and (12). By adding (11) multiplied by ℎ to (13) multiplied by ℎ , 2 , 2 is eliminated to yield: which can be further transformed into: Thus, the PMCs (11) and (13) imply (10). In fact, not only do they imply the rate condition, but it is clear that the PMCs represent more restrictive constraints than rate conditions. Finally, eliminating from the PMCs through the combination of (11) multiplied by ℎ 1 , with (13) multiplied by ℎ , yields: (27) which can be rearranged into: Once again, the PMCs for the decoding of 2 and at and 1 imply their rate condition counterparts. Note that the necessary channel condition that appears from eqs. (26) and (27) is the same as in the case of HD-SIC in the second half-time slot: Also, the combinations of (15) with (13), and (17) with (11), while eliminating , give the following conditions: These inequalities yield two other necessary, but not sufficient, channel conditions for the application of full SIC to the system: Repeating the same procedure for the second decoding order delivers the same results: 1) the PMCs encompass the rate conditions, 2) the same necessary four channel conditions (6), (8), (28), and (29) are obtained.
Therefore, in the FD-SIC scenario, the system checks the validity of eqs. (6), (8), (28) and (29) prior to solving the PA problem for each decoding order. If the channel conditions are not valid or no solution is obtained for (2), the FD-SIC algorithm reverts to the FD-NoSIC procedure described in section III-A. As a conclusion for this section, Problem (2) is now only equipped with the PMC set corresponding to the decoding order (i.e. eqs. (11), (13), (15) and (17), or eqs. (22) to (25)), in addition to constraints eqs. (2a) to (2d). This reduces the number of combinations of active/inactive constraints from 2 12 − 1 to 2 8 − 1 which is still considerable. The aim of the next section is to workaround the need of a full search over the corresponding 255 cases for determining the optimal PA. This is done by efficiently determining the meaningful constraint combinations, based on the geometrical interpretation of the FD-SIC PA problem. Considerable complexity reductions arise from this approach as shown next.

VIII. SOLUTION FOR FD-SIC OPTIMAL PA
The proposed geometrical resolution of the FD-SIC D2D rate maximization PA problem is presented in detail for the first decoding order. First, the geometrical representation of the solution space satisfying the PMC and power limit constraints is provided. Then, a procedure is elaborated leading to the reduction of the search space to the minimum required. Afterwards, the optimization is conducted on the resulting reduced search space. At last, a quick summary of the optimal PA procedure is presented including the required changes to obtain the optimal PA for the second decoding order.

A. 3D Solution Space Representation
The four PMCs that must be satisfied in eq. (11), (13), (15) and (17) are written in the following form: In the 3D space of axes , , representing variables 1 , 2 and respectively, we introduce the planes PL 1 , PL 2 , PL 3 and PL 4 whose equation is given by the PMCs 1, 2, 3 and 4 when the conditions are met with equality. In the following, we refer to PL as the plane derived from, or equivalently, corresponding to, or simply, as the plane of . Each PMC restricts the search space either to the half space below its corresponding plane like for 1 and 3 , or to the half space above its corresponding plane as for 2 and 4 . On the other hand, the transmit power limits restrict the search space to the region within the parallelepiped defined by the sides = 1, , = 0, = 2, , = 0, = , , = , . To have a non-empty search space (i.e. FD-SIC is feasible), the pentahedron defined by the space region above PL 2 and PL 4 and below PL 1 and PL 3 must be non-empty, and it must have a common region with the parallelepiped.
• Non-empty pentahedron: The pentahedron is non-empty if the intersection lines of PL 1 with PL 2 ( 1,2 ) and PL 4 ( 1,4 ) are below the intersection line of PL 1 with PL 3 ( 3 ), as shown in Fig. 4. Let ì be the direction vector of 1,2 ; we must have ( ì) ( ì) < ℎ , 1 /ℎ , , which yields: leading to the following channel condition: Doing the same for 1,4 with respect to 3 , we get the following channel condition: Note that (30) and (31) are more restrictive than the necessary conditions of eq. (28) and (29), which is normal since they turn them into sufficient channel conditions. • Pentahedron ∩ parallelepiped: For the pentahedron to have a non-empty intersection with the parallelepiped, it is sufficient to make sure that the intersection line of PL 1 with PL 3 ( 3 ) intersects the plane = , within the 1, and 2, limits. These conditions on the , coordinates of PL 3 ∩ PL 1 ∩ , yield the constraints: Conditions (30), (31) and (32) form the necessary and sufficient constraints for the existence of a solution to the FD-SIC PA problem according to the first decoding order.

B. Search Space Reduction
We prove in this section that the optimal solution lies on the intersection line of PL 2 , PL 4 or the lower side of the parallelepiped , with one of the outer sides of the parallelepiped 1 , 2 , (cf. Fig. 4), respectively defined by: Proof. The D2D rate is given by: For any couple ( 1 , 2 ), and ∀ > 1, the throughput of ( 1 , 2 ) is greater than 2 ( 1 , 2 ) since: Therefore, given an initial triplet ( 1 , 2 , ), a higher throughput-achieving triplet can be obtained by simply multiplying the components by a factor larger than 1. The higher the , the higher the throughput, meaning that should be increased until reaching the boundaries of the region, which can be either 1, , 2, or , .
Moreover, it is clear that the D2D rate is independent of . This means that when moving on a vertical line in the solution space, 2 is constant and only affects the CU rate. To keep the CU rate as close as possible to , , we select the smallest value from the range of admissible values for a given ( 1 , 2 ) couple. Since every point in the solution space must be on top of PL 2 and PL 4 , the minimum allowed value of is given by forcing the equality either on 2 or on 4 , according to the one that delivers the higher minimum value of for the considered ( 1 , 2 ) couple. As a conclusion, the optimal solution lies on the intersection segment of one of the outer sides of the parallelepiped 1 , 2 , or , with one of the planes PL 2 , PL 4 , or , resulting in a total of eight possible combinations (eight segments). Given the shape of the solution space, some of the combinations are mutually exclusive. The aim of the next section is to determine which subset of segments should be accounted for in the power optimization process, depending on the channel conditions of the D2D-CU couple.

C. Selection of the Useful Intersections
As can be seen from Fig. 4, some of the eight intersections can be discarded. For example, the intersection of with 1 and 2 is not relevant, since the value of is decided by 2 and 4 , whose planes are on top of near the sides 1 and 2 . Fig. 5 shows the projection on the ( 1 , 2 ) plane of the partition of the space into two vertical regions where 4 encompasses 2 for region 1, and 2 encompasses 4 for region 2. The plane separating the two regions is the vertical plane passing through the straight line PL 4 ∩ PL 2 . Therefore, for the case of Figs. 4 and 5, the D2D rate optimization is to be conducted over the segment 1 ∪ 4 which is included in 1 , over the segment 4 5 included in , and on the segment 5 2 over 2 . By doing so, the optimization over the segments 1 , 3 , 3 2 and 2 2 is avoided. Therefore, the first step in reducing the number of intersections to be considered lies in determining which of 4 and 2 encompasses the other, and for which region of the space. To that end, a schematic of PL 2 and PL 4 is presented in Figs. 6a and 6b, showing their intersection with the planes  1) Interplay between 2 and 4 : Depending on the angles Ω, , and , four cases are identified to determine the interplay between 2 and 4 : 1) Ω > , > : 2 encompasses 4 ( 2 ⇒ 4 ) over all the positive ( 1 , 2 ) plane. 2) Ω < , < : 4 encompasses 2 ( 4 ⇒ 2 ) over all the positive ( 1 , 2 ) plane. 3) Ω < , > : 4 encompasses 2 in region 1 and 2 encompasses 4 in region 2, (cf. Fig. 5). 4) Ω > , < : 4 encompasses 2 in region 2 and 2 encompasses 4 in region 1. Before proceeding, it must be noted that even for cases 3) and 4), a PMC may encompass the other on the entire search space if the actual search space is included either in region 1 or 2. This is depicted in the examples of Fig. 7 which take back the conditions of Fig. 5 with some modifications. In Fig. 7a, PL 1 is such that 1,4 is at the right side of ( 1,4 is in region 2), then the search space is included in region 2 and only 2 needs to be accounted for. The other scenario is represented in Fig. 7b where PL 3 ∩ PL 2 is at the left side of (in region 1), hence 4 encompasses 2 over the entirety of the search space. The first scenario occurs when is on top of PL 1 , and the second one occurs when is on top of PL 3 . The explicit channel conditions enabling each scenario are derived in detail in Appendix A, where it is shown that two simple tests are required to determine if the search space is included in any region.
As a conclusion, by comparing Ω to and to , and following the discussion in Appendix A, the number of intersections to be considered is reduced by selecting the appropriate PMC between 2 and 4 in the corresponding space region. For the sake of clarity, we introduce 2,4 as the efficient combination of 2 and 4 , given by: , elsewhere.
2) Selection of the useful parallelepiped sides: With 2,4 at hand, the next step is to reduce the unnecessary sides of the parallelepiped, which do not intersect with PL 2,4 , or that do intersect with PL 2,4 but not inside the range allowed between PL 1 and PL 3 . To that end, we study 1 and 3 which do not affect the intersection segments (of PL 2,4 with the parallelepiped sides) as such, but rather the end points of these intersection segments. A typical example is given in Fig. 7a where 1 sets the end point 1 from the side 1 , and 3 sets the end point 2 from the side 2 .
Let 1 regroup the intersection lines 1,2 and 1,4 such that 1 = PL 2,4 ∩ PL 1 , and let 3 be the intersection line of PL 3 with PL 2,4 (cf. Fig. 8). Each of 1 and 3 may intercept the sides 1 or 2 or , yielding a total of nine potential combinations. Since each side is a rectangular surface within the infinite plane S of equation = , , then 1 and 3 can intercept only one side of the parallelepiped ( aside) for a given channel configuration. Let and be the intersection points of 1 and 3 with , we have: To determine which sides are hit by 1 (resp. 3 ), i.e. to determine if we have 1 , 2 or (resp. 1 , 2 or ), we consider the points (resp. ), intersections of 1 (resp. ( 1 ∩ S 1 ) ( 1 ∩ S 1 ) Then, two tests are needed to determine which of 1 , 2 or occurs for the given channel states (cf. Algorithm 1). Note that and have positive coordinates as shown in Appendix B-C. The same tests are replicated for . From the nine possibilities, only six combinations are actually viable Algorithm 1: 1 intersection with the parallalepiped input : 1, , 2, , , , , , ℎ , , ℎ , 1 , 2 , ℎ 1, , ℎ 2 , , ℎ 1 , , ℎ 2 , Result: = , keep end end because the pairs ( , 1 ), ( 2 , ) and ( 2 , 1 ) cannot be achieved without violating (30) or (31) as can be seen in Fig. 8. Indeed, the three cases shown in Fig. 8 lead to empty search spaces. The six viable pairs are given in Table II with the correspondence between the pairs and the parallelepiped sides hosting the useful intersection segments. Note that if PL 1 and PL 3 intercept PL 2,4 at the same side, then the search space can be reduced to a single segment as it is the case for the first, the fourth and the fifth rows in Table II. For the second and third rows, two segments are involved in the D2D rate optimization. Finally, in the case where 1 intercepts 1 and 3 intercepts 2 (as in Fig. 5), the segment 4 5 belonging to is to be included in the D2D optimization process -in  addition to the segments in 1 and 2 -if and only if the value of obtained from PL 2,4 at 1 = 1, and 2 = 2, is greater than , .
3) Segments endpoints: Having determined the relevant intersection segments (a maximum of three segments) for the D2D rate optimization using PMCs 1 and 3, we detail hereafter how the endpoints of every segment are determined for each side of the parallelepiped. For the sake of clarity, let 1 , 2 , 3 , 4 , 5 be the edges of the parallelepiped (cf. Fig. 4) given by: Also, let the three families of points , , and be the intersections of PL 2 , PL 4 and PL 2,4 with : Examples of such points can be seen in Fig. 7b for 4 , 5 , 3 and 2 . Note that points are only used to designate the points or depending on whether we are in region 1 or 2. We can now efficiently designate the segment endpoints on each side. a) Side 2 : The optimization over 2 translates into an optimization over 1 , since 2 is equal to 2, . It is clear that the minimal value of 1 is bound to 3 . In Fig. 5 for example, the minimal value of 1 is obtained for the point 2 , intersection of PL 3 with PL 2,4 . However, if 4 were the one encompassing 2 in region 2, the considered segment over 2 would have been 3 2 ∪ 2 2 ; then, the minimum 1 value would have been given by the point 2 = PL 3 ∩ 2 ∩ . Thus, we can generalize by stating that the minimal value of 1 in the segment over 2 is given by: Regarding the maximum value of 1 , it can be due to the intersection of 2 ∩ PL 2,4 with either (like for 5 in Fig.  5), 1 (like for 3 ), or with PL 1 (for the case of 2 and 2 , in the fifth row of table II). Also, the maximum 1 value may be simply set by 2 , the intersection of 2 ∩ PL 1 with . The maximum value of 1 is given by: Note that for the side 2 , 3 is involved in the minimum 1 value, and 1 in the maximum value.
b) Side 1 : Regarding the side 1 , 3 is now involved in the maximum 2 value, while 1 settles the minimum 2 value; its expression is given by: The maximum value of 2 depends on which plane intercepts first PL 2,4 among the three candidates: PL 3 , (cf. Fig. 5), or 2 . We have: In the example of Fig. 5, the intersection segment starts at 1 and ends at 4 passing by . Although 1 ∪ 4 is a different segment from 1 4 , their projections over the ( 1 , 2 ) plane are identical, thus we are only interested in segment ends over both sides 1 and 2 (the projection of 3 2 ∪ 2 2 is the same as that of the segment 3 2 ).

c) Side
: Unlike for the other sides, none of 1 or 2 is fixed, but 2 can be expressed in terms of 1 ; therefore, we evaluate the position of the endpoints of the segments on in terms of maximum 1 and minimum 1 . When = PL 2 ∩ PL 4 does not intercept as it is the case for Fig. 5 for example, the intersection of PL 2,4 with yields a unique segment ( 4 5 in case of Fig. 5). The endpoint corresponding to the minimum value of 1 is due to the intersection of PL 2,4 with either 2 or PL 3 .
The maximum value of 1 is due to the intersection of PL 2,4 with either PL 1 or 1 .
If resides on , then the intersection segment of PL 2,4 with is broken into two segments. In that case, if we let and be the points given by (33) and (34) respectively, then the optimization over has to be conducted separately over from the side of region 1, and over from the side of region 2. In this case, corresponds to the max 1 point in and to the min 1 point in . Assuming the conditions of the last row in Table II, this is the only case where 4 segments in total have to be checked to find the optimal D2D throughput achieving point. The coordinates of all the points mentioned in this section are provided in Appendix B-C.

D. D2D Throughput Optimization
At last, given the segments locations and endpoints, the analytical power optimization can be conducted. The mathematical formulation varies according to the side the segment is included in.
1) Side 1 : The optimization variable is 1 and the problem formulation is the following: * Taking the derivative of ( 1 ) = 2 ( 1 , 2, ) with respect to 1 , we get: The sign of / 1 is equal to the sign of the following second-degree polynomial of 1 : It is shown in Appendix B-A that independently of the sign of the polynomial's discriminant, it is sufficient to test which of min 1 or max 1 delivers the best throughput and then select the corresponding segment endpoint. The coordinates of the endpoint form the optimal triplet ( * 1 , * 2 , * ) maximizing the D2D throughput over the side 1 . The endpoint coordinates are given in Appendix B-C.
2) Side 2 : Following the same reasoning as for 1 (with the only difference that the optimization variable is now 2 instead of 1 , and 1 = 1, ), the same conclusion is reached, i.e. the maximum D2D throughput is delivered by the points corresponding either to min 2 or to max 2 .
3) Side : Since is accounted for in the maximum and minimum values of 1 for each intersection segment, the optimization can thus be conducted over each segment independently.
The D2D throughput maximization problem over the intersection segment of PL 2 with can be written as follows: * , we get: .
Since 1 ∈ U, it can be easily verified that both denominators are positive; therefore, only the numerator is needed to evaluate the sign of / 1 : After some simplifications and re-arrangements, the sign of / 1 can be written as the sign of a second-degree polynomial of 1 of the form 2 1 + 1 + with: Given the root 1 = (− − √ 2 − 4 )/2 of the polynomial, we show in Appendix B-B that * 1 is either given by min 1 , max 1 , or 1 (when it is included in the interval [min 1 , max 1 ]), according to the value delivering the highest throughput. Regarding the optimization over the intersection segment of PL 4 with , the same steps are followed to determine the optimal value of 1 : we start by writing the expression of ( 1 ) by replacing 2 in ( 1 , 2 ) with ( , ℎ 2 , − 1 ℎ )/ 2 ). Then, the study of the sign of / 1 turns into the study of the sign of another seconddegree polynomial 2 1 + 1 + with: Also, following the different channel conditions concerning sgn( 1 2 − ℎ 2 ), and considering all the possible relative positions between max 1 , min 1 , and 1 , the same result as previously is obtained, which can be cast as: * As a conclusion, the optimization over the sides 1 and 2 resides in selecting the corresponding endpoint achieving the highest throughput. On the side , a maximum of three additional points ( , 1 , 1 ) may need to be considered to get the highest D2D throughput.

E. Summary of the PA Procedure and Extension to the Second Decoding Order
In this section, the geometrical representation of the FD-SIC PA problem allowing for a drastic reduction of the search space size was described. It was shown that the initial search volume in Sec. VIII-A can be reduced to a set of intersection segments (Sec. VIII-B) from which a subset is selected (VIII-C). These segments search spaces are then further reduced to become a finite set of points (Secs. VIII-C3, VIII-D). In the worst case scenario, the original PA problem, which had 2 12 − 1 variants, is converted into the search for the maximum throughput of a list of seven elements: two elements from 1 , two from 2 and three additional elements from ( 4 is a common endpoint to 1 and , and 5 is common to 2 and ). The global PA procedure to determine the optimal D2D rate for the first decoding order of FD-SIC is summarized in algorithm 2.
Algorithm 2: Optimal PA procedure for FD-SIC input : 1, , 2, , , , , , ℎ , , ℎ , 1 , 2 , Test Ω, , , and build 2,4 ; Execute Algorithm 1 to determine and ; Follow Table II  Regarding the resolution for the second decoding order, the PA procedure itself is unchanged, but the changes in 1 and 3 lead to some modifications. Here is the list: • Modification in the expressions of 1 and 3 : • The necessary and sufficient conditions (30), (31) and (32) become: • Concerning section VIII-C3, the roles of 1 and 3 are interchanged concerning the settlement of the segment endpoints. • The three non-occuring ( , ) pairs of section VIII-C2 become: ( 1 , 2 ), ( 1 , ), ( , 2 ). Sections VIII-C1 and VIII-D are kept unchanged because building 2,4 is independent of 1 and 3 , and given the endpoints of the segments subset, the optimization of section VIII-D is not affected by the change in 1 and 3 .

IX. CHANNEL ALLOCATION
In this section, the procedure for optimal channel allocation to D2D devices is conducted. Recalling that the D2D system is underlaying a pre-established CU network, D2D channel allocation is equivalently referred to as D2D-CU pairing.
Having determined the analytical PA solutions for all the transmission scenarios, their resolution cost is a constanttime operation. Therefore, filling the D2D rate tables s.t. the constraints of (1) are verified This assignment problem is efficiently solved by the Kuhn-Munkres (KM) algorithm [33], also called the Hungarian method, with a complexity of ( 2 ) [34]. Note that the global resource allocation complexity is now dominated by that of the channel assignment after the important PA complexity reduction. The KM can be directly applied in our study to yield the optimal channel assignment by rewriting the problem as a minimization of the opposite objective function (−R 2 ). As a conclusion, the optimal PA procedures allowed for an efficient filling of the rate tables which are then fed to the KM solver. This delivers the global optimal solution of the joint channel and power allocation problem formulated in section II-A.
X. NUMERICAL RESULTS In our simulation setup, the BS is positioned at the center of a hexagonal cell with an outermost radius of 300 m. The D2D users and the CU users are randomly located within the cell. The distance between the D2D users of every pair is below a maximum value . The propagation model includes large-scale fading with a path-loss exponent = 3.76, and an 8 dB zero mean lognormal shadowing. The maximum transmit power of the devices and CU is 24 dBm. The system bandwidth is 20 MHz, divided into = 64 channels, leading to a UL bandwidth of = 312.5 kHz, with a noise power of −119 dBm. The minimum required rate , is the same for all the CU users, and the SI cancellation factor is the same for all D2D pairs, its value being varied between −130 and −80 dB. The results are averaged over 1000 different realizations of the devices and CU positions. Unless specified otherwise, , is set to 1.5 Mbps, is set to 20 CUs, = 100 m, and = 5 D2D pairs. Figure 9 presents . On the other hand, the decrease of reduces the interference terms in the D2D throughput expression, which translates into a higher achieved throughput.
As expected, when comparing the performance for different required CU rates between Figs. 9a and 9b, the increase of from 1.5 Mbps to 3 Mbps decreases the achieved D2D throughput for all proposed methods. However, the percentage gain in the performance of SIC procedures with respect to NoSIC increases from 41 % to 86 % for the HD case, and from 33 % to 70 % for the FD case (for = −130 dB). The reason behind this gain increase is that NoSIC algorithms are highly affected by the value of (≥ , ) since they suffer from its interference, which is not the case of SIC techniques. In fact, even though the total number of FD-SIC enabled D2D-CU pairs decreases with also be observed from another perspective in Fig. 10, where the total D2D throughput is presented as a function of the CU required rate. In the conditions of Fig. 10, the gap between FD-NoSIC and HD-NoSIC is large enough so that no intersection occurs between FD-SIC and HD-SIC. However, it can still be observed that the gap between FD-SIC and HD-SIC reduces as the CU required rate increases.
In Fig. 11, the variation of the total D2D throughput is presented as a function of the D2D maximum user distance . The increase of leads to a significant decrease in the performance of all proposed methods since ℎ , the channel gain of the direct link between 1 and 2 , is reduced on average. However, this increase of is accompanied by a greater percentage increase in performance due to mutual SIC for FD and HD transmission scenarios, with respect to No-SIC scenarios. Indeed, FD-SIC achieves a D2D throughput 128 % higher than FD-NoSIC for = 100 m, compared to the 81 % increase achieved for = 20 m. This is due to having more FD-SIC enabled D2D-CU pairs when distancing the D2D users further apart from one another, since an average of 1.96 pairs successfully apply FD-SIC for = 20 m as opposed to 3.33 pairs for = 100 m. The reason behind this increase is the decrease in ℎ which relaxes the sufficient conditions (30) and (31), thereby enabling more FD-SIC cases.   12 presents the evolution of the D2D throughput as a function of the number of CU users in the cell. Although the channel properties of the D2D users (i.e. ℎ , ℎ 1 , and ℎ 2 , ) are unchanged, the total D2D throughput of all techniques benefits from the additional diversity provided by the greater number of CU users. This also favors the FD-SIC enabled pairs, as their average number grows from 1.98 for = 20 to 2.46 for = 50. We can therefore conclude that the important performance gain achieved by SIC methods, with respect to No-SIC methods, can be obtained without requiring the implementation of SIC at all D2D and CU receivers. Indeed, generally only 2 or 3 triplets need to perform SIC which is enough to boost the D2D system capacity, while the others can settle for the simple classical No-SIC receivers. Therefore, the additional complexity is localised at the level of the users performing SIC for which the major throughput increase is worth the incurred SIC complexity. Finally, the  total and average throughput variations are presented in Fig.  13 as a function of the number of D2D pairs in the system, for a fixed value of = 50. In Fig. 13a, the average throughput per D2D pair is shown to slightly decrease with the increasing number of D2D pairs. In a sense, this is the dual of the behavior observed in Fig. 12, since the ratio / decreases with and thus the system diversity -in terms of the average number of possible CU channel choices for every D2D pair to be collocated on -decreases, thus reducing the achievable throughput per D2D pair. Nonetheless, the total throughput follows a quasi linear progression with the number of D2D pairs because the additional D2D pairs are allocated on orthogonal channels, therefore each D2D pair can be associated more or less to an additional D2D rate unit. Figs. 12 and 13 indicate that, for a fixed number of 2 users or CU users, the effect of the proportion / on the average D2D throughput per D2D pair is rather limited. The most dominant factors remain the distance between D2D users, the SI cancellation capabilities of the receivers (for FD-SIC), and the required CU rate.

XI. CONCLUSION
In this paper, the use of NOMA with mutual SIC was proposed for the first time between cellular users and FD-D2D devices underlaying the cellular channels. The necessary and sufficient conditions for applying FD-SIC were derived and a highly efficient PA procedure was elaborated to solve, in constant time operation, the throughput maximization problem of significant original complexity. The optimal, yet simple, PA resolution allowed for achieving global optimal resource allocation by conveniently combining the Kuhn-Munkres channel assignment with the proposed PA methods. The results show important performance gains obtained by applying SIC in D2D underlay systems in both HD and FD transmission schemes, promoting thereby the use of mutual SIC NOMA for D2D systems whenever possible. When applying mutual SIC, the comparison between HD and FD transmission scenarios showed that FD-SIC is more efficient for average to high SI cancellation capabilities, moderate CU rate requirements and significant D2D distances, while HD-SIC performs better especially at low SI cancellation capabilities. Future work directions of the study could be to adapt the PA procedure for D2D underlay to downlink cellular systems which present different sets of challenges and networking paradigms. Moreover, the integration of our work in UAV-assisted NOMA networks for IoT could also be considered, as in the context of [35].
ACKNOWLEDGMENT This work has been funded with support from IMT Atlantique and the Lebanese University Research Support Program.
• If ℎ 2 > 1 2 , is increasing with 1 ⇒ Set * 1 to max 1 • If ℎ 2 < 1 2 , is decreasing with 1 ⇒ Set * 1 to min 1 However, if Δ > 0, then we have the two solutions • If ℎ 2 > 1 2 ⇒ 1 < 2 and 1 < 0. But not much can be said about the sign of 2 and how it compares to min 1 and max 1 . However, we note that the right side of the variation table (where 1 > 1 ) is similar to the variation table in Fig.  14. Therefore, we conclude that: * 1 = arg max[ (min 1 ), (max 1 )].
• if ℎ 2 < 1 2 ⇒ 2 < 1 and Since 1 is a local maximum, ( 1 ) > ( 1 ), ∀ 1 > 2 . Then, the only values of 1 which might give a better throughput than 1 are those at the left of 2 . We can distinguish the following three cases: if max 1 < 1 , set * 1 to max 1 . if 1 < min 1 , set * 1 to min 1 . if 1 ∈ [min 1 , max 1 ], then: * if min 1 > 2 , set * 1 to 1 . * if min 1 < 2 , set * 1 = arg max[ (min 1 ), ( 1 )]. To sum up, in the optimization over the intersection segment of PL 2 with , all the possible channel conditions lead at some point to choosing * 1 from the values min 1 , max 1 , and 1 (when it is included in the interval U) according to the one delivering the highest throughput.