# Utility Maximization for Load Optimization in Cellular/WLAN Interworking Network Based on Generalized Benders Decomposition

Authors: Fanqin Zhou; Lei Feng; Peng Yu; Wenjing Li; Luoming Meng
Journal: Mobile Information Systems
 (2016)

---

## Abstract

Load steering is widely accepted as a key SON function in cellular/WLAN interworking network. To investigate load optimizing from a perspective of system utilization maximization more than just offloading to improve APs’ usage, a utility maximization (UTMAX) optimization model and an ASRAO algorithm based on generalized Benders Decomposition are proposed in this paper. UTMAX is to maximize the sum of logarithmic utility functions of user data rate by jointly optimizing user association and resource allocation. To maintain the flexibility of resource allocation, a parameterβ is added to the utility function, where smaller β means more resources can be allocated to edge users. As a result, it reflects a tradeoff between improvements in user throughput fairness and system total throughput. UTMAX turns out to be a mixed integer nonlinear programming, which is intractable intuitively. So ASRAO is proposed to solve it optimally and effectively, and an optional phase for expediting ASRAO is proposed by using relaxation and approximation techniques, which reduces nearly 10% iterations and time needed by normal ASRAO from simulation results. The results also show UTMAX’s good effects on improving WLAN usage and edge user throughput.


---

## Body

## 1. Introduction

With the popularization of intelligent mobile terminals and enrichment of Internet services, it arouses surging demands for mobile Internet access in mobile users. So higher requirements on access capacity and data rate are put forward to commercially operated mobile networks. As spectrum band authorized to a mobile network is so scarce, academic and industrial research groups are continuously working on more spectrum-efficient radio communication technologies. Besides endeavors in this direction, taking advantages of available extra frequency bands seems much more economical and practical, which leads to a common deed among network operators of deploying WLAN in public areas, because WiFi works on the open ISM frequency band and wireless traffic can be offloaded from cellular network (CN) to WLAN.A typical traffic offloading case in cellular/WLAN coexisted scenario is illustrated in Figure1. In Figure 1(a), users (or UEs without distinction hereafter) are served dominantly by cellular base stations (BSs), while APs are underutilized. As service areas of CN and public WLAN overlap, to perform offloading, we can simply adjust user association and optimize resource allocation accordingly. The optimized scenario is expected to be as in Figure 1(b), where some users handover from BSs to nearby APs.Figure 1
Offloading in cellular/WLAN interworking network.
(a) (b)Initially, public access points (APs) in WLAN were deployed mainly in traffic hot-spots where many mobile users gathered, resulting in high usage of APs. However, with public WLAN gradually being a scale large enough to cover almost all urban areas, the problem is increasingly highlighted that most public WLAN APs are underutilized. The strict power limit of ISM devices is a root cause, while a lack of load distribution optimization approaches in present cellular/WLAN coexisting networks also plays a role.So load steering between cellular and WLAN was an important research case in SEMAFOUR, a project in EU Framework Program 7 (EU-FP7) for developing multi-RAT/multilayer self-optimizing network (SON) function and integrated SON management system [1]. Its resulting approaches pay more attention to applicability in practice, while a theoretical performance bound is not well investigated. A typical approach proposed is to make AP serve users as many as possible, which is sure to improve AP’s usage to the most but will not produce the highest system utilization or overall users’ satisfaction. Therefore, offloading is deemed as a basic purpose, but load optimizing is a further goal in cellular/WLAN interworking networks.Moreover, providing better satisfactory data rate for users with the same infrastructure means more economical utilizing of investment, which is of more significance to commercially operated networks. To this end, optimizing user association and resource allocation to achieve ultimate overall users’ satisfaction in cellular/WLAN interworking network are to be investigated in this paper.Many prior works related to load optimizing in cellular/WLAN interworking networks have been surveyed in [2, 3]. In this paper, we roughly group the related approaches into three categories, according to whether limited, fair, or sufficient information exchange between CN and WLAN is required.(1)
UE-decision based approaches are widely used where information exchange is difficult, and statistical methods are used to estimate proper handover trigger criteria for individual UE on aspects, such as the average amount of data that could be transferred through CN or WLAN [4], the holding and continuing of an existing traffic data stream [5], the hysteresis for handover between CN and WLAN [6], the load status in a BS or AP [7], and the partial packet success [8]. These approaches focus on designing practical rules for individual UE to make proper handover decision without assists from networks.(2)
Network-assisted UE-decision based approaches become possible with more information attainable about a target network through information exchange, such as the achievable data rate in a target eNB or AP [9, 10], UE’s velocity, and received signal strength of target AP [11]. These works use different methods, such as TOPSIS [10] and fuzzy logical system [11] to process multidimensional information obtained from network to help UE make better handover choice.(3)
Centralized optimization approaches are applicable when sufficient information is achievable. They usually focus on providing theoretical performance bounds for coexisting networks on various optimization goals, such as balancing resource utilization [12], minimizing power consumption [13, 14], maximizing quality of experience [15], and maximizing multilink aggregation [16]. These approaches may not be intended for application in actual networks, but they provide theoretical bounds for target optimization goals, which is significant for the design of practical network optimizing methods.The existence of various types of approaches is partially due to the different requirements in actual practices of network optimizing and theoretical analyses in academic researches and largely due to the development of network technologies. In the early stage of cellular/WLAN coexisting network, the lack of information exchange between noninterworking CN and WLAN caused great difficulty in controlling user access or handover in a centralized way. Therefore, works at that time focused on UE-decision based approaches.To improve the performance of coexisting networks, information exchange between networks is a prerequisite. For this purpose, 3GPP proposed solutions for interworking between 3GPP network and WLAN in [17], and a media independent handover (MIH) standard was developed by IEEE in [18]. With more information attainable, especially that about load status in a target BS or AP (we use AN to denote a BS or AP hereafter), individual UE can even estimate its achievable data rate in a target BS or AP, so better handover decision can be made. On these bases, network-assisted UE-decision based became popular.As multiple users may contend for access chances and resources, UE-decision based schemes neither can guarantee a resulting decision always to be the optimal nor can tell where a performance bound goes from a whole network perspective. Therefore, centralized optimization approaches are widely investigated to derive such a theoretical performance bound. These approaches are usually of high computational complexity; however, deriving such a performance bound will give guidance and direction for more practical approaches with better tradeoff between optimality and complexity.Therefore, this paper is to derive a centralized logarithmic utility maximization model to improve overall users’ satisfaction by jointly optimizing user association and resource allocation. The key contributions are summarized as follows.(i)
A utility maximization (UTMAX) optimization model is formulated to depict the users’ satisfaction degree maximization by jointly planning association selection and resource allocation (ASRA). Instead of directly using a logarithmic function of data rate to model user’s satisfaction of a data rate, which produces identical resource allocation in a same BS and thus lack flexibly in resource allocation, we add a control parameterβ to it. With smaller β value, the users with poor channel conditions are to be assigned more resources and vice versa. Therefore, a tradeoff between user throughput fairness and system throughput gain can be achieved by selecting proper β value.(ii)
The derived UTMAX has both binary and continuous variables for user association and resource allocation, respectively, which is a mixed integer nonlinear programming (MINLP) and seems intractable intuitively. To solve it optimally and effectively, we devise an ASRAO algorithm on the basis of Bender Decomposition (BD) framework [19, 20]. We also propose optional procedures for ASRAO to expedite the computation of UTMAX using relaxation and approximation techniques [21], which dramatically reduces the number of iterations needed before convergence, especially when a UTMAX instance is in large scale. In addition, the convergence of our algorithms is validated both theoretically and experimentally in this paper.(iii)
Our work contributes to SON. The UTMAX and its corresponding ASRAO can serve directly in the scheme-planning module, which forms a self-optimizing closed loop together with monitor, evaluation, and execution modules, to generate the optimal ASRA scheme through a centralized way. It also provides a performance upper bound for load optimization on purposes of a system preferred tradeoff between overall throughput and end user satisfaction. So It can guide the development of novel suboptimal but less computationally complex algorithms by evaluating the performance gaps between them.The rest of the paper are arranged as follows. In Section2 we formulate the user data rate and propose the UTMAX optimization model. Section 3 analyzes a special property of UTMAX, and, to exploit this property, we propose the ASRAO algorithm based on generalized Benders Decomposition (GBD) to solve it. Optional procedures for ASRAO to expedite the computing of UTMAX are also given in this section. In Section 4, we evaluate the performance of (A-)ASRAO and effect of UTMAX through extensive simulations. The paper will be concluded in Section 5.


## 2. System Model and Problem Formulation

A cellular/WLAN interworking network is composed ofJ BSs and K APs, serving I UEs. We denote by J, K, and I the sets of all BSs, APs, and UEs, respectively. UE is served uniquely by a BS or AP, from which it gets largest received power. The association between UE i and AN l is indicated by binary variable xil, which equals 1 if UE i is served by AN or 0 otherwise as in (1). A matrix X=(xil)i∈I,l∈J∪K represents associations between all UEs and ANs. UE is uniquely served by one AN, namely, a BS or AP, which forms the constraint in (2).(1)xil=1,UEiisservedbyANl,0,otherwise,(2)∑l∈J∪Kxil=1.


### 2.1. Channel and Resource Model of CN

The access rate of a UE in CN is determined by two factors, the channel condition and the allocated resource.ΛJ=λij, i∈I, j∈J represents the resource allocation scheme in CN, in which λij denotes the ratio of resource allocated to UE i from BS j. To assure a basic quality of service (QoS), a minimal ratio of resource, λε, should be allocated to UE i from BS j if xij=1, and λij=0 if xij=0, which can be summarized in (3). The channel condition between UE i and BS j can be indicated by SINRij, as in (4), in which Prij is the received signal power of UE i from BS j, and α is a coefficient of cochannel interference. According to Shannon Formula, the data rate of UE i served by BS j, rij can be expressed in (5), wherein Rij=Wjlog(1+SINRij).(3)λεxij≤λij≤xij,(4)SINRij=Prijα·∑i′∈B/iPri′j+σ2,(5)rij=λijWjlog21+SINRij=λijWjγij=λijRij.


### 2.2. Channel and Resource Model of WLAN

