Combining Time-Flexible GEO Satellite Payload With Precoding: The Cluster Hopping Approach

High throughput geostationary orbit (GEO) satellite systems are characterized by a multi-beam wide coverage. However, developing efficient resource management mechanisms to meet the heterogeneous user traffic demands remains an open challenge for satellite operators. Furthermore, the spectrum shortage and the ever-increasing demands claim for more aggressive frequency reuse. In this article, we combine the time-flexible payload capabilities known as beam hopping (BH) with precoding techniques in order to satisfy user traffic requests in areas of high demand (i.e. hot-spot areas). The proposed framework considers a flexible beam-cluster hopping where adjacent beams can be activated if needed, forming clusters with various shapes and sizes. In this context, we present three strategies to design the beam illumination patterns. First, a max-min user demand fairness satisfaction problem; second, a penalty-based optimization is considered to penalize the occurrence of adjacent beams in an attempt to avoid precoding whenever possible. Third, seeking a low-complexity design, we propose a queuing-based approach to solve the problem in a time-slot-by-time-slot basis trying to provide service to users based on the requested demands. The three methods are discussed in detail and evaluated via numerical simulations, confirming their effectiveness versus benchmark schemes and identifying the pros and cons of each proposed design.

Combining Time-Flexible GEO Satellite Payload With Precoding: The Cluster Hopping Approach Vaibhav Kumar Gupta , Member, IEEE, Vu Nguyen Ha , Member, IEEE, Eva Lagunas , Senior Member, IEEE, Hayder Al-Hraishawi , Senior Member, IEEE, Lin Chen , Graduate Student Member, IEEE, and Symeon Chatzinotas , Fellow, IEEE Abstract-High throughput geostationary orbit (GEO) satellite systems are characterized by a multi-beam wide coverage.However, developing efficient resource management mechanisms to meet the heterogeneous user traffic demands remains an open challenge for satellite operators.Furthermore, the spectrum shortage and the ever-increasing demands claim for more aggressive frequency reuse.In this article, we combine the time-flexible payload capabilities known as beam hopping (BH) with precoding techniques in order to satisfy user traffic requests in areas of high demand (i.e.hot-spot areas).The proposed framework considers a flexible beam-cluster hopping where adjacent beams can be activated if needed, forming clusters with various shapes and sizes.In this context, we present three strategies to design the beam illumination patterns.First, a max-min user demand fairness satisfaction problem; second, a penalty-based optimization is considered to penalize the occurrence of adjacent beams in an attempt to avoid precoding whenever possible.Third, seeking a low-complexity design, we propose a queuing-based approach to solve the problem in a time-slot-by-time-slot basis trying to provide service to users based on the requested demands.The three methods are discussed in detail and evaluated via numerical simulations, confirming their effectiveness versus benchmark schemes and identifying the pros and cons of each proposed design.Index Terms-Very high throughput satellite system, beam hopping, cluster hopping, precoding, illumination pattern, demand satisfaction.

I. INTRODUCTION
T HE continuous and significant growth in broadband sub- scribers and bandwidth-hungry applications have shown the integration of terrestrial networks with the satellite broadband systems very appealing [1].The main focus of 5G and beyond systems is to ensure uninterrupted and high-speed connectivity across the globe to in-motion terminals on aircraft, ships, and vehicles.However, due to the fact that high-speed and broadcast-enabled satellite links will supplement existing terrestrial connectivity, this goal cannot be fulfilled solely by terrestrial networks [2].Indeed, the satellite resources are exceptionally scant and very expensive, and thus, these resources must be utilized effectively [3].Further, satellite traffic changes with time as well as with coverage areas [4].Hence, equal resource allocation across the whole coverage area all the time results in inefficient resource utilization, e.g., satellite power [5].
Onboard resource management in conventional satellite communication systems employs constrained pre-design methods [6].On one hand, this can simplify the satellite payload design.However, such limited flexibility in practice barely accounts for heterogeneous traffic and changing demands.Mismatches between supplied capacity and requested traffic in practical situations lead to inefficient resource utilization.For next-generation very high throughput satellites (VHTS), it is crucial to develop enhanced resource management techniques that take advantage of multi-dimension flexibilities to solve this prevalent problem [7].
The space industry is on the verge of a technological breakthrough to extend satellite system performance beyond current levels [8].Many commercial players are working hard to launch the flexible VHTS system that is able to offer up to 10 3 Gbps data rates even to the remote and unserved areas [9].These high data rates can only be achieved by evolving the multi-beam satellite architecture which uses a smaller beam size.In this context, we are currently going through a very interesting and exciting period where the main goal is two-fold: (i) to drastically reduce the operational cost of the satellite system by increasing the satellite capacity while improving the system spectral efficiency, and (ii) to deliver satellite services in a more flexible manner by allocating the satellite resources based on the per-beam traffic demands.
Beam Hopping (BH) and linear precoding are two disruptive technologies developed to achieve the above-mentioned two goals for the satellite systems [10], [11], [12], [13], [14].While BH enables flexibility to adapt the offered capacity to the heterogeneous demand, precoding strives to improve the spectral efficiency [13], [14].Conventional BH exploits the time-domain flexibility to provide all available spectrum to a selected set of beams as long as they are not adjacent to each other [15].With conventional BH, all the available satellite resources are employed to provide service to a certain subset of beams, which is active till the requested demand in each beam is satisfied.In conventional BH illumination patterns, the active spot beams are designed such that they do not interfere to each other.Hence, in conventional BH, there are strict constraints in the possible illumination patterns due to interference avoidance and uneven traffic distribution [13].
The introduction of precoding in BH systems is required to deal with co-channel interference when a set of adjacent beams is illuminated at the same time and operating at the same frequency band to satisfy high beam demand in hot-spot scenario [14].Such a set of adjacent beams is known as cluster and this technique is known as "Cluster Hopping" [16].The design of the optimal illumination pattern, i.e. deciding which clusters to illuminate together and for how long, is the key design challenge of Cluster Hopping (CH).A specific arrangement of all the clusters is called a snapshot.The main advantage of the proposed cluster-based hopping is the resource allocation flexibility within a cluster, meaning that when a group of adjacent beams is simultaneously illuminated, one can implement a beam-free approach where all the resources of the cluster can be shared across the cluster's coverage area.However, the performance of CH strongly depends on the design of the cluster size and shape.In the case of overlapping clusters, there will be a very large number of possible clusters resulting in a huge search space for the proposed problem and rendering a high-complexity solution.Furthermore, it is difficult to regulate the supplied capacity.The shape of the clusters has an impact on the precoding performance.Furthermore, compact clusters are preferred as they will allow sharing resources from neighbouring beams.Finally, the channel state information (CSI) estimation procedure (which is fundamental for the precoding part) needs to be adapted to the time-varying nature of the CH technique.

