Ngoc Lam

Dinh

Ibtissam Labriji
email: ibtissam.labriji@renault.com

Mickael Maman
email: mickael.maman@cea.fr

Emilio Calvanese Strinati
email: emilio.calvanese-strinati@cea.fr

Toward URLLC with Proactive HARQ Adaptation

In this work, we propose a dynamic decision maker algorithm to improve the proactive HARQ protocol for beyond 5G networks. Based on Lyapunov stochastic optimization, our adaptation control framework dynamically selects the number of proactive retransmissions for intermittent URLLC traffic scenarios under time-varying channel conditions without requiring any prior knowledge associated with this stochastic process. It then better exploits the trade-off between Radio Access Network (RAN) latency, reliability and resource efficiency, which is still limited in its realization on current HARQ designs. We then evaluate the performance of several HARQ strategies and show that our proposal further improves latency over the reactive regime without affecting the resource efficiency such as fixed proactive retransmission while maintaining target reliability.

I. INTRODUCTION

The 5G and beyond network enables the exploitation of new emerging use cases, such as Ultra-Reliable Low Latency Communication (URLLC) which is applicable for mission critical services, remote surgery and wireless factory control. However, ensuring URLLC in wireless communication is a challenging task as it requires the simultaneous enhancement of low latency and high reliability (e.g., industrial control demands a packet error rate of 10 -5 with a latency less than 1 ms). Under such stringent requirement, the retransmission regime in Hybrid Automatic Repeat reQuest (HARQ) is inadequate since each Round Trip Time (RTT) could take several time slots [START_REF] Kim | An implementation of downlink asynchronous HARQ for LTE TDD system[END_REF] for the corrupted packet to be retransmitted again and thus the latency requirement is no longer met if many retransmissions (RTXs) are required.

In order to tackle with the problem of long RTTs for mission critical applications, a K-repetition regime [START_REF] Le | An Overview of Physical Layer Design for Ultra-Reliable Low-Latency Communications in 3GPP Releases 15, 16, and 17[END_REF] and a proactive regime with early termination [START_REF] Liu | Analyzing Grant-Free Access for URLLC Service[END_REF] have been proposed. By performing such fixed repetition of a packet multiple times, regardless of the RTX success or failure, one can opportunistically decode the packet at the receiver in a shorter time at the expense of inefficient resource utilization [START_REF] Jacobsen | System Level Analysis of Uplink Grant-Free Transmission for URLLC[END_REF]. To cope with the over-estimation of packet repetition, Le et al. [START_REF] Le | Optimal reserved resources to ensure the repetitions in Ultra-Reliable Low-Latency Communication Uplink Grant-free transmission[END_REF] proposed a reserved resource scheme to guarantee that each HARQ process can perform an adequate number of repetitions under a configured period. However, in a time-varying radio channel, the lack of an adaptive packet RTX strategy evidently degrades the system performance since the number of RTXs is applied imprecisely and both resource utilization and latency are negatively impacted. In the current literature, the adaptation of the HARQ strategy can be achieved in several ways such as adapting the modulation and coding scheme [START_REF] Pfletschinger | Adaptive HARQ With Non-Binary Repetition Coding[END_REF], the transmission power [START_REF] Jabi | Energy Efficiency of Adaptive HARQ[END_REF] and the maximum number of RTXs [START_REF] Mukhtar | Content-aware and occupancy-based hybrid ARQ for video transmission[END_REF]. By optimally tuning these parameters under dynamic channel conditions, the performance of adaptive HARQ can be substantially improved.

Considering the benefits of proactive HARQ and the adaptation strategy to dynamic channel conditions, in this paper, we proposed an adaptation control algorithm to improve proactive HARQ. In addition, we evaluate our design in an intermittent traffic scenario that is close to the use cases of URLLC in 5G. To the best of our knowledge, there is still limited amount of research that considers such dynamics in the results. In [START_REF] Han | Clarq: A dynamic arq solution for ultra-high closed-loop reliability[END_REF], they proposed a Closed-Loop ARQ protocol that dynamically re-allocates the remaining resource between uplink and downlink slots upon the result of last uplink transmission. In contrast to the current state-of-theart, our proposition will robustly select the number of RTXs sufficient to maintain a good level of resource efficiency and ensure a shorter packet delivery time. To design such an optimal control action, knowledge of stochastic processes such as traffic behaviour and instantaneous channel quality is required. However the acquisition of such information is usually difficult and their behaviours are unpredictable, which is an obstacle to design a feasible optimization algorithm. Fortunately, the framework of so-called Lyapunov stochastic optimization does not require the prior probabilities associated with these processes [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF]. Furthermore, given the compatibility of the time-slot system in the 5G New Radio (NR), the Lyapunov optimization framework becomes a promising candidate to deal with our problem. The optimal performance of the balance between shorter packet delivery and better resource usage is then achieved by fine-tuning the V parameter of the Lyapunov framework.

Based on rich information perceived prior to the control action to be taken, such as the current state of the queue backlog, our objective is to minimize the long-term average of total resource allocation needed for each packet, subject to the stability of the queue and the reduction of the risk of bad decisions. To formalize our problem, we presented the dual objectives of resource allocation and queue stability in drift-plus-penalty problem, and then introduced a virtual queue whose evolution is associated with the risk of taking an inadequate proactive retransmission.

The remainder of the paper is organized as follows: Section II presents a detailed scenario setup where many physical (PHY) and Medium Access Control (MAC) layers assumption are described in our work and is followed by problem formulation. Then, Section III details our designed algorithm. In Section IV, the simulation results demonstrate the performance of the system and Section V concludes the paper.

II. SYSTEM MODEL

A. Scenario Description

The network contains 1 gNb and 1 UE. During the simulation, the UE node moves around the gNB and varies its distance to it. Downlink (DL) traffic is assumed to be generated by a remote host located near the gNb, in which the variable arrival rate and packet length follow an exponential distribution. The application packets are then queued in the transmission buffer Q 1 (t) of the Radio Link Control (RLC) layer. After completing the scheduler operation at MAC layer, the gNb prepares a transport block (TB) whose data is extracted from Q 1 (t) and sends it over the air. At the same time, the scheduler keeps a copy of this TB and associates it with an identifier to construct an identical HARQ process [START_REF] Gpp | Medium Access Control (MAC) protocol specification[END_REF], which will give the receiver information about the processed data. In order to demonstrate the state of ongoing HARQ processes that are not yet decoded at the UE side, we defined a Q 2 (t) that contains them. Due to the dynamic nature of not only the traffic model but also the channel behaviour, the lengths of Q 1 (t) and Q 2 (t) are influenced and can be considered as random variables. The state of Q 1 (t) and Q 2 (t) demonstrated a twostage queuing system whose length should be minimized.

TB = N OF DM × BW 8 × SCS(N um) × M (m) × CR(m) (1)
where SCS(N um) is the sub-carrier spacing at the numerology N um, M (m), CR(m) are the modulation order and code rate of the selected MCS m, respectively. These values are given in Table 5.1.3.1 of 3GPP TS 38.214 [START_REF] Gpp | Physical layer procedures for data[END_REF]. Concerning the improvement of the diversity order of the communication, U tx = M BS × N BS and S rx = M U E × N U E Uniform Linear Array (ULA) antennas are installed at the gNB and UE, respectively. The transmission power at the gNb is set to P tx dBm. In order to evaluate the performance of the system in a highly changing environment, we modeled the Indoor Factory-sparse clutter scenario (InF-SL) according to the 3GPP technical report [START_REF] Gpp | Study on channel model for frequencies from 0.5 to 100 GHz[END_REF]. This model clarifies the channel behaviours based on antenna modeling, slow-fading factors such as shadowing, propagation loss and fast-fading factors due to multi path. Initially, the propagation loss is calculated based on whether the path is in visibility or not (LOS or NLOS) and can be derived:

P L LOS (d) = 31.84+21.50×log 10 (d)+19.00×log 10 (f c ) (2) P L N LOS (d) = 33+25.50×log 10 (d)+20.00×log 10 (f c ) (3)
where d is the instantaneous 3D distance between the UE and gNb in meter and f c is the center frequency in GHz.

According to the channel statistics, the probability of LOS communication P r LOS between the transmitter and receiver which is based on the distance d and the scenario parameter k = -d clutter ln(1-r) can be derived as follows.

P r LOS = e -d k (4)
where d clutter = 10 is the typical clutter size in meter and r = 0.3 is the clutter density. In this study, we consider the UE mobility with respect to the gNb, so the shadowing effect plays an important role in the channel status and it has been modelled as in [START_REF] Itu | The prediction of the time and the spatial profile for broadband and mobile services using UHF and SHF bands[END_REF].

To better understand the benefits of good decision making for proactive HARQ in improving Radio Access Network (RAN) latency and resource efficiency, Figure 2 shows the schemes of (A) classical reactive HARQ procedure, (B) fixed repetitions of 3 RTXs and (C) dynamic redundancy applied to the number of RTXs. In Figure 2A, a delay L 12 is introduced to demonstrate the TB preparation time from the gNb scheduler to the antenna. Then, a feedback will be encoded together with an Uplink Control Information (UCI) message and sent back to the gNB after T f b = K 1 slot(s), thus illustrating the processing time at the UE. In 5G NR standard, this processing time reflects a delay between the reception of the UL grant in the DL (PDCCH) and the transmission of the corresponding UL data (Physical Uplink Shared Channel (PUSCH)) [START_REF] Patriciello | The Impact of NR Scheduling Timings on End-to-End Delay for Uplink Traffic[END_REF]. Afterwards, the gNB has the information about the corrupted HARQ process on the UE side and decides to retransmit the erroneous TB after L 12 slots. This process continues 5 times until the corrupted TB is successfully decoded by the UE. By doing this, the resource are perfectly utilized, but the latency could be unacceptable for URLLC communications. Figure 2B shows the fixed repetition of corrupted TB with a fixed redundancy level R = 3. This means that the first cluster which will occupy 3 consecutive slots, will be reserved for retransmissions. If it fails, the second cluster with the same redundancy level will be allocated. As shown in Figure 2B, latency is significantly improved at the cost of resource efficiency (RTX6 is useless when the TB has been decoded after 5 retransmissions). By dynamically selecting the proactive redundancy level for each cluster, Figure 2C shows better performance in terms of reduced latency and improved resource efficiency. In this design, to reduce the control overhead due to multiple feedbacks to the transmitter, we grouped their feedbacks into a single feedback that represents the current proactive retransmission status.