Due to the sparse deployment of public WLAN APs and the limited transmit power, for simplification, the interferences between WLAN APs will not be taken into consideration. Reference [9] gives the way to evaluate the service rate of UE i∈I accessing public WLAN AP k∈K, as in(6)rik=CWLAN·MCSWLANPrikmax⁡MCSWLAN·i∈I,k∈K,whereinCWLAN denotes the average data rate evaluated by AP when the modulation and coding scheme (MCS) with the highest spectrum efficiency is used and Prik is the received power strength of UE i from WLAN AP k. MCSWLAN(·) is corresponding to each possible Prik in MCS table. Taking resource allocation into consideration, we have CWLAN=CWLAN_FULLλik, where CWLAN_FULL is the available data rate if the observing UE occupies resources exclusively in AP. As interference to AP is not taken into consideration and only a limited number of UEs are served in each AP, we assume CWLAN_FULL/max{MCSWLAN(·)}=1. Therefore, rik=λikRik.AP uses Round Robin (RR) for packet data scheduling, which from long time average perspective can be viewed as UEs in AP equally sharing channel resources (7). Network Operators commonly limit the data rate of UEs in WLAN APs to avoid resources being occupied dominantly by user requesting big volume of data with good channel condition. For simplification, in this paper we assume the ratio of resources occupied by UE in AP has an upper bound λm≤1 as in (8). In addition, UEs contend to access AP via shared wireless channel, and too many UEs would severely deteriorate the throughput performance, so the number of UEs served simultaneously by AP has an upper bound Nm as in (10), which amounts to setting a lower bound to resources allocated to UE. Hence, we have the following constraints on resource allocation in AP.(7)λik=λkavrxik,i∈Ik,k∈K,(8)λik≤λmxik,(9)∑i∈Iλik≤1,(10)∑i∈Ixik≤Nm.As UE is served uniquely by a BS or AP, the overall data rate of UEi can be formulated as ri=∑i‍λilRil. It should be noted that Ril is not a variable to be optimized in our model that means more accurate approaches for full band rate estimation can be used, which is not a key point in this paper and we will pay no more attention to it. After formulating the data rate of UE in CN and WLAN, we take a utility function perspective to investigate a utility maximization problem for the data rate ri of each UE i in the system, to get the optimal association and resource allocation.


### 2.3. Utility Maximization Problem

For utility maximization, utility functionUi(·) is often chosen to be continuously differentiable, monotonically increasing, and strictly concave. UE i gets utility Ui(ri) when its receiving rate is ri. The overall system utility can then be represented as ∑i‍Ui(ri). Maximize it with some practical constraints, we can get the general utility maximization problem. Ui(·) should be chosen carefully. If Ui(·) is a linear function, max∑i‍Ui(ri) turns out to be a throughput maximization and it produces a trivial solution, where each BS allocates resources only to its user with best channel quality. It is not a satisfactory solution for a multiple user system. So we seek a utility function that would naturally encourage fairness resource allocation among users. To this end, a logarithmic function is widely used. It is concave and has diminishing returns. This property approximates the fact that a well-served user has low priority in getting more resources, so it encourages load balancing and improves overall users’ satisfaction.However, logarithmic utility function, like in [15], always results in equal resource allocation among users in a same BS. Because ∑i‍logRiλi=∑i‍logRi+∑i‍logλi, the utility of individual BS has nothing to do with Ri, and ∑i‍logλi=log∏i‍λi≤nlog∑i‍λi/n always hold. It reaches maximum only when resources are all shared equally, λi1=λi2, ∀i1,i2∈Il. This lacks flexibility in resource allocation, so we slightly adjust this logarithm utility function and add a BS related parameter βl to ri, to make resource allocation scheme represent certain system preference, like caring more about fairness or demanding higher system throughput. If bl=0, it means an identical resource allocation. If bi>0, users with better channel condition get more resources and vice versa. The logarithmic utility function then becomes log(∑l‍Rilλil+βlxil). As xil is binary variable and ∑l‍xil=1, log(∑l‍Rilλil+βlxil) can be transformed into ∑l‍xillog(∑l‍Rilλil+βl). And our utility maximization problem is formed as P1 in P1, where for simplicity we assume in CN λm=1, and in WLAN λε=0.(11)P1:maxΛ,XΦΛ,X,s.t.ΦΛ,X=∑i∑lxillog⁡∑lRilλil+βl,∑lxil=1,∑iλil≤1,λεxil≤λil≤λmxil,λik=λkavrxik,∑i∈Ixik≤Nm,xil∈0,1,i∈I,k∈K,l∈J∪K.P1 has both binary and continuous variables and nonlinear part in goal function. These properties make it an MINLP, which is intractable intuitively. However, this problem has a special structure, utilizing which we propose ASRAO algorithm based on GBD to solve it effectively. The problem analyses and ASRAO algorithm will be presented in next section.


## 2.1. Channel and Resource Model of CN

The access rate of a UE in CN is determined by two factors, the channel condition and the allocated resource.ΛJ=λij, i∈I, j∈J represents the resource allocation scheme in CN, in which λij denotes the ratio of resource allocated to UE i from BS j. To assure a basic quality of service (QoS), a minimal ratio of resource, λε, should be allocated to UE i from BS j if xij=1, and λij=0 if xij=0, which can be summarized in (3). The channel condition between UE i and BS j can be indicated by SINRij, as in (4), in which Prij is the received signal power of UE i from BS j, and α is a coefficient of cochannel interference. According to Shannon Formula, the data rate of UE i served by BS j, rij can be expressed in (5), wherein Rij=Wjlog(1+SINRij).(3)λεxij≤λij≤xij,(4)SINRij=Prijα·∑i′∈B/iPri′j+σ2,(5)rij=λijWjlog21+SINRij=λijWjγij=λijRij.


## 2.2. Channel and Resource Model of WLAN

Due to the sparse deployment of public WLAN APs and the limited transmit power, for simplification, the interferences between WLAN APs will not be taken into consideration. Reference [9] gives the way to evaluate the service rate of UE i∈I accessing public WLAN AP k∈K, as in(6)rik=CWLAN·MCSWLANPrikmax⁡MCSWLAN·i∈I,k∈K,whereinCWLAN denotes the average data rate evaluated by AP when the modulation and coding scheme (MCS) with the highest spectrum efficiency is used and Prik is the received power strength of UE i from WLAN AP k. MCSWLAN(·) is corresponding to each possible Prik in MCS table. Taking resource allocation into consideration, we have CWLAN=CWLAN_FULLλik, where CWLAN_FULL is the available data rate if the observing UE occupies resources exclusively in AP. As interference to AP is not taken into consideration and only a limited number of UEs are served in each AP, we assume CWLAN_FULL/max{MCSWLAN(·)}=1. Therefore, rik=λikRik.AP uses Round Robin (RR) for packet data scheduling, which from long time average perspective can be viewed as UEs in AP equally sharing channel resources (7). Network Operators commonly limit the data rate of UEs in WLAN APs to avoid resources being occupied dominantly by user requesting big volume of data with good channel condition. For simplification, in this paper we assume the ratio of resources occupied by UE in AP has an upper bound λm≤1 as in (8). In addition, UEs contend to access AP via shared wireless channel, and too many UEs would severely deteriorate the throughput performance, so the number of UEs served simultaneously by AP has an upper bound Nm as in (10), which amounts to setting a lower bound to resources allocated to UE. Hence, we have the following constraints on resource allocation in AP.(7)λik=λkavrxik,i∈Ik,k∈K,(8)λik≤λmxik,(9)∑i∈Iλik≤1,(10)∑i∈Ixik≤Nm.As UE is served uniquely by a BS or AP, the overall data rate of UEi can be formulated as ri=∑i‍λilRil. It should be noted that Ril is not a variable to be optimized in our model that means more accurate approaches for full band rate estimation can be used, which is not a key point in this paper and we will pay no more attention to it. After formulating the data rate of UE in CN and WLAN, we take a utility function perspective to investigate a utility maximization problem for the data rate ri of each UE i in the system, to get the optimal association and resource allocation.


## 2.3. Utility Maximization Problem

For utility maximization, utility functionUi(·) is often chosen to be continuously differentiable, monotonically increasing, and strictly concave. UE i gets utility Ui(ri) when its receiving rate is ri. The overall system utility can then be represented as ∑i‍Ui(ri). Maximize it with some practical constraints, we can get the general utility maximization problem. Ui(·) should be chosen carefully. If Ui(·) is a linear function, max∑i‍Ui(ri) turns out to be a throughput maximization and it produces a trivial solution, where each BS allocates resources only to its user with best channel quality. It is not a satisfactory solution for a multiple user system. So we seek a utility function that would naturally encourage fairness resource allocation among users. To this end, a logarithmic function is widely used. It is concave and has diminishing returns. This property approximates the fact that a well-served user has low priority in getting more resources, so it encourages load balancing and improves overall users’ satisfaction.However, logarithmic utility function, like in [15], always results in equal resource allocation among users in a same BS. Because ∑i‍logRiλi=∑i‍logRi+∑i‍logλi, the utility of individual BS has nothing to do with Ri, and ∑i‍logλi=log∏i‍λi≤nlog∑i‍λi/n always hold. It reaches maximum only when resources are all shared equally, λi1=λi2, ∀i1,i2∈Il. This lacks flexibility in resource allocation, so we slightly adjust this logarithm utility function and add a BS related parameter βl to ri, to make resource allocation scheme represent certain system preference, like caring more about fairness or demanding higher system throughput. If bl=0, it means an identical resource allocation. If bi>0, users with better channel condition get more resources and vice versa. The logarithmic utility function then becomes log(∑l‍Rilλil+βlxil). As xil is binary variable and ∑l‍xil=1, log(∑l‍Rilλil+βlxil) can be transformed into ∑l‍xillog(∑l‍Rilλil+βl). And our utility maximization problem is formed as P1 in P1, where for simplicity we assume in CN λm=1, and in WLAN λε=0.(11)P1:maxΛ,XΦΛ,X,s.t.ΦΛ,X=∑i∑lxillog⁡∑lRilλil+βl,∑lxil=1,∑iλil≤1,λεxil≤λil≤λmxil,λik=λkavrxik,∑i∈Ixik≤Nm,xil∈0,1,i∈I,k∈K,l∈J∪K.P1 has both binary and continuous variables and nonlinear part in goal function. These properties make it an MINLP, which is intractable intuitively. However, this problem has a special structure, utilizing which we propose ASRAO algorithm based on GBD to solve it effectively. The problem analyses and ASRAO algorithm will be presented in next section.


## 3. ASRAO Algorithm

### 3.1. Problem Analyses