A. Related Works
We now briefly present the existing literature on BH and precoding techniques for multi-beam satellite systems.In the conventional BH technique, all the available resources are allocated to a subset of beams which are illuminated in each time slot to serve the users in such a way that any two activated beams are far apart so that there is no interference.The performance of BH is hampered by the limitation that adjacent beams can not be activated simultaneously using the same spectrum resource [13].The main task in BH technique is to efficiently select a subset of beams to be illuminated in each time slot which can be accomplished by considering a binary variable for each beam and setting it to one if the beam is illuminated.This helps to formulate the problem into a mixed integer non-linear programming problem (MINLP) which is very hard to solve.A supplied capacity maximization problem in a multi-beam satellite system was considered in [17] and iterative sub-optimal heuristic algorithms were proposed to solve the problem.Similarly, an efficient beam capacity utilization problem via beam hopping technique was considered in [18], and a solution was proposed based on optimization and genetic algorithms.
Furthermore, a resource allocation problem in a beamhopping satellite system in which weights are assigned to the users to categorize them was studied in [19].An optimization framework was formulated to maximize the weights of the served users and a Cuckoo search-based algorithm was proposed.An optimal capacity allocation problem in gateway satellite systems where gateways offer the capacity to satisfy beam demands was investigated in [20].A fair scheme based on Monge arrays to maximize the minimum demand satisfaction ratio was proposed.A resource allocation problem in low Earth orbit (LEO) beam-hopping satellite systems that share the spectrum of the geostationary orbit (GEO) satellite was considered in [21].The problem was decomposed into three sub-problems relating to the optimal selection of frequency, cell, and transmit power for LEO satellite beams.In [22], a resource allocation problem was formulated with the objective of minimizing the total number of carriers and power allocated to each satellite beam.The problem was decomposed into two sub-problems; first to select the optimal number of carriers required for each beam and then power allocation for selected carriers to meet the beam demand.
The fact that the binary selection variable for beam activation results in an outsized search space is undoubtedly the main issue in early investigations.When there is a large number of beams, the issue gets exceptionally worse.Therefore, Machine Learning (ML) tools are being widely used to solve complex optimization problems.Authors in [23] proposed a deep learning-based procedure, to identify the optimal beam illumination pattern for satellite systems, which results in low computation complexity.In [24], a multi-agents Deep Reinforcement Learning (DRL) based method was proposed to optimally solve the joint beam illumination pattern identification and bandwidth allocation problem in beam hopping satellite systems.Similarly, a double loop learning-based multi-action selection policy was proposed in [25] for beam hopping in DVB-S2X broadband satellite system with multiple purposes.A DRL-based genetic algorithm was proposed in [26] to optimize the beam hopping decision for the DVB-S2X GEO satellite system.
The fundamental disadvantage of conventional BH is the need for spatial separation across active beams, which restricts its capacity to adjust to any demand distribution.It is evident that in certain scenarios where more than one adjacent beam has high demand, i.e., hot-spot scenario, it is required to simultaneously activate a group of more than one adjacent beam, i.e., a cluster of beams.This was the driving force behind the exploratory research conducted in [27], where precoding was initially integrated with BH to reduce the inter-beam interference caused by adjacent beams that are illuminated at the same time.This framework, in which a subset of adjacent beams forming a cluster is activated at a time, was first analyzed by the European Space Agency (ESA) [28].Based on this study, the concept of Cluster Hopping (CH) in which the overall multibeam coverage area is divided into a number of mutually exclusive beam clusters that are illuminated for just long enough to satisfy the demand in each cluster was first proposed in [16].Precoding is employed within a cluster to mitigate inter-beam interference, while inactive clusters are physically positioned between active clusters to reduce inter-cluster interference.
Additionally, an efficient precoded cluster hopping-based illumination pattern design method was proposed in [16], [29] to maximize the minimum supply-capacity-demand ratio.However, fixed shape and size clusters were assumed to reduce the search space and so as complexity.Although, employing precoding techniques within clusters enhance the computational complexity of the system.Therefore, the authors in [30] have formulated an interference penalty minimization framework with demand satisfaction constraints and proposed a low complexity dynamic beam illumination mechanism with selective precoding used only for the active beams which create strong inter-beam interference to each other.Authors in [31] investigated different kinds of precoding for high throughput satellite hot spot scenarios.However, the 7-beam cluster was considered in which a high-demand beam is surrounded by six zero-demand beams, and an average throughput maximization problem was formulated with the assumption that interference at the receiver is treated as noise only.