In this work, a TB is successfully decoded by the receiver if the Signal-to-Noise-plus-Interference-Ratio (SINR) reaches the target SINR. We consider an incremental redundancy HARQ, thus each RTX contains different coded bits than the previous one. According to [START_REF] Lagen | New Radio Physical Layer Abstraction for System-Level Simulations of 5G Networks[END_REF], a singular value of SINR (r) after r RTXs is computed quantitatively based on a set of received Resource Blocks (RBs) ω, their SINR (r) x (∀x ∈ ω) and the previous RTX, i.e. SINR (r-1) tbn as follows:

SINR (r) tbn = -β × ln( 1 |ω| × x∈ω e - SINR (r) x +SINR (r-1) tbn β ) (5)
β is a constant depending on the MCS as specified in [START_REF] Lagen | New Radio Physical Layer Abstraction for System-Level Simulations of 5G Networks[END_REF].

B. Problem Formulation

In this work, we limit the number of decisions into c max clusters of proactive RTXs. Each cluster c j ∈ {c 0 ...c max } can include r cj RTXs. c 0 corresponds to the initial TB transmission. The decision maker we designed will dynamically choose the size of each cluster c j (i.e, r cj ) to reduce both RAN latency and resource waste. With respect to the resources allocated for proactive RTXs of a T B n , the decision maker selects an element-wise positive resource allocation vector (r n,c0 , r n,c1 , ..., r n,cmax ) that satisfies the following condition:

1 ≤ j r n,cj ≤ R max (6)
where R max is the maximum number of RTXs for a T B n .

In the case where a TB is not decoded at receiver after R max RTXs, the TB that contains many application packets is considered a loss. r min ≤ r n,cj ≤ r max constrains the number of proactive RTXs at cluster c j not to exceed a value r max . We can define an objective function f obj that is akin to the average resource allocation to be provided for each TB as follows:

f obj = lim N →∞ 1 N N -1 n=0 ( j r n,cj ) (7)
The decision of how to optimally select r cj is based on various factors, such as the current status of the Q 1 (t), Q 2 (t), the current cluster index c j and the current aggregated SINR.

The first decision at time t is to select the queue Q i (t) ∀i ∈ {1, 2} to proactively (re)transmit the TB over r cj consecutive time slots. Then, a series of actions a ∈ {0, 1, 2, ...} is made at the corresponding action slot t a in which the proactive RTXs span r a time slots. The time interval (t a , t a+1 ) forms the a -th time frame denoted by F a . In order to control at frame F a , which queue among Q i (F ) will be served, we introduced the control variable α a where α a = 1 means serving Q 1 (F a ) and α a = 0 means serving Q 2 (F a ). Knowing that HARQ processes have a higher priority, α a = 0 when Q 2 (t) > 0. Then, the queuing dynamic will be given as follows:

Q 1 (F a+1 ) = max{Q 1 (F a ) -α a .T B rc 0 a , 0} + A 1 (F a ) (8) Q 2 (F a+1 ) = max{Q 2 (F a )-(1 -α a ).1 T B rc j a .T B rc j a , 0} + A 2 (F a ) (9)
where Q i (F a+1 ) are the backlogs of the queue i at time frame F a+1 . The value A 1 (F a ) = Fa t=1 A 1 (t) represents the total amount of high layer packets that accumulated in Q 1 during frame F a . During this frame, an amount of T B rc j a that corresponds to the proactive RTX of T B in r cj slots, will be served. The indicator function 1 T B rc j a , in Equation 9, is equal to 1 if after the r cj proactive RTXs, T B is successful and is 0, otherwise. If the first transmission of T B rc 0 a is a failure, A 2 (F a ) = T B a will be added as a HARQ process to Q 2 , otherwise A 2 (F a ) = 0 as the HARQ process of T B rc 0 a will be removed from Q 2 . According to Little's law, the average delay is related to the queue length. Therefore, the first constraint of our problem is both to minimize the average long-term queue length Q 1 (F a ) + Q 2 (F a ) and to make their average rate stable. By definition, a stochastic process Q(t) is stable at the average rate if:

lim t→∞ E{Q(t)} t = 0 (10)
In the following, we present an additional constraint that is associated with the average long term risk. The visible risk in our problem is when the decision maker at cluster c max chooses to serve the HARQ process in Q 2 by insufficient allocation of r n,cmax to recover the corrupted TB:

ζ n = P[(SINR ( max j=0 rn,c j ) tbn ≤ SINR t ) | SINR ( max-1 j=0 rn,c j ) tbn ] (11)
where SINR ( max j=0 rn,c j ) tbn is the SINR of T B n after c max clusters of RTXs and SINR t is the target SINR.

So, the long-term average risk across all TBs is:

ζ = lim N →∞ 1 N N -1 n=0 ζ n ( 12
)
The constraint is then to guarantee the expectation of the long-term risk ζ below a predefined value ζ o . For that purpose, we introduce a virtual queue Z(t) which is incremented each time slot

t by ζ -ζ o each time ζ exceeds ζ o . Z(t + 1) = max{Z(t) + ζ -ζ o , 0} (13)
Hence, the constraint of satisfying the long-term risk becomes a stability constraint on the average rate of the queue Z(t).

lim t→∞ E{Z(t)} t = 0 (14)
To summarize, our optimization problem P 1 is to minimize the objective function f obj subject to several constraints:

minimize f obj (P 1 ) s.t. lim t→∞ E{Q i (t)} t = 0, ∀i ∈ {1, 2} (C 1,2 ) lim t→∞ E{Z(t)} t = 0, (C 3 ) 1 ≤ cmax j r n,cj ≤ R max , (C 4 ) r min ≤ r n,cj ≤ r max (C 5 )
III. DYNAMIC DECISION MAKER ALGORITHM In this section, we propose a dynamic decision maker algorithm based on Lyapunov's optimization tools to solve the optimization problem P 1 . To simplify the design of the sequential decision maker, we assume that (i) the first decision is the same as the reactive decision (i.e. only one slot is reserved for the TB transmission) and (ii) proactive RTXs in the same cluster share the same HARQ feedback.

First, we define the current state in the slot t as Θ(t) = (Q 1 (t), Q 2 (t), Z(t)) and the Lyapunov function as follows:

L(Θ(t)) ≜ 1 2 Q 2 1 (t) + Q 2 2 (t) + Z 2 (t) (15)
Next, we define the one-slot conditional Lyapunov drift ∆(Θ(t)) representing the expected change of the Lyapunov function over a slot as follows:

∆(Θ(t)) = E{L(Θ(t + 1)) -L(Θ(t)) | Θ(t)} (16)
By minimizing both ∆(Θ(t)) and f obj , we can solve the problem P 1 because the queues are stable in terms of average rate and the objective function is minimized . However, according to [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF], there is a performance-delay trade-off between these dual objective optimizations that can be parameterized by a constant V . By setting a large positive value to V , the control algorithm will favor minimizing the objective function f obj over the stability of the average rate queues. Our ultimate objective now is to minimize the following Lyapunov-driftplus-penalty function:

g(t) = ∆(Θ(t)) + V.E{f obj | Θ(t)} (17)
Its upper bound, γ(t), can be derived as follows for any action, any possible value of Θ(t) and any parameter V > 0 [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF]:

γ(t) =V.E{f obj | Θ(t)} + E{Z(t).(ζ -ζ o ) | Θ(t)} + B + 2 i=1 Q i (t).E{A i (t) -b i (t) | Θ(t)} ( 18
)
where B is a constant that satisfies:

B ≥ 1 2 2 i=1 E{A 2 i (t) -b 2 i (t) | Θ(t)} + 1 2 E{(ζ -ζ o ) | Θ(t)} - 2 i=1 E{A i (t). min{Q i (t), b i (t)} | Θ(t)} (19)
and b i (t) is the quantity of T B n at time slot t that the queue i can process.

b i (t) =    T B rn,1 n if i = 1 T B rn,c j n .1 T B rn,c j n if i = 2 0 otherwise (20)
Under Slater's condition [START_REF] Neely | A Simple Convergence Time Analysis of Drift-Plus-Penalty for Stochastic Optimization and Convex Programs[END_REF], the drift-plus-penalty algorithm provides a time-averaged O(1/ϵ) of the expected size on all actual queues. lim sup

t→∞ 1 t . t-1 τ =0 Q 1 (τ ) + Q 2 (τ ) ≤ B + C + V [f obj,max -f obj,min ] ϵ (21)
where C, ϵ > 0 are constant values.

According to Equations ( 17) and ( 21), by increasing the value of V minimizes f obj at a cost of higher average queue length Q 1 (t) + Q 2 (t) and vice versa.

Through the opportunistic minimization framework of a conditional expectation [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF], by minimizing γ(t), the upper bound of the dual objective optimization, we can guarantee that the optimization problem P 1 will be satisfied. Therefore, the design of our algorithm will be based on the control action a at the decisive time t a and will choose the control action that exhaustively minimizes function γ(t) as follows:

Algorithm 1 Control algorithm 1: Observe time slot t 2: if t ̸ = t a then 3: t = t + 1 4: else 5:
Observe the concatenated queue Θ(t)

6:
Choose optimal action : a = arg min a {γ(t)} 7: end if 8: Observe the outcomes of taken action 9: Update the queue Θ(t + 1)

Since the optimal action is chosen based on the exhaustive search minimizing the designed function γ(t) and the action space is bounded from r min to r max , the computation complexity is relatively low and does not lead to an increase in the transmission processing time.

IV. RESULTS AND DISCUSSIONS

In this section, we compare the performance of several HARQ strategies (i.e. our proactive adaptation algorithm, a fixed proactive and a reactive) between a gNB and single mobile user. Initially, the UE is placed at a distance d 0 from the gNb and moves away from it with constant speed and direction (v x , v y ). In this work, packets are generated in distinct ON and OFF periods that follow the Internet Protocol (IP) traffic model [START_REF] Marvi | On the Use of ON/OFF Traffic Models for Spatio-Temporal Analysis of Wireless Networks[END_REF]. The average duration of the ON and OFF periods are t ON and t OF F , respectively. In the ON state, packets of variable size are generated with an arrival rate of λ ON and fill Q 1 (t). The performance is evaluated in terms of RAN latency, reliability outage, packet loss and resource efficiency. We define the resource efficiency as the ratio of the number of radio resources required for a TB to be received by the receiver to the number of radio resources allocated by the scheduler. We also define the RAN latency as the times between the arrival of IP packets in the RLC layer of the gNb and their arrival in the IP layer at the UE side. In the scheduling process, K 1 and L 12 [START_REF] Gpp | Medium Access Control (MAC) protocol specification[END_REF] are modelled to illustrate the feedback processing time and data preparation time at the UE and gNB, respectively. For simplicity, we assumed that the core network latency and propagation delay are negligible. Table I summarizes all parameters including communication band, transmit power P tx , target Block Error Rate (BLER) ϵ t , maximum number of proactive cluster c max and proactive retransmission per cluster r max and risk threshold ζ o .