Assume association indicatorX is fixed, then each UE i is served by a specific AN l, and each AN l has a fixed UE collection Il={i∣xil=1}. UE i contributes log(Rilλil+βl) to utility of AN l. Its value is determined by λil, which is only related to AN l. So resource allocation in each AN can be solved independently with given user association X˙. Therefore, Φ(Λ,X˙) can be identically transformed to(12)ΦΛ,X˙=∑i∈I,l≡argxil=1logRilλil+βl=∑l∈J∪K∑i∈IllogRilλil+βl=∑l∈J∪KΦlλl,x˙l,where Φl(λl,x˙l)=∑i∈Il‍logRilλil+βl is the utility of AN l. So we get a subproblem of P1 on condition of X˙.(13)maxΛ∑i∈IllogRilλil+bls.t.∑iλil≤1,λεx˙il≤λil≤λmx˙il,λik=λkavrx˙ik.This is a classical problem in information theory and a water filling (WF) principle is specially developed to solve it, which is derived from Karush-Kuhn-Tucker (KKT) optimal conditions. It provides direction solutions to problems with the form like (13). To enable its applicability in P1, we extend WF to a multi-AN case and denote its results as MWF, short for multi-AN water filling. With X fixed as X˙ in P1, we have a MWF applicable problem as in (14)maxΛΦΛ,X˙,s.t.ΦΛ,X˙=∑i∑lx˙illog∑lRilλil+βl,λεx˙il≤λil≤λmx˙il,∑iλil≤1,λik=λkavrxik.According to generalized duality theory, the dual function ofΦ(Λ,X˙) is as in (15)LΛ,ν,ω,γ,X˙=∑i∑lx˙illog∑lRilλil+βl-∑lνl∑iλil-1-∑i∑lωilλil-x˙il+∑i∑lγilλil-εx˙il-∑i∑kμikλik-x˙ikλkavr,where ν, ω, γ are nonnegative relaxed variables and μ is nonzero. The local optimal is reached, when the following KKT conditions are satisfied.(16)∂ΦΛ∂λij=Rij∑jRijλij+∑lx˙ijbj-νj-ωij+γij=0,∂ΦΛ∂λik=Rik∑kRikλik+∑lx˙ikbk-νk-ωik+γik-μik=0,∂ΦΛ∂λkavr=μikx˙ik=0,ωilλil-x˙il=0,νl∑iλil-1=0,γilλil-εx˙il=0,μikλik-x˙ikλkavr=0,νl,ωil,γil∈0,+∞,μik≠0.The solutions to (16) are as follows. MWF results for CN is as in (17)νj∗=∑ix¯˙ijρ¯+∑ibjx¯˙ijRij-1,λij∗=maxx˙ij1νj∗-βjRij,0+λεx˙ij,γij=-RijRij′λij′+βj′-νj∗,RijRij′λij′+βj′≤νj∗,0,else,ωij=RijRij′λij′+βj′-νj∗+γij.The following (18) is MWF results for WLAN.(18)λkavr=minλm,1∑ix˙ik,∑ix˙ik>0,0∑ix˙ik=0,λik∗=λkavrx˙ik,νk∗=minRikβk+∑ix˙ik,λm∑ix˙ik≥1,0,λm∑ix˙ik<1,μik=-RikRik′λik′+βk′-νk∗,λm∑ix˙ik>1,1,λm∑ix˙ik≤1,γik=-RikRik′λik′+βk′-νk∗,RikRik′λik′+βk′<νk∗,0,else,ωik=RikRik′λik′+βk′-νk∗+γik+μik.Using above MWF results, we can directly calculate optimalΛ with a given X˙. To utilize MWF in solving problem P1, an optimal X should be retrieved. However, as X is a matrix of binary variable entries and Λ is closely related with X, it is difficult to find such an optimal. In this paper, instead of directly solving an optimal X, we devise ASRAO algorithm to solve it iteratively and effectively. ASRAO is originated from BD, which provides a framework for addressing mixed integer programming problems by decomposing it into two smaller subproblems and solving them iteratively. The details are in the sequel.


### 3.2. Generalized Benders Decomposition

BD is originally proposed to solve mixed integer linear programming (MILP). Instead of solving all variables and constraints simultaneously, it decomposes an MILP into two subproblems, a master problem (MP) and a slave problem (SP), and the original MILP can be solved by solving MP and SP iteratively. In BD frameworks, MP is an MILP problem which consists of all integer variables in the original problem, while SP is a Linear Programming problem with all the continuous variables from the original problem. In each iteration, BD utilizes an extreme point or extreme ray derived from the dual of SP to generate an optimality cut or feasibility cut to trim the feasible domain of MP, and MP will eventually reach the same optima as the original MILP.GBD extends the BD approach to a more general class of problems by adopting nonlinear duality theory, and some nonlinear problems are brought into range. When solving MINLP, GBD follows the framework of BD, but SP is derived as NLP. To generate the feasibility cuts and optimality cuts, GBD utilized KKT conditions to calculate an extreme point or extreme ray of the Lagrange dual functions of NLP SP. A GBD applicable problem has a general form as in(19)maxx,yFx,y,s.t.Gx,y≥0,Hx,y=0,x∈0,1m×1,y∈D⊆R+n×1.It should meet following situations: (a) for fixedx, (19) separates into a number of independent optimization problems; (b) for fixed x, (19) assumes a well-known special structure that efficient solution procedures are available; (c) (19) may be not concave program in x and y jointly, but fixing x renders it so in y. According to the analyses, UTMAX obviously satisfies these conditions.The initial MP of (19) in GBD procedures is as follows:(20)maxx,ηη,s.t.Gx≥0,Hx=0,η∈R,x∈E⊆0,1m×1,where G(x)≥0 and H(x)=0 are constraints from (19) that only related with x. Solving a solution (x(t),η) from MP, SP can be derived by simply fixing x in (19) as x(t),(21)maxyFxt,y,s.t.Gxt,y≥0,Hxt,y=0,y∈D⊆R+n×1,where x(t) is the solution to MP in tth iteration. The Lagrange function to problem (14) is L(y,ν,υ)=F(x(t),y)+νG(x(t),y)+υH(x(t),y). If optimal solution y∗ to SP exists, it can be derived by solving KKT point (y(u),ν(u),μ(u)), and due to the strictly convex property of SP, L(y(u),ν(u),μ(u))=max(y)F(x(t),y). After getting (y(u),ν(u),μ(u)), an optimality cut Γ(x,y(u),ν(u),υ(u))≥η can be generated and added to MP, which actually constructs an upper bound for MP. If the optimal solution to SP does not exist, a feasibility cut should be added to cut off the infeasible x(v)=x(t). Specifically, the feasible cut is ∑i∈Yv‍xi-∑i∈Nv‍xi≤Yv-1, derived from x≠x(v), wherein Yv is the set of indexes indicating the nonzero elements in x(v) and Nv is that indicating the zero ones.Therefore, the MP intth iteration is as (22)maxx,ηη,s.t.Γx,yu,νu,υu≥η,∀u=1,2,…,t1,∑i,l∈Yv-∑i,l∈Nv≤Yv-1,∀ν=1,2,…,t2,Gx≥0,Hx=0,x∈0,1m×1,where t1+t2=t. Based on above problem decomposition framework, we propose the following ASRAO algorithm to solve P1.


### 3.3. ASRAO Algorithm

As is depicted in BD framework, problem P1 can be decomposed into MP and SP as in (26) and (24), separately. The initial MP to be solved is(23)maxXη,s.t.∑l∈J∪Kxil=1,∑i∈Ixik≤Nm,xil∈0,1,i∈I,l∈J∪K.In the tth iteration, denote the solution to MP as (X(t),η(t)), and the SP has the form(24)maxΛΦΛ,Xt,s.t.Λ,Xt=∑i∑lxiltlog∑lRilλil+βl,λεxilt≤λil≤λmxilt,λik=λikavrxik,∑i∈Iλij≤1,λij∈R+.If optimal solution exists in problem (24), the KKT point is (Λ(u),ν(u),ω(u),γ(u),μ(u))=(ΛJ(t),ν(t),ω(t),γ(t),μ(t)), and an optimal cut Γ(X,Λ(u),ν(u),ω(u),γ(u),μ(u))≥η should be added to MP, wherein(25)ΓX,Λu,νu,ωu,γu,μu=∑i∑lxillog∑lRilλilu+βl-∑lνlu∑iλilu-1-∑i∑lωiluλilu-λmxil+∑i∑lγiluλilu-λεxil-∑i∑kuμikuλiku-xikλkavru.Due to strong duality of SP in (24), we have Φ(Λ(u),X(t))=Γ(X(t),Λ(u),ν(u),ω(u),γ(u),μ(u)). If an optimal solution does not exist, a feasibility cut ∑(i,l)∈Yv‍xil-∑i,l∈Nv‍xil≤Yv-1 is added to MP to cut off the infeasible X(v)=X(t), where Yv={(i,l)∣xil=1,i∈I,l∈J∪K} and Nv={(i,l)∣xil=0,i∈I,l∈J∪K}. And the derived MP in tth iteration is as (26) and t1+t2=t.(26)maxXJ,XK,ηη,s.t.ΓX,Λu,νu,ωu,γu,μu≥η,∀u=1,…,t1,∑i,l∈Yvxil-∑i,l∈Nvxil≤Yv-1,∀v=1,…,t2,∑l∈I∪Jxil=1,∑i∈Ixik≤Nm,xil∈0,1,i∈I,l∈I∪J.Denote the lower bound and upper bound of problem P1 asLB(t) and UB(t), respectively. They can be derived from problem (24) and (26), which will be proved later in Lemma 1. The iteration procedure terminates when the gap between UB(t) and LB(t) becomes zero.The pseudocodes of the above procedures are as shown in Algorithm1.Algorithm 1:The ASRAO algorithm.
(1)Initialize:  UB(0)=+∞, LB(0)=-∞, Let t=0, u=0, v=0. (2)repeat (3)  Lett=t+1 and solve (26) to obtain current optimal solution (Xt,ηt),  UB(t)=η(t) (4)if (24) with X(t) is bounded: u=u+1, and solve it with MWF to get a KKT point (Λ(u),ν(u),ω(u),γ(u),μ(u))= (Λ(t),ν(t),ω(t),γ(t),μ(t)). The lower bound is set to LB(t)=max{max1≤s≤t{Φ(Λ(s),X(s))},LB(0)} and addΓ(X,Λu,νu,ωu,γu,μu)≥η to (26) (5)elseif (24) with Xt is infeasible: v=v+1,   X(v)=X(t)UB(t)=UB(t-1), and add ∑(i,l)∈Yv‍xil-∑(i,l)∈Nv‍xil≤|Yv|-1 to (24). (6)else as (24) is unbounded, P1 is unbounded, ASRAO stops and return an indication of infeasible problem. (7)until  UB(t)-LB(t)≤0 (8) Return(X∗,Λ∗)=(Xt,Λt) as the optimal solution to problem P1.Lemma 1.
Intth iteration, denote the values of goal functions in problem (24) and (26) as L(t) and U(t), and then LBt=max0≤s≤t{Lt}=max0≤s≤t{Φ(Λs,Xs)} and UB(t)=Lt=η(t).Proof.
First, we proveUB(t)=U(t). Due to the strictly convex property, P1 has the strong duality, by which we have P1 which is equivalent to(27)maxXminν,ω,γ,μ⁡supΛ⁡ΓXJ,XK,ΛJ,ν,s.t.∑l∈I∪Jxil=1,∑i∈Ixik≤Nm,xil∈0,1,i∈I,l∈I∪J.Therefore, problem P1 and (27) have the same optimal solution and objective value of goal function. Denote the optimal solution of P1 as (X∗,Λ∗) and the corresponding goal function value as Φ(Λ∗,X∗)=(∑i‍∑l‍(xil∗log(∑l‍Rilλil∗+βl)))=M∗. In tth iteration, the optimal solution to problem (26) is (X(t),η(t)). Since (26) is a relaxation of (27), it follows that Φ(Λ∗,X∗)≤U(t)=η(t). Therefore, U(t) is a upper bound of the objective function in P1.
Next, we proveLB(t)=max0≤s≤t{L(s)} and L(s)=Φ(Λs,Xs) is an upper bound of P1. Whether L(t) is bounded is decided by (24). If L(s)=-∞ (∀0≤s≤t), LB(t)=max1≤s≤t{Φ(Λ(s),X(s))}=-∞, which is obviously the lower bound of P1, so we focus on the case in which problem (24) is bounded. If (24) is bounded, LB(t)>-∞. We denoted s_=argmax1≤s≤t{Φ(Λ(s),X(s))}, and thereby LB(t)=L(s_)=Φ(Λ(s_),X(s_)). If LB(t) is not the lower bound of P1, LB(t)>M∗. Due to the strong duality of problem (24), L(s_)=Φ(Λ(s_),X(s_))>M∗=Φ(Λ∗,X∗). This means (Λ(s_),X(s_)) generates a larger objective function value of P1 than that (Λ∗,X∗) does, which is in contradiction to the fact that (Λ∗,X∗) is the optimal solution of P1. Therefore, LB(t)=max0≥s≥t{L(s)} is a lower bound of P1.Theorem 2.
The ASRAO algorithm converges to a global optimal solution to P1 with finite number of iterations.Proof.
As is indicated in [19, 20], there always exists a pointed convex polyhedral cone C that makes problem (24) and (26) equivalent to(28)maxX,ηη,s.t.X,η∈C,∑l∈J∪Kxil=1,∑i∈Ixik≤Nm,xil∈0,1,∀i∈I,∀k∈K,∀l∈J∪K,
whereinC can be expressed as a convex hull of finitely many half-lines. In each iteration of ASRAO the set of solution X,η to problem (26) is shrunk by introducing a new extreme half-line being either an optimality cut Γ(X,Λu,νu,ωu,γu,μu)≥η or a feasibility cut ∑i,l∈Yv‍xil-∑i,l∈Nv‍xil≤Yv-1. Because C is a convex hull of finally many extreme half-lines and the new introduced extreme half-line is different from the preceding ones, the complete set of constraints determining the C can be obtained within a finite number of iterations. This means the optimal UE-AN associations X∗ and solution to problem (24) can be obtained within finite iterations by solving(29)ΦΛu,Xt=ΓXt,Λu,νu,ωu,γu,μu.Therefore, the optimal UE-AN associations X∗ and solution to problem P1 can also be obtained within a finite number of iterations. As is analyzed above, the sequence {UB(t)} is nonincreasing, and after the set C is determined by constraints in (29), UB(k)=M∗, where M∗ denotes the optimal objective function value. The sequence {LB(t)} is nondecreasing and satisfies LB(t)=M∗ after X∗ is found through solving (29). Therefore, we can claim that UB(t)-LB(t)=0 will guarantee that the ASRAO converges to the optimal solution to problem P1 within a finite number of iterations.