B. Contributions
Satellite data traffic has a spatiotemporal nature that is highly heterogeneous and variable.Indeed, the traffic demand for each beam may vary significantly due to geographical and temporal differences, especially when considering mobile user terminals.The fixed cluster size limits the capability of adapting to the highly heterogeneous traffic demand distributions.In particular, fixed cluster size leads to under utilization of the satellite resources where demand is comparatively low (e.g., in the middle of an ocean) and poor demand satisfaction in hot-spot scenarios (e.g., surrounding of an international airport with multiple high demand mobile users) where demand is very high.Therefore, to reduce the complexity, compact-shaped, non-overlapping, and equal size clusters result in sub-optimal resource utilization [29].In fact, cluster size and shape should be as close as possible to the expected demand clusters in order to reach the best demand-matching performance.Furthermore, illumination of flexible sized clusters was considered in [30] to minimize the inter-beam interference penalty but it does not guarantee the efficient utilization of satellite resources to best satisfy the beam demands.This motivates us to combine the promising frameworks developed in [29], [30].Hence, in this article, we design a more flexible CH illumination pattern allowing clusters of different shapes and sizes and employing precoding within a cluster to improve the resource utilization and supply-demand matching as compared to that was achieved with the pre-defined clustering used in [29].
The main contributions of this article are summarized as follows: r A max-min fair satellite resource allocation framework for per beam supply-demand matching problem is investigated.We propose a solution called flexible cluster hopping (FCH) in which the design goal has been extended from illuminating a few clusters of fixed shape and size to jointly illuminating a number of beams satisfying the beam demand constraints and clustering the illuminated beams for precoded transmission to mitigate the inter-beam interference.
r We relax the formulated problem by simplifying the non-linear constraints and reformulate the non-linear and non-convex optimization problem to make it tractable.
r Formulate the beam hopping problem as a penalty mini- mization problem while ensuring beam demand satisfaction.The problem is converted into a convex optimization problem using Semi-Definite Programming (SDP) relaxation and Successive Convex Approximation (SCA) based methods.
r Propose a low complexity queuing-based method that is a throughput optimal beam hopping technique with the objective of minimizing the remaining demand of each beam in each time slot.
r The performance of the proposed techniques is investigated herein in terms of the beam supply-demand matching profile, system efficiency, and other key performance indicators (KPIs).Simulation results including performance comparisons are provided to demonstrate the validity and gains of the proposed method over the benchmark CH scheme.The remainder of this article is organized as follows.Section II introduces the system, channel and precoding models.A maxmin fair beam supply-demand matching problem is formulated and a flexible cluster hopping scheme is proposed in Section III.The problem is reformulated as a penalty minimization framework and two solution approaches namely: SDP based and SCA based are proposed to convexify the problem in Section IV.Next, we reformulate the problem as maximization of the sum weighted throughput of all the illuminated beams in each time slot and two queuing-based methods are proposed in Section V.The complexity analysis of the proposed schemes is provided in Section VI.Section VII reports on detailed numerical results and a concluding section ends the article.

II. SYSTEM MODEL
We consider a downlink scenario of a high-throughput multibeam GEO satellite system where the satellite can generate N beams for serving ground users within a time window T H , which is also known as a hopping window.The hopping window is segmented in T equal time slots of length Δ T , which is known as the dwelling time of the hopping system.Herein, an illumination pattern is created over a particular hopping window and is then repeatedly used over time.The users in every beam should be served at their data demand.This demand is represented as an amount of data bits that should be transmitted within the hopping window T H .In this work, beam-level frameworks are considered; hence, one user per beam scenario is assumed. 1 In such scenario, the data demand raising within beam n ∈ {1, . .., N} can be represented by the aggregated data transmission demand of only one virtual user in this beam, which is denoted as D n 2 (bits to be transmitted in T times slots).The scenarios where there is more than one user located within the coverage of each beam will be considered in future works dealing with user-level strategies.Regarding the hardware limitation or practical implementation requirements, we also assume that the number of illuminated beams in time slot t cannot be larger than K3 (K ≤ N ).Also, all the illuminated beams operate over the same spectrum B employing full frequency reuse.It is assumed that multiple users are multiplexed in a time division multiplexing (TDMA) manner.The group of adjacent and simultaneously illuminated beams is called a cluster and each permuted arrangement of the illuminated clusters is called "a snapshot".The size of the clusters has an impact on the complexity and flexibility of the CH operation.When the clusters are small, it approaches the conventional BH solution where precoding does not provide any gain.On the other hand, if the size of the clusters is very high, we will have much lower flexibility in terms of supplying capacity, as we will have very less snapshots.Therefore, in this work, we consider efficient and flexible cluster size such that each beam demand must satisfy.Additionally, it is assumed that clusters are non-overlapping and compact in shape.

A. Channel Model
Suppose y n [t] denotes the received signal at user n in time slot t and the corresponding received signal vector at the users in cluster i in time slot t is represented by , where c i is the size of cluster i at time slot t.Therefore, the received signal vector at all the users in time slot t, i.e., T can be represented as follows: represents the channel matrix, which is assumed to be perfectly known at the transmitter, W [t] ∈ C N ×N represents the corresponding precoding factors for all the beams, w l is the precoding factor designed for the cluster l, i.e., the lth column of W [t] which corresponds to the precoding factor of the beams in cluster l, L t denotes the number of clusters formed in time slot t, and N ∈ C N ×1 denotes the zero-mean additive Gaussian noise.In particular, we assume T I N , where σ T = √ κT Rx B; κ denotes the Boltzmann constant and T Rx is the clear sky noise temperature of the receiver.
The downlink channel gain coefficient from lth satellite beam to kth users is modeled as, where G k R is the receiver antenna gain at the kth user, G k,l is the gain from lth beam to the kth user, φ k,l is a phase component resulting from the beam-pattern, d k is the distance of the kth user from the satellite, and λ denotes the wavelength.