Figure 3 illustrates how the objective function and resource efficiency behave as a function of the parameter V. For small value of V, our algorithm tends to minimize the expected changes in the Lyapunov function, i.e L(θ(t)) rather than the number of allocated resources i.e. f obj . By doing so, a large number of radio resources are generously provided to each TB for proactive RTXs, so the average resource allocation is high and the resource efficiency is low. When the value of V increases, the focus is on minimizing the objective function and fewer resources are allocated. At a certain value of V (i.e. around 50), we found that a minimum value of the objective function is reached (i.e. around 3.2). This is because at a high value of V, the algorithm favors allocating smaller number of parallel RTXs on each cluster and thus more proactive RTX clusters are required for each TB. Therefore, in our future evaluations, we will set V to 60. 4 shows that the reactive HARQ scheme performs the worst compared to the other solutions due to the high RTT cost associated with triggering many RTXs. By enabling two parallel RTXs per received NACK, the latency can be improved at the cost of decreasing resource efficiency to around 0.82 compared to 1 in the reactive case. The latency can be further improved with 3 proactive RTXs per cluster, but the resource efficiency drops significantly to 0.69 when we redundantly provide radio resources for a TB to be successfully decoded. Concerning the performance of our algorithm with the adaptation of proactive HARQ, the performance of latency is slightly better than the others cases. As the channel condition and traffic behaviour varies dynamically, our dynamic decision maker will wisely select different redundancy level according to the current queues status and instantaneous channel condition. As the result, we can enhance the latency while keeping a good level of resource efficiency at around 0.8. To show the temporal variation in MAC layer delay of the different HARQ schemes for intermittent URLLC traffic scenarios under time-varying channel conditions, Figure 5 compares the latency between the departure of a TB at the sender's MAC layer and its successful arrival at the receiver's MAC layer. The reactive HARQ scheme experiences more peaks as it needs more time to successfully decode a packet when critical errors occur. When fixed proactive HARQ schemes are applied, 2-parallel strategy reduces the peak delay and 3-parallel strategy further improves the delay when the radio resources allocated in parallel help decode the packet faster. Our proposed algorithm also reduces the number and amplitudes of peaks when it dynamically selects an optimal level of proactive redundancy for each cluster, thus achieving a better trade-off between reserving more radio resources to recover the corrupted packet faster and less to maintain a minimum objective function.

Figure 6 shows the latency as a function of V value for an outage level of 0.9 and 0.95, respectively. Certain values of V make the outage latency minimal because they allow a better balance between the number of resources allocated for each packet and the stability of the queues. As the result, the effect of overloaded buffer is better managed and leads to a good latency outage. On the other hand, by setting V too high or too low, this balance is no longer efficiently controlled and it causes a queuing effect due to the fact that too many clusters are used (large V) or too many retransmissions in each cluster are set (small V).

In addition, various performance metrics such as application packet loss (APP loss), resource efficiency as well as average delay and its standard deviation for different HARQ schemes are discussed. First, the APP loss is guaranteed to be between 0.8 % and 1 % because the maximum number of retransmissions (i.e. R max = 10) is applied in all cases. Second, we find the largest average delay and its standard deviation in reactive HARQ whose values are around 10.77 ms and 6.57 ms, respectively. These coupled values are enhanced to roughly 8.14 ms, 4.91 ms and 8.1 ms, 4.5 ms in fixed 2-parallel and 3 parallel HARQ schemes at the expense of resource efficiency which decreases to 0.82 and 0.69, respectively. By dynamically deciding the number of proactive retransmissions at V=60, the average delay and jitter are slightly better than other HARQs which are about 7.2 ms and 4 ms while maintaining a good resource utilization level at about 0.8.

V. CONCLUSIONS

In this study, we demonstrate the application of dynamic decision maker framework that enables a novel proactive HARQ design to cope with a time-varying channel and intermittent traffic source rate. Based on Lyapunov stochastic optimization tool, a mathematical framework is proposed to understand the performance-delay trade-off by minimizing the objective function of the total resource allocation and the total queue length that is parameterized by a V value. Our results reveal that an appropriate selection of V enables the dynamic selection of proactive retransmission to overcome the defect of a long RTT in the reactive scheme while maintaining a good level of resource efficiency that is considered a drawback of the fixed proactive schemes.

In future work, this study will be extended to dynamic resource scheduling. The resource efficiency-latency-reliability trade-off will be achieved by the number, timing and intensity of the decisions. Moreover, we will implement our optimization for real experimentation with OpenAirInterface.

Fig. 1 .

 1 Fig. 1. System model To calculate the TB size, the TDMA-based scheduler takes into account the number of OFDM symbols allocated to user in a single time slot (1 ≤ N OF DM ≤ 12), as well as the choice of the Modulation and Coding Scheme (MCS) m, the numerology N um, the bandwidth BW , as follows:

Fig. 2 .

 2 Fig. 2. Different DL HARQ procedures: (A) classical, (B) fixed R=3 parallel RTXs and (C) Dynamic and proactive RTXs.

Fig. 3 .Figure 4

 34 Fig. 3. Objective function and resource efficiency as a function of V Figure 4 compares the Complementary Distribution Function (CDF) of latency for different HARQ schemes: (i) reactive HARQ scheme (i.e. R = 1) with a maximum number of 10 RTXs, (ii) 5 RTX clusters where each cluster contains 2 proactive RTXs (2-2-2-2-2 clusters), (iii) 4 RTX clusters where the first three clusters contain 3 proactive RTXs and the last cluster contains a single RTX (3-3-3-1 clusters) and (iv) our proposition with a maximum of 5 clusters.Figure4shows that the reactive HARQ scheme performs the worst compared to the other solutions due to the high RTT cost associated with triggering many RTXs. By enabling two parallel RTXs per received NACK, the latency can be improved at the cost of decreasing resource efficiency to around 0.82 compared to 1 in the reactive case. The latency can be further improved with 3 proactive RTXs per cluster, but the resource efficiency

Figure

  Figure 4 compares the Complementary Distribution Function (CDF) of latency for different HARQ schemes: (i) reactive HARQ scheme (i.e. R = 1) with a maximum number of 10 RTXs, (ii) 5 RTX clusters where each cluster contains 2 proactive RTXs (2-2-2-2-2 clusters), (iii) 4 RTX clusters where the first three clusters contain 3 proactive RTXs and the last cluster contains a single RTX (3-3-3-1 clusters) and (iv) our proposition with a maximum of 5 clusters.Figure4shows that the reactive HARQ scheme performs the worst compared to the other solutions due to the high RTT cost associated with triggering many RTXs. By enabling two parallel RTXs per received NACK, the latency can be improved at the cost of decreasing resource efficiency to around 0.82 compared to 1 in the reactive case. The latency can be further improved with 3 proactive RTXs per cluster, but the resource efficiency

Fig. 4 .

 4 Fig. 4. CDF of latency for reactive RTX, fixed 2-parallel RTX, fixed 3-parallel RTXs and our dynamic adaptation algorithm with V=60

Fig. 5 .Fig. 6 .

 56 Fig. 5. Evolution of MAC delay between HARQ schemes

ACKNOWLEDGMENT

This work was supported by the European Union H2020 / Taiwan Project 5G CONNI [19] under grant n.861459.


Ngoc Lam

Mickael Dinh

Emilio Calvanese Maman

Strinati

Proactive Resource Scheduling for 5G and Beyond Ultra-Reliable Low Latency Communications

Effective resource use in Ultra-Reliable and Low-Latency Communications (URLLC) is one of the main challenges for 5G and beyond systems. In this paper, we propose a novel scheduling methodology (combining reactive and proactive resource allocation strategies) specifically devised for URLLC services. Our ultimate objective is to characterize the level of proactivity required to cope with various scenarios. Specifically, we propose to operate at the scheduling level, addressing the trade-off between reliability, latency and resource efficiency. We offer an evaluation of the proposed methodology in the case of the well-known Hybrid Automatic Repeat reQuest (HARQ) protocol in which the proactive strategy allows a number of parallel retransmissions instead of the "send-wait-react" mode. To this end, we propose some deviations from the HARQ procedure and benchmark the performance in terms of latency, reliability outage and resource efficiency as a function of the level of proactivity. Afterwards, we highlight the critical importance of proactive adaptation in dynamic scenarios (i.e. with changing traffic rates and channel conditions).

I. INTRODUCTION

5G and beyond networks enable the exploitation of new emerging use cases, such as Ultra-Reliable Low Latency Communication (URLLC) which plays a pivotal role in the connectivity of highly demanding services and applications in factory automation, autonomous driving, etc [START_REF] Li | 5G URLLC: Design Challenges and System Concepts[END_REF]. In addition to presenting the key performance indicators related to URLLC, namely reliability and latency, it should be note that resource usage is also not a simple matter, as it is still unclear how to efficiently utilize resources while maintaining reliable, low latency communications [START_REF] She | Radio Resource Management for Ultra-Reliable and Low-Latency Communications[END_REF]. Given the requirements for delay and reliability, two approaches can activate available resources for the communication. On the one hand, reactive strategies activate additional resources on demand, which allows for efficient resource utilization, but significantly increases latency, as the demand for additional resources is not instantaneous. On the other hand, the proactive approach is proposed to systematically apply additional resources to stretch the latency below the deadline and it usually activates for the worst case with a margin. As a result, this approach implies a high cost in terms of resource utilization, especially when worst-case impairments are very rare. Our goal is then to define the level of proactivity to apply, to cope with various scenarios, at the scheduling level, taking into account the trade-off between reliability, latency and resource efficiency. This study will help in the design of early decision maker, as patented in [START_REF] Maman | Procede et dispositif d'orchestration de l'execution de mecanismes dans un reseau sans fil[END_REF], dynamically adapting reactive-proactive modes required in dynamic scenario.