### 3.4. Accelerated ASRAO Algorithm

Problem (26) is an MILP, which is generally calculated by branch and bound approach. Thus, the computation of problem (24) dominates the computation complexity of the whole BD process. In order to reduce the cost of solving problem (26), in this paper we propose to accelerate ASRAO algorithm by relaxing integer constraints in (26) into continuous constraints in intermediate iterations. Specifically, intth iteration, the MP can be relaxed into a LP problem, as is in (30)maxXJ,XK,ηη,s.t.ΓX,Λu,νu,ωu,γu,μu≥η,∀u=1,…,t1,∑i,l∈Yvxil-∑i,l∈Nvxil≤Yv-1,∀v=1,…,t2,∑l∈I∪Jxil=1,∑i∈Ixik≤Nm,xil∈0,1,i∈I,k∈K,l∈I∪J,where xil is a relaxed continuous variable ranged in [0,1]. As problem (30) is a LP problem, it can be solved utilizing standard algorithms, simplex method for instance. Denote the solution derived from problem (30) as (X^,η^), in which the entries in X^ are possibly not integers. Theorem 3 addresses the fact that relaxation of problem (26) will not exclude the optimal solution to P1.Theorem 3.
Denote(X^,η^) as an arbitrary feasible solution to the relaxed problem. The optimality cuts and feasibility cuts generated with X^ will not exclude the optimal solution (Λ∗,X∗) from problem P1.Proof.
If problem (24) is bounded, denote (Λ^,ν^) as the optimal solution to problem (24) with X^ and an optimality cut Γ(X,Λ^J,ν^,ω^,γ^,μ^)≥η is generated. Otherwise, a feasibility cut ∑i,l∈Yv‍xil-∑i,l∈Nv‍xil≤Yv-1 is generated. To prove Theorem 3, we show that, in either cases, optimal solution (Λ∗,X∗) to problem P1 will not violate the newly introduce constraints. Hence, it will not be excluded from the feasible set by the generated cuts.
Denote the optimal solution to problem (24) as (Λ∗,ν∗,ω∗,γ∗,μ∗) with X∗; the corresponding optimal value of the objective function is M∗=η∗=Φ(Λ∗,X∗). In the case where problem (24) is bounded with X being X^, suppose X∗,η∗ violates Γ(X,Λ^,ν^,ω^,γ^,μ^)≥η, that is Γ(X∗,Λ^,ν^,ω^,γ^,μ^)<η∗, due to the completeness of the constraints in (26), it means a smaller η′=Γ(X∗,Λ^,ν^,ω^,γ^,μ^) is found, which violates the fact that M∗=η∗ is the optimal solution to problem P1. Hence, the optimality cut will not be violated. In the case where the problem (24) is unbounded, suppose X∗ violates the feasibility cut ∑i,l∈Yv‍xil-∑i,l∈Nv‍xil≤Yv-1, which means X∗ makes problem (24) unbound, and therefore, X∗ will not be the optimal solution to P2. This conflicts the fact that (X∗,η∗) is the optimal solution to P1. So the feasibility cut will not be violated.
Theorem3 means we can safely replace (26) by (30) in our ASRAO algorithm without concerning about losing the optimal solution of P1. However, it should be noted that X^ in the optimal solution to problem (30) can be used to calculate the lower bound of P1, but when problem (24) is solved with X^, the upper bound derived from the objective function value of problem (24) may not be a valid upper bound for P1, since entries in X^ are not integers and thereby infeasible to constraints in P1. One way to cope with the problem is to round X^ to the nearest matrix X~ with binary entries. However, as the summation of elements in each row of X equals 1, a direct approximation would produce an all-zero matrix. Hence, in view of the special property of X that the entries in the matrix are binary variables and only one element in each row equals 1, we use a technique inspired by feasible pump method proposed in [21] to get X~.
To select a 1 element inxi, the ith row of X^, we divide a range [0,1] into J+K slots, and each element in xi successively occupies a slot, whose length equals to the value of element xil. Then a randomizer uniformly generates a pointer valued within [0,1]. This element with the slot where the pointer locates is selected, and all other elements are 0 elements. To control randomness, a parameter T is used to decide only 1/T rows are allowed to perform randomly selection and other rows just select their largest elements as 1 element.
Thereby, fromX^ we get an approximate matrix X~ which is possibly feasible to problem (26) and P1, but with very low probability X~K fails the constraint ∑i∈I‍xik≤Nm, and therefore, X~K should be checked. And to avoid the fact that X~ stays at a fixed point, novelty of the solution should be checked. If no feasible solution or solution stays unchanged, we revert back to the normal ASRAO algorithm, and the current solution calculated from problem (30) serves as an initial X~ to the following normal ASRAO process.
The pseudocodes of accelerated ASRAO algorithm are as shown in Algorithm2.Algorithm 2:The A-ASRAO algorithm.
(1)Initialize:  UB(0)=+∞, LB(0)=-∞, Let t=0, u=0, v=0. PHASE-I (2)repeat (3)     Lett=t+1 and solve (30) to obtain current optimal solution (X^(t),η^(t)), UB(t)=η^(t) (4)if  X^(t)=X^(t-1) or feasible X~(t) can not be found: t=t-1, break (5)else: get approximative matrix X~(t), and X(t)=X~(t) (6)if (24) with Xt is bounded: u=u+1, and solve it with Xt to obtain KKT point Λu,νu,ωu,γu,μu =(Λ(t),ν(t),ω(u),γ(u),μ(t)). The lower bound is set to LB(t)=max{max1≤s≤t{Φ(Λ(s),X(s))},LB(0)}. AddΓ(X,Λu,νu,ωu,γu,μu)≥η to (26) and (30) (7)else: v=v+1,   X(v)=X(t), LB(t)=LB(t-1), and add ∑(i,l)∈Yv‍xil-∑(i,l)∈Nv‍xil≤|Yv|-1 to (26) and (30). (8)until  UB(t)-LB(t)≤0 PHASE-II (9)   ASRAO step (2) to (7). (10) Return(X∗,Λ∗)=(Xt,Λt) as the optimal solution to problem P1.


## 3.1. Problem Analyses