B. Precoding Model
This work mainly focuses on developing flexible beamhopping techniques where satellite beams are clustered and illuminated in clusters of different sizes to meet user demands and improve system efficiency.When a cluster has two or more illuminated beams simultaneously, precoding is used to mitigate inter-beam interference.The precoding factor design tasks are lessened by simply employing some benchmark techniques such as MMSE-based precoding [14].It is noteworthy that precoding is performed at the gateway to ensure reduced complexity of the satellite's on-board processing effort and user equipment.In particular, for the cluster i having at least two or more beams, the precoding factor is defined as follows: where H i represents the channel matrix of all the users in cluster i and α stands for a predefined regularization factor.Then, w i is determined based on ŵi by normalizing every column vector of the matrix to meet the per-beam transmit power constraints P n .For all the beams which are not activated in a time slot, the precoding factor must be zero.Similarly, when there is only a single beam in a cluster, its precoding factor is simply √ P n .We tend to use dynamic beam scheduling and selective precoding to match demands with supplied capacity.This extended design goal considering these constraints can be stated in the following general problem. s.t. where ) is the objective function which will be defined thoroughly in certain scenarios with specified objective designing goals.For instance, it can be the achievable rate, the transmission power, the fairness ratio, the transmission time, the operation cost, or system revenue.Herein, the first constraint indicates that the number of illuminated beams at each time slot cannot be more than K due to payload architecture limitations.In the second constraint, Δ T is the time duration of one time-slot, and f (.) is a rate function in terms of Signal-to-Interference-plus-Noise Ratio (SINR).SINR experienced by beam n, i.e., SINR n (z [t], W [t]) is expressed as in (6) shown at the bottom of this page, where C n (z [t]) denotes the cluster consisting of beam n which can be figured out based on z [t].With clustering transmission, the precoding can exploit the spatial dimension of the channels from beams and users within clusters to maximize the received power or to mitigate inter-beam interference.Without a precoding design, the clustering problem is still challenging.This can be illustrated by the number of feasible snapshots turning out to be exceptionally large.Now, selecting K beams and grouping them into different clusters is equivalent to choosing K elements of z[t] and setting them to ones while remaining others to zeros.Note that the clustering solution can be permuted without affecting the results.This yields a total of ( (NK)! K!(NK−K)!)/K! possible snapshots.For example, assume N = 67 and K = 17 beams can be simultaneously activated in a time slot.The total number of reliable snapshots is 15 which is an extremely substantial number just for only one time slot.This number will increase exponentially if the whole window with T time slots is considered.Therefore, efficient and flexible clustering algorithms are needed which will be presented in the following sections.

III. FAIR SUPPLY DEMAND MATCHING FRAMEWORK USING FLEXIBLE CLUSTER HOPPING
The first natural step is to find a solution where adjacent beams are simultaneously activated (and precoded) only if needed.First, we propose a flexible cluster hopping (FCH) technique, which is a simple extension of the CH solution used in [29].Then, we start with the formulation of a fair demand-satisfaction problem making use of the max-min of the ratios (supplied capacity/demand).The achievable data transmission of beam n during all T time slots can be expressed as Next, the normalized amount of achievable data transmission of beam n can be calculated as the ratio between the amount of data transmission and the data demand as R n D n .Then, the extension of the CH design problem is stated as in problem P 0 : where we have removed the constraint that ensured the activation of non-adjacent beams.Note that, L = I T ×T ⊗ 1 1×N and the maximum number of simultaneously active beams in each time slot instance over the hopping window is equal, i.e., k = [K, . .., K] T and x = [x 1 [1], . . ., x N [1], . . ., x 1 [T ], . . ., x N [T ]] T in which the variables x n [t] 's are the binary variable which are used to define the beam illumination solution, i.e., x While in the existing CH technique, the potential snapshot candidates are limited to the snapshots that do not contain adjacent active beams or fixed size clusters, herein we do not impose such constraint thus making the algorithm more "flexible" in selecting active beams and forming clusters whenever needed.Obviously, to make sure this works, we need to add a final step where adjacent beams are precoded together to mitigate inter-beam interference.

A. Proposed Solution
The formulated problem P 0 is a mixed integer, non-convex, and non-linear programming (MINLP) mainly due to the presence of a binary variable (x) and the interference term in the SINR expression.By introducing the variable η, this problem is equivalent to the following one max x,η η (P 1 ) In other words, problem P 0 can be solved by defining the maximum value of η such that the following problem is feasible.
Both the P 1 and P 2 problems are still MINLPs.As a result, no globally optimal solution to these problems can be obtained using the well-known convex optimization techniques [32].
Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
Regarding the interference, the interference coming from adjacent beams can be mitigated by the use of precoding, while the interference coming from non-adjacent beams is generally low and can be assumed negligible.Therefore, the first simplification step is to re-formulate P 1 4 using Signal-to-Noise Ratio (SNR) instead of SINR, where, γ n [t] is the SNR of beam n in time slot t.We further simplify P1 by considering the demand vector d = [D 1 , . . ., D N ] T and writing the offered capacity vector as follows where ] T , denotes the element-wise product, and R = [r 1 , . . ., r N ] T , is a vector that contains the interference-free offered data rate per beam (which is independent of the time-instance as it ignores the interference term), i.e, Note that P2 is equivalent to P 0 iff interference is negligible.Hence, the solution of P2 is optimal only when SINR ≈ SNR.
Although P2 still contains a binary variable, some state-of-theart solvers, e.g., MOSEK,5 can be applied to obtain a solution.

IV. PENALTY MINIMIZATION FRAMEWORK USING FLEXIBLE CLUSTER HOPPING
Herein we present a different approach where the objective is generalized to a minimization of a given penalty function subject to demand constraint.This approach is cast by the following problem: where Ω ∈ R (NT×NT) stands for the penalty matrix whose element represents a penalty cost when the corresponding two beams are activated in a time slot t and other notation have a similar interpretation as in Section III.In problem P 3 , again f (.) stands for the rate-SINR function which can be the Shanon capacity function or the rate function used in DVB standard, and , where p n [t] is the transmission power of beam n and h j n [t] stands for the channel gain from beam j to the user n.Since we consider the beam-level design, the SINR are estimated based on the channels from the satellite beams to the users located at the center of the beams.In the user-level design, the achievable data rate of one beam can be defined as the total rate achieved by all users located inside that beam's coverage area.Once the channel gain and the transmission power are assumed unchanged during the window time, the suffix [t] can be omitted.
As can be seen, problem P 3 is formulated in a general form where the penalty matrix can be flexibly determined according to the design goals, e.g., inter-beam interference, and precoding computation cost.The problem P 3 is again non-convex, non-linear optimization problem due to demand satisfaction constraint as in the problem P 0 .