In order to highlight the importance of the above-mentioned trade-off for improving URLLC, we investigate the direct application of reactive and proactive strategies to the well-known Hybrid Automatic Repeat reQuest (HARQ) protocol. In a dynamic and rapidly changing environment, the implementation of HARQ protocol is a simple and effective method to ensure communication resiliency against deep channel fading that occurs unexpectedly during transmission and provides us with a good level of efficient resource utilization [START_REF] Shi | Analysis of HARQ-IR Over Time-Correlated Rayleigh Fading Channels[END_REF]. In the usual reactive HARQ scheme, retransmission (RTX) is triggered as soon as the sender receives the negative acknowledgement (NACK) from the receiver. Thus all allocated resources are used efficiently at the cost of a long Round Trip Time (RTT) trigger if many RTXs are required. In the literature, adaptation of the HARQ strategy is usually achieved by adapting the modulation and coding scheme [START_REF] Pfletschinger | Adaptive HARQ With Non-Binary Repetition Coding[END_REF], the transmission power [START_REF] Jabi | Energy Efficiency of Adaptive HARQ[END_REF] and the maximum number of RTXs [START_REF] Mukhtar | Content-aware and occupancy-based hybrid ARQ for video transmission[END_REF] but rarely at the RTX scheduler level. A K-repetition scheme [START_REF] Le | An Overview of Physical Layer Design for Ultra-Reliable Low-Latency Communications in 3GPP Releases 15, 16, and 17[END_REF] and a proactive scheme with early termination [START_REF] Liu | Analyzing Grant-Free Access for URLLC Service[END_REF] have been proposed, allowing for a number of redundant retransmissions upon receipt of the acknowledgment by the sender. By doing so, one can opportunistically decode the packet at the receiver in a shorter time at the expense of inefficient resource usage [START_REF] Jacobsen | System Level Analysis of Uplink Grant-Free Transmission for URLLC[END_REF]. However, the adaptation of HARQ strategies at the scheduling level in dynamic scenarios is limited in current research.

The contributions of this paper are as follows: (1) We show how proactive resource scheduling can be included in the regular HARQ procedure, validating the implementation before a real experiment. [START_REF] She | Radio Resource Management for Ultra-Reliable and Low-Latency Communications[END_REF] We consider end-to-end performance by developing a system-level simulator based on NS-3 [START_REF] Patriciello | An E2E simulator for 5G NR networks[END_REF] applying to the 5G New Radio (NR). This simulator handles several HARQ processes and measures the latency between the transmitter Radio Link Control (RLC) layer and the receiver RLC layer assuming that the transmission buffer size is infinite. We therefore consider both the queuing latency at the scheduler (due to reactive/proactive approaches) and (re)transmission latency (i.e. PHY/MAC). (3) We highlight the crucial importance of proactive HARQ adaptation for enhancing the reliability, latency and resource efficiency tradeoff in dynamic scenarios. We consider several levels of proactivity, different traffic rates and channel conditions. The results of this study pave the way for the design of an online decision maker framework enabling to dynamically choose the appropriate level of proactivity based on the network state [START_REF] Dinh | Toward URLLC with Proactive HARQ Adaptation[END_REF].

The remainder of the paper is organized as follows: Section II presents a detailed scenario setup where PHY and MAC layer are fully described. Then, Section III formulates the problem. In Section IV, the simulation results demonstrate the performance of the system and Section V concludes the paper.

II. SYSTEM MODEL A. Scenario Description

The network contains 1 gNb and 1 UE at a fixed distance d. It is assumed that a downlink (DL) traffic is generated by a remote host located near the gNb with an arrival rate of λ and a packet size of L Bytes. Given the slot duration t s , the amount of traffic that fills the buffer at each slot n is A(n) = λ.L.t s . In order to improve the diversity order of the communication, U tx = M BS × N BS antennas and S rx = M U E × N U E antennas are installed at the gNB and UE, respectively.

The performance of the system is evaluated for an Indoor Factory-sparse clutter scenario (InF-SL) modeled according to the 3GPP technical report [START_REF]Study on channel model for frequencies from 0.5 to 100 GHz[END_REF]. This model clarifies the channel behaviours based on antenna modeling, slow-fading factors such as shadowing, propagation loss and fast-fading factors due to multi path. Initially, the propagation loss is calculated based on whether the path is in visibility or not (LOS or NLOS) and can be derived as follows:

P L LOS (d) = 31.84+21.50×log 10 (d)+19.00×log 10 (f c ) (1) P L N LOS (d) = 33+25.50×log 10 (d)+20.00×log 10 (f c ) (2)
According to the channel statistics, the probability of LOS communication P r LOS between the transmitter and receiver can be derived as follows:

P r LOS = e d×ln(1-r) d clutter (3)
where d is the 2D distance, d clutter = 10 is the typical clutter size in meter and r = 0.3 is the clutter density which provides the statistical knowledge of the Indoor factory environment.

In this study, the position between the gNB and the UE is fixed, so the shadowing effect is negligible. Concerning the fast fading model, the sub-channel between the gNB antenna u and the UE antenna s is modeled by channel impulse response H u,s (t, τ ). Then, the MIMO channel matrix at transmission slot n e.g., H n is formed by the superposition of single subchannels [H u,s (t, τ )] which were described in [START_REF]Study on channel model for frequencies from 0.5 to 100 GHz[END_REF]. In this paper, we also consider the scenario where the channel fading is too fast to apply instantaneous channel quality feedback to the transmitter, resulting in the transmitter has little knowledge of the channel statistics or only having the outdated delay versions of the feedback.

B. MAC configurations

In our study, HARQ -Incremental Redundancy (IR) is assumed. Each RTX contains different coded bits than the previous RTX. The quality of data transmission over different Resource Blocks (RBs) is modeled using Effective Exponential SINR Mapping (EESM) method in which a single effective signal-to-noise ratio SINR ef f is combined from individual SINR received from individual RBs. In HARQ-IR, the SINR r ef f after r RTXs is derived as follows [START_REF] Lagen | New Radio Physical Layer Abstraction for System-Level Simulations of 5G Networks[END_REF]:

SINR r ef f = -β × log( 1 |ω| × x∈ω e - SINRx,r +SINR r-1 ef f β ) (4)
where SINR r-1 ef f is the effective SINR after the previous RTX, i.e., r -1 retransmissions, SINR x,r is the SINR experienced by the x-th RB in the r-th RTX, and ω is the set of RBs. The value of β depends on the MCS selection.

Afterwards, the BLER of a MAC-layer transport block (TBLER) as a function of the radio channel, the link-to-system mapping (L2SM) method was exploited for the system level simulation. In doing so, a SINR-BLER lookup table can be used for the calculation of the TBLER with respect to the MCS selection m associated with LDPC coding that corresponds to multiple lifting sizes following the 3GPP specifications in [START_REF]Multiplexing and channel coding[END_REF].

In terms of MAC procedure, the communication between the gNB and the UE is done on a slot basis, where the slot duration t s is preconfigured by selecting the numerology (N um). The TDMA-based scheduler is capable of reserving a maximum of 12 OFDM data symbols in a slot, depending on the state of the transmission buffer queue. Based on the number of OFDM symbols allocated to users in a single slot (1 ≤ N OF DM ≤ 12), as well as the choice of MCS m, the numerology N um, the bandwidth BW , the Transport Block Size TBS(m) is defined as follows:

TBS(m) = N OF DM × BW 2 N um × SCS(N um) × M (m) 8 × CR(m) (5)
where SCS(N um) is the sub-carrier spacing at the numerology N um, M (m), CR(m) are the modulation order and code rate of the selected MCS m, respectively. These values are given in Table 5.1.3.1 of 3GPP TS 38.214 [START_REF]Physical layer procedures for data[END_REF].

In order to prepare for data transmission, a delay of a slots is configured at the scheduler to put the data over-the-air. At the receiver, the packet is considered corrupted if the TBLER is greater than a reliability constraint given by a target error t i.e, the actual SINR r ef f received after RTX r is less than a target SINR t . In this case, the NACK will be returned to the sender after a delay of b slots. Basically, depending on the channel statistics and the selected MCS m, the probability of a packet being successfully received after r RTXs is:

P s (m, r) = P(SINR r ef f ≥ SINR t ) (6)
Concerning the length of transmission buffer queue, the accumulated data in the transmission buffer queue until time slot n + 1 is determined as follows:

Q(n + 1) = Q(n) + A(n) - n+1 i=0 1 × TBS(m) (7)
where Q(n) and A(n) are the queue state and the arrival event at time slot n, respectively. For the departure event, a quantity of transport block size TBS(m) will be served to current slot when it received an ACK feedback.

For the queue stability condition, the departure rate (D) must be greater than or equal to the arrival rate (A) to avoid buffer overflow. Assuming that out of N time slots, there are N 1 slot(s) (N 1 ≤ N ) that do not have an ongoing NACK feedback, the queue stability condition can be expressed as follows:

E(A) E(D) = lim N →∞ N × L × λ × t s N 1 × TBS(m) ≤ 1 (8)
Equation 8 shows that either N 1 or TBS(m) must be large enough to prevent an adverse queue delay that is greater than the target delay. This equation also shows the trade-off between packet latency and reliability when running RTX. By triggering HARQ protocol in the current slot, we can improve the communication reliability at the cost of increasing the queue delay caused by the departure process to zero.

III. PROBLEM FORMULATION In this section, we make a tradeoff between end-to-end latency and resource usage when the R-parallel RTXs scheme is applied. By opportunistically triggering retransmissions before their feedback, the corrupted packet can be recovered faster than the classical approach i.e, the over-the-air latency L c of the transmitted packet is shortened. On the other hand, the protocol has to reserve unnecessary slots to retransmit the already decoded message, which results in low resource utilization. We formulate the optimization problem to find an optimal R-parallel RTXs scheme that reduces the end-to-end latency L e2e and guarantees good resource usage, as follows:

min R E( L e2e η R,ef f ) = min R E( L c (R) + L q (Q) η R,ef f ) s.t. R ∈ [0, 1, ..., R max ] E(A) E(D) ≤ 1 (9)
We defined the resource efficiency η R.ef f as the ratio of resources required to accomplish packet transmission to those provided. The end-to-end delay L e2e of the generated packets consists of two main parts, namely the queuing delay at the transmission buffer L q and the over-the-air delay L c which includes the retransmission process.

Figure 1 depicts the DL HARQ scheme in which reactive single retransmission (R = 1) and proactive parallel retransmission (R = 3) are compared. The details of how the corrupted packet is resent in these 2 different cases are revealed below. The classical HARQ procedure corresponds to a reactive retransmission procedure (i.e. R = 1). Figure 1A depicts DL data flow from an application server, which is close to the gNB, and the UL feedback from the UE side. The procedure of reactive HARQ can be divided into several steps:

• Step 1: Data packets are generated periodically at the application layer with a size of L bytes and are queued in the gNb transmission buffer. • Step 4: The UE feedback is acquired by gNb. A NACK feedback triggers a retransmission of the corrupted TB (RV=1) that was previously recorded with the sender. From there, steps 2 and 3 can be repeated until the maximum number of allowed RTXs is reached or the packet is successfully received by the receiver (e.g., RV=2 as shown in Figure 1A). In contrast, an ACK gives the scheduler the opportunity to collect new TBs from the queue and repeat step 1.

Based on the statistic of the probability that the packet is successfully decoded after being retransmitted r times with a given MCS m e.g P s (m, r), the average latency L c,r of reactive HARQ is given as follows:

L c,r = P s (m, 0) × a × t s + Kmax r=1 P s (m, r) × r × T × t s ( 10
)
where T = a + b + 1 refers to the RTT of the retransmission.

Since the receiver requires retransmission of the error packet after receiving the NACK, resource efficiency is always guaranteed, η R.ef f = 1 (i.e., the number of resources requested to accomplish the retransmission is equal to those allocated).

B. R-Parallel HARQ (R > 1)

R-Parallel HARQ deviates from the last step of regular reactive HARQ as follows:

• Step 4': Once NACK is acquired at the sender, instead of systematically retransmitting the single corrupted TB, multiple TBs is scheduled consecutively in adjacent transmission slots. The level of redundancy R is determined by the decision maker taking into account that a higher value of R improves latency in exchange for resource efficiency and vice versa. This approach is illustrated in Figure 1B where R = 3 refers to the amount of resources the scheduler gives to the sender. In this case, instead of processing a single transmission, the receiver must process a stream of RTXs (RV=1,2,3) to generate their corresponding feedbacks. Again, these feedbacks must undergo a constant delay of b slots before reaching the gNb's scheduler. Since the first feedback is a NACK and the others are ACKs, a wasted RTX (RV=4) will be triggered before the scheduler realized that the packet has already been decoded at RV=2. Therefore, RV=3 and RV=4 are considered wasted resources. The resource efficiency can be derived as follows:

η R.ef f = Kmax r=1 r × P s (m, r) min{ Kmax r=1 r × P s (m, r) + R, K max } (11)
The underlying over-the-air latency in this HARQ scheme is based on many factors that take into account the redundancy level R chosen by the gNb scheduler and the over-the-air queue caused by the scheduler operation. In general, the sender opportunistically retransmits the corrupted packet R times in R consecutive slots and thus postpones the scheduling of others packets. In this case, the over-the-air latency is given by:

L c,p = P s (m, 0) × a × t s + E Kmax r=1 r+2 R ×T +(r+2)modR ×P s (m, r)×t s ( 12
)
where a is the nearest integer less than or equal to a and mod is the modulo operator. This formulation consists of two parts: the average over-the-air latency of the initial transmission and the average latency of the RTXs. In the example shown in Figure 1B where a = 2, b = 4 (then T = 7) and R = 3, the distance in TTIs between the initial transmission (RV=0) and second retransmission (RV=2) is 8 slots.

IV. RESULTS AND DISCUSSION In this section, we evaluated the performance in terms of latency, reliability outage and resource efficiency, for different traffic rates and channel conditions. The simulation results are averaged over 10 random seeds. Latency is measured from the time the packet arrives in the RLC buffer of the gNB to the time it leaves the RLC buffer at the UE side, assuming that core network latency, signal processing latency and propagation delay are negligible. Each packet generated by a remote server has a fixed length of 60 Bytes. The arrival rates λ are 100 and 800 for low and high rate data generation, respectively. In terms of PHY/MAC configurations, the bandwidth of DL communications is 50 MHz in the sub-6GHz band (i.e. f c = 3.5GHz). The uniform linear array antenna designs are operated on the gNB and UE side, in which U tx = 8 × 4 and S rx = 4 × 4 antenna elements are equipped respectively. Since the gNb scheduler does not know the instantaneous link quality of the communication, a fixed transmission power P tx = 10 dBm is distributed equally among the transmitted antenna elements. Besides, a fixed MCS m = 5 with its corresponding spectral efficiency η S.ef f = 0.7402 and a code rate CR = 0.3701 are chosen according to [START_REF]Study on channel model for frequencies from 0.5 to 100 GHz[END_REF]. The numerology N um is set to 1 for the slot configuration. The TDMA-based scheduler takes a = 2 slots to prepare the data stored in the transmission buffer and push it to the air. Once the received payload is processed by the receiver, a delay b = 4 slots is required to send the feedback (ACK/NACK) to the transmitter. The TBLER target threshold is set to t = 1E -4. If this reliability requirement is not met before the maximum number of RTXs (i.e. K max = 10) is reached, the packet is considered lost. Finally, the performance is evaluated in a bad/good channel where the UE is placed far/near the gNb at the distance d = 150m and d = 100m, respectively. Table I summarizes the simulation parameters. Figure 2 shows the Complementary Distribution Function (CDF) of the latency for different HARQ schemes: classical (i.e. R = 1) and proactive (i.e. R > 1). In this case, a low application traffic rate is associated with a bad channel condition to perform DL communications between the gNB and the UE. Compared to the reactive approach (R = 1), a significant latency gain is achieved when a higher redundancy level R is set. By increasing R, the retransmission process will be accelerated by consecutively scheduling the radio resources for the retransmission and prior to the ACK/NACK reception. As a result, the delivery of the corrupted packets from the source to the destination is both reliably guaranteed and considerably faster than the classical method. To evaluate the impact of channel conditions and arrival traffic rates, Figure 3 plots the latency CDFs for a redundancy level R = 3. Under very good channel condition, a low traffic regime performs best compared to the other scenarios, as shown in particular by the probability of successful transmission on the first communication. When packets are generated by the application at very high rates, and a packet is undecodable at the receiver, the aggregation of queuing effects occurs at the transmit buffer and the radio at the same time, resulting in degraded system performance. In this regard, a poor channel condition or a high traffic rate scenario will severely prolong the overall packet transmission delay. Figure 4 shows the gain in latency and jitter as a function of the redundancy level R for different channel conditions and traffic rates. Latency and jitter are improved compared to the classical HARQ (R = 1) regardless of the channel condition and traffic generation behaviour. In low traffic mode, the increased redundancy further improves performance. In contrast, in high traffic mode, performance is optimized at a particular redundancy level before it becomes saturated and degrades. This is because queue aggregation causes buffer latency to increase and moderates the latency gain. Figure 5 shows the outage probability over latency (i.e. the latency required to reach an outage from 0.9 to 0.99999 with different R). It demonstrates how the latency gap is optimized to reach a more critical outage from a more relaxed outage by increasing R. In general, the cost to reach to critical outages in the R-parallel RTX scenarios is lower than the classical method (R = 1). The benefit of significantly reduced outage is clearly demonstrated in the low traffic regime where latency outages continue to decrease as R increases. For high traffic, the latency required is less expensive than reactive HARQ (R = 1) and starts increasing after an optimal value of R. Regarding resource efficiency, we observe a significant decrease when switching from a reactive (R = 1) to a proactive (R > 1) strategy. In particular for R = 2, the resource efficiency is reduced to about 0.64 and 0.59 for the bad and good channel, respectively. This couple of values continues to decrease to (0.5, 0.43), (0.4, 0.36) and (0.32, 0.26) for R = 3, R = 4 and R = 5 accordingly. When the channel condition is good, proactively enabling additional RTXs often results in redundancy of allocated resources. Thus, we are more likely to cause unnecessary resource allocation and efficiency reduction.

V. CONCLUSION In this study, we focused on different scheduling strategies of HARQ retransmissions for URLLC in a system level sim-ulator. We proposed some deviations of HARQ procedure and evaluated the tradeoff between reliability, latency and resource efficiency by comparing the performance of reactive HARQ and proactive HARQ as a function of traffic source rates and channel conditions. This studies clarified the correlation between traffic source rate and PHY/MAC configurations to serve the data over the air in a highly dynamic channel and highlighted the aforementioned tradeoff by analyzing the application of different redundancy levels to retransmit corrupted packets in a classical or R-parallel manner. This study has been extended to dynamic redundancy adaptation based on the current network state to aim for optimal system performance and we are currently implementing it for real experimentation with OpenAirInterface.

Fig. 1 .

 1 Fig. 1. Schema of (A) classical DL HARQ scheme and (B) R=3-parallel RTXs DL HARQ scheme

Fig. 2 .

 2 Fig. 2. CDF of latency for different redundancy levels R

Fig. 3 .

 3 Fig. 3. CDF of latency for different channel conditions and traffic rates.

Fig. 4 .

 4 Fig. 4. Latency as a function of the redundancy levels, channel conditions and traffic source rates

Fig. 5 .

 5 Fig. 5. Outage probability on latency as a function of the redundancy levels, channel conditions and traffic source rates

ACKNOWLEDGMENT

This work was supported by the European Union H2020 / Taiwan Project 5G CONNI [17] under grant n.861459.


Mattia Merluzzi
email: mattia.merluzzi@cea.fr

Serge Bories
email: serge.bories@cea.fr

Emilio Calvanese Strinati
email: emilio.calvanese-strinati@cea.fr

Energy-Efficient Dynamic Edge Computing with Electromagnetic Field Exposure Constraints

Keywords: Green communications, EMF, Edge Computing

We present a dynamic resource allocation strategy for energy-efficient and Electromagnetic Field (EMF) exposure aware computation offloading at the wireless network edge. The goal is to maximize the overall system sum-rate of offloaded data, under stability (i.e. finite end-to-end delay), EMF exposure and system power constraints. The latter comprises end devices for uplink transmission and a Mobile Edge Host (MEH) for computation. Our proposed method, based on Lyapunov stochastic optimization, is able to achieve this goal with theoretical guarantees on asymptotic optimality, without any prior knowledge of wireless channel statistics. Although a complex long-term optimization problem is formulated, a per-slot optimization based on instantaneous realizations is derived. Moreover, the solution of the instantaneous problem is provided with closed form expressions and fast iterative procedures. Besides the theoretical analysis, numerical results assess the performance of the proposed strategy in striking the best trade-off between offloading sum-rate, power consumption, EMF exposure, and E2E delay. To the best of our knowledge, this is the first work addressing the problem of energy and exposure aware computation offloading.

I. INTRODUCTION

The sixth generation of mobile communication networks (6G) targets a radical long-term transformation of wireless systems, with a pervasive deployment of computing resources at the edge, as well as a myriad of new connections spanning from conventional users to Internet of Things (IoT) and vertical sectors such as autonomous vehicles and Industry 4.0 [START_REF] Strinati | 6G networks: Beyond Shannon towards semantic and goal-oriented communications[END_REF]. In this new perspective, human, physical and digital worlds will be embedded into the same ecosystem [START_REF] Uusitalo | Hexa-X The European 6G flagship project[END_REF], making network management increasingly complex. In this context, extremely large amounts of capillary data will be wirelessly transferred between intelligent edge devices to operate data distillation, computation offloading services, training and inference of distributed (and/or federated) Artificial Intelligence (AI) and Machine Learning (ML) models, etc. In this direction, Multi-access Edge Computing (MEC) [START_REF] Pham | A Survey of Multi-Access Edge Computing in 5G and Beyond: Fundamentals, Technology Integration, and State-of-the-Art[END_REF] supports the network evolution toward 6G [START_REF] Strinati | 6G: The next frontier: From holographic messaging to artificial intelligence using subterahertz and visible light communication[END_REF], by bringing data processing at the network edge, thus close to the end consumers. Among several services, MEC networks allow resource poor end devices to offload computational tasks to nearby processing units. This comes with the challenge of a communicationcomputation co-design of wireless networks [START_REF] Merluzzi | Dynamic Computation Offloading in Multi-Access Edge Computing via Ultra-Reliable and Low-Latency Communications[END_REF], since computation offloading services involve both transmission and computation delays and energy consumption. Several works addressed the problem of striking the best trade-off between energy consumption, service delay and/or application level This work has been partly funded by the European Commission through the H2020 project Hexa-X (Grant Agreement no. 101015956).