Assume association indicatorX is fixed, then each UE i is served by a specific AN l, and each AN l has a fixed UE collection Il={i∣xil=1}. UE i contributes log(Rilλil+βl) to utility of AN l. Its value is determined by λil, which is only related to AN l. So resource allocation in each AN can be solved independently with given user association X˙. Therefore, Φ(Λ,X˙) can be identically transformed to(12)ΦΛ,X˙=∑i∈I,l≡argxil=1logRilλil+βl=∑l∈J∪K∑i∈IllogRilλil+βl=∑l∈J∪KΦlλl,x˙l,where Φl(λl,x˙l)=∑i∈Il‍logRilλil+βl is the utility of AN l. So we get a subproblem of P1 on condition of X˙.(13)maxΛ∑i∈IllogRilλil+bls.t.∑iλil≤1,λεx˙il≤λil≤λmx˙il,λik=λkavrx˙ik.This is a classical problem in information theory and a water filling (WF) principle is specially developed to solve it, which is derived from Karush-Kuhn-Tucker (KKT) optimal conditions. It provides direction solutions to problems with the form like (13). To enable its applicability in P1, we extend WF to a multi-AN case and denote its results as MWF, short for multi-AN water filling. With X fixed as X˙ in P1, we have a MWF applicable problem as in (14)maxΛΦΛ,X˙,s.t.ΦΛ,X˙=∑i∑lx˙illog∑lRilλil+βl,λεx˙il≤λil≤λmx˙il,∑iλil≤1,λik=λkavrxik.According to generalized duality theory, the dual function ofΦ(Λ,X˙) is as in (15)LΛ,ν,ω,γ,X˙=∑i∑lx˙illog∑lRilλil+βl-∑lνl∑iλil-1-∑i∑lωilλil-x˙il+∑i∑lγilλil-εx˙il-∑i∑kμikλik-x˙ikλkavr,where ν, ω, γ are nonnegative relaxed variables and μ is nonzero. The local optimal is reached, when the following KKT conditions are satisfied.(16)∂ΦΛ∂λij=Rij∑jRijλij+∑lx˙ijbj-νj-ωij+γij=0,∂ΦΛ∂λik=Rik∑kRikλik+∑lx˙ikbk-νk-ωik+γik-μik=0,∂ΦΛ∂λkavr=μikx˙ik=0,ωilλil-x˙il=0,νl∑iλil-1=0,γilλil-εx˙il=0,μikλik-x˙ikλkavr=0,νl,ωil,γil∈0,+∞,μik≠0.The solutions to (16) are as follows. MWF results for CN is as in (17)νj∗=∑ix¯˙ijρ¯+∑ibjx¯˙ijRij-1,λij∗=maxx˙ij1νj∗-βjRij,0+λεx˙ij,γij=-RijRij′λij′+βj′-νj∗,RijRij′λij′+βj′≤νj∗,0,else,ωij=RijRij′λij′+βj′-νj∗+γij.The following (18) is MWF results for WLAN.(18)λkavr=minλm,1∑ix˙ik,∑ix˙ik>0,0∑ix˙ik=0,λik∗=λkavrx˙ik,νk∗=minRikβk+∑ix˙ik,λm∑ix˙ik≥1,0,λm∑ix˙ik<1,μik=-RikRik′λik′+βk′-νk∗,λm∑ix˙ik>1,1,λm∑ix˙ik≤1,γik=-RikRik′λik′+βk′-νk∗,RikRik′λik′+βk′<νk∗,0,else,ωik=RikRik′λik′+βk′-νk∗+γik+μik.Using above MWF results, we can directly calculate optimalΛ with a given X˙. To utilize MWF in solving problem P1, an optimal X should be retrieved. However, as X is a matrix of binary variable entries and Λ is closely related with X, it is difficult to find such an optimal. In this paper, instead of directly solving an optimal X, we devise ASRAO algorithm to solve it iteratively and effectively. ASRAO is originated from BD, which provides a framework for addressing mixed integer programming problems by decomposing it into two smaller subproblems and solving them iteratively. The details are in the sequel.


## 3.2. Generalized Benders Decomposition

BD is originally proposed to solve mixed integer linear programming (MILP). Instead of solving all variables and constraints simultaneously, it decomposes an MILP into two subproblems, a master problem (MP) and a slave problem (SP), and the original MILP can be solved by solving MP and SP iteratively. In BD frameworks, MP is an MILP problem which consists of all integer variables in the original problem, while SP is a Linear Programming problem with all the continuous variables from the original problem. In each iteration, BD utilizes an extreme point or extreme ray derived from the dual of SP to generate an optimality cut or feasibility cut to trim the feasible domain of MP, and MP will eventually reach the same optima as the original MILP.GBD extends the BD approach to a more general class of problems by adopting nonlinear duality theory, and some nonlinear problems are brought into range. When solving MINLP, GBD follows the framework of BD, but SP is derived as NLP. To generate the feasibility cuts and optimality cuts, GBD utilized KKT conditions to calculate an extreme point or extreme ray of the Lagrange dual functions of NLP SP. A GBD applicable problem has a general form as in(19)maxx,yFx,y,s.t.Gx,y≥0,Hx,y=0,x∈0,1m×1,y∈D⊆R+n×1.It should meet following situations: (a) for fixedx, (19) separates into a number of independent optimization problems; (b) for fixed x, (19) assumes a well-known special structure that efficient solution procedures are available; (c) (19) may be not concave program in x and y jointly, but fixing x renders it so in y. According to the analyses, UTMAX obviously satisfies these conditions.The initial MP of (19) in GBD procedures is as follows:(20)maxx,ηη,s.t.Gx≥0,Hx=0,η∈R,x∈E⊆0,1m×1,where G(x)≥0 and H(x)=0 are constraints from (19) that only related with x. Solving a solution (x(t),η) from MP, SP can be derived by simply fixing x in (19) as x(t),(21)maxyFxt,y,s.t.Gxt,y≥0,Hxt,y=0,y∈D⊆R+n×1,where x(t) is the solution to MP in tth iteration. The Lagrange function to problem (14) is L(y,ν,υ)=F(x(t),y)+νG(x(t),y)+υH(x(t),y). If optimal solution y∗ to SP exists, it can be derived by solving KKT point (y(u),ν(u),μ(u)), and due to the strictly convex property of SP, L(y(u),ν(u),μ(u))=max(y)F(x(t),y). After getting (y(u),ν(u),μ(u)), an optimality cut Γ(x,y(u),ν(u),υ(u))≥η can be generated and added to MP, which actually constructs an upper bound for MP. If the optimal solution to SP does not exist, a feasibility cut should be added to cut off the infeasible x(v)=x(t). Specifically, the feasible cut is ∑i∈Yv‍xi-∑i∈Nv‍xi≤Yv-1, derived from x≠x(v), wherein Yv is the set of indexes indicating the nonzero elements in x(v) and Nv is that indicating the zero ones.Therefore, the MP intth iteration is as (22)maxx,ηη,s.t.Γx,yu,νu,υu≥η,∀u=1,2,…,t1,∑i,l∈Yv-∑i,l∈Nv≤Yv-1,∀ν=1,2,…,t2,Gx≥0,Hx=0,x∈0,1m×1,where t1+t2=t. Based on above problem decomposition framework, we propose the following ASRAO algorithm to solve P1.


## 3.3. ASRAO Algorithm

As is depicted in BD framework, problem P1 can be decomposed into MP and SP as in (26) and (24), separately. The initial MP to be solved is(23)maxXη,s.t.∑l∈J∪Kxil=1,∑i∈Ixik≤Nm,xil∈0,1,i∈I,l∈J∪K.In the tth iteration, denote the solution to MP as (X(t),η(t)), and the SP has the form(24)maxΛΦΛ,Xt,s.t.Λ,Xt=∑i∑lxiltlog∑lRilλil+βl,λεxilt≤λil≤λmxilt,λik=λikavrxik,∑i∈Iλij≤1,λij∈R+.If optimal solution exists in problem (24), the KKT point is (Λ(u),ν(u),ω(u),γ(u),μ(u))=(ΛJ(t),ν(t),ω(t),γ(t),μ(t)), and an optimal cut Γ(X,Λ(u),ν(u),ω(u),γ(u),μ(u))≥η should be added to MP, wherein(25)ΓX,Λu,νu,ωu,γu,μu=∑i∑lxillog∑lRilλilu+βl-∑lνlu∑iλilu-1-∑i∑lωiluλilu-λmxil+∑i∑lγiluλilu-λεxil-∑i∑kuμikuλiku-xikλkavru.Due to strong duality of SP in (24), we have Φ(Λ(u),X(t))=Γ(X(t),Λ(u),ν(u),ω(u),γ(u),μ(u)). If an optimal solution does not exist, a feasibility cut ∑(i,l)∈Yv‍xil-∑i,l∈Nv‍xil≤Yv-1 is added to MP to cut off the infeasible X(v)=X(t), where Yv={(i,l)∣xil=1,i∈I,l∈J∪K} and Nv={(i,l)∣xil=0,i∈I,l∈J∪K}. And the derived MP in tth iteration is as (26) and t1+t2=t.(26)maxXJ,XK,ηη,s.t.ΓX,Λu,νu,ωu,γu,μu≥η,∀u=1,…,t1,∑i,l∈Yvxil-∑i,l∈Nvxil≤Yv-1,∀v=1,…,t2,∑l∈I∪Jxil=1,∑i∈Ixik≤Nm,xil∈0,1,i∈I,l∈I∪J.Denote the lower bound and upper bound of problem P1 asLB(t) and UB(t), respectively. They can be derived from problem (24) and (26), which will be proved later in Lemma 1. The iteration procedure terminates when the gap between UB(t) and LB(t) becomes zero.The pseudocodes of the above procedures are as shown in Algorithm1.Algorithm 1:The ASRAO algorithm.
(1)Initialize:  UB(0)=+∞, LB(0)=-∞, Let t=0, u=0, v=0. (2)repeat (3)  Lett=t+1 and solve (26) to obtain current optimal solution (Xt,ηt),  UB(t)=η(t) (4)if (24) with X(t) is bounded: u=u+1, and solve it with MWF to get a KKT point (Λ(u),ν(u),ω(u),γ(u),μ(u))= (Λ(t),ν(t),ω(t),γ(t),μ(t)). The lower bound is set to LB(t)=max{max1≤s≤t{Φ(Λ(s),X(s))},LB(0)} and addΓ(X,Λu,νu,ωu,γu,μu)≥η to (26) (5)elseif (24) with Xt is infeasible: v=v+1,   X(v)=X(t)UB(t)=UB(t-1), and add ∑(i,l)∈Yv‍xil-∑(i,l)∈Nv‍xil≤|Yv|-1 to (24). (6)else as (24) is unbounded, P1 is unbounded, ASRAO stops and return an indication of infeasible problem. (7)until  UB(t)-LB(t)≤0 (8) Return(X∗,Λ∗)=(Xt,Λt) as the optimal solution to problem P1.Lemma 1.
Intth iteration, denote the values of goal functions in problem (24) and (26) as L(t) and U(t), and then LBt=max0≤s≤t{Lt}=max0≤s≤t{Φ(Λs,Xs)} and UB(t)=Lt=η(t).Proof.
First, we proveUB(t)=U(t). Due to the strictly convex property, P1 has the strong duality, by which we have P1 which is equivalent to(27)maxXminν,ω,γ,μ⁡supΛ⁡ΓXJ,XK,ΛJ,ν,s.t.∑l∈I∪Jxil=1,∑i∈Ixik≤Nm,xil∈0,1,i∈I,l∈I∪J.Therefore, problem P1 and (27) have the same optimal solution and objective value of goal function. Denote the optimal solution of P1 as (X∗,Λ∗) and the corresponding goal function value as Φ(Λ∗,X∗)=(∑i‍∑l‍(xil∗log(∑l‍Rilλil∗+βl)))=M∗. In tth iteration, the optimal solution to problem (26) is (X(t),η(t)). Since (26) is a relaxation of (27), it follows that Φ(Λ∗,X∗)≤U(t)=η(t). Therefore, U(t) is a upper bound of the objective function in P1.
Next, we proveLB(t)=max0≤s≤t{L(s)} and L(s)=Φ(Λs,Xs) is an upper bound of P1. Whether L(t) is bounded is decided by (24). If L(s)=-∞ (∀0≤s≤t), LB(t)=max1≤s≤t{Φ(Λ(s),X(s))}=-∞, which is obviously the lower bound of P1, so we focus on the case in which problem (24) is bounded. If (24) is bounded, LB(t)>-∞. We denoted s_=argmax1≤s≤t{Φ(Λ(s),X(s))}, and thereby LB(t)=L(s_)=Φ(Λ(s_),X(s_)). If LB(t) is not the lower bound of P1, LB(t)>M∗. Due to the strong duality of problem (24), L(s_)=Φ(Λ(s_),X(s_))>M∗=Φ(Λ∗,X∗). This means (Λ(s_),X(s_)) generates a larger objective function value of P1 than that (Λ∗,X∗) does, which is in contradiction to the fact that (Λ∗,X∗) is the optimal solution of P1. Therefore, LB(t)=max0≥s≥t{L(s)} is a lower bound of P1.Theorem 2.
The ASRAO algorithm converges to a global optimal solution to P1 with finite number of iterations.Proof.
As is indicated in [19, 20], there always exists a pointed convex polyhedral cone C that makes problem (24) and (26) equivalent to(28)maxX,ηη,s.t.X,η∈C,∑l∈J∪Kxil=1,∑i∈Ixik≤Nm,xil∈0,1,∀i∈I,∀k∈K,∀l∈J∪K,
whereinC can be expressed as a convex hull of finitely many half-lines. In each iteration of ASRAO the set of solution X,η to problem (26) is shrunk by introducing a new extreme half-line being either an optimality cut Γ(X,Λu,νu,ωu,γu,μu)≥η or a feasibility cut ∑i,l∈Yv‍xil-∑i,l∈Nv‍xil≤Yv-1. Because C is a convex hull of finally many extreme half-lines and the new introduced extreme half-line is different from the preceding ones, the complete set of constraints determining the C can be obtained within a finite number of iterations. This means the optimal UE-AN associations X∗ and solution to problem (24) can be obtained within finite iterations by solving(29)ΦΛu,Xt=ΓXt,Λu,νu,ωu,γu,μu.Therefore, the optimal UE-AN associations X∗ and solution to problem P1 can also be obtained within a finite number of iterations. As is analyzed above, the sequence {UB(t)} is nonincreasing, and after the set C is determined by constraints in (29), UB(k)=M∗, where M∗ denotes the optimal objective function value. The sequence {LB(t)} is nondecreasing and satisfies LB(t)=M∗ after X∗ is found through solving (29). Therefore, we can claim that UB(t)-LB(t)=0 will guarantee that the ASRAO converges to the optimal solution to problem P1 within a finite number of iterations.