A. Convexification of the Problem P 3
The problem P 3 has a non-convex constraint related to demand satisfaction: The non-convexity comes in part from the interference in the denominator of the SINR.The strong interference from the adjacent beams can be eliminated while the interference from non-adjacent beams is exceptionally low thanks to the narrow beam pattern design of satellite systems.Therefore, in an interference-free model, the SINR can be approximated as SNR and expressed as follows: where α ≥ 1 is an adjustable scale factor that can be used for further implementation (iteration algorithms if necessary).Then, the data-demand constraint can be transferred to where d = [D 1 , . . ., D N ] T is the demand vector and When the channel and power transmission are assumed to be unchanged during T time slots, we can simplify and express C as 1 1×T ⊗ diag(r 1 , . . ., r N ), where r n stands for the achieved data rate of the user in beam n within one time slot.The approximated problem can be stated as Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
This is a binary quadratic programming, which can be solved by several methods or optimization solvers such as Gurobi, Mosek.
Other solution approaches are also introduced in the following subsections.
Next, we investigate the optimal solution of problem P3 .To do so, we adopt a new method in which the demand D n is translated into the minimal number of time slots that beam n must activate in order to match its demand over the hopping window.Let r n (in bps) be the average achievable rate of the user in beam n in the interference-free model, i.e., f (γ n ) = r n .The achievable rate of the user in beam n over the hopping window can be written as R n = T t=1 (x n [t]r n Δ T ).Now, let us consider the following Lemma [11].
Lemma 1: Constraint (C2) in problem P3 can be reformulated as where dn = D n Δ T r n , and .represents the ceiling operator.
Since, T t=1 x n [t] is an integer, the right hand side of (8) can be replaced by dn = D n Δ T r n .This completes the proof.In constraint ( Ĉ2), dn can be considered as the minimum number of time slots in each of which beam n is illuminated.Next, the following theorem [11] relates the optimal solution of problem P3 and dn .
Theorem 1: Let {x n [t]} be the optimal solution of problem P3 , then T t=1 x n [t] = dn , ∀n, if r n 's are estimated accurately, i.e., constraints ( Ĉ2) hold for all beams.
Proof: We use the contradiction method to prove this theorem.Suppose that there exists at least one beam that the corresponding constraint ( Ĉ2) does not hold.Let n be the beam for which T t=1 x n [t] > dn .Selecting any time slot t such that x n ,t = 1, we generate the new solution of P3 , {x n [t]}, such that xn [t] = x n [t] ∀(n, t) = (n , t ) and xn [t ] = 0.It is easy to observe that xn [t] 's satisfies constraints ( Ĉ2).Moreover, this new solution satisfies constraints (C1) while resulting in the lower objective function.It follows by a contradiction since {x n [t]} is the optimal solution.Therefore, constraint ( Ĉ2) holds for all beams with any optimal solution of problem P3 .
More precisely, under the optimal solution of problem P3 , constraint (C2) holds with equality for all beams, i.e., the supplied capacity and demand are exactly matched for all beams.Therefore, the solutions obtained from the proposed algorithms are compared with the optimal solution in terms of the average beam demand satisfaction in Section VII where 100% average beam demand satisfaction denotes the optimal solution.

B. Proposed Solution 1) Semidefinite Programming (SDP) Relaxation Method:
Problem P3 can be transformed into the semidefinite form as follows: where X = xx T .This problem can be solved with advanced optimization toolboxes supporting mixed-integer models like CVX.