performance (e.g. AI/ML training/inference accuracy) [START_REF] Mao | Stochastic joint radio and computational resource management for multi-user mobileedge computing systems[END_REF]- [START_REF] Merluzzi | Wireless Edge Machine Learning: Resource Allocation and Trade-Offs[END_REF]. However, these works do not consider Electromagnetic Field (EMF) human exposure, which will increase due to continuous communications between the myriad of sensors and wireless Access Points (APs) and a drastic uplink traffic growth [START_REF] Oueis | Uplink traffic in future mobile networks: Pulling the alarm[END_REF]. Therefore, sustainable operations of future wireless networks impose to jointly consider: i) Performance; ii) Energy efficiency; iii) EMF exposure awareness. While all these aspects have been previously disjointly investigated, a holistic view is lacking and needed towards future networks deployment and operation phases. Indeed, while traditional network optimization takes into account Key Performance Indicators (KPIs) such as data rate and latency, 6G networks will be designed to also strengthen less tangible indicators, termed as Key Value Indicators (KVIs) [START_REF] Uusitalo | Hexa-X The European 6G flagship project[END_REF], such as the EMF exposure, with the usual monitoring of new recommendations [START_REF]ICNIRP guidelines for limiting exposure to electromagnetic[END_REF], fundamental to run operators' networks with the best performance while not exceeding exposure limits, and also for a smooth public acceptance of new technologies. Related works. Previous works focused on EMF-aware networking and hardware perspectives [START_REF] Zhao | EMF Exposure Study Concerning mmWave Phased Array in Mobile Devices for 5G Communication[END_REF]- [START_REF] Chiaraviglio | 5G Network Planning under Service and EMF Constraints: Formulation and Solutions[END_REF]. The authors of [START_REF] Chiaraviglio | 5G Network Planning under Service and EMF Constraints: Formulation and Solutions[END_REF] present the problem of network planning, considering a weighted function of gNB installation costs and 5G service coverage as objective, while [START_REF] Galán-Jiménez | Measuring the impact of ICNIRP vs. stricter-than-ICNIRP exposure limits on QoS and EMF from cellular networks[END_REF] analyzes the impact of more or less strict regulations on network performance, showing how stricter limitations degrade network performance. This straightforward yet important result motivates us to investigate trade-offs between network power consumption, EMF exposure, and performance of MEC-enabled services. In [START_REF] De Domenico | EMF-Aware Cell Selection in Heterogeneous Cellular Networks[END_REF], the authors study a load balancing scheme in heterogeneous networks to achieve target data rate and EMF exposure. The work has been carried out within the project LEXNET [START_REF] Tesanovic | The LEXNET Project: Wireless Networks and EMF: Paving the Way for Low-EMF Networks of the Future[END_REF]. Our contribution. We focus on dynamic computation offloading, with end devices continuously generating data to be uploaded to a Mobile Edge Host (MEH) and processed [START_REF] Merluzzi | Dynamic Computation Offloading in Multi-Access Edge Computing via Ultra-Reliable and Low-Latency Communications[END_REF]- [START_REF] Merluzzi | Wireless Edge Machine Learning: Resource Allocation and Trade-Offs[END_REF]. Exploiting Lyapunov stochastic optimization, we design an online algorithm able to jointly allocate radio and computing resources, achieving the maximum system sum-rate under stability constraints, as well as system power constraints, and exposure constraints over intended areas. Our method does not require any prior knowledge of wireless channel statistics, but only the solution of simple per-slot problems that enjoy closed form expressions or fast iterative procedures. Also, it comes with theoretical guarantees on sum-rate asymptotic optimality, system stability, as well as power and EMF exposure guarantees. To the best of our knowledge, this is the first work with a joint view of offloading performance, system power consumption, and EMF exposure. We consider a dynamic scenario, with time organized in slots t = 1, 2, . . . of duration τ . In our setting, computation offloading entails two phases: i) Uplink data transmission; ii) Computation. Thus, we design a 2-hops queueing system, to model uplink transmission and computation. Let us define the long-term average of a random variable X as

X = lim T →∞ 1 T T t=1 E {X t } , (1)
where the expectation is generally taken with respect to random context parameters (e.g. time-varying wireless channels).

The above notation (X) will be used throughout the paper to denote the long-term average of all involved variables. In each time slot t, every device k = 1, . . . K stores A k,t new bits into an uplink buffer Q u k,t . As it will be clarified later on, A k,t will be optimized online to achieve the maximum sumrate of offloaded data. At the same time, device k possibly transmits data, buffered in previous slots, with data rate R u k,t . Thus, the uplink queue evolves as follows:

Q u k,t+1 = max 0, Q u k,t -τ R u k,t + A k,t , ∀k (2)
The data transmitted by the end devices join a remote queue Q m k,t , k = 1, . . . , K, which is then drained by computing the tasks. Therefore, denoting by J k the number of bits computed by each CPU cycle, the queue evolution can be written as:

Q m k,t+1 = max 0, Q m k,t -τ f k,t J k + min(τ R u k,t , Q u k,t ), (3)
where f k,t is the CPU frequency dedicated to device k during slot t, which comes from the (optimized) MEH scheduling.

A. Average E2E delay

The simple queueing system described above is very useful to characterize the first two KPIs of our study: the average sum-rate and the E2E delay of computation offloading, measured as the time elapsed from the generation of new arrivals, until their computation at the MEH. Specifically, let us recall that a queueing system is stationary if all queues are stable, i.e. their long-term average is bounded. In particular, denoting by Q u k , Q m k the long-term average of the involved queues (cf. (1)), the stability of the system is formalized as follows:

Q u k < ∞, ∀k; Q m k < ∞, ∀k. (4)
Moreover, if the system is stable, thanks to Little's law [START_REF] Little | A Proof for the Queuing Formula: L = λW[END_REF], the average E2E delay is finite and can be written as follows:

D k = τ (Q u k + Q m k )/A k , (5)
where A k is the average number of new arrivals A k,t (cf.

(2)), computed as in [START_REF] Strinati | 6G networks: Beyond Shannon towards semantic and goal-oriented communications[END_REF], to be optimized as we will clarify.

B. Average data rate

Since a queueing system is stable if the average departure rate of each queue is greater than the average arrival rate, by designing a policy that stabilizes the system and, at the same time, maximize the average arrival rate, we guarantee the maximum system sum-rate [START_REF] Neely | Fairness and optimal stochastic control for heterogeneous networks[END_REF], [START_REF] Lakshminarayana | Throughput maximization with channel acquisition in energy harvesting systems[END_REF]. Then, in this paper, our aim is to maximize the average sum-rate K k=1 A k , to achieve the maximum sum-rate with finite E2E delay, and with long-term power and EMF exposure constraints as described later on. In particular, intuitively speaking, the idea is to admit, into the uplink queue, the maximum number of bits that the system is able to support (i.e. without instability), by limiting them through a fictitious control valve that chokes the input data rate to match the system capacity [START_REF] Neely | Fairness and optimal stochastic control for heterogeneous networks[END_REF], as depicted in Fig. 1. For this purpose, to achieve stability, a joint resource allocation policy for uplink transmission and MEH computation will also be designed, since, a priori, the system rate can be limited by the uplink, the processing capacity, or the EMF exposure limits to be guaranteed. The challenge is to devise a joint communication and computation resource management policy able to solve this problem without a priori knowledge of channel statistics.

C. System power consumption

In this paper, we consider the devices and the MEH as sources of power consumption, with the aim of keeping their average power expenditure under predefined thresholds. Given the uplink data rate R u k,t a time t, the instantaneous power consumption of user k can be written by inverting the well-known Shannon formula as

p u k,t = N 0 B u k,t h u k,t exp R u k,t ln(2) B u k,t -1 , where B u k,t is the uplink bandwidth assigned to user k, h u k,t
is the time-varying uplink channel power gain, and N 0 is the noise power spectral density at the receiver. For simplicity, we assume a frequency division multiple access, without focusing on spectrum allocation. However, our method can be easily extended to incorporate bandwidth assignment. Finally, the MEH power consumption depends on the CPU core frequency f s,t , which we assume to be dynamically scalable to save power [START_REF] Merluzzi | Dynamic Computation Offloading in Multi-Access Edge Computing via Ultra-Reliable and Low-Latency Communications[END_REF]- [START_REF] Liu | Dynamic Task Offloading and Resource Allocation for Ultra-Reliable Low-Latency Edge Computing[END_REF]. Thus, the power consumption during time slot t can be written as p c,t = κf 3 s,t , where κ is the effective switched capacitance of the processor [START_REF] Burd | Processor design for portable systems[END_REF]. Note that, in this case, the constraint on the CPU should satisfy K k=1 f k,t ≤ f s,t , where f k,t is defined in (3), i.e. the sum overall all users cannot exceed f s,t . Our aim, besides system stability, is to ensure that the long-term average of p u k,t , ∀k and p c,t do not exceed predefined thresholds. We summarize the constraints as follows (cf. (1)):

p u k ≤ p th k , ∀k, p c ≤ p th c , (6)
where p th k and p th c are set a priori, based on the specific needs of devices (e.g. battery) and the MEC operator (e.g. OPEX).

D. The EMF exposure as a Key Value Indicator

As already mentioned, none of the previous works on dynamic computation offloading considered the EMF exposure. Inspired by [START_REF] Chiaraviglio | 5G Network Planning under Service and EMF Constraints: Formulation and Solutions[END_REF], where only communication aspects are tackled, let us consider a generic area populated by machines/sensors and humans. The area is divided in pixels as in the upper part of Fig. 1 and, in each pixel, we aim to guarantee a level of exposure below a predefined threshold. Moreover, the threshold can vary across space, due to specific needs or sensitive areas (e.g. with children), where the exposure has to be kept lower (e.g. over the red areas in Fig. 1). To this aim, following ICNIRP recommendations [START_REF]ICNIRP guidelines for limiting exposure to electromagnetic[END_REF], we use the time-average incident power density, a measure that allows us to write the EMF exposure in closed form, given the channel power gains between devices and pixels. Since we consider the uplink exposure, let us write the total instantaneous power density at time t over pixel i as [21]