## 3.4. Accelerated ASRAO Algorithm

Problem (26) is an MILP, which is generally calculated by branch and bound approach. Thus, the computation of problem (24) dominates the computation complexity of the whole BD process. In order to reduce the cost of solving problem (26), in this paper we propose to accelerate ASRAO algorithm by relaxing integer constraints in (26) into continuous constraints in intermediate iterations. Specifically, intth iteration, the MP can be relaxed into a LP problem, as is in (30)maxXJ,XK,ηη,s.t.ΓX,Λu,νu,ωu,γu,μu≥η,∀u=1,…,t1,∑i,l∈Yvxil-∑i,l∈Nvxil≤Yv-1,∀v=1,…,t2,∑l∈I∪Jxil=1,∑i∈Ixik≤Nm,xil∈0,1,i∈I,k∈K,l∈I∪J,where xil is a relaxed continuous variable ranged in [0,1]. As problem (30) is a LP problem, it can be solved utilizing standard algorithms, simplex method for instance. Denote the solution derived from problem (30) as (X^,η^), in which the entries in X^ are possibly not integers. Theorem 3 addresses the fact that relaxation of problem (26) will not exclude the optimal solution to P1.Theorem 3.
Denote(X^,η^) as an arbitrary feasible solution to the relaxed problem. The optimality cuts and feasibility cuts generated with X^ will not exclude the optimal solution (Λ∗,X∗) from problem P1.Proof.
If problem (24) is bounded, denote (Λ^,ν^) as the optimal solution to problem (24) with X^ and an optimality cut Γ(X,Λ^J,ν^,ω^,γ^,μ^)≥η is generated. Otherwise, a feasibility cut ∑i,l∈Yv‍xil-∑i,l∈Nv‍xil≤Yv-1 is generated. To prove Theorem 3, we show that, in either cases, optimal solution (Λ∗,X∗) to problem P1 will not violate the newly introduce constraints. Hence, it will not be excluded from the feasible set by the generated cuts.
Denote the optimal solution to problem (24) as (Λ∗,ν∗,ω∗,γ∗,μ∗) with X∗; the corresponding optimal value of the objective function is M∗=η∗=Φ(Λ∗,X∗). In the case where problem (24) is bounded with X being X^, suppose X∗,η∗ violates Γ(X,Λ^,ν^,ω^,γ^,μ^)≥η, that is Γ(X∗,Λ^,ν^,ω^,γ^,μ^)<η∗, due to the completeness of the constraints in (26), it means a smaller η′=Γ(X∗,Λ^,ν^,ω^,γ^,μ^) is found, which violates the fact that M∗=η∗ is the optimal solution to problem P1. Hence, the optimality cut will not be violated. In the case where the problem (24) is unbounded, suppose X∗ violates the feasibility cut ∑i,l∈Yv‍xil-∑i,l∈Nv‍xil≤Yv-1, which means X∗ makes problem (24) unbound, and therefore, X∗ will not be the optimal solution to P2. This conflicts the fact that (X∗,η∗) is the optimal solution to P1. So the feasibility cut will not be violated.
Theorem3 means we can safely replace (26) by (30) in our ASRAO algorithm without concerning about losing the optimal solution of P1. However, it should be noted that X^ in the optimal solution to problem (30) can be used to calculate the lower bound of P1, but when problem (24) is solved with X^, the upper bound derived from the objective function value of problem (24) may not be a valid upper bound for P1, since entries in X^ are not integers and thereby infeasible to constraints in P1. One way to cope with the problem is to round X^ to the nearest matrix X~ with binary entries. However, as the summation of elements in each row of X equals 1, a direct approximation would produce an all-zero matrix. Hence, in view of the special property of X that the entries in the matrix are binary variables and only one element in each row equals 1, we use a technique inspired by feasible pump method proposed in [21] to get X~.
To select a 1 element inxi, the ith row of X^, we divide a range [0,1] into J+K slots, and each element in xi successively occupies a slot, whose length equals to the value of element xil. Then a randomizer uniformly generates a pointer valued within [0,1]. This element with the slot where the pointer locates is selected, and all other elements are 0 elements. To control randomness, a parameter T is used to decide only 1/T rows are allowed to perform randomly selection and other rows just select their largest elements as 1 element.
Thereby, fromX^ we get an approximate matrix X~ which is possibly feasible to problem (26) and P1, but with very low probability X~K fails the constraint ∑i∈I‍xik≤Nm, and therefore, X~K should be checked. And to avoid the fact that X~ stays at a fixed point, novelty of the solution should be checked. If no feasible solution or solution stays unchanged, we revert back to the normal ASRAO algorithm, and the current solution calculated from problem (30) serves as an initial X~ to the following normal ASRAO process.
The pseudocodes of accelerated ASRAO algorithm are as shown in Algorithm2.Algorithm 2:The A-ASRAO algorithm.
(1)Initialize:  UB(0)=+∞, LB(0)=-∞, Let t=0, u=0, v=0. PHASE-I (2)repeat (3)     Lett=t+1 and solve (30) to obtain current optimal solution (X^(t),η^(t)), UB(t)=η^(t) (4)if  X^(t)=X^(t-1) or feasible X~(t) can not be found: t=t-1, break (5)else: get approximative matrix X~(t), and X(t)=X~(t) (6)if (24) with Xt is bounded: u=u+1, and solve it with Xt to obtain KKT point Λu,νu,ωu,γu,μu =(Λ(t),ν(t),ω(u),γ(u),μ(t)). The lower bound is set to LB(t)=max{max1≤s≤t{Φ(Λ(s),X(s))},LB(0)}. AddΓ(X,Λu,νu,ωu,γu,μu)≥η to (26) and (30) (7)else: v=v+1,   X(v)=X(t), LB(t)=LB(t-1), and add ∑(i,l)∈Yv‍xil-∑(i,l)∈Nv‍xil≤|Yv|-1 to (26) and (30). (8)until  UB(t)-LB(t)≤0 PHASE-II (9)   ASRAO step (2) to (7). (10) Return(X∗,Λ∗)=(Xt,Λt) as the optimal solution to problem P1.


## 4. Performance Evaluation

In this section, we first introduce network scenarios and parameter settings used in our simulations. After this, the performance of ASRAO is evaluated in scenarios with different AP and UE settings. The effectiveness of the accelerative technique is also discussed here. Then, we analyze effects the value ofβ may cause to our UTMAX model by evaluating the optimized network performances on metrics, such as system throughput and access fairness. In the last part of this section, our UTMAX model is evaluated under different AP densities to validate its effectiveness and compatibility.


### 4.1. Scenario and Parameter Settings

The simulations are carried out in a scenario consists of 2 fixed BSs and 6 or 12 APs. The distance between the two BSs is about 300 m. The number of UEs varies from 10 to 20 by a step size of 2. In order to reduce random errors, each point of the simulation results is averaged over 100 different network topologies among which the positions of APs and UEs differ but their amounts stay unchanged. According to [15], we simplify the propagation loss as 37+37.6logdkm+21log(2.6/2)dB in CN, and 34+37.6logdkm+21log(2.4/2)dB in WLAN. To model the shadowing effects in CN a log-normal random variable with zero mean and 8 db variation is introduce in CN, while in WLAN shadowing effects are neglected. Moreover, a received power maximization (PRMAX) user association scheme, which is commonly used in present commercial networks, serves as a performance benchmark. Another scheme maximally uses WLAN by getting AP to serve each possible user in its range. We denote it as APMAX. Our UTMAX optimization model with different β values for different system preferences, such as system throughput and user data rate fairness, is compared with PRMAX and APMAX to validate its performance. The main parameters and their values are as listed in Table 1, from [7, 15].Table 1
Simulation parameter.
Items (unit) Values Transmitting power of BSs (dBm) 46 Transmitting power of WiFi AP (dBm) 20 Distance between BSs (m) 300 LTE bandwidth (MHz) 10 WiFi bandwidth (MHz) 20 eNB antenna height from eNB (m) 40 AP antenna height from AP (m) 5 LTE operating frequency (MHz) 2600 WiFi operating frequency (MHz) 2400 Number of eNode BSs and APs 2 : 6, 2 : 12 Thermal noise (dBm/Hz) −174 Shadowing for LTE (dB) 8 LTE cochannel interference coefficientα 1 Number of UEs 10 : 2 : 20 Min resource for UE in eNBλε 0.01 Max resource for UE in APλm 0.3 Fairness coefficientβ −0.05, 0, 1, 5