2) Successive Convex Approximation (SCA) Based Solution:
The main challenging issue for solving problem P3 comes from the binary constraint.To overcome this issue, we apply feasible point pursuit successive convex approximation (FPP-SCA), which not only linearizes the non-convex parts of the problem as conventional SCA does, but also adds slack variables to sustain feasibility, and a penalty to the objective function to ensure the slacks are sparingly used.The procedures are given as in the following.
First, we can rewrite binary constraint x n [t] ∈ {0, 1} as Using a tolerance factor ξ n,t (ξ n,t ≥ 0), these conditions can be further expressed into two following inequalities convex, the other is not.By using the first-order Taylor series to linearize this concave constraint, the constraint can be converted into where stands for the value of x n [t] achieved in the iteration k of the applied iterative solution procedure.Here, it is worth mentioned that the starting points {x Then problem P3 can be transferred into the following convex problem, which can be solved optimally in each iteration of the iterative solution procedure.
Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
where δ is the penalty parameter.Usually, the higher value it is set to, the more closed to zeros or ones the solution can be achieved.Hence, δ should be set sufficiently large in order to obtain a good/acceptable solution.However, if the value of δ is set too high, the δ ξ component will be the main term of the objective function, and the solution process may only focus on minimizing this term.This may lead to an inefficient solution of x.Basically, we can start with a small value of δ, and then δ is increased step by step if the continuous solution of x is still not close to zeros or ones.

V. QUEUING-BASED TIME-SLOT BY TIME-SLOT SOLUTION
As can be observed, the complexity of defining the BH for the whole window time is scaled up heavily with the numbers of beams and time-slots.The optimization-relating task becomes extremely challenging when these numbers are sufficiently large [11].Considering the remaining data as the queue length for every user, this section aims to exploit the queue-lengthbased resource allocation method presented in [33], [34], [35] to develop two effective low-complexity BH approaches.In particular, these approaches aim to decompose the window-based BH optimization problem into several sub-problems each of which is stated for BH design in a specific time slot based on the remaining data at users.Generally, to satisfy users' demand on delivering all their data amount during the time window, the beams should be illuminated in each time slot to maximize the weighted sum rate of all users.Herein, the weights can be defined carefully according to the design perspective.
To begin with, one can express the remaining data demand of beam n at time slot t as Let W (q n [t]) represent the weight factor relating to the queue length corresponding to beam n in time-slot t, which is a definable function of q n [t].Utilizing these weights, i.e., W (q n [t])'s, assuming the interference free model as in (7), the BH design task corresponding to time slot t can be cast as sum-weightedrate maximization problem as follows [33], [34], [35]: Regarding the fairness and the satisfaction rate among all beams, there are two weighting approaches will be introduced for setting the weight factor W (q n [t]) for every beam in every time slot.

A. Linear Weighting Approach
Thanks to (9), one can observe that the data demand of beam n is fully transmitted before time slot t if and only if q n [t] = 0. Normally, fairness among all users with different demands can be cast by focusing on minimizing the window time within that all beams can achieve their demand.This problem can be stated as In this scheme, the higher priority is set to the users having the longer queue length.This method aims to reduce the demand for the high-load beam as much as possible.In particular, the weight factor of each beam in each time slot t in this scheme is determined as [33] W (q Following the result given in [33], determining the weights as the linear form of queue lengths can achieve fairness where the remaining demand data among all beams are kept equal even though each beam has different demand.However, there exist some crucial scenarios in that no feasible solution for problem (MIN_BH) can be defined, e.g., the system cannot transmit all data demand for all beams within T time slot.In such cases, satisfying as many users at their demand data within T time slots as possible can be considered as a suitable decision, which will be presented in the following section.

B. Hyperbolic Weighting Approach
As can be observed, by setting the weight factor as in (11), the queueing-based framework aims to select the beam demanding high traffic for illumination pattern in every time lot.When the scheme is overloaded, this weight setting leads to the fact that the beams with extremely high traffic demand will be activated in all-time lots.In such cases, the beam-illumination fairness should be violated and the satisfaction rates related to beams having low capacity demands are very low.To deal with such problems occurring in the overload schemes, we consider another weighted setting as follows.
By defining the weight factor as in (12), the framework will take higher priority for users who have shorter queue lengths in every time slot.Then, the number of users satisfying their traffic demands can be increased.Hereafter, we name the queuingbased BH method using the weight factor as in (11) and in (12) as Linear Weight Queuing-based (LWQ) and Hyperbolic Weight Queuing-based (HWQ) methods, respectively.

VI. COMPLEXITY ANALYSIS OF THE PROPOSED SCHEMES
In this section, we evaluate the computational complexity of the proposed schemes.
1) FCH Scheme: The complexity of the FCH scheme depends on the number of snapshots in each time slot.Let N ss be the number of snapshots for the highest illumination ratio, then the complexity of the FCH scheme can be estimated as Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
O(T N ss ) >> O(T N), where T and N are the numbers of time slots and beams, respectively.
2) SDP Based Scheme: As can be observed, the number of variables in problem (P _SDP) are T N. Therefore, the time complexity of the SDP scheme is O([T N] 4.5 ) [36].
3) SCA Based Scheme: Generally, (P _SCA) is a convex quadratic problem with T N number of variables.Hence, this problem can be solved in polynomial time with the complexity of O([T N] 3 ) [32].Then, the complexity of the SCA-based scheme can be estimated as O([T N] 3 ) × I SCA , where I SCA denotes the number of iterations required to converge which is around 10-20 due to our implementation.
4) Queuing Based Scheme: scheme aims at solving T sub-problem corresponding to the obtain a beam illumination pattern in each time slot as follows: The problem is equivalent to selecting K beams that maximize the objective function value in each time slot.Therefore, the complexities of two queuing-based schemes can be estimated as O(T Nk).

VII. NUMERICAL RESULTS AND DISCUSSION
In this section, we present extensive numerical simulation results to test the proposed algorithms.For this purpose, we consider a simulation scenario consisting of 67-beam GEO satellite beam pattern, i.e., N = 67, covering Europe and part of America and Africa which is provided by ESA [28].Super-user per-beam abstraction is considered, i.e., each beam has a single representative user who aggregates the total demand of randomly distributed users per beam.Suppose K is the maximum number of simultaneously active beams in each time slot, then illumination ratio Q is defined as K N .Detailed performance comparison of the proposed schemes with the benchmark CH technique 6 [29] which uses fixed size, pre-defined clusters is demonstrated.We consider two clustering options, Clustering Option-1 and Clustering Option-2 corresponding to the cluster sizes of 4 beams/cluster and 6 beams/cluster, respectively, for the CH technique similar to that was used in [29].Hereafter, the CH technique with clustering Option-1 and Option-2 is also abbreviated as CH 1 and CH 2, respectively.Also, we consider three illumination ratio values, i.e., 1  4 , 1 6 and 1 8 .The beam configuration for the three illumination ratios is provided in Table I.
The per-beam power depends on the illumination ratio (or on the number of active beams per snapshot) and can be expressed as, P beam [dB] = 10 • log 10 P T Number of active beams 6 Since the CH technique outperforms the conventional BH scheme as already shown in [29], we have selected the best up-to-date illumination pattern design scheme, i.e., CH scheme, as a benchmark scheme.Remark 1: Attenuation due to rainfall can be significant at Kaband.Rain attenuation is already considered in actual operating GEO VHTS at the radio link design level.When GEO VHTS started moving towards higher frequency bands, the adaptive code modulation (ACM) in the DVB air interface standard was introduced for rainfall fade mitigation.
The following KPIs are used to measure the performance of the proposed algorithms: r System Unmet Capacity, defined as where d i and s i denote the demand and the supplied capacity of the i-th beam, respectively.
r System Unused Capacity, defined as r The system supplied capacity is the total average supplied capacity of considering all beams: Note that the supplied capacity C supplied = C unused + C met , where C met represents the amount of system demand that is satisfied.Average beam demand satisfaction, which is defined as, Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.Sometimes it is interesting to evaluate the minimum of the beam demand satisfaction percentage, i.e.BDS min = min(BDS 1 , . .., BDS N ), since it gives us an indication of how well the worst users on the system are addressed.
r System efficiency, i.e., the fraction of the supplied capacity utilized to satisfy the system demand is defined as,