P i,tot d,t = K k=1 4π λ 2 p u k,t h i k,t , i ∈ P (7)
with P the set of all pixels, h i k,t is the instantaneous channel power gain between device k and pixel i, and λ is the wavelength. As already mentioned, our aim is to ensure that the long-term average of P i,tot d,t does not exceed a threshold:

P tot,i d ≤ P i,th d , ∀i ∈ P (8)
As an example, ICNIRP recommendation for communications at 2-300 GHz is P i,th d = 10 W/m 2 for the general public, with P i,tot d,t averaged over 30 minutes. However, different values are used in national regulations or in specific conditions. The goal of this paper is to propose a method able to guarantee a predefined constraint P i,th d , rather than focusing on the values of P i,th d over different scenarios. Then, we keep P i,th d generic on purpose to specify it in the numerical results.

III. PROBLEM FORMULATION AND SOLUTION

Having defined and formalized all the KPIs and KVIs, let us now formulate the long-term optimization problem, aimed at maximizing the total average arrival rate, subject to stability constraints of all the queues involved in the offloading procedure, as well as exposure and power constraints. In this way, the long-term maximum sum-rate of the MEC system is attained [START_REF] Neely | Fairness and optimal stochastic control for heterogeneous networks[END_REF], [START_REF] Lakshminarayana | Throughput maximization with channel acquisition in energy harvesting systems[END_REF]. The problem, is formulated as follows:

max {Ψt}t K k=1 A k (9)
subject to (4), ( 6), ( 8), (a) 0

≤ A k,t ≤ A max k , ∀k, t; (b) 0 ≤ p u k,t ≤ p max k , ∀k, t (c) f s,t ∈ F; (d) f k,t ≥ 0, ∀k, t (e) K k=1 f k,t ≤ f s,t , ∀t,
where

Ψ t = [{A k,t } K k=1 , {p u k,t } K k=1 , f s,t , {f k,t } K k=1 ].
Besides the long-term constraints already discussed in Section II, the instantaneous constraints in [START_REF] Oueis | Uplink traffic in future mobile networks: Pulling the alarm[END_REF] have the following meaning: (a) The new arrivals to the uplink queues are non negative and not greater a maximum value, to limit congestion [START_REF] Neely | Fairness and optimal stochastic control for heterogeneous networks[END_REF], [START_REF] Lakshminarayana | Throughput maximization with channel acquisition in energy harvesting systems[END_REF]; (b) The device uplink transmit power is non negative and lower than a maximum value; (c) The MEH CPU cycle frequency is chosen from a discrete set F; (d) The CPU frequency assigned to each user is non negative; (e) The sum of the CPU frequencies assigned to each user does not exceed the selected MEH CPU cycle frequency. In principle, problem [START_REF] Oueis | Uplink traffic in future mobile networks: Pulling the alarm[END_REF] would require the knowledge of the uplink channel statistics to be solve optimally, thus making it very challenging. Moreover, even in case of perfect knowledge of the statistics, the problem is difficult to solve due to its non-convexity and the discrete nature of f s,t over a longterm horizon. However, we propose a simple solution based on Lyapunov stochastic optimization, inspired by [START_REF] Neely | Fairness and optimal stochastic control for heterogeneous networks[END_REF], [START_REF] Lakshminarayana | Throughput maximization with channel acquisition in energy harvesting systems[END_REF], [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF]. We will show how ( 9) is first reduced to a pure stability problem that can be solved in a per-slot fashion, without requiring knowledge of the channel statistics. Furthermore, the per-slot problem is decoupled into radio and computation resource allocation problems that admit closed form solutions or fast iterative algorithms. Based on the theoretical foundations in [START_REF] Neely | Fairness and optimal stochastic control for heterogeneous networks[END_REF], [START_REF] Lakshminarayana | Throughput maximization with channel acquisition in energy harvesting systems[END_REF], [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF], this solution is asymptotically optimal with respect to the original problem (9) (i.e. the joint communication and computation co-design), with a single hyperparameter tuned to explore the trade-off between system sum-rate, E2E delay, power consumption, and EMF exposure. The reason for using Lyapunov tools to solve ( 9) is twofold: i) It comes with theoretical guarantees on system stability, constraint violations, and asymptotic optimality; ii) Being a framework that guarantees long-term averages, it naturally fits the definition and ICNIRP recommendations of EMF exposure, in terms of average incident power density over time, as well as device and MEH power constraints.

A. The Lyapunov approach

Let us now introduce the tools needed to solve the complex problem in [START_REF] Oueis | Uplink traffic in future mobile networks: Pulling the alarm[END_REF], with low complexity. Our first aim is to transform (9) into a pure stability problem. Therefore, we define three virtual queues used to guarantee constraints [START_REF] Mao | Stochastic joint radio and computational resource management for multi-user mobileedge computing systems[END_REF] and [START_REF] Merluzzi | Wireless Edge Machine Learning: Resource Allocation and Trade-Offs[END_REF]. In particular, for each devices k, we define a virtual queue Y k,t evolving as follows (cf. ( 6)):

Y k,t+1 = max 0, Y k,t + y p u k,t -p th k , ∀k, (10)
where y > 0 is a step size. The role of the virtual queue is simple and intuitive: since Y k,t grows if the power constraint is instantaneously violated, and decreases otherwise, thus capturing the behavior of the system in terms of constraint violations. Similarly, for the MEH power consumption (cf. ( 6)), we define a virtual queue H t evolving as

H t+1 = max 0, H t + h p c,t -p th c , (11)
with h > 0. Finally, for constraint (8) we define, for each pixel, a virtual queue Z i,t , i ∈ P evolving as

Z i,t+1 = max(0, Z i,t + z (P i,tot d,t -P i,th d )), (12)
with z > 0. From [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF], we know that the mean rate stability 1 of a virtual queue guarantees to meet the corresponding constraint. For instance, the mean rate stability of Z i,t , ∀i ∈ P would ensure the EMF exposure limit. Therefore, our problem is to guarantee physical and virtual queues' stability (cf. ( 2), ( 3), ( 4), ( 10), ( 11), ( 12)). To this aim, letting

Θ t = [{Q u k,t } k , {Q m k,t } k , {Y k,t } k , H t , {Z i,t } i ],
we define the Lyapunov function as follows [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF]:

L(Θ t ) = 1 2 K k=1 [(Q u k,t ) 2 + (Q m k,t ) 2 + Y 2 k,t ] + H 2 t + i∈P Z 2 i,t , (13)
which is a measure of the overall congestion state of the system in terms of both physical and virtual queues. Based on (13), we can define the drift-plus-penalty (DPP) function:

∆ p,t = E L(Θ t+1 ) -L(Θ t ) -V K k=1 A k,t Θ t . ( 14
)
The DPP in ( 14) represents the conditional expected change of the Lyapunov function over one slot, with a penalty factor, weighted with parameter V , assigning more or less priority to the objective function of ( 9) with respect to queue backlogs. Interestingly, if the DPP is upper bounded by a finite constant for each time slot, the stability of the physical and virtual queues is guaranteed, thus guaranteeing constraints ( 6), ( 8), and finite delay (cf. ( 5)). Then, to push the network towards low congestion states, we proceed by minimizing a suitable upper bound of ( 14) [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF], which in this case reads as

∆ p,t ≤ζ + E k k=1 Q u k,t A k,t -τ R u k,t + Q m k,t τ R u k,t -τ f k,t J k + y Y k,t (p u k,t -p th k ) + i∈P z Z i,t P i,tot d,t -P i,th d + h H t (p c,t -p th c ) -V K k=1 A k,t Θ t , (15)
where ζ is a positive constant omitted due to the lack of space, as well as the derivations leading to [START_REF] De Domenico | EMF-Aware Cell Selection in Heterogeneous Cellular Networks[END_REF], which directly follow from [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF], using the following two inequalities: For a generic physical queue evolving as

Q t+1 = max(0, Q t - b) + A, we can write Q 2 t+1 -Q 2 t ≤ A 2 + b 2 + 2Q(A -b); for a generic virtual queue G t+1 = max(0, G t + g t -g), we can write G 2 t+1 -G 2 t ≤ (g max t -g) 2 + 2G t (g t -g
), where g max t is the instantaneous finite maximum value attainable by g t , bounded by hypothesis. In particular, instantaneous power consumption and EMF exposure are clearly bounded by a finite constant. As in [START_REF] Neely | Fairness and optimal stochastic control for heterogeneous networks[END_REF], [START_REF] Lakshminarayana | Throughput maximization with channel acquisition in energy harvesting systems[END_REF], [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF] and based on stochastic optimization arguments, minimizing [START_REF] De Domenico | EMF-Aware Cell Selection in Heterogeneous Cellular Networks[END_REF] in each time slot theoretically leads to system stability under the assumption of i.i.d channel realization over consecutive time slots. Moreover, the optimal solution of (9) (i.e. the maximum system sum-rate) is asymptotically achieved as the Lyapunov parameter V increases, with the cost of increased average delay and convergence time [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF]. In particular, the strategy that minimizes [START_REF] De Domenico | EMF-Aware Cell Selection in Heterogeneous Cellular Networks[END_REF] in each time slot achieve a [O(1/V ), O(V )] performance-delay trade-off, i.e. the distance between the value obtained by the online strategy and the optimal value of ( 9) decreases as O(1/V ), with a cost of O(V ) increased upper bound on queue backlogs (i.e. delay and convergence time). More details follow from the theoretical results in [START_REF] Neely | Fairness and optimal stochastic control for heterogeneous networks[END_REF], [START_REF] Lakshminarayana | Throughput maximization with channel acquisition in energy harvesting systems[END_REF], [START_REF] Neely | Stochastic Network Optimization with Application to Communication and Queueing Systems[END_REF], and are omitted due to the lack of space. Now, let us notice that, once the problem is cast in a per slot fashion, the flow control (i.e. data arrivals), uplink radio resource allocation, and CPU scheduling at the MEH, are decoupled into three independent sub-problems, which we present in the following, along with their respective solutions. Due to the lack of space, the derivations of the solutions will be summarized, but mainly come from simple mathematical manipulations involving classical convex optimization arguments that we will briefly recall [START_REF] Boyd | Convex Optimization[END_REF]. Note that the solutions of the following problems require the instantaneous knowledge of wireless channels, as well as physical and virtual queues' states. We assume them to be centrally solved by the MEH, although distributed implementations that go beyond the scope of this paper can be envisioned.

B. Flow control sub-problem

The flow control problem involves the data arrivals (i.e. A k,t , ∀k) and is formulated as follows at time slot t:

min {A k,t } ∀k K k=1 Q u k,t -V A k,t s.t. 0 ≤ A k,t ≤ A max k , (16)
Problem ( 16) is linear and admits the straightforward simple closed-form solution

A * k,t = A max k • 1{Q u k,t ≤ V },

C. Uplink Radio Resource Allocation sub-problem

The uplink resource allocation problem aims at optimizing the transmit power of each device, and is formulated as

min {p u k,t } ∀k - K k=1 Q u k,t -Q m k,t τ B u k,t log 2 1 + h u k,t p u k,t N 0 B u k,t + i∈P z Z i,t K k=1 4π λ 2 p u k,t h i k,t + K k=1 y Y k,t p u k,t subject to 0 ≤ p u k,t ≤ p max k , ∀k; (17)
Assuming the bandwidth to be assigned, the above problem can be decoupled among different users. Once decoupled, for a generic user k let us notice that, if [START_REF] Little | A Proof for the Queuing Formula: L = λW[END_REF] are monotone non-decreasing functions of p u k,t . Therefore, the optimal solution for the devices having

Q u k,t ≤ Q m k,t , all terms in
Q u k,t ≤ Q m k,t is p u, * k,t = 0.
For all other devices, the problem is strictly convex. Therefore, the optimal solution can be found in closed-form by solving the Karush-Kuhn-Tucker (KKT) conditions [START_REF] Boyd | Convex Optimization[END_REF] (omitted due to the lack of space), and is:

p u, * k,t = σ u k,t - N 0 B u k,t h u k,t p max k 0 , ∀k : Q u k,t > Q m k,t , where σ u k,t = τ Q u k,t -Q m k,t B u k,t ln(2) i∈P 4π λ 2 z Z i,t h i k,t + y Y k,t
.

Besides its theoretical mathematical derivation, this solution is also intuitive: if the local queue of device k grows compared to its computation queue, device k increases the transmit power to drain it; at the same time, if the virtual queues Z i,t and/or Y k,t grow, the transmit power is reduced to meet the long-term power and exposure constraints. Finally, the transmit power increases as h u k,t increases, to exploit good radio channel conditions when transmitting to efficiently Algorithm 1 Energy efficient EMF-aware comp. offloading 2), ( 3), ( 10), [START_REF] Zhao | EMF Exposure Study Concerning mmWave Phased Array in Mobile Devices for 5G Communication[END_REF], [START_REF] Thors | Exposure to rf emf from array antennas in 5g mobile communication equipment[END_REF]. drain local queues. Let us notice that, even though for simplicity we assumed the bandwidth to be assigned a priori, for all devices transmitting at time t, the objective function of ( 17) is convex also on B u k,t , making it suitable to optimize the spectrum allocation under the same optimization framework.

Input: {Y k,0 } K k=1 , {Z i,0 } i , H 0 , {Q u k,0 } K k=1 , {Q m k,0 } K k=1 . In each time slot, observe Q u k,t , Q m k,t , Z i,t , Y u k,t ,

D. Computation Resource Allocation sub-problem

For the CPU scheduling at the MEH, from [START_REF] De Domenico | EMF-Aware Cell Selection in Heterogeneous Cellular Networks[END_REF], the subproblem is cast as follows:

min fs,t,{f k,t } k h H t κf 3 s,t - K k=1 Q m k,t τ f k,t J k (18) subject to (a) f s,t ∈ F; (b) f k,t ≤ min f s,t , Q m k,t J k τ (c) f k,t ≥ 0, ∀k, t (d) K k=1 f k,t ≤ f s,t , ∀t,
where, for efficiency purposes, we added constraint (b) to ensure that each user cannot be assigned a clock frequency higher than the one necessary to empty the queue. Problem ( 18) is a mixed integer non-linear program. However, For a fixed f s t ∈ F, it is linear in {f k,t } K k=1 , and the optimal frequencies can be found using a simple procedure (requiring at most K steps) that iteratively assigns largest portions of f s t to the devices with largest values of Q m k,t , as also proved in [START_REF] Merluzzi | Dynamic Computation Offloading in Multi-Access Edge Computing via Ultra-Reliable and Low-Latency Communications[END_REF], [START_REF] Merluzzi | Wireless Edge Machine Learning: Resource Allocation and Trade-Offs[END_REF]. In the case a device's queue is completely drained, the remaining devices are allocated the available frequency until no frequency is left or all devices are served. Then, letting {f * k,t (f s,t )} be the optimal frequencies assigned to the devices for each f s t ∈ F, the optimal f s t is

f * s,t = arg min fs,t∈F - K k=1 Q m k,t τ J k f * k,t (f s,t ) + h H t κ(f s,t ) 3 ,
which involves |F| evaluations of the objective function in [START_REF] Neely | Fairness and optimal stochastic control for heterogeneous networks[END_REF]. Finally, the overall procedure for Energy-efficient EMFaware computation offloading is summarized in Algorithm 1.

IV. NUMERICAL RESULTS

For the numerical results, we consider an AP, located at the center of a squared area of side 15 m, operating at f 0 = 3.5 GHz, with a total available bandwidth of 20 MHz, equally shared among all devices. Around the AP, in the same area, we consider 100 sensors uniformly randomly distributed in space, and transmitting with a maximum power p max k = 20 dBm. Also, we assume a noise power spectral density N 0 = -174 dBm/Hz. The path-loss model is taken from [START_REF] Gpp | Study on channel model for frequencies from 0.5 to 100 GHz[END_REF], considering a dense factory plant. A Rayleigh fading with unit variance and coherence time equal to the slot duration τ = 5 ms is assumed, while A max k = τ R u,max k (cf. constraint (a) of ( 9)). At the MEH side, the CPU clock frequency can be selected from the discrete set F = {0, 0.1, 0.2, . . . , 1}×f max s , with f max s = 4.5 GHz. The effective switched capacitance of the processor is set to κ = 10 -27 . We use conversion factor J k = 10 -1 , ∀k (cf. ( 3)). The area around the AP is divided into 1 m 2 pixels. All results are averaged for 2000 slots and over 100 different realizations of devices' positions. As numerical results, we show the trade-off between sumrate, EMF exposure, power consumption, and the E2E delay of the offloading service. We present three cases: i) Unconstrained: No limits on power and EMF are imposed, i.e. no virtual queues are used (cf. ( 12), ( 10), [START_REF] Zhao | EMF Exposure Study Concerning mmWave Phased Array in Mobile Devices for 5G Communication[END_REF]), but only physical queues (cf. ( 2), (3)) for system stability. Therefore, this is the solution achieving the highest possible sum-rate, and represents our benchmark for the second setting; ii) Constrained: Constraints on EMF exposure over the space and power consumption (including sensors and MEH) are imposed, thanks to physical and virtual queues stability; iii) Constrained (no MEH power constraint): Constraints on EMF exposure over the space and power consumption (including only sensors) are imposed. Let us show the results in Fig. 2, first with focus on the unconstrained setting, represented by the blue curve in all plots. In particular, Fig. 2a shows the EMF exposure as a function of the average sum-rate, obtained by increasing the Lyapunov parameter V (cf. ( 14)) from left to right. As it is also clear from Fig. 2c, we tuned V to achieve similar E2E delay in the three cases, with a maximum value around 50 ms. We can notice how, by increasing the average sum-rate, the EMF exposure increases as expected in the unconstrained setting, up to a maximum value around 90 mW/m 2 , corresponding to the maximum average sum-rate (around 450 Mb/s). At the same time, in Fig. 2b, we show the average device and MEH power consumption as a function of the average sum-rate. Again, by increasing the sum-rate, the average power consumption of both devices and MEH show an increasing behaviour towards its maximum value (100 mW and 90 W, respectively), corresponding again to the maximum average sum-rate on the abscissa. Finally, in Fig. 2c, we show the average E2E delay (cf. ( 5)) as a function of the average sum-rate, showing its monotone increasing behaviour as expected by the theory ([O(1/V ),(V )] trade-off). Now that the unconstrained results have been described, let us focus on the constrained case, represented by the red curve in all plots. In particular, we impose: i) An EMF exposure threshold P i,th d = 40 mW/m 2 , represented by the black dashed horizontal line in Fig. 2a; ii) An MEH average power threshold p th c = 45 W, and a device average power threshold p th k = 50 mW, represented by the black dashed lines in Fig. 2b. First of all, we can notice how the online method is able to guarantee all the required constraints, in terms of exposure and power consumption, as also guaranteed by the theoretical arguments in Section III-A. Moreover, EMF and sensor power consumption values are strictly lower than the thresholds, while the MEH power consumption approaches the desired threshold. This result suggests that the system sum-rate (around 350 Mbps in this case) is limited by the MEH in this setting. To validate this conclusion, we run the same simulation, removing the MEH power constraint, thus achieving the result represented by the green curve in all plots. As expected, since the system is limited by the MEH capacity in the previous case, by removing this constraint, the method is able to guarantee sensor power consumption and EMF constraints, without degrading the performance in terms of sum-rate and E2E delay, as visible in Fig. 2c. Therefore, we can conclude that our method is able to strike the best trade-off between average sum-rate, power consumption, EMF exposure and delay. Also, depending on the limiting capacity (computation or wireless), the method can dramatically reduce individual power consumption and/or EMF exposure, with negligible impact on network performance. More results involving different combinations of constraints are omitted due to the lack of space, but are foreseen for future investigations.

V. CONCLUSIONS In this paper, we proposed an online algorithm to jointly optimize radio and computing resources for computation offloading purposes, to maximize the average sum-rate of offloaded data under EMF exposure and network power constraints, involving multiple end devices and an MEH. Our method does not require any prior knowledge of wireless channel statistics, and comes with theoretical guarantees on stability, exposure and power constraints, as well as on the asymptotic optimality of the solution in terms of average sumrate. Besides the theoretical analysis, numerical results show the effectiveness of our method in striking the best trade-off between sum-rate, power consumption, EMF exposure, and service E2E delay. Future work includes multi-cell investigations with interference, directive communications, downlink communications and AP power consumption.

Fig. 1 :

 1 Fig. 1: Scenario and arrivals II. SYSTEM MODEL AND KEY PERFORMANCE/VALUES

  ∀k where 1{•} denotes the indicator function.

Fig. 2 :

 2 Fig. 2: Trade-off between average sum-rate, EMF exposure, network power consumption and E2E delay

  ∀k, i, H t : S1. Observe h u k,t , ∀k, and find A * k,t and p u, * k,t as in Sections III-B and III-C, respectively; S2. Find f * s (t) and {f * k (t)} k as in Section III-D; S3. Update all queues as in (

For a virtual queue Gt, it is defined as lim T →∞ E{G T }/T = 0