### 4.2. Algorithm Performance Analyses

In this subsection, we validate the convergence of ASRAO and A-ASRAO in simulations and analyze the effects the value of a randomness control parameterT produces on the convergence of the latter. The simulations are carried out in Matlab under the settings mentioned above, and our algorithms are implemented with YALMIP [22], an optimization toolbox, which is able to solve our UTMAX model directly, but at unsatisfying speed. We denote it by B&B, because an MINLP solver based on branch and bound is used. It works differently in each detailed iteration with (A-)ASRAO, so only the computation time spent by these approaches is compared.Figure2 illustrates the convergent processes of ASRAO and A-ASRAO in a UTMAX instance where 20 UEs exist. In essence, ASRAO performs fixed path searching, whose interim solutions are decided by their preceding solutions. ASRAO is sure to converge, but it lacks the ability to try better searching paths. While A-ASRAO adopts a random searching strategy, with wider searching directions, it shows more steady descending upper bound. From Figure 2, it is obvious that A-ASRAO converges much faster in the later part of the iterations than ASRAO as is depicted by lines with hollow marks. As these upper bounds reach the same optimal objective function value and all lower bounds get to optima much earlier than upper bounds, A-ASRAO will converge much faster than ASRAO.Figure 2
The converge process of ASRAO and A-ASRAO withT=(1,3,9).In addition, A randomness control parameterT is set in A-ASRAO to control its behavior, and A-ASRAO with T is denoted by A-ASRAO(T). With different values of T, A-ASRAO has different converging rates. This is probably caused by random generation of approximation matrix X~. Smaller T means more rows in X~ perform random 1-element selection, while large T will make less rows perform random selection. Extremely small or large values of T will make too many or few rows perform random selection, but neither of them is good for finding ideal initial solution for Phase II in A-ASRAO. So A-ASRAO(3) performs better than A-ASRAO(1) and A-ASRAO(9).In Figures3 and 4, respectively, the number of iterations performed and computation time spent by ASRAO and A-ASRAO are plotted against the number of UEs ranging from 10 to 20. Each point is obtained by averaging over 100 different network topologies. It shows that ASRAO needs a little fewer iterations than A-ASRAO when UE number is less than 14. However, when UE number is larger, A-ASRAO reduces the number of iterations by up to 15% compared with ASRAO in average.Figure 3
The number of iterations versus the number of UEs.Figure 4
Computational time versus the number of UEs.The reason for the different convergent properties lies in the scale of a UTMAX instance. A-ASRAO works in two phases and in Phase II A-ASRAO works the same as ASRAO, which implies that A-ASRAO in Phase I works actually to obtain an initial solution to Phase II in A-ASRAO. That is to say, if a better solution, which produces smaller gap between the upper bound and lower bound, found by the end of Phase I in A-ASRAO than that by ASRAO within the same iterations, A-ASRAO will have much larger probability to converge within fewer iterations than ASRAO. Therefore, the dominant factor relating to the different convergent rates of ASRAO and A-ASRAO is whether a better solution would be found by the end of Phase I in A-ASRAO.When a UTMAX instance is small-scale problem, its optimal solution can be derived within a few iterations even through exhaustive search, not alone ASRAO. While A-ASRAO tries a more steady but random searching, the finding of optimal solution may thus be postponed. So it hardly performs better than ASRAO. When a UTMAX instance is of large scale, a large solution space means a bolder searching practice may produce better searching path. So A-ASRAO(3) performs better than ASRAO and A-ASRAO(9). Nevertheless, too strong randomness goes against the finding of optimal solution, like in A-ASRAO(1) where all rows inX~ perform random 1-element selection. We speculate that is because a too wild searching will waste many iterations wandering around rather than heading towards a converging direction. That is why A-ASRAO(3) outperforms A-ASRAO(1). Moreover, in Figure 4 it shows that A-ASRAO spends less time computing optimal solution than ASRAO, and the performance gap is even larger than iteration numbers. This benefits from relaxing the MILP master problem into LP in Phase I of A-ASRAO, and LP can be solved much faster than MILP. Figure 4 also shows that ASRAO computes UTMAX much quicker than a common B&B method. Above results prove the effectiveness of (A-)ASRAO in solving UTMAX, and in most cases A-ASRAO(3) performs better than others.


### 4.3. Effect ofβ to UTMAX

As is mentioned in Section2, β indicates UTMAX’s preference on throughput fairness among UEs served in a BS by setting effects on resource allocation. As UE’s attainable data rate is changed, user associations are inevitably affected, which is reflected in Figure 5. And as a consequence, the performance of the entire interworking network is influenced. To analyze these effects, the performance of optimized network is evaluated on throughput produced and the number of UEs served in WiFi tier and cellular tier in Figure 6 and throughput gains of UEs in Figure 7 through APMAX and UTMAX with β in {-0.05,0,1,5}. PRMAX is a default user association scheme and adopts identical resource allocation in each BS and AP.Figure 5
Performance at AN granularity of PRMAX, APMAX, and UTMAX(β).Figure 6
Performance at network granularity of PRMAX, APMAX, and UTMAX(β).Figure 7
Performance at UE granularity of PRMAX, APMAX, and UTMAX(β).Figure5 plots the throughput and number of UEs served in each BS and AP with different approaches marked by different colors and symbols. It shows that, from a AN perspective, both APMAX and UTMAX do well in improving AP’s utilization, and APMAX performs even better, as it tends to maximize AP’s usage. The improvements of WLAN’s offloading performance by different approaches are much more obvious in Figure 6.Figure6 plot the throughput and number of UEs served in cellular tier and WiFi tier. APMAX gets most UEs served in WLAN, while, as resources are fully occupied, each UE in WLAN will get less resources and the throughput of WLAN will not be improved greatly compared with UTMAX. However, due to lack of careful UE selection from CN to WLAN in APMAX, UE amount decreases greatly in CN and the left UEs are probably with poor wireless channel condition or too far from a BS, which causes the poor CN throughput performance in APMAX. As for UTMAX, with different β value, it performs differently. Larger β value will get BS to allocate more resources to UEs with higher unit band data rate, which probably makes the contribution to the throughput growing in CN. However, as UEs with poor channel condition in CN get much less resources, they turn to WLAN to compensate the throughput loss. That would cause a minor UE amount increase in WLAN, while, as AP’s resources are equally shared, new coming UEs with poor channel condition share resources with indigenous WLAN UEs, causing the decrease in WLAN’s throughput. Nevertheless, the whole system throughput shows a continuous growth with the increase of β. It seems that larger β performs better on total throughput and throughput fairness between WLAN and CN, whereas it is on the expense of throughput fairness among UEs.Figure7 depicts UE throughput gain versus probability P(r<α) for various approaches versus PRMAX, where P(r<α) represents the ratio of UE whose throughput is less than α in all UE. The UEs’ throughput gains are quite large at 10% ratio point by all approaches, among which UTMAX(-ε) brings the largest throughput gain. With the growing of β, UE throughput gain decreases. However, UTMAX always has higher throughput gain over APMAX. With the growth of probability P(r<α) UE throughput gain decreases, because PRMAX also has relatively larger UE throughput at larger ratio points. At the meanwhile, UTMAX with larger β value has larger UE throughput gain, as larger β makes UTMAX allocate more resources to UEs with better channel conditions.


### 4.4. Effects of AP Density to UTMAX

The performances of UTMAX with different values ofβ are evaluated on larger AP density, to validate whether the same conclusion holds with different AP settings. The results are in Figures 8 and 9.Figure 8
Performance at network granularity of PRMAX, APMAX, and UTMAX(β).Figure 9
Performance at UE granularity of PRMAX, APMAX, and UTMAX(β).In Figure8, compared with Figure 6 the most obvious changes are the WLAN tier throughput increase for PRMAX and minimization of performance gap between UTMAXs with each of two different values of β. Due to the increase of AP density, more UEs access WLAN and the average distance from a UE to an AP is shorten. This is a direct cause to WLAN tier throughput increase for PRMAX and APMAX. As for the latter difference, it is possibly because of much fewer UEs in each BS or AP, and β probably has weaker effect on resource allocation among fewer UEs in a serving node. Figure 9 is plotted by adding curves of UE throughput gain where 12 APs exist. The additional UE throughput gain is largely contributed by newly added APs. A similar performance relation is shown between UTMAX with larger AP density, while the performance gap is narrowed and the curves meet and overlap earlier than that with smaller AP density.From these results, the effect ofβ to UTMAX is further checked that a smaller β makes UTMAX perform better on improving throughput of the part of UEs with poor channel conditions, while a larger β does better on improving overall system throughput.


## 4.1. Scenario and Parameter Settings

The simulations are carried out in a scenario consists of 2 fixed BSs and 6 or 12 APs. The distance between the two BSs is about 300 m. The number of UEs varies from 10 to 20 by a step size of 2. In order to reduce random errors, each point of the simulation results is averaged over 100 different network topologies among which the positions of APs and UEs differ but their amounts stay unchanged. According to [15], we simplify the propagation loss as 37+37.6logdkm+21log(2.6/2)dB in CN, and 34+37.6logdkm+21log(2.4/2)dB in WLAN. To model the shadowing effects in CN a log-normal random variable with zero mean and 8 db variation is introduce in CN, while in WLAN shadowing effects are neglected. Moreover, a received power maximization (PRMAX) user association scheme, which is commonly used in present commercial networks, serves as a performance benchmark. Another scheme maximally uses WLAN by getting AP to serve each possible user in its range. We denote it as APMAX. Our UTMAX optimization model with different β values for different system preferences, such as system throughput and user data rate fairness, is compared with PRMAX and APMAX to validate its performance. The main parameters and their values are as listed in Table 1, from [7, 15].Table 1
Simulation parameter.
Items (unit) Values Transmitting power of BSs (dBm) 46 Transmitting power of WiFi AP (dBm) 20 Distance between BSs (m) 300 LTE bandwidth (MHz) 10 WiFi bandwidth (MHz) 20 eNB antenna height from eNB (m) 40 AP antenna height from AP (m) 5 LTE operating frequency (MHz) 2600 WiFi operating frequency (MHz) 2400 Number of eNode BSs and APs 2 : 6, 2 : 12 Thermal noise (dBm/Hz) −174 Shadowing for LTE (dB) 8 LTE cochannel interference coefficientα 1 Number of UEs 10 : 2 : 20 Min resource for UE in eNBλε 0.01 Max resource for UE in APλm 0.3 Fairness coefficientβ −0.05, 0, 1, 5