A. Performance Evaluation of the Proposed and the Benchmark Schemes for Pre-Fixed Demand Instances Generated by Traffic Emulator
In this set of experiments, pre-fixed traffic demand distributions are generated using Satellite Traffic Simulator in [37].As a part of the traffic simulator, the three credible input data sets, i.e., population, aeronautical, and maritime, are analyzed and processed simultaneously to simulate the spatial-temporal traffic demand and distribution during 24 hours of the day.To account for the time variability of the traffic, we have considered the traffic demand instances at three different time stamps, i.e., 1:00 am, 8:00 am, and 12:00 pm as shown in Fig. 1 which correspond to low, medium, and high demand scenarios, respectively.The total system demand of the three traffic instances is summarized in Table III.Fig. 2. Comparison of per beam supply-demand matching using various techniques corresponding to the demand instance at 12:00 pm (i.e., high demand) and illumination ratio 1  4 .
Fig. 3. Comparison of per beam supply-demand matching using various techniques corresponding to the demand instance at 1:00 am (i.e., low demand) and illumination ratio 1 6 .
Figs. 2 and 3 show the demand-supply matching profile of each beam corresponding to demand instance at 12:00 pm (i.e., high demand) and 1:00 am (i.e., low demand), respectively, for different considered techniques.In addition, high (i.e., 1/4) and low (i.e., 1/6) illumination ratios are used corresponding to the high and low demand cases, respectively.It can be seen that there is a large deviation in the beam supply-demand matching profile obtained from the benchmark techniques CH 1 and CH 2 as compared to the proposed techniques.In particular,  SCA based technique outperforms all other techniques in terms of the beam supply-demand matching metric.
Figs. 4 and 5 report on the CDF of the per slot transmit power utilized to satisfy beam demands corresponding to 1:00 am and 12:00 pm, respectively, for all the considered techniques.The benchmark CH scheme activates a fixed number of beams per time slot, and hence, the transmitted power in each slot is also fixed, i.e., independent of the clustering option as well as the illumination ratio.As observed in Figs. 4 and 5, the CH scheme always performs the worst whilst the proposed FCH scheme outperforms for low as well as for high demand instances.With the proposed schemes, the algorithm decides how many beams per time slot to activate, typically resulting in transmit power savings as fewer beams are activated whenever possible.
Tables IV and V summarize the performance of the various schemes in terms of all the KPIs for low and high demand cases, respectively.From both tables, it is apparent that the proposed schemes perform better than the benchmark CH scheme in terms of the average beam demand satisfaction, and system efficiency.Similarly, the benchmark CH scheme performs the worst in terms of the total of unmet and unused capacity for both the low and high-demand cases.
The results given in Table VI represent the Matlab implementation running time complexity of the considered schemes which is estimated by calling the "Tic-Toc Matlab function" beforeafter the CVX on PC with Intel Core i7-1085H CPU @ 2.7 GHz for the demand instance case at 12:00 pm (i.e., high demand).Note that we set-up a limit of 120 seconds to conclude the task.This is because CVX uses an iterative approximation method   to obtain the optimal solution and it might get stuck trying to achieve the target accuracy.As observed, the implementation time complexities of the queuing methods (i.e., LWQ and HWQ) are extremely less as compared to that of the other schemes and the SCA scheme has the highest time complexity due to its iterative behavior.Focusing on the queue-based schemes, we identified that the LWQ method works better for low-loaded scenarios and the HWQ method works better for highly-loaded scenarios.In general, assuming that the satellite resources have been well planned for the expected demand, the LWQ approach would be the preferred option.

B. Performance Evaluation of the Proposed and the Benchmark Schemes Using Monte Carlo Simulations for Randomly Generated Demand Instances
Herein, we present the numerical results corresponding to Monte-Carlo simulations for the proposed and benchmark techniques for 100 instances of the three types of traffic demands which are randomly generated from different distributions.Essentially, the beams are classified into three categories, i.e. high-demand beams (hot), medium-demand beams (warm), and low-demand beams (cold) based on the demand limit values shown in Table VII.Once beams are classified, their demand instance is generated uniformly randomly in the limited range.The following distributions are considered: (i) Clustered Demand 1, (ii) Clustered Demand 2, (iii) Clustered Demand 3.
Clustered demand 1 consists of 1 big cluster of "Hot" beams at the center.Clustered Demand 2 has 3 randomly located clusters each consisting of 6-7 "Hot" beams.Finally, Clustered demand 3 has randomly located multiple small-sized clusters which consist of high-demand "Hot" beams.Herein, we want to investigate the effect on the cluster sizes obtained from FCH algorithm with different demand instances.Table VIII depicts the variation of the average number of clusters of different sizes for two different kinds of randomly generated demand instances, i.e., Clustered Demand 2 and Pre-fixed 1:00 am, which corresponds to high and low demand scenarios with a total demand of 26.24 Gbps and 24 Gbps, respectively.It can be clearly seen from Table VIII that the number of large-size clusters are less and increases when demand increases.This also shows the superiority of the proposed FCH algorithm since less number of large-size clusters results in the low computation cost required for precoding.Moreover, Clustered Demand 2 consists of a large number of beams with high demand (i.e., warm and hot beams) as compared to the Pre-fixed 1:00 am demand.Therefore, the FCH algorithm results in an illumination pattern in which large-size clusters are generated to match the high demand of those warm and hot beams.As illustrated in Table VIII, there exist clusters with sizes 14, 15, 16, and 17 in the case of Clustered Demand 2 as compared to the case of Pre-fixed 1:00 am demand in which the maximum cluster size is 13.Hence, the proposed FCH algorithm adapts supplied capacity to the beam demand distribution and simultaneously activates the required number of multiple beams (i.e., cluster) to match the heterogeneous demands.
Fig. 6 describes the CDF variation of the unmet capacity, minimum beam demand satisfaction (i.e., BDS min ) and average beam demand satisfaction (i.e., BDS avg ) using different considered schemes for the Clustered Demand 2 and illumination ratios 1/4 and 1/6.As observed from the Fig. 6 for the given illumination ratio, the proposed FCH scheme outperforms and than LWQ scheme performs second best among the remaining schemes in terms of all the unmet capacity, BDS min and BDS avg for Clustered Demand 2. The same experiments are repeated also for Clustered Demand 1 and Clustered Demand 3. Figs.7 and 8 report the CDF of the unmet capacity and average beam demand satisfaction rate BDS avg for Clustered Demand 1 and 3, respectively.It can be confirmed that the proposed FCH and queuing-based schemes outperform the benchmark CH scheme.Next, Fig. 9 reports the system efficiency for Clustered Demand 2 using all the considered schemes and it can be inferred that the FCH scheme with an illumination ratio of 1/6 gives the best efficiency.Also, the FCH scheme with an illumination ration of 1/4 results in poor efficiency because it provides excess supplied capacity compared to the system demand.
After evaluating the proposed techniques with very different demand instances, we can confirm that the proposed techniques are able to provide very good demand satisfaction with any type of input demand.The FCH and the queuing-based schemes achieve good demand matching with a slightly better match for FCH scheme.However, the computational time of the queuingbased scheme and its generally good performance make it very attractive for practical satellite systems.