## 4.2. Algorithm Performance Analyses

In this subsection, we validate the convergence of ASRAO and A-ASRAO in simulations and analyze the effects the value of a randomness control parameterT produces on the convergence of the latter. The simulations are carried out in Matlab under the settings mentioned above, and our algorithms are implemented with YALMIP [22], an optimization toolbox, which is able to solve our UTMAX model directly, but at unsatisfying speed. We denote it by B&B, because an MINLP solver based on branch and bound is used. It works differently in each detailed iteration with (A-)ASRAO, so only the computation time spent by these approaches is compared.Figure2 illustrates the convergent processes of ASRAO and A-ASRAO in a UTMAX instance where 20 UEs exist. In essence, ASRAO performs fixed path searching, whose interim solutions are decided by their preceding solutions. ASRAO is sure to converge, but it lacks the ability to try better searching paths. While A-ASRAO adopts a random searching strategy, with wider searching directions, it shows more steady descending upper bound. From Figure 2, it is obvious that A-ASRAO converges much faster in the later part of the iterations than ASRAO as is depicted by lines with hollow marks. As these upper bounds reach the same optimal objective function value and all lower bounds get to optima much earlier than upper bounds, A-ASRAO will converge much faster than ASRAO.Figure 2
The converge process of ASRAO and A-ASRAO withT=(1,3,9).In addition, A randomness control parameterT is set in A-ASRAO to control its behavior, and A-ASRAO with T is denoted by A-ASRAO(T). With different values of T, A-ASRAO has different converging rates. This is probably caused by random generation of approximation matrix X~. Smaller T means more rows in X~ perform random 1-element selection, while large T will make less rows perform random selection. Extremely small or large values of T will make too many or few rows perform random selection, but neither of them is good for finding ideal initial solution for Phase II in A-ASRAO. So A-ASRAO(3) performs better than A-ASRAO(1) and A-ASRAO(9).In Figures3 and 4, respectively, the number of iterations performed and computation time spent by ASRAO and A-ASRAO are plotted against the number of UEs ranging from 10 to 20. Each point is obtained by averaging over 100 different network topologies. It shows that ASRAO needs a little fewer iterations than A-ASRAO when UE number is less than 14. However, when UE number is larger, A-ASRAO reduces the number of iterations by up to 15% compared with ASRAO in average.Figure 3
The number of iterations versus the number of UEs.Figure 4
Computational time versus the number of UEs.The reason for the different convergent properties lies in the scale of a UTMAX instance. A-ASRAO works in two phases and in Phase II A-ASRAO works the same as ASRAO, which implies that A-ASRAO in Phase I works actually to obtain an initial solution to Phase II in A-ASRAO. That is to say, if a better solution, which produces smaller gap between the upper bound and lower bound, found by the end of Phase I in A-ASRAO than that by ASRAO within the same iterations, A-ASRAO will have much larger probability to converge within fewer iterations than ASRAO. Therefore, the dominant factor relating to the different convergent rates of ASRAO and A-ASRAO is whether a better solution would be found by the end of Phase I in A-ASRAO.When a UTMAX instance is small-scale problem, its optimal solution can be derived within a few iterations even through exhaustive search, not alone ASRAO. While A-ASRAO tries a more steady but random searching, the finding of optimal solution may thus be postponed. So it hardly performs better than ASRAO. When a UTMAX instance is of large scale, a large solution space means a bolder searching practice may produce better searching path. So A-ASRAO(3) performs better than ASRAO and A-ASRAO(9). Nevertheless, too strong randomness goes against the finding of optimal solution, like in A-ASRAO(1) where all rows inX~ perform random 1-element selection. We speculate that is because a too wild searching will waste many iterations wandering around rather than heading towards a converging direction. That is why A-ASRAO(3) outperforms A-ASRAO(1). Moreover, in Figure 4 it shows that A-ASRAO spends less time computing optimal solution than ASRAO, and the performance gap is even larger than iteration numbers. This benefits from relaxing the MILP master problem into LP in Phase I of A-ASRAO, and LP can be solved much faster than MILP. Figure 4 also shows that ASRAO computes UTMAX much quicker than a common B&B method. Above results prove the effectiveness of (A-)ASRAO in solving UTMAX, and in most cases A-ASRAO(3) performs better than others.


## 4.3. Effect ofβ to UTMAX

As is mentioned in Section2, β indicates UTMAX’s preference on throughput fairness among UEs served in a BS by setting effects on resource allocation. As UE’s attainable data rate is changed, user associations are inevitably affected, which is reflected in Figure 5. And as a consequence, the performance of the entire interworking network is influenced. To analyze these effects, the performance of optimized network is evaluated on throughput produced and the number of UEs served in WiFi tier and cellular tier in Figure 6 and throughput gains of UEs in Figure 7 through APMAX and UTMAX with β in {-0.05,0,1,5}. PRMAX is a default user association scheme and adopts identical resource allocation in each BS and AP.Figure 5
Performance at AN granularity of PRMAX, APMAX, and UTMAX(β).Figure 6
Performance at network granularity of PRMAX, APMAX, and UTMAX(β).Figure 7
Performance at UE granularity of PRMAX, APMAX, and UTMAX(β).Figure5 plots the throughput and number of UEs served in each BS and AP with different approaches marked by different colors and symbols. It shows that, from a AN perspective, both APMAX and UTMAX do well in improving AP’s utilization, and APMAX performs even better, as it tends to maximize AP’s usage. The improvements of WLAN’s offloading performance by different approaches are much more obvious in Figure 6.Figure6 plot the throughput and number of UEs served in cellular tier and WiFi tier. APMAX gets most UEs served in WLAN, while, as resources are fully occupied, each UE in WLAN will get less resources and the throughput of WLAN will not be improved greatly compared with UTMAX. However, due to lack of careful UE selection from CN to WLAN in APMAX, UE amount decreases greatly in CN and the left UEs are probably with poor wireless channel condition or too far from a BS, which causes the poor CN throughput performance in APMAX. As for UTMAX, with different β value, it performs differently. Larger β value will get BS to allocate more resources to UEs with higher unit band data rate, which probably makes the contribution to the throughput growing in CN. However, as UEs with poor channel condition in CN get much less resources, they turn to WLAN to compensate the throughput loss. That would cause a minor UE amount increase in WLAN, while, as AP’s resources are equally shared, new coming UEs with poor channel condition share resources with indigenous WLAN UEs, causing the decrease in WLAN’s throughput. Nevertheless, the whole system throughput shows a continuous growth with the increase of β. It seems that larger β performs better on total throughput and throughput fairness between WLAN and CN, whereas it is on the expense of throughput fairness among UEs.Figure7 depicts UE throughput gain versus probability P(r<α) for various approaches versus PRMAX, where P(r<α) represents the ratio of UE whose throughput is less than α in all UE. The UEs’ throughput gains are quite large at 10% ratio point by all approaches, among which UTMAX(-ε) brings the largest throughput gain. With the growing of β, UE throughput gain decreases. However, UTMAX always has higher throughput gain over APMAX. With the growth of probability P(r<α) UE throughput gain decreases, because PRMAX also has relatively larger UE throughput at larger ratio points. At the meanwhile, UTMAX with larger β value has larger UE throughput gain, as larger β makes UTMAX allocate more resources to UEs with better channel conditions.


## 4.4. Effects of AP Density to UTMAX

The performances of UTMAX with different values ofβ are evaluated on larger AP density, to validate whether the same conclusion holds with different AP settings. The results are in Figures 8 and 9.Figure 8
Performance at network granularity of PRMAX, APMAX, and UTMAX(β).Figure 9
Performance at UE granularity of PRMAX, APMAX, and UTMAX(β).In Figure8, compared with Figure 6 the most obvious changes are the WLAN tier throughput increase for PRMAX and minimization of performance gap between UTMAXs with each of two different values of β. Due to the increase of AP density, more UEs access WLAN and the average distance from a UE to an AP is shorten. This is a direct cause to WLAN tier throughput increase for PRMAX and APMAX. As for the latter difference, it is possibly because of much fewer UEs in each BS or AP, and β probably has weaker effect on resource allocation among fewer UEs in a serving node. Figure 9 is plotted by adding curves of UE throughput gain where 12 APs exist. The additional UE throughput gain is largely contributed by newly added APs. A similar performance relation is shown between UTMAX with larger AP density, while the performance gap is narrowed and the curves meet and overlap earlier than that with smaller AP density.From these results, the effect ofβ to UTMAX is further checked that a smaller β makes UTMAX perform better on improving throughput of the part of UEs with poor channel conditions, while a larger β does better on improving overall system throughput.


## 5. Conclusion

To solve a key research case in multi-RAT SON about alleviating the biased load distribution in cellular/WLAN interworking network, we propose to optimize access load from users satisfaction perspective, and a UTMAX optimization model with sum of logarithmic utilities is derived. This is an MINLP and seems hard to solve intuitively. In order to solve it optimally and efficiently, we devise an ASRAO algorithm on the basis of general Benders Decomposition and prove its convergence in finite iterations to optimal. To reduce the computation complexity in ASRAO, especially when dealing with large-scale problem instance, we propose to accelerate the ASRAO algorithm with a relaxation and approximation technique inspired by feasible pump. Simulation results validate the convergence and effectiveness of both ASRAO and A-ASRAO and prove the effects of UTMAX optimization model on improving throughput fairness among users or the system throughput, and in either case, the utilization of WLAN is improved obviously especially on the number of users and throughput in it.As shown in simulation results, when dealing with large-scale UTMAX instance, A-ASRAO has distinct advantage over ASRAO on convergent speed, which implies that A-ASRAO is a proper algorithm to solve UTMAX. The resource allocation among users in BS is affected byβ, whereas simulation results imply it, and also functions as an indicator for preference on throughput fairness among users or total system throughput, which is a tradeoff or preference from viewpoints of network operators. However, due to computation complexity limitation, it is hard to perform extensive simulations on much larger UTMAX instances. So our future work includes a heuristic approach to produce solution close to that of ASRAO but with much lower computational complexity, and therefore effects of UTMAX can be further checked with larger UE density.



---

*Source: 6929575-2016-11-06.xml*