1) Effect of Hopping Window Size:
Herein we analyze the impact of the parameter T , i.e., the number of time slots in the hopping window on the performance of the various proposed schemes.Fig. 10 presents the BDS avg and unmet capacity versus the various values of the number of time slots in the hopping window.As can be observed, the higher value of T results in the best demand-matching numbers.In particular, when the hopping window size increases, the proposed schemes have more flexibility in assigning the active time-slots per beam resulting in better performance in terms of unmet capacity and average beam satisfaction rate.In practical systems, however, we need to consider a trade-off between window length (periodicity of the illumination pattern) and per-beam latency of activation.The longer the beam-hopping window, the more flexibility in Fig. 10.Effect of the parameter T , i.e., the number of time slots in the hopping window on the performance of the proposed schemes with illumination ratio 1/6 and for the Clustered Demand 2. matching the demand but a longer window may have a negative impact on latency and revisiting times.From the results, it seems that, after some point, increasing the window size does not have a significant impact on the performance; which evidences that an optimal point may be found by balancing the different aspects of the final users' quality of experience.
2) Effect of the Number of Simultaneously Illuminated Beams K in a Time Slot: Herein, we investigate the performance of the proposed algorithms with the maximum number of illuminated beams, K, in a time slot.Table IX summarizes the performance of the various algorithms in terms of Unmet and Unused capacity for Clustered Demand 2 case.It can be seen that the FCH algorithm requires less number of illuminated/active beams in each time slot to best satisfy the user's demand as compared to the case of using the SCA, LWQ, and HWQ algorithms.Furthermore, it can be deduced that FCH algorithm outperforms for low value of K and the LWQ algorithm outperforms for the high value of K.

VIII. CONCLUSION
In this article, highly flexible beam illumination pattern designs were proposed for efficient utilization of resources in multi-beam satellite systems which allow clusters of different shapes and sizes depending on the user demands.A maxmin fair beam supply-demand framework was developed with the goal of reducing interference using precoding by ensuring the maximum number of illuminated beams in each time slot.The interference minimization problem was also formulated as a penalty minimization problem, and two solution approaches, i.e., SDP-based and SCA based, were proposed without the need for complex precoding designs.Additionally, two queuingbased approaches, LWQ and HWQ, were proposed to maximize the sum weighted rate in each time slot.
The performance of the proposed schemes has been evaluated via extensive numerical simulations.The experimental results have revealed that the proposed schemes outperform the relevant benchmark CH scheme in terms of beam supply-demand matching.Finally, the proposed schemes are practically important as they address the challenge of efficiently designing beam illumination patterns while considering the heterogeneous user demands and inter-beam interference, which is crucial in satellite communication systems where interference can severely impact the quality of service.

n
[t]} 's can be chosen freely.

r
Beam demand satisfaction (BDS) rate of beam b, which is defined as, BDS b (%) = 100 × min s b d b , 1 .

Fig. 4 .
Fig. 4. Comparison of the CDF of the transmit power per slot to satisfy beam demands corresponding to 1:00 am demand instance for various techniques.

Fig. 5 .
Fig. 5. Comparison of the CDF of the transmit power per slot to satisfy beam demands corresponding to 12:00 pm demand instance for various techniques.

Fig. 6 .
Fig. 6.Comparison of the benchmark CH scheme with the proposed schemes in terms of the CDF of unmet capacity, minimum beam demand satisfaction (BDS min ) and average beam demand satisfaction (BDS avg ) obtained from Monte Carlo simulations corresponding to the Clustered Demand 2.

Fig. 7 .
Fig. 7. Comparison of the benchmark CH scheme with the proposed schemes in terms of the CDF of unmet capacity and average beam demand satisfaction rate (BDS avg ) obtained from Monte Carlo simulations corresponding to the Clustered Demand 1.

Fig. 8 .
Fig. 8.Comparison of the benchmark CH scheme with the proposed schemes in terms of the CDF of unmet capacity and average beam demand satisfaction rate (BDS avg ) obtained from Monte Carlo simulations corresponding to the Clustered Demand 3.

Fig. 9 .
Fig. 9. Comparison of the benchmark CH technique with the proposed techniques in terms of the CDF of the system efficiency obtained from Monte Carlo simulations corresponding to the Clustered Demand 2.

TABLE I NUMBER
OF CLUSTERS AND NUMBER OF ACTIVE BEAMS PER SNAPSHOT FOR DIFFERENT CLUSTERING OPTIONS − Total loss[dB].

TABLE IV COMPARISON
OF VARIOUS TECHNIQUES IN TERMS OF DIFFERENT KPIS FOR THE DEMAND INSTANCE 1:00 AM

TABLE V COMPARISON
OF VARIOUS TECHNIQUES IN TERMS OF DIFFERENT KPIS FOR THE DEMAND INSTANCE 12:00 PM

TABLE VI RUN
TIME COMPLEXITY OF DIFFERENT SCHEMES FOR DIFFERENT ILLUMINATION RATIOS AND FOR HIGH DEMAND INSTANCE AT 12:00 PM

TABLE IX COMPARISON
OF VARIOUS ALGORITHMS IN TERMS OF UNMET AND UNUSED CAPACITY FOR CLUSTERED DEMAND 2 FOR DIFFERENT VALUES OF K