INTERFERENCE MINIMIZATION FOR DEVICE-TO-DEVICE (D2D)
COMMUNICATIONS
PEASH RANJAN SAHA
Bachelor of Science, University of Chittagong, 2015
A Thesis
Submitted to the School of Graduate Studies
of the University of Lethbridge
in Partial Fulfillment of the
Requirements for the Degree
MASTER OF SCIENCE
Department of Mathematics and Computer Science
University of Lethbridge
LETHBRIDGE, ALBERTA, CANADA
©c PEASH RANJAN SAHA, 2018
INTERFERENCE MINIMIZATION FOR DEVICE-TO-DEVICE (D2D)
COMMUNICATIONS
PEASH RANJAN SAHA
Date of Defence: December 07, 2018
Dr. Daya Gaur
Supervisor Professor Ph.D.
Dr. Shahadat Hossain
Committee Member Professor Ph.D.
Dr. Robert Benkoczi
Committee Member Associate Professor Ph.D.
Dr. Howard Cheng
Chair, Thesis Examination Com- Associate Professor Ph.D.
mittee
Dedication
To my parents & grandmother
iii
Abstract
In a cellular network, a central base station manages the cellular users. Direct device-to-
device (D2D) communication within short range without going through the base station
can improve the spectral efficiency of a cellular network. Thus D2D communication under-
laying cellular networks can play a crucial role in the fifth generation (5G) network. D2D
communication will also help the evolution of inter-device location-based applications ev-
erywhere such as emergency social services and advertising networks. D2D pair shares the
spectrum resources with the cellular user. As a result, a significant amount of interference
is generated in the system. We study the allocation of the resources from the cellular users
to the D2D pairs such that the total interference is minimized while guaranteeing a mini-
mum target sum rate. We propose a two-phase combinatorial algorithm which computes an
allocation subject to the sum rate constraint. For the case when the interference generated is
uniform irrespective of which D2D pair is communicating, the algorithm finds an optimal
solution in polynomial time. We also evaluate the algorithm empirically both on synthetic
and random data.
iv
Acknowledgments
It is my pleasure to acknowledge my thesis supervisor Dr. Daya Gaur. He guided me to-
wards research with his positive approach from the beginning. I also want to thank my thesis
committee members Dr. Shahadat Hossain and Dr. Robert Benkoczi for their thoughts and
suggestions.
I am very much thankful to Dr. Salimur Choudhury for his valuable suggestions. I would
also like to thank all the members of the Optimization Research Group for their helping
hand towards me every time.
I am grateful to my parents, wife, and parents-in-law for their continuous support and en-
couragement for the last two years while I was studying in Canada.
Finally, I would like to thank School of Graduate Studies, University of Lethbridge, Dr.
Daya Gaur and the external funding agency for the financial assistantship and giving me
the opportunity to learn and enrich my skill.
v
Contents
Contents vi
List of Tables viii
List of Figures ix
1 Introduction 1
1.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2
1.2 Contributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.3 Organization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
2 Preliminaries and Related Concepts 5
2.1 Definitions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
2.2 Maximum weighted bipartite matching . . . . . . . . . . . . . . . . . . . . 8
2.3 Integer linear program (ILP) . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.4 Integrality gap . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
2.5 NP-Completeness . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
3 System Model and Related work 13
3.1 SINR model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
3.1.1 Channel model . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
3.1.2 Network model . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
3.2 Related Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
3.2.1 Local search based solution . . . . . . . . . . . . . . . . . . . . . 15
3.2.2 Auction based solution . . . . . . . . . . . . . . . . . . . . . . . . 17
3.2.3 Other strategies . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
3.3 Problem Formulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
4 Two-phase Approach 23
4.1 The Characterization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
4.2 First phase . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
4.3 Second Phase . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
4.3.1 Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
4.4 Proof of Correctness . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
4.5 Complexity and Efficiency . . . . . . . . . . . . . . . . . . . . . . . . . . 34
4.6 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
vi
CONTENTS
5 Implementation & Performance Comparison 38
5.1 CPLEX Implementation . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
5.2 Evaluation and discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
5.2.1 Performance on the real network data . . . . . . . . . . . . . . . . 40
5.2.2 Performance on random data . . . . . . . . . . . . . . . . . . . . . 42
5.3 Improvement in Phase 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
6 An ILP for Phase 2 46
6.1 Formulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
6.2 Implementation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
6.3 Evaluation and discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
6.4 CORS conference presentation . . . . . . . . . . . . . . . . . . . . . . . . 50
7 Conclusion and Future Work 51
7.1 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
7.2 Future work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
Bibliography 53
A Table for different δ values 57
vii
List of Tables
5.1 Performance on data from network simulator . . . . . . . . . . . . . . . . 41
5.2 Performance on random data for δ = 0.4 . . . . . . . . . . . . . . . . . . . 42
5.3 Performance on random data for δ = 0.1 . . . . . . . . . . . . . . . . . . . 43
5.4 Performance on random data for δ = 0.9 . . . . . . . . . . . . . . . . . . . 43
6.1 Time needed for the second phase integer linear program . . . . . . . . . . 49
A.1 Performance on random data for δ = 0.2 . . . . . . . . . . . . . . . . . . . 57
A.2 Performance on random data for δ = 0.3 . . . . . . . . . . . . . . . . . . . 57
A.3 Performance on random data for δ = 0.5 . . . . . . . . . . . . . . . . . . . 58
A.4 Performance on random data for δ = 0.6 . . . . . . . . . . . . . . . . . . . 58
A.5 Performance on random data for δ = 0.7 . . . . . . . . . . . . . . . . . . . 58
A.6 Performance on random data for δ = 0.8 . . . . . . . . . . . . . . . . . . . 59
viii
List of Figures
2.1 Vertex-induced subgraph . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
2.2 Alternating Path . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
2.3 Maximum Weighted Matching . . . . . . . . . . . . . . . . . . . . . . . . 8
3.1 Scenario of network model with downlink resources . . . . . . . . . . . . . 15
4.1 Decrementing path . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
4.2 Special triple . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
4.3 Instance . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
5.1 Percentage of time in each phase . . . . . . . . . . . . . . . . . . . . . . . 44
5.2 Frequency of instances . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
6.1 Decision matrix for an edge . . . . . . . . . . . . . . . . . . . . . . . . . . 49
ix
Chapter 1
Introduction
Device-to-device (D2D) communication in an underlaying cellular network is becoming
increasingly common. D2D pairs can communicate directly without going through the
base station. This short range communicating feature makes D2D communication energy
efficient [1]. D2D communication also enables transmission of data at a very high rate.
D2D communication, therefore, increases the overall throughput in an LTE based cellular
network. Research in both academia and industry shows that D2D communication will play
a vital role in the fifth generation (5G) network and beyond [2].
In the recent past, several algorithms have been proposed to solve the interference mini-
mization problem in the presence of D2D communicating pairs in a cellular network. Local
search based methods have been proposed in [3], [4], and [5]. A few auction based methods
have been proposed in [6] ,[7] and [8]. Other notable work in minimizing interference in
D2D communications are in [9, 10].
A significant amount of interference is generated in the system when D2D pairs try to
share resources with the cellular users. Thus minimization of interference in D2D commu-
nications underlaying mobile networks is an important research question being studied in
recent years [5]. We study the problem of assigning resources from cellular users to D2D
pairs such that the total interference generated in the system is minimum and a minimum
target sum rate is also satisfied. We propose a two-phase resource allocation algorithm
where one D2D pair can share resources from at most one cellular user and vice-versa. Our
algorithm works for the general case where the interference generated by any pair is arbi-
1
1.1. MOTIVATION
trary. The algorithm generates an optimal solution when the interference from all the pairs
is uniform. We use a two-phase based approach. We use the maximum weighted bipar-
tite matching algorithm to generate a feasible assignment in the first phase. The total sum
rate in the first phase satisfies the minimum target level; otherwise, no solution is possible.
However, the interference generated by the allocation in the first phase can be arbitrarily
large. We show the existence of a particular structure which we call decrementing structure
(defined in chapter 4) for the uniform interference case, if the solution generated in the first
phase is not optimal. This interference decreasing structure is then used to improve the
solution. The main issue that arises is whether the interference decreasing structures can
be identified efficiently. We show that this is indeed the case. The second phase of our
algorithm determines the decrementing structures efficiently to improve the solution while
maintaining the same sum rate from the first phase. The output of our algorithm for the case
of uniform interference is an assignment of cellular users to D2D pairs where the total in-
terference is minimum, and the target sum rate is also satisfied. In general, the interference
is only locally optimal and may not be the minimum possible.
1.1 Motivation
In recent years, the number of smart handheld devices has exploded. Proximity-based
application and services are becoming increasingly popular. Location-based advertisements
and communication between emergency vehicles are such applications. The current cellular
network is four generation old and follows the traditional mechanism for transmission. The
infrastructure is also limited and cannot support the demand for the increased data rates
with small delays for the billions of connections [11]. The need for increasing the network
capacity is thus unavoidable.
Device-to-Device (D2D) communication is a significant option. It can help cope with
the increase in the demand of the spectrum to provide a better experience for the end users.
The interference generated in the system while a D2D pair shares the resources from the
2
1.3. ORGANIZATION
cellular users is unavoidable. As mentioned earlier, interference minimization has been
extensively studied in the recent past. We are not aware of any work on the uniform in-
terference scenario. We give a combinatorial two-phase based algorithm to minimize the
interference which works well in general and is optimal for the uniform interference sce-
nario.
1.2 Contributions
The resource assignment problem is NP-complete, shown by a simple reduction from
the Knapsack problem when the interference is arbitrary. The complexity is unknown for
the case where the interference is uniform. In this thesis, we study the case of uniform in-
terference and propose an exact algorithm to generate a resource assignment with minimum
interference. We show the resource allocation can be computed in polynomial time for this
scenario.
Ours is a two-phase an optimal combinatorial algorithm. We model the first phase as
the maximum weighted bipartite matching problem. We prove the existence of a decre-
menting structure if the solution of the first phase is not optimal. The second phase uses the
decrementing structures to improve the solution. Our algorithm also works well for the case
of nonuniform interference. We tested our algorithm on network simulator data as well as
random data. We analyze the performance of our algorithm with the solution of an integer
linear program (ILP) model for the problem.
We also model the second phase of our algorithm as an ILP. We presented some ex-
perimental results at the 60th Canadian Operational Research Society (CORS) conference,
2018, Halifax, Nova Scotia, Canada based on the ILP model for the second phase [12].
1.3 Organization
There are seven chapters including this chapter. We start by discussing the related topics
in chapter 2. The chapter is mostly concerned with maximum weight matching in bipartite
3
1.3. ORGANIZATION
graphs. We explain the system model and the network environment in chapter 3. We also
discuss some of the previous work using a different strategy in this chapter. We also give
some of the significant work in the broad research area of interference minimization for
D2D communication.
We describe our two-phase approach in chapter 4. We prove the correctness and analyze
the complexity of the algorithm in this chapter.
In chapter 5, we discuss the implementation technique. We compare our algorithm
performance with an ILP formulation of the interference minimization problem. In chapter
6, we describe an ILP model for the second phase and some experiments. Results were
presented at the CORS conference.
We conclude in chapter 7 with a summary of our findings and the limitations. We also
discuss some directions for future research.
4
Chapter 2
Preliminaries and Related Concepts
In this chapter, we discuss the preliminaries.
2.1 Definitions
Consider a graph G = (V,E), where V is the set of vertices and E is the set of edges.
We define the following notions using graph G.
Bipartite Graph
G is called a bipartite graph, if V can be partitioned into two set C and D, such that
V = (C∪D) and no two vertices from C or D have an edge between them.
Complete bipartite graph
if every vertex of C is adjacent to every vertex of D, then G is a complete bipartite graph.
The number of edges in G is (|C| ∗ |D|), where |C| is the cardinality of the set C and |D| is
the cardinality of the set D. A complete bipartite graph is denoted by K|C|,|D|.
Bipartite matching
A matching M in G is a subset of edges (M ⊆ E), such that no two edges in M share an
endpoint. The size of M is the number of edges in M and represent by |M|.
5
2.1. DEFINITIONS
Maximum matching
A matching M with the maximum number of edges is called maximum matching.
Perfect Matching
A matching M is perfect if every vertex in the graph G is an endpoint of an edge e ∈M.
Vertex-Induced Subgraph
A graph G′ = (V ′,E ′) is a vertex-induced subgraph of G, if it contains the vertex set
V ′ ⊆ V and an edge set E ′ ⊆ E, such that all the edges in E ′ have both the endpoints in
V ′. A vertex-induced subgraph of a complete bipartite graph is also a complete bipartite
graph. We denote it as a subgraph induced by V ′. Figure 2.1 show the subgraph induced by
V ′ = {2,3,5,6}.
1 4
2 5 2 5
3 6 3 6
(a) Complete bipartite graph G (b) Induced subgraph by {2,3,5,6}
Figure 2.1: Vertex-induced subgraph
Vertex-Induced Subgraph Matching
A matching in V ′ is called an induced matching. We define it as matching in the sub-
graph induced by V ′. The maximum size of a induced matching in V ′ is min(|C′|, |D′|)
where C′ and D′ are two set of vertex in V ′ and V ′ = (C′∪D′).
6
2.1. DEFINITIONS
Alternating path
Vertices in G which are not an endpoint of any edge e ∈M, are the unmatched vertices.
Alternating path P with respect to a matching M is a connected path with edges that alternate
in M, and not in M.
c1 d1
c2 d2
c3 d3
Figure 2.2: Alternating Path
Example of an alternating path is given in Figure 2.2. The path (c3→ d1→ c1→ d2→
c2 → d3) alternate between matched and unmatched edges. Here solid line represents a
matched edge and dashed line represents an unmatched edge. No two consecutive edges
are in M, or are not in M.
Augmenting path
Augmenting path P is an alternating path which starts with an unmatched vertex and
ends with an unmatched vertex. The sequence of edges in P alternate between unmatched
and matched. The number of edges in an augmenting path is odd. The number of un-
matched edges is one more than the number of matched edges.
The path (c3→ d1→ c1→ d2→ c2→ d3) in Figure 2.2 is an augmenting path. There
are three edges (c3,d1),(c1,d2) and (c2,d3) which are not in M and two edges (c1,d1) and
(c2,d2) are in M which alternates.
Weighted bipartite graph
If every edge e∈E of G has a positive weight we associated with it, then G is a weighted
bipartite graph.
7
2.2. MAXIMUM WEIGHTED BIPARTITE MATCHING
Path loss
Path loss is the attenuation of the transmission signal as it propagates from the transmit-
ting antenna to the receiving antenna. It can be described as a function of the propagation
distance.
Channel gain
The space between the transmitting antenna and the receiving antenna is considered as
a channel. Channel gain is defined by the difference of the transmitter output power and the
power received at the receiver end. Path loss is a component to consider while calculating
the channel gain. Some other factors are the hardware impairments, multipath fading, and
imperfect connections.
2.2 Maximum weighted bipartite matching
The matching M in G, for which the total weight of all of the edges in M = ∑e∈M we
is maximum, is called a maximum weighted matching. The size of M may or may not be
maximum in a maximum weighted matching.
C D C D
7 4 7 4
8 8
7 5 7 5
3 3
(a) (b)
      (a)  bipartite graph with edge weight        (b) Maximum weighted matching
Figure 2.3: Maximum Weighted Matching
Figure 2.3(a) is a weighted bipartite graph. The total weight of the matching in (b) is
20, and it is of maximum weight.
The Hungarian algorithm proposed in [13] computes a maximum weighted bipartite
8
2.2. MAXIMUM WEIGHTED BIPARTITE MATCHING
matching, and it is one of the earliest combinatorial optimization algorithms. It is also
known as the Kuhn-Munkres algorithm. The algorithm can be described using an augment-
ing path.
Hungarian algorithm using augmenting paths
The algorithm manipulates the weights of a bipartite graph to compute a maximum
weight matching. It searches for the augmenting paths until there exists any and computes
a perfect matching with the maximum weight.
Vertex labeling of a graph is a mapping between the set of vertices and a set of labels.
An integer value is assigned to every vertex which represents the label of the vertex. Let us
define the valid labeling of the vertices in a graph. Labeling is valid if the following relation
is true for any two vertices c ∈C and d ∈ D and an edge (c,d).
lc + ld ≥ w(c,d) ∀c ∈C,∀d ∈ D (2.1)
Where lc(Id) is the label of the vertex c(d), and w(c,d) is the weight of the edge (c,d).
The steps are given below:
• The algorithm starts with an empty matching and performs valid labeling on one
set of vertices by assigning the maximum weight among all the edges incident on a
vertex. The other set of vertices are labeled by 0.
• It then traverses the edges according to the non-increasing order of the weight. Two
vertices which are the endpoint of an edge are assigned to that edge if the vertices are
not an endpoint of any other matched edges.
• If the matching at this stage is perfect, then it is the final solution with maximum
weight and size. Otherwise, the algorithm searches for the augmenting path. Every
unmatched edge is not considered while searching for the augmenting paths. Only
the unmatched edges where the label of both endpoints of the edges is nonzero are
9
2.3. INTEGER LINEAR PROGRAM (ILP)
considered. if no augmenting path exists at the current stage, the algorithm improves
the labeling by splitting the weight of an edge into both sets of vertices such that
more unmatched edges can be considered while searching for the augmenting paths.
It then repeats the steps again.
• Whenever the algorithm finds an augmenting path, it switches the matching edges
with the unmatched edges of that augmenting path and repeats the previous step.
The algorithm stops when it generates a perfect matching. If necessary , dummy vertices
are used so that |C| = |D|. All the edges associated with the dummy vertices have zero
weight. The necessity to search for an augmenting path is followed from the Hopcroft-
Karp maximum matching theorem [14].
2.3 Integer linear program (ILP)
An Integer linear program (ILP) is a mathematical optimization program in which some
or all of the variables are restricted to be integers. The goal of an ILP is to optimize (max-
imize or minimize) a given objective function while satisfying some constraints. ILPs can
be classified as:
1. Pure ILP: All of the variables are integral.
2. Binary ILP: All of the variables are 0 or 1.
3. Mixed ILP: Some of the variables are restricted to be integral, and others are not.
The standard form of an integer linear program is :
minimize cT x
sub ject to Ax≥ b
x≥ 0
x ∈ (Z+∪0)
10
2.4. INTEGRALITY GAP
Here c and b are column vectors and A is a matrix. A,c,b,x can be represent as follows
[15].
 c1
 
  x1
    c2 

 

x2 

 .   . 
a11 a12 ... a    1n 1
  
b1  .   .  a21 a22 ... a2n 1   
b2
 

   
   . . 
 . . . . A =  
 .     
, b =  , c = c  and x =  x 
 . . . .  
   n  n  . 

 0 xn +1 . . . .   .  

 
  
 .   . am1 a  m2 ... amn 1 bm    

 . .    .   . 
0 xn+m
The last constraints restricts the decision variables to the set of integers.
2.4 Integrality gap
ILPs are known to be NP-hard [16]. The restrictions on the decision variables of ILPs
can be relaxed by assigning fractional values which are called linear programming relax-
ation (LP relaxation) of ILPs. Linear programming relaxation is a widely used technique to
design approximation algorithms for optimization problems. The solution to the LP relax-
ation is called a fractional solution. Integrality gap is the maximum ratio between the value
of the optimal solution when all variables are restricted to be an integer and the optimal
solution to the LP relaxation. For a minimization program p, if ip(I) is the optimal integer
solution and fp(I) is the optimal fractional solution for instance I, then the integrality gap
11
2.5. NP-COMPLETENESS
can be written as,
i (I)
IG p=
fp(I)
For a maximization problem, IG fp(I)= i (I) . The integrality gap is always at least 1.p
2.5 NP-Completeness
A problem whose solution is either a yes or a no is called a decision problem. Class P
is the set of decision problems which can be solved in polynomial time. The complexity
class NP (Nondeterministically polynomial) contains all the decision problems which can
be solved in polynomial time using a non-deterministic Turing machine. The correctness of
a solution of the problems in NP can be verified in polynomial time. NP-completeness is the
theory of determining the hardness of a problem using reduction. A problem x is reducible
to another problem y if some modification of input definition of x can be represented as an
input instance of problem y. If two problems are reducible to each other, then there exists
an algorithm to solve one problem which uses the subroutine of the algorithm to solve
another problem. NP-complete problems are the hardest problems in the class NP. If any
problem is polynomial time reducible to an NP-complete problem, then the problem also
belongs to the NP-complete problem class. The satisfiability problem is an example of an
NP-complete problem. Further information on NP-completeness can be found in [16].
12
Chapter 3
System Model and Related work
Modeling of interference in wireless networks is a very active area of research. Several
methods have been studied for the interference modeling. One unifying approach has been
given in [17]. They analyzed two well-known approaches given by the standard interference
function and general interference function. They established a bridge between both the
functions. A model for D2D interference based on stochastic geometry is presented in [18].
It takes into account the ability of BS to measure before selecting and configuring the D2D
communications. The most commonly used interference model in wireless communication
is the signal-to-interference-plus-noise ratio (SINR) model. Here we describe the SINR
model and some related work using this model. We use this model, as SINR models are
considered as the minimum level of detail needed to study interference [19].
3.1 SINR model
SINR is a mathematical technique to determine channel quality in wireless networks. It
defines the ratio between the gain of a channel based on a certain signal and the sum of the
interference and the external noise. The SINR is expressed in decibels (dB).
The external background noise and the interference generated by the other links are
unavoidable in wireless communication. During transmission, the distance between the
transmitter and the receiver is an essential factor that determines the channel gain at the
receiver. In a propagation model, a path loss function is defined based on the distance
between transmitter and receiver. Multipath fading has been introduced in the propagation
13
3.1. SINR MODEL
model. The path loss function from the propagation model in [20] is used to define the
SINR model. We describe the SINR channel model and the network model next.
3.1.1 Channel model
We consider an Urban Micro System as [3]. This also follows the Rayleigh fading path
loss model. The path loss (PL) equation is
PL = 36.7log10(d)+22.7+26log10( fc)[21] (3.1)
Where fc is the frequency measured in GHz and d is the distance between a D2D transmitter
t and receiver r. The channel gain between D2D transmitter and receiver is
t,r
Gt,r
−PL
= 10 10 [21] (3.2)
PLt,r is the distance dependent path loss between t and r.
3.1.2 Network model
We consider a cellular network with n cellular users and m D2D pairs where n ≥ m.
Realistically, n>>m. C is the set of cellular user, and D is the set of D2D pairs. We assume
that a single cellular user can share the radio spectrum with at most one D2D pair, and a
D2D pair can accept resource from up to one cellular user. The network has one central
base station which handles the allocation of resources to D2D pairs from cellular users [22].
There are two types of radio spectrum resources (uplink and downlink) is accessible for a
cellular user. We consider downlink resources in this thesis. The base station uses downlink
resources for signal transmission to the cellular users [21].
Two types of intra-cell interference are generated in two directions in this scenario when
a D2D pair share the resources with a cellular user. No inter-cell interference exists as the
downlink resources use an orthogonal channel. Victims of the interference are cellular
14
3.2. RELATED WORK
Cellular
Interfering  user
link
D2D
transmitter cellular link
Direct 
link
Base station
D2D
receiver Interfering link
Figure 3.1: Scenario of network model with downlink resources
user due to the D2D transmitter and the D2D receiver from the base station [3]. Figure
3.1 represents a scenario for the network model. The interfering link exists from D2D
transmitter to the cellular user and from the base station to the D2D receiver.
3.2 Related Work
For the last few years, many algorithmic solutions have been proposed for the research
allocation problems in D2D communication. Much research has been done to gain maxi-
mum throughput in D2D communication underlaying cellular networks while minimizing
the interference. The proposed solutions are mainly based on local search and auctions. In-
terference is modeled using the SINR model. We also discuss some other solution methods
in the literature.
3.2.1 Local search based solution
Local search is a frequently used method to solve computationally hard optimization
problems [23]. A local search method iteratively improves the initial feasible solution using
local changes. Local search based reassignment approach was first proposed by Islam et al.
15
3.2. RELATED WORK
in [3]. They proposed an algorithm to improve the system sum rate given an initial feasible
solution. System sum rate is the total sum rate generated by the assignment. SINR can be
computed for the victims of the interference (B and D2D receiver). They also set a lower
bound on the SINR at both B and the D2D receiver for each assignment. The higher the
SINR is, the less is the interference. Their approach is to restrict the interference for each
assignment at a target level. Their algorithm starts with a feasible one to one assignment.
It then searches for a pair of the cellular users such that swapping the assignment of their
corresponding D2D pairs improves the system sum rate while maintaining the target for
interference. They find that their algorithm performs better than the greedy heuristic based
resource allocation algorithm presented in [22].
The method due to Islam et al. [3] was extended by Hussain et al. [24]. They con-
sider three different approaches (one-to-one, one-to-many, many-to-many) for a D2D pair
to share resources from cellular users. They proposed a solution for the one-to-one ap-
proach to maximize the system sum rate using the weighted bipartite matching algorithm.
The difficult part is to maintain a quality of service requirement for the interference. Some
sharing can decrease the system sum rate. Their algorithm avoids such sharing and gen-
erates an assignment with the maximum total sum rate. They also proposed a resource
allocation algorithm for the other two approaches. The total interference generated can be
arbitrarily large in both of the algorithm proposed in [24].
Hassan et al. [5] also proposed an algorithm to minimize the interference based on the
local search method in two stages. The first stage is a combinatorial stage, which computes
a feasible allocation by solving a weighted bipartite matching problem. Interferences are
the weight on the edges. Sharing of resources which can decrease the sum rate are dis-
carded. Their algorithm stops if the total system sum rate meets the target. Otherwise, they
compute a new matching with the maximum total sum rate using the weighted bipartite
matching algorithm. In the second stage, local moves are used to decrease the interfer-
ence. Their algorithm tries to swap the assignment of a pair of a cellular user (for a pair
16
3.2. RELATED WORK
of D2D pair), if the swapping decreases the interference while maintaining the target sum
rate. This continues until no such improving pair exists. They implemented their algorithm
in a network simulator and compared the performance with an auction-based fair resource
allocation algorithm given in [6]. They find that the solution generated by their algorithm
is either optimal or very near to optimal in all the cases.
3.2.2 Auction based solution
A few auction based algorithms have been proposed to minimize the total interference
while satisfying the target sum rate. Islam et al. [6] have proposed a two-stage algorithm
for resource allocation. They ensure a fair resource allocation such that every D2D pair is
assigned to at least one cellular user. Instead of starting with a solution that meets the target
sum rate but has high interference, they begin with a solution with minimum interference
and try to modify it using an auction to increase the sum rate. In the first stage, they
assign resources to the selected D2D pairs using the minimum weighted bipartite matching
algorithm. If the target sum rate is satisfied, they return the solution. Otherwise, in the
second stage, D2D pairs bid according to a fixed bidding strategy. They try to increase
the sum rate by bidding for resources from an unallocated cellular user. They consider a
maximum increment to the sum rate and a maximum decrement to the interference in the
new bid. Their algorithm swaps the cellular user out if the bid for an unallocated cellular
user improves the solution. This algorithm may miss a feasible solution in some cases due
to a lack of synchronization between a D2D pair bid priority and the allocation of a cellular
user. A D2D pair selected based on minimum interference which might not be the best
choice for a D2D pair. Examples can be given where the algorithm misses an assignment
which increases the overall throughput of the system. It was shown by Hassan et al. [5] that
the solution returned by the two-stage algorithm of Islam et al. [6] can be trapped in local
optima. Hassan et al. [5] also showed the performance ratio of Islam et al. [6] is unbounded
in the worst case.
17
3.2. RELATED WORK
3.2.3 Other strategies
Now we mention some solutions based on different strategies not directly related to our
work, but they are essential in the broader general context of resource allocation in cellular
networks.
Feng et al. [25] have considered an approach based on admissibility. They proposed
a resource allocation algorithm to maximize the network throughput with guarantees to
satisfy the minimum quality-of-service. Admission control determines the admissibility of
a D2D pair to obtain power for transmission. The parameter to determine the admissibility
of a D2D pair is the distance from the base station. Computing a solution in their scheme
is a three-stage process. In the first stage, they perform admission control on D2D pairs
based on the proposed minimum distance metric. A D2D pair is considered admissible
if it does not violate the sum rate quality requirement. In this way, they also impose a
constraint on the access rate for a D2D pair. They allocate power to each admissible D2D
pair and the possible cellular users in the second stage. In the third stage, an allocation
is generated using the maximum weighted bipartite matching algorithm to assign suitable
cellular users to each admissible D2D pair. The proposed scheme significantly improves
the performance in terms of network throughput. However, some admissible D2D pairs can
have a significant amount of interference which can reduce the system capacity.
Hussain et al. [21] observed that at times sharing of resources between a cellular user
and a D2D pair can decrease the sum rate. This imposes the constraints of which cellular
users can be matched with which D2D pairs. They use this observation to rule out a match-
ing that decreases the overall sum rate. In effect when any cellular user shares the resource
with any D2D pair, the sum rate of the system is at least the base sum rate of the cellular
user without any resource sharing. They first generate an assignment matrix by satisfying
the quality constraints where each entry in the matrix represents a sum rate between the
associated cellular user and the D2D pair in the row and column. They replace the entry
with base sum rate of the cellular user if it decreases the system sum rate. In the next step,
18
3.3. PROBLEM FORMULATION
the weighted bipartite matching algorithm is used to find an allocation. The interference
in such a solution can be high as the algorithm did not seek to minimize it. Previously
mentioned techniques can be used to decrease interference though.
Janis et al. [26] discuss a solution based on local awareness between the cellular users
and D2D pairs. Doppler et al. [27] proposed a mechanism for D2D communication ses-
sion setup and management in LTE advanced network in an area with limited interference
to increase the overall throughput. A reverse iterative combinatorial auction-based radio
resource management scheme has been proposed in Song L. et al. [28]. They proposed a
non-monotonic descending price auction algorithm. Janis et al. [29] proposed a power con-
trol mechanism for the D2D pairs while sharing resources to ensure minimum performance
quality for the cellular users. Chia-Hao Yu et al. [30] also applied a power control method
to put constraints on the SINR degradation on the channel.
3.3 Problem Formulation
We use the SINR model and describe an integer linear program (ILP) for the interference
minimization problem as defined in [5]. Let us define the different types of channel’s link
and their signal gain in the network whenever a cellular user c ∈ C is sharing resources
with a D2D pair d ∈ D. Two transmission links exist in c and d to transmit the signal.
One between B and c and another between the transmitter and the receiver of d. As the
transmitter of d creates interference on c, there exists an interference channel link between
the transmitter of d and c. Similarly, another interference channel link exists between B and
the receiver of d. The corresponding channel gains are defined as follows:
Gbc = channel gain of the link between B and c.
Gtr = channel gain of the link between the transmitter and the receiver of d.
Gbd = channel gain of the interference link between B and the receiver of d.
Gdc = channel gain of the interference link between the transmitter of d and c.
19
3.3. PROBLEM FORMULATION
Let Pc and Pd are the transmission power of the cellular user c and the D2D pair d
respectively. The SINR at c while B is transmitting signal to c is represented by γ(c,d) (as c
is affected by d). γ(c,d) can be calculated by the following formula.
PcGbc
γ(c,d) = (3.3)T +PdGdc
c bc
Where T is the thermal noise. SINR at c can be written as γ = P Gc T if c does not shares
resources with d. Similarly SINR at D2D receiver is,
PdGtr
γ(d,B) = T +PB bd
(3.4)
G
Here PB is the transmission power of the base station B. The total interference generated
in the system due to sharing between c and d is given by i d dc(c,d) = P G +PBGbd . The sum
rate denoted by s(c,d) and can be calculated using Shannon’s formula [5] as follows:
s(c,d) = Bnlog2(1+ γ(c,d))+Bnlog2(1+ γ(d,B)) (3.5)
Where Bn is the channel bandwidth. if c does not share resources, then the base sum rate of
c can be written as sc = Bnlog2(1+ γc) [31].
We now describe an integer linear programming formulation for the interference mini-
mization problem. Let us define a decision variable x(c,d), which is 1, when a cellular user
c shares resources with a D2D pair d, otherwise it is 0. Total system sum rate S is given by
the following equation:
S = ∑ ∑ x(c,d)s(c,d)+ ∑ (sc(1− ∑ x(c,d))) ∀c ∈C,∀d ∈ D (3.6)
c∈C d∈D c∈C d∈D
The total system sum rate is calculated by adding sum rate for each cellular users that shares
20
3.3. PROBLEM FORMULATION
resources with any D2D pair. If any c does not share resources with any d, then
(1− ∑ x(c,d)) = 1
d∈D
. Let us assume the target sum rate is P. The interference generated in the system due to
resource sharing between cellular user c and D2D pair d is i(c,d), as earlier. The optimization
problem can be written as:
min ∑ ∑ x(c,d)i(c,d) (3.7)
c∈C d∈D
S≥ P (3.8)
∑ x(c,d) ≤ 1 ∀d ∈ D (3.9)
c∈C
∑ x(c,d) ≤ 1 ∀c ∈C (3.10)
d∈D
x(c,d) ∈ {0,1} ∀c ∈C,∀d ∈ D (3.11)
The objective function is to minimize the total interference. The first constraint ensures
that the target sum rate is satisfied. The second and third constraints impose the matching
constraint. This ensures that each D2D pair shares resources with at most one cellular
user and vice-versa. This ILP is NP-hard to solve. We propose a two-phase combinatorial
algorithm to find a solution to this ILP. For the case when the interference is uniform, we
prove that our algorithm ( algorithm 1) finds an optimal solution in polynomial time. We
also evaluate the performance of our combinatorial algorithm empirically on synthetic and
random data. The synthetic data is generated by the network simulator using different
parameters such as cell radius, maximum D2D pair distance, D2D transmit power.
Next, we show that the integrality gap of the local solution generated by Hassan et al.
[5] is unbounded.
21
3.3. PROBLEM FORMULATION
Consider a complete bipartite graph G = (V,E), where V = 2n, n cellular users and
n D2D pairs. Every edge e ∈ E has sum rate se and interference ie. Let assume se1 = 1
and ie1 = 1 and every other edges has value (
1
n+1 ,0) and the target sum rate R is 1. The
optimal solution for this instance is (1,1) as it selects only e1. However, the algorithm
proposed in [5] selects (n−1) edges with sum rate 1n+1 and
2
n+1 fraction of e1 such that the
total sum rate of the solution is (n−1)∗ ( 1 2n+1)+(n+1)∗1 = 1 and the total interference is
(n−1)∗0+( 2n+1)∗1 =
2
n+1 . So, the integrality gap is
1
2 =
n+1
2 , which is unbounded.
n+1
22
Chapter 4
Two-phase Approach
This chapter is similar in text to the manuscript [32]. The characterization and the two-
phase solution in the thesis are same as in the manuscript [32].
We study this optimization problem and propose a two-phase solution. We analyze
the features of the original graph with respect to the optimal solution of the interference
minimization problem (ILP in 3.7 - 3.11). We show the existence of a particular feature in
the graph for the uniform interference scenario which is the primary focus of our approach.
Sharing of any resources generates the same amount of interference irrespective of which
D2D pair is communicating in the uniform interference scenario. In the first phase, we
use the maximum weighted bipartite matching algorithm. The combinatorial algorithm
(Algorithm 1) in the second phase is used to improve the solution if the solution of the first
phase is not optimal. Our algorithm works well in general, and It generates the optimal
solution when the interferences are uniform. Let us assume that the interference generated
by any sharing of resources is uniform for the time being. We will show the correctness of
our algorithm under the uniformity assumption. We then discuss how the algorithm can be
used for the general case at the end of this chapter. In the next few subsections, we describe
the approach behind our algorithm, the first phase, the algorithm in the second phase, the
motivation of the algorithm design, and its complexity, and we prove the correctness.
23
4.1. THE CHARACTERIZATION
4.1 The Characterization
First, we characterize the optimal solution of the interference minimization problem
(ILP) when the interference is uniform. We analyze how the matched and the unmatched
edges relate to each other to form some specific structures if the solution is not optimal.
We establish the existence of certain structures that we call a decrementing structure if
the solution is not optimal. The total interference can be minimized using the decre-
menting structures. Let M be the current solution with maximum weight which is not
an optimal solution in terms of interference. M is a matching with total sum rate S(M) =
∑e∈M se ≥ (P− (sc(1−∑d∈D x(c,d)))), and total interference I(M) = ∑e∈M ie = k|M|, for
some k. Removing any edges from the matching violates the target sum rate constraints as
we assume the target sum rate satisfies in equality and all the edges in the matching has
nonzero sum rate. Some edges might exists where one endpoint of the edges is in V (M)
and other endpoint is not in V (M) and those edges has more sum rate than some edges in
M. Adding any of those edges in the matching can only remove one edge from the current
matching, otherwise the edge would have been added to the matching before. To describe
the decrementing structure we need some notation. Given a solution M, we call a path
Z = e1,e2, . . . ,ek alternating with respect to M, if the edges alternate in M and E \M, i.e.,
e1,e3, . . . ,ek ∈ M, and e2,e4, . . . ,ek−1 6∈ M. If an alternating path is such that ek ∈ M and
s(e1)+ s(e3)+ ....+ s(ek) ≤ s(e2)+ s(e4)+ ....+ s(ek−1), then the path is called a decre-
menting path. The number of edges in Z is odd and the number of matched edges is one
more than the number of unmatched edges. Let assume N(e1,e3,e5, ....,ek) is the number
of edges in the sequence e1,e3,e5, ....,ek, then N(e1,e3,e5, ....,ek) = N(e2,e4, ...,ek−1)+ 1
in Z. The decrementing path is an augmenting path starting at a matched edge. An example
of a decrementing path of length five is given in Figure 4.1.
Here solid line represents the matched edges and dashed line represents the unmatched
edges. Sum rates are given with each edge and the interferences are q for every edges. The
decrementing path is ordered as e1 = (c1,d1),e2 = (c2,d1),e3 = (c2,d2),e4 = (c3,d2),e5 =
24
4.1. THE CHARACTERIZATION
c 11 d1
2
c2 d5 2
10
c3 d6 3
Figure 4.1: Decrementing path
(c3,d3). The matched edges are e1,e3,e5 and e2,e4 are the unmatched edges. Furthermore,
s(e1)+ s(e3)+ s(e5) = s(e2)+ s(e4) = 12 and N(e1,e3,e5) = N(e2,e4)+1 in Figure 4.1.
Now we show the existence of the decrementing path whenever the solution given by
the matching M can be improved.
Lemma 4.1. If M is not an optimal solution, then there exists a decrementing path with
respect to M.
Proof. Since, M is not optimal there is another (optimal) solution M′ such that S(M′) ≥
(P− (sc(1−∑ ′d∈D x(c,d)))) and i(M ) < i(M). Let us examine the symmetric difference of
edges in M′ and M. These are edges that are either only in M or only in M′ but not in both.
The symmetric difference consists of even cycles or alternating paths. As the interferences
are uniform, the number of edges in M′ is smaller than the number of edges in M. Therefore,
one of the paths in the symmetric difference has more edges from M then from M′. This
path starts and ends with edges in M, and it is alternating. Therefore it is a decrementing
path.
If we have a decrementing path Z with respect to M, then (Z⊕M) = {Z \M}∪{M \Z}
is a new solution which satisfies the target sum rate constraint and has lesser interference
(number of edges) than M. Therefore the approach of our algorithm is: start with a matching
which has maximum sum rate and satisfies the target sum rate constraint, and while there
is a decrementing path, find one and update the current solution. Improving the solution
iteratively with existing decrementing path leads the solution to be optimal. We show an
efficient way to find a decrementing path if one exists.
25
4.1. THE CHARACTERIZATION
Lemma 4.2. Suppose Z = e1,e2,e3, . . . ,ek is a decrementing path with respect to M where
e1,e3 ∈ M and e2 ∈6 M. Let V (M) be the set of vertices covered by edges in M. Let,
V ′ = V (M \ e1)\V (e2∧ e3). If M′ is a maximum sum rate matching in the bipartite graph
induced by V ′ (Note that the endpoints of any edge in M′ are in V ′), then M′ ∪{e2} is a
matching in V such that S(M′∪ e2)≥ S(M).
Proof. Edge e2 is not an incident on any vertex in V ′, therefore, M′ ∪ e2 is a matching.
By Lemma 4.1, M1 = {Z \M}∪{M \Z} is a matching with sum rate S(M1) ≥ S(M) and
|M ′ ′1| < |M|. M1 \ {e2} is a matching induced by V . V has 2 ∗ (|M| − 2) + 1 vertices,
therefore a maximum matching in V ′ has size (|M|−2) and its sum rate is at least the sum
rate of M1 \{e2}. Therefore, for a maximum matching M′ induced by V ′ , M′∪ e2 satisfies
the sum rate constraint in V and its cardinality is at most (|M| − 1). At some point on Z,
there exists such triple of edges e1,e2,e3.
The efficiency of any algorithm which updates the current solution can be determined
by the technique to find e1,e2,e3, which are the starting triple of edges in Z. The other part
of Z is connected to e3 and can be determined by solving a subproblem once e1,e2,e3 is
identified. Let us call an update using a decrementing path Z as one step of the algorithm.
Since we do not know Z, we do not know e1,e2,e3. However, one can iterate through all
the possible choices for the first three edges of Z. The running time can be improved in
practice by discarding some choices for the first three edges. We observe any three edges of
a decrementing structure has some special property. We define this three edges as a special
triple, and this triple can be the starting point of searching for the decrementing path.
A special triple with respect to a matching M is a triple of three edges (e1 =(v1,v2),e2 =
(v ,v ),e = (v ,v )), where e ,e ∈ M and e ∈6 M, s(e ) ≥ s(e1)+s(e3)2 3 3 3 4 1 3 2 2 2 . Example of a
special triple is shown in Figure 4.2.
The special triple shown in Figure 4.2 is a part of the decrementing path of Figure 4.1.
There are two such triples of edges as (c1,d1),(c2,d1),(c2,d2) and (c2,d2),(c3,d2),(c3,d3)
in the decrementing path of Figure 4.1, but only the triple (c2,d2),(c3,d2),(c3,d3) is spe-
26
4.2. FIRST PHASE
c2 d5 2
10
c3 d6 3
Figure 4.2: Special triple
cial. This special triple can be defined as the starting triple such that e1 = (c2,d2),e2 =
(c3,d2),e3 = (c3,d3) for that decrementing path. Next, we show that every decrementing
path has at least one special triple.
Lemma 4.3. If there exists a decrementing path Z with respect to M, then there exists a
special triple in Z with respect to M.
Proof. Suppose, there is no special triple in Z. For all triples, (ei,ei+1,ei+2) on the decre-
menting path Z, where i+1 is even,
(s(ei)+ s(es e i+2
))
( i+1)< 2
Summing the previous equation over all even i+1 we get,
s1 +2∗ s3 +2∗ s5 + . . .+2∗ sk−1 + sk+1 > 2∗ s2 +2∗ s4 + . . .+2∗ sk
s1 s
+ s3 + s
k+1
5 + . . .+ s2 k−1
+ > s2 + s4 + . . .+ s2 k
This implies that the sum rate of the edges on Z but not in M is strictly less than the sum
rate of the edges in Z and M, a contradiction.
Now, instead of trying out every triple of edges defined by some two edges in the match-
ing to determine the decrementing path, only special triples are considered. We use this
approach to develop the solution. The first phase of our solution generates a matching with
maximum total sum rate. In the second phase, our algorithm searches for the decrementing
paths using only special triples. We describe the first phase and the second phase next.
27
4.2. FIRST PHASE
4.2 First phase
We model the first phase as the maximum weighted bipartite matching problem as in
[21]. For the rest of the thesis, we define it as the maximum sum rate matching problem
as the sum rate corresponds to the weight on the edges. We use the same decision variable
x(c,d) for the edge between a cellular user c ∈ C and a D2D pair d ∈ D as we used in the
ILP for the interference minimization problem. x(c,d) = 1, if c shares resources with d,
otherwise 0. The integer programming formulation for the first phase is
max ∑ ∑ x(c,d)s(c,d) (4.1)
c∈C d∈D
∑ x(c,d) ≤ 1 ∀d ∈ D (4.2)
c∈C
∑ x(c,d) ≤ 1 ∀c ∈C (4.3)
d∈D
x(c,d) ∈ {0,1} ∀c ∈C,∀d ∈ D (4.4)
Here the objective is to find a matching with maximum total sum rate where s(c,d) is the
sum rate of the edge between c and d. The constraints in (3.9) and (3.10) ensures the output
of that program is a matching M. We assume the graph is complete and each cellular user is
sharing resources with at most one D2D pair. The maximum weighted bipartite matching
algorithm always returns the solution with the maximum number of edge [33]. If the total
sum rate of all edges in M is less than the target sum rate P, then there is no solution to the
ILP for the interference minimization problem. Otherwise, the target sum rate is achievable.
A slack value α is defined by the difference between the total sum rate of the edges in M
and the target sum rate P. The equation to calculate α is
α = ∑ se−P (4.5)
e∈M
Although the total sum rate in the matching of the first phase satisfies the target, the total
28
4.3. SECOND PHASE
interference may not be optimal and arbitrarily large. The slack α measures to what extent
we can relax the sum rate while minimizing interference in the second phase. Without loss
of generality, we assume the sum rate cannot be decreased after the first phase. In the second
phase, we set α = 0. We also assume that the interference over all edges are uniform. We
discuss the use of the algorithm in general where the interferences are arbitrary at the end
of the chapter. The number of edges is less in the final matching than the number of edges
in current matching M if the interference can be minimized for the uniform interference.
Now the issue is whether the number of edges can be decreased while keeping the same
sum rate in the second phase.
4.3 Second Phase
We propose a combinatorial algorithm for the second phase to improve the solution.
The output of the first phase is matching with the maximum total sum rate. We describe an
integer linear program for the second phase in chapter 6. We have also proved the existence
of a decrementing path if the solution of the first phase is not optimal. Any decrementing
path has one less unmatched edges than matched edges where the total sum rate of the
unmatched edges is equal to the total sum rate of the matched edges. Each decrementing
path has at least one special triple (Lemma 4.3). We propose an iterative algorithm to de-
termine the decrementing paths using the special triples. In our algorithm, we search for
the special triples on one set of vertices say C. Special triples are ”identical” from both
sides of the graph. The unmatched edge of a special triple connects two vertices from two
sets. Any special triple can be determined by searching the edges from both of the vertices
endpoints of the unmatched edge. So, searching in one set is sufficient. The existence of a
special triple does not guarantee the existence of a decrementing path. It is possible for the
algorithm to find a special triple, but there may be no decrementing path associated with
it. We aim to update the current matching M only if the special triple is a part of a decre-
menting path. To accomplish this, whenever the algorithm finds a special triple with edges
29
4.3. SECOND PHASE
e1,e2,e3, it generates a temporary maximum sum rate matching M′ induced by vertex set
V ′ =V (M \e1)\V (e2∧e3). Since both e1 and e3 of the special triple is a candidate as e1 of
the decrementing path, we swap the position of e1 and e3 and generates the temporary max-
imum sum rate matching M′ at most two times. The algorithm updates current matching
M by M′ only if S(M′∪ e2)≥ P, and start searching again from the beginning to determine
another decrementing path if there exists any.
In a special triple, edges e1,e3 are in M and e2 is not in M. A special triple has four
vertices, and all the four vertices are matched and the endpoint of at least one of the edges
e1,e2,e3 of the decrementing path. Let us assume M′ is the improved matching where
I(M′)< I(M). The edge e ∈M′2 as e2 ∈/ M. It is not necessary to consider the two vertices
which are the endpoint of e2 as they can not be connected to any other edge. The rest of the
decrementing path can either be connected to the vertex which is an endpoint of e3, but not
an endpoint of e2 or to the vertex which is an endpoint of e1, but not an endpoint of e2. To
find the connection with rest of the decrementing path, we discard three vertices which are
endpoints of e1 and e2 and generates a maximum sum rate matching in the induced subgraph
given by V ′=V (M\e1)\V (e2∧e3). The other vertex is kept alive as it can be the connector
between the special triple and rest of the decrementing path. If s(M′∪ e2) < P, then there
is no decrementing path through that vertex. We then swap e3 and e1 to check whether the
connection can be established through the other side of the special triple. The algorithm
executes the previous steps again. If the special triple is not a part of a decrementing path,
then the algorithm continues searching for the next special triple. The algorithm returns the
solution whenever there exists no more special triple. We describe the algorithm below.
30
4.3. SECOND PHASE
4.3.1 Algorithm
Algorithm 1 Two-phase algorithm
1: procedure ALGORITHM(resource allocation from C(c1,c2, ...,cn) to D(d1,d2, ...,dm))
2: Phase 1 : M = maximum sum rate matching in input graph G.
3: Phase-2 :
4: start :
5: order the vertices of C from 1 to n in G.
6: loop :
7: for i = 1 to n do
8: for j = (n+1) to (n+m) do
9: if ((s(i, j))≥
s(i,ip)+s( j, jp)
2 ) then
10: e1 = (i, jp),e2 = (i, j),e3 = ( j, ip)
11: New matching :
12: Let M′ = maximum sum rate matching induced by V ′ =V (M \ e1)\V (e2∧ e3).
13: if S(M′∪ e2)≥ P then
14: M = (M′∪ e2)
15: goto start
16: else if (e3! = (i, jp)), then
17: e1 = ( j, ip) and e3 = (i, jp).
18: goto New matching
19: end if
20: end if
21: end for
22: end for
23: end procedure
The algorithm starts with the given sets of cellular users C and D2D pairs set D, which
we denote as original input graph G. The corresponding sum rates give the weights on the
edges in G. The output of the first phase is a maximum sum rate matching M on G. The
second phase starts at line 4, based on the matching M from the first phase. Algorithm 1
first orders the vertices of C in line 5 and starts executing the loop over vertices of C from
line 6. In the algorithm, numbers 1 to n denote vertices of C and (n+1) to (n+m) are for
vertices of D. The pair vertex of a vertex i is ip, such that edge (i, ip) ∈M. For each vertex,
it checks all edges incident on it and searches for a special triple in lines 7−9. Once it finds
a special triple, it initializes the edge number of that triple (e1,e2,e3) in line 10. It then
generates a temporary matching M′ induced by vertex set V ′ = V (M \{e1,e2}) in line 12.
31
4.4. PROOF OF CORRECTNESS
If no decrementing path exists connected to one side of the special triple, then we swap the
assignment of edge number between e1 and e3 in line 17. The algorithm repeats the steps
from line 11 in this case. If any of the two initialization leads to a decrementing path, then
the condition S(M′∪ e2)≥ P at line 13 is satisfied. The algorithm then updates the current
matching as M = (M′∪ e2) in line 14. If no decrementing path exists associated with that
special triple, it continues the loop over the next vertex to search for another special triple.
The algorithm returns the solution at line 23 when there are no more special triples.
4.4 Proof of Correctness
We show that our two-phase algorithm always returns the matching as the final out-
put with the minimum number of edges while satisfying the target sum rate. A maximum
sum rate matching M is returned in the first phase to ensure the target sum rate is achiev-
able. There exists a decrementing path if M is not optimal (Lemma 4.1). The matching M
can be updated using a decrementing path. The second phase of the algorithm repeatedly
updates M using the decrementing paths. It first finds the special triples associated with
decrementing paths and generates a new matching. Next, we show the Algorithm 1 always
successfully determines the decrementing paths until there exist any and updates the current
matching.
Lemma 4.4. If there exists a decrementing path with respect to the current matching M,
then the algorithm 1 always updates M.
Proof. Let us assume, there exists a decrementing path Z where set of matched edges X =
{(x1,x2, ....|xi) ∈ (M ∩ Z)} and set of unmatched edges Y = {(y1,y2, ...|y j) ∈ (Z −M)}.
Based on the definition of Z, the following relation is true.
(|M \X ∪Y |= (|M|−1))∩ (S(M−X +Y )≥ S(M)) (4.6)
Here M \X ∪Y is a new matching and S(M \X ∪Y ) is the total sum rate of all of the
32
4.4. PROOF OF CORRECTNESS
edges in that matching. It is possible to find a matching M′ with size at most (|M| − 2)
without using an edge e ∈Y , and after adding e to M′, the matching M \X ∪Y exists, and it
satisfies the target sum rate.
In a special triple, there are two edges from X and one edge from Y . When the algorithm
finds the maximum sum rate matching with V ′=V (M\e1)\V (e2∧e3) in line 12, it discards
the two matched edges (i, ip),( j, jp) from M and the two vertices from C. The total number
of vertices in V ′ is 2∗ (|M|−2)+1 where one side has |M|−2 vertices and other side has
|M|−2+1 vertices. The maximum size of the matching M′ generated in this step can be at
most |M|−2 and all the edges are selected in a way such that total sum rate is maximum.
In the next step, the algorithm adds e2 to M′, and S(M′∪ e2)≥ P∩ (|M′∪ e2| ≤ (|M|−1))
holds. The algorithm then updates M.
The algorithm always successfully determines the decrementing paths and updates the
current matching. The matching size decreases by at least one in each update. Algorithm 1
repeats this process iteratively as long as any decrementing path exists.
Theorem 4.5. Algorithm 1 always returns the optimal solution.
Proof. The algorithm finds a matching M as an output of the first phase. After that, only a
decrementing path with respect to M can decrease the number of edges in M while keeping
the sum rate unchanged. Our algorithm always updates M and reduce its size if there exists
a decrementing path (Lemma 4.4). It searches for the special triples until there exists any
and uses the decrementing paths associated with those triples to update M. In this way,
the algorithm decreases the number of edges of M by at least one in every update. The
iterative process continues until there exists any decrementing path. Finally, the algorithm
returns the matching M which has the minimum number of edges while satisfying the target
sum rate. M is of minimum cardinality if there exists no decrementing path with respect to
M.
We showed the correctness of our two-phase combinatorial algorithm. Next, we discuss
33
4.5. COMPLEXITY AND EFFICIENCY
the complexity of the update in the second phase and how the use of special triples affects
the running time in practice.
4.5 Complexity and Efficiency
First, we address the complexity of the update in phase 2. Let us call an update to the
current solution using a decrementing path Z, as one step of the algorithm. Since we do
not know Z, we do not know e1,e2,e3. However, (one)can iterate through all the possible
choices for the first three edges of Z. There are O( |V |2 ) choices for e1,e3 and two choices
for e2 in each step. Each step of the algorithm reduces the interference by at least one;
therefore the worst case number of steps is O(|V |). Each step requires computation of a
maximum weighted matching which takes O(|V |3) steps [34]. Therefore the total running
time is O(|V |6). Instead of trying out every possible choice for e1,e2,e3, the Algorithm 1
tries only the special triples. It generates the maximum sum rate matching using a special
triple to determine the decrementing path. Since we do not know the frequency of the
special triple in the real network scenario, we can not determine the exact complexity for
the updating technique of the Algorithm 1. In practice, the frequency of the special triples
are below 10% and solving the subproblem model using only the special triples improves
the running time noticeably.
In the general case when the interferences are arbitrary, the decrementing structures
may either be alternating paths or alternating cycles. The number of unmatched edges
is equal to the number of matched edges in an alternating cycle, and the total number of
edges is even. Algorithm 1 only searches for alternating paths. Therefore, it might miss
the alternating cycles which can also be used to minimize the interference. In such cases,
Algorithm 1 might not return an optimal solution. The interference minimization problem
can be proved as NP-complete by a reduction from the Knapsack problem. Algorithm 1
finds an optimal solution for the uniform interference scenario in polynomial time, but we
should not expect the optimal solution for the general case. However, in some preliminary
34
4.6. DISCUSSION
experiments for the general case with this solution, Algorithm 1 does well when it comes to
the quality of the solution. In Lemma 4.3, we showed the existence of at least one special
triple in a decrementing path. This is also true for the general case. The technique to
determine the decrementing paths using special triples in the second phase is still useful
to decrease the interference. One can, of course, obtain an exact algorithm by searching
for decrementing structure (path, cycle both) using another integer linear program. It is
an interesting open problem to design a combinatorial algorithm for finding decrementing
structures that are paths and cycles.
4.6 Discussion
In this section, we discuss the motivation behind the design technique of Algorithm
1. We show examples where determining the decrementing path in the second phase is
necessary. We use the maximum weighted bipartite matching algorithm in the first phase.
The goal of using the maximum weighted bipartite matching algorithm is to ensure the
target sum rate can be satisfied. The case where the solution to the first phase is not optimal
frequently occurs in a network.
The maximum weighted bipartite matching algorithm always returns the solution with
the maximum number of edge. If there are two different solution possible where one has
more edges than the other, then the maximum weighted bipartite matching algorithm al-
ways returns the solution with more number of edges. The Hungarian algorithm is a well
known maximum weight matching algorithm and always tries to find a perfect matching
with maximum weight [33]. Let us assume that M is the matching returned by the maxi-
mum weighted bipartite matching algorithm, and S(M) is the total sum rate of all the edges
in M. We assume that an edge exists between every cellular user and every D2D pair. A
perfect matching M′ have edges with zero-sum rate and edges with non-zero sum rate can
exist, where S(M′) = S(M). The solution with the less number of edges can still generate
a perfect matching by adding the edges with zero-sum rate. M′ has less number of edges if
35
4.6. DISCUSSION
we discard edges with the zero-sum rate. It still has a total sum rate S(M′). The maximum
weighted bipartite matching algorithm might fail to select the matching M′ and symmetric
difference of M, and M′ is a decrementing path.
c 21 d1
3
c2 d2 2
3
c3 d2 3
Figure 4.3: Instance
For example, if we execute the maximum weighted bipartite matching algorithm on the
graph in Figure 4.3, it always returns a matching with three edges of sum rate 2 each. The
same amount of sum rate can be generated with two edges with sum rate 3 by discarding
edges with the zero-sum rate. The edges (c1,d3),(c2,d1),(c2,d3),(c3,d2), do not exist in
the graph are assumed to have the zero-sum rate. We verified such example instances using
the Hungarian algorithm [33].
We also consider the complete bipartite graph scenario in which all D2D pairs are eli-
gible to share resources from any cellular user. In a real network scenario, a D2D pair does
not need to pull resources from all the cellular users. As a result, edges with the zero-sum
rate frequently appear in the network.
The distance between the cellular users and the D2D pairs and the frequent change of
channel state determines the transmission power of devices [35]. As a result, sum rates
vary between a cellular user and different D2D pair connected to it. Variation in the sum
rate on the edges and the existence of edges with the zero-sum rate increases the chances
of missing the optimal combination in the first phase of Algorithm 1. In those cases, the
maximum weighted bipartite matching algorithm might fail to generate the optimal solution
in terms of minimizing the number of edges.
For the uniform case, in the decrementing path, less number of unmatched edges have
36
4.6. DISCUSSION
the same total sum rate as a larger sized matching. This means at least one of the unmatched
edges has a larger sum rate than one or more of the matched edges. Otherwise, the total
sum rate of the unmatched edges can not sum up to total sum rate of the matched edges.
We define a special triple (e1,e2,e3) in a way such that s(e ) ≥ s(e1)+s(e3)2 2 , where e2 is an
unmatched edge. s(e2) may be more than the sum rate on both of the matched edges. We
aim to find a triple first such that we can get an unmatched edge which has at least the sum
rate as of one matched edge and then the algorithm searches for the rest of the decrementing
path. In this way, the chances are higher for the Algorithm 1 to determine the decrementing
paths.
37
Chapter 5
Implementation & Performance
Comparison
We explain the implementation of our two-phase algorithm for the uniform interference sce-
nario in this chapter. We discuss the simulation environment and the different datasets, on
which we perform the experiments. We show the performance comparison of our algorithm
with the optimal solution of the ILP (shown in 3.7 - 3.11) for the interference minimization
problem.
5.1 CPLEX Implementation
We implement the simulation process in IBM ILOG CPLEX 12.8 [36]. IBM ILOG
CPLEX optimizer can solve very larger integer linear programs defined using the opti-
mization programming language (OPL). A large variety of mathematical descriptions is
supported in OPL to specify an integer linear program. CPLEX model specified using OPL
defined as a model in this chapter.
The model where the graph is defined with the initialization of the parameters is called
the primary model. The primary model also calls the other models which are called sec-
ondary models and use the results returned by those models. The graph is complete as
we assume every D2D device can share resources from every cellular user. We define a
secondary model to create an instance on the graph. We generate an incident matrix for
the complete graph instance in that model. Let us assume that the graph has a total of 2n
vertices. Vertices from 1 to n are assigned to the first set (cellular users), and vertices from
38
5.1. CPLEX IMPLEMENTATION
n+ 1 to 2n are assigned to the second set (D2D pairs). The total number of edges is n2,
each vertex is the endpoint of n edges. One incident matrix is generated to keep track of the
edges incident on a vertex.
We define an array to store the sum rate assigned to each of the edges in the model.
The incident matrix and the sum rate array is then stored into an external data file so that
other models can use them. Another secondary model is defined to solve the integer linear
program of the first phase of Algorithm 1. The primary model calls the secondary model
of the first phase after a graph instance is created. The first phase secondary model uses
the incident matrix stored in the external file. The primary model then uses the matching
computed in the first phase in order to implement the second phase. The primary model put
a flag on the matched and the unmatched vertices and edges using the first phase matching
and the incident matrix. The purpose of using the flag of an edge is to determine whether
the edge is matched while solving the second phase. The primary model uses the external
file to retrieve the graph instance information. It then implements the iterative process to
determine the decrementing paths using special triples. The primary model also updates
the incident matrix whenever there is an update of the matching.
Every time a special triple is found, a subproblem need to solve to generate a temporary
maximum sum rate matching. The primary model uses the model defined for the first phase
to solve the subproblem since the subproblem can be solved using the first phase integer
linear program. Whenever the current solution is updated, the primary model also updates
the associated matched and unmatched flags of the vertices and edges. Algorithm 1 repeats
the steps until it completes traversing the last vertex. Every time Algorithm 1 determines a
decrementing path and updates the matching, it restarts from the first vertex again to search
for another decrementing path. When there is no more decrementing path, then Algorithm
1 stops. The matching returned from the primary model after the iterative process is the
final output of our algorithm.
We also define another secondary model to implement the integer linear program (ILP)
39
5.2. EVALUATION AND DISCUSSION
of the interference minimization problem. The primary model calls that model after execu-
tion of Algorithm 1. The graph instance stored in the external file is used in this model to
generate a matching for this problem. The primary model compares the result returned by
this model with the solution of Algorithm 1.
5.2 Evaluation and discussion
We perform multiple simulations for graphs with vertices within the range of 4 to 500.
We consider only even numbers in this range and split the vertices into two equal sized sets.
We perform two different experiments. The first experiment uses the sum rates generated
by the network simulator. We use the parameter as in [6] and [31] to generate the sum
rate for all of the edges. In the second experiment, we use random data to assign the sum
rate to the edges. A random value is generated within the range [0− 50] using a uniform
distribution as the sum rate for each of the edges. We set 50 as an upper bound of the
random value as the sum rate values generated by the network simulator is less than 50 in
most cases. We assign zero-sum rate using a threshold point δ in this case. We generate
another random value β within the range [0−1] for each edge. We fix a value for threshold
point δ and assign the zero-sum rate to the edge if the β value of that edge is less than the
threshold point δ. We compare the output of our algorithm with the optimal solution to the
ILP in both of the cases. Next, we note some results and the performance of our algorithm
on real network data and random data.
5.2.1 Performance on the real network data
We verified Algorithm 1 with the sum rates generated by network simulator using the
parameters as in [6] and in [31]. Whenever the interferences are uniform, and the first
phase solution is not optimal, Algorithm 1 successfully improves the solution in the second
phase. Results show that our algorithm always produces the same interference as the ILP
solution after the second phase. We show some notable instances and the comparison of
40
5.2. EVALUATION AND DISCUSSION
our algorithm with ILP in Table 5.1.
Table 5.1: Performance on data from network simulator
I I Time: Time: Time:
I Time:
| | | | Time: after after phase 2 phase 2 Algo-V E in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
10 25 4 0.0106 5 4 0.0062 0.0150 0.0013 0.0303
28 196 12 0.0446 14 12 0.0123 0.0570 0.0382 0.1599
64 1024 28 0.5358 32 28 0.0667 2.8637 0.0774 6.2092
90 2025 42 2.1359 45 42 0.1746 9.0720 0.4826 17.8165
150 5625 71 15.0942 75 71 0.7025 104.2081 3.7 204.8778
174 7569 84 24.2191 87 84 0.9544 135.633 4.0689 226.055
220 12100 104 31.121 110 104 1.4437 137.8935 4.1672 277.428
310 24025 149 71.414 155 149 2.8666 273.7927 8.2742 429.096
Notes: All times are given in second. |V |= number of vertex, |E|= number of edge.
I = number of edge in matching.
ILP refers to the integer linear program of (3.7) - (3.11).
We show the total time needed for Algorithm 1 and the ILP in Table 5.1. We break down
the total time for Algorithm 1 into phase 1 time and phase 2 time. Phase 2 time is further
subdivided into the time spent on computing matchings that were parts of some decrement-
ing paths and those that were not part of any decrementing path. We define the temporary
maximum sum rate matching model using a special triple as a subproblem. The algorithm
tries out every special triple as a possible choice to determine the decrementing paths. How-
ever, every special triple is not a part of a decrementing path. Algorithm 1 also solves the
subproblem model in phase 2 for cases where the solution cannot be improved using a spe-
cial triple. We divide the total time needed in phase 2 into the time when a subproblem
result was used to decrease the interference (update M) and the time the subproblem result
was not used to update the current solution.
As shown in Table 5.1, the runtime for Algorithm 1 is significant than the time it takes to
solve the ILP. In phase 2, Algorithm 1 searches for the special triples and solves a subprob-
lem model using those triples. Solving the subproblem model several times to determine
a decrementing path is the reason for the increase in time. Most of the time in phase 2 is
41
5.2. EVALUATION AND DISCUSSION
spent on solving maximum weighted matching problems for triples which are not part of
the decrementing paths. If special triples could be identified with certainty, then Algorithm
1 would be faster by a factor of 10. This motivates the development of a faster algorithm to
detect special triples that are part of some decrementing path. Another approach to speed
up Algorithm 1 would be to try to develop a primal-dual algorithm [37] for the uniform
case.
5.2.2 Performance on random data
We also tested our algorithm on random sum rate data. We assign a random value
generated uniformly at random in the range [0− 50] as the sum rate on each of the edges.
Algorithm 1 produces the same interference as the solution to the ILP does for this case
also. We define a new threshold value δ which is used to set the zero-sum rate on some of
the edges. We generate another random value β within the range of [0− 1] for each sum
rate. We assign zero-sum rate on an edge if the β value of the edges is less than δ. Below ee
show Tables with different δ values where Algorithm 1 is compared with the ILP. Results
for δ = 0.4 is shown in Table 5.2.
Table 5.2: Performance on random data for δ = 0.4
I I Time: Time: Time:
I Time:
| | | | Time: after after phase 2 phase 2 Algo-V E in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
16 64 7 0.0166 8 7 0.0096 0.0209 0.0042 0.0427
40 400 18 0.1019 20 18 0.0224 0.2871 0.0209 0.5684
74 1369 35 0.9082 37 35 0.0969 4.3370 0.1162 8.4565
120 3600 57 6.2196 60 57 0.3679 26.7939 1.0241 59.0941
168 7056 80 10.3026 84 80 0.6262 69.4826 1.7574 115.4051
218 11881 106 13.0327 109 106 1.0062 91.1816 4.4604 148.868
266 17689 129 16.6180 133 129 1.4109 134.7855 5.9573 212.7633
340 28900 164 56.9295 170 164 2.9534 303.9525 8.2212 374.0301
Notes: Same as Table 5.1
Algorithm 1 always improves the solution to optimal value when the phase 1 solution
42
5.2. EVALUATION AND DISCUSSION
is not optimal. However, the time needed by Algorithm 1 is significant for this case also
because of the subproblems to determine a decrementing path.
Results for δ = 0.1 and δ = 0.9 is shown in the Tables 5.3 and 5.4 respectively.
Table 5.3: Performance on random data for δ = 0.1
I I Time: Time: Time:
I Time:
| | | | Time: after after phase 2 phase 2 Algo-V E in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
22 121 10 0.141 11 10 0.0097 0.3946 0.0175 0.531
60 900 29 1.266 30 29 0.1044 2.9871 0.0312 3.777
92 2116 45 2.0419 46 45 0.1868 9.394 0.6712 13.4204
136 4624 66 5.040 68 66 0.3706 25.0732 1.4059 32.5067
168 7056 82 7.4088 84 82 0.4995 30.6145 1.7056 43.7356
270 18225 131 16.4429 135 131 1.6174 103.8161 2.3254 129.1139
Notes: Same as Table 5.1
Table 5.4: Performance on random data for δ = 0.9
I I Time: Time: Time:
I Time:
| | | | Time: after after phase 2 phase 2 Algo-V E in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
14 49 6 0.0402 7 6 0.0192 0.0641 0.0067 0.11
26 169 11 0.216 13 11 0.0662 1.1516 0.2871 1.494
68 1156 31 0.9445 34 31 0.4528 2.5653 0.2686 3.6067
114 3249 55 1.2840 57 55 0.7536 8.3754 0.6829 13.959
150 5625 70 3.899 75 70 0.7799 17.3131 1.3069 23.4999
204 10404 95 5.5068 100 95 1.3121 40.060 4.93 61.631
260 16900 124 20.256 130 124 1.5428 109.551 5.3254 152.155
324 26244 156 56.2371 162 156 2.6873 309.2314 7.4656 362.3571
Notes: Same as Table 5.1
The number of edges zero-sum rate varies based on the value of δ. More edges have
zero-sum rate if the value of δ is high. This also increases the chances of phase 1 returns
the non-optimal solution. Algorithm 1 produces the same interference as the ILP model for
each case.
43
5.3. IMPROVEMENT IN PHASE 2
Based on the results shown in the Tables, the total time needed by Algorithm 1 is much
more than the total time needed to solve the ILP of the interference minimization problem.
As we mentioned earlier, this is due to the subproblem models solved while searching for
decrementing paths. The percentage of total time spent in phase 1 and the percentage to
update the solution in phase 2 are shown in Figure 5.1. The percentages shown in Figure
5.1 are based on the average of the results presented in Table 5.1, Table (5.2 - 5.4), and
Appendix A.
search time in phase 2
92.88%
2.07%
phase 1 5.04%
update time in phase 2
Figure 5.1: Percentage of time in each phase
From the results in figure 5.1, we see that more than 90 percent of the time is spent
to determine the decrementing paths in the second phase. The time needed for the first
phase is only around 2 percentage of the total time. We can minimize the time needed by
the second phase if the special triples which are the part of the decrementing paths can be
identified without searching over all the special triples. Every time Algorithm 1 solves the
subproblem model to check for the existence of a decrementing path. The number of times
the subproblem model is solved can be decreased using more efficient techniques, possibly.
5.3 Improvement in Phase 2
We analyze the experimental results and observe that as the number of vertices and
edges increase in the graph, the solution of phase 1 is further away from the optimal so-
lution. As a consequence, there is more improvement in phase 2 as the graph increase in
44
5.3. IMPROVEMENT IN PHASE 2
size. The chances that phase 1 misses the optimal solution also increases as the number of
edges with zero-sum rate increases. The graph with more vertices is more likely to have
more edges with the zero-sum rate in both of the real network and random data scenario.
We show the improvement in phase 2 in Figure 5.2 for both types of data. We only consider
the random data when δ = 0.4.
10
Network data
Random data
6
4
2
0
4 100 200 300 500
Number of vertex
Figure 5.2: Frequency of instances
Phase 2 improvement is more significant in a real network as shown in Figure 5.2. The
frequency of the solution instance in a real network is notable as compared to random data
for different δ values as shown in Tables (5.3) - (5.4) and Appendix A. The result also
establishes the significance of determining the decrementing structures in the real network
scenario.
45
Improvement in phase 2
Chapter 6
An ILP for Phase 2
We describe an integer linear programming formulation for the second phase for the uni-
form interference scenario in this chapter. We explain the implementation and discuss some
significant experimental results.
The target of the interference minimization problem is to minimize the total interference
while the D2D pairs shares resources with the cellular users. A minimum target sum rate is
also needed to be satisfied. The first phase gives a maximum weighted bipartite matching
such that the target sum rate is satisfied. A solution does not exist, if the total sum rate
returned in the first phase is less than the target sum rate. The goal of the second phase
is to minimize the interference while maintaining the target sum rate. We assume that
the sum rate cannot be decreased further after the first phase. In the uniform interference
scenario, the goal of the second phase is to find a matching with less number of edges while
maintaining the same sum rate.
The objective of the second phase is to maximize the difference between the removed
matched edges and added unmatched edges if the solution can be improved after the first
phase while keeping the same sum rate. We model the goal of the second phase as an
integer linear program. We assume the interferences are uniform.
6.1 Formulation
Let assume the bipartite graph is G = (V,E) and se is the sum rate for an edge e. Con-
sider an indicator variable ce for the edge e. If an edge e is in the matching, then ce = 1,
46
6.1. FORMULATION
otherwise 0. Now, we introduce a new decision variable ye for an edge e in the second
phase which is 1 if the edge e is selected to either add into the matching or remove from the
matching, otherwise 0. The optimization problem for the second phase can be written as,
max ∑ ceye−∑ (1− ce)ye (6.1)
e∈E e∈E
∑ yese ≥ ∑ yese (6.2)
e∈/M e∈M
∑ ceye ≥ ∑ (1− ce)ye ∀v ∈V (6.3)
e∈δ(v) e∈δ(v)
∑ ye(1− ce)≤ 1 ∀v ∈V (6.4)
e∈δ(v)
ye ∈ {0,1} ∀e ∈ E (6.5)
Here the objective is to maximize the difference between the number of edges selected
to be removed from the matching and the number of edges selected to be added into the
matching. Constraint (6.2) ensures that the total sum rate of unmatched edges added is
at least the total sum rate of matched edges that are removed. If slack α > 0 after the
first phase, then the constraint in (6.2) can be written as ∑e∈/M yese ≥ ∑e∈M yese−α. The
equation (6.3) and (6.4) are the matching constraints. Constraint (6.3) ensures that if an
edge is added which is incident on a vertex v, then any matched edge incident on the vertex
v must be removed. Constraint (6.4) says that the number of edges added to every vertex
cannot be more than one. The last constraint says that the decision variable ye is either 0
or 1. The output of the program is the maximum decrease in interference that can be made
after the first phase. The solution of the first phase can be improved using a decrementing
path (Lemma 4.1). Algorithm 1 searches for a decrementing path to improve the solution.
For every special triple, it solves a subproblem to determine a decrementing path. The
output of the above ILP in (6.1) - (6.5) can be used in each iteration of Algorithm 1. In this
47
6.2. IMPLEMENTATION
way, subproblem for the special triples which are not a part of a decrementing path can be
avoided in Algorithm 1.
The output of the ILP mentioned in equation (6.1) - (6.5) can also be used to determine
whether the first phase solution is optimal. If the ILP returns a solution with value more
than 0, then the solution of the first phase is not optimal.
6.2 Implementation
We implement the integer linear program for the second phase using IBM ILOG CPLEX
12.8 [36]. The specification of the model for phase 1 and phase 2 is similar to the imple-
mentation of Algorithm 1. The primary model call the secondary models to create the graph
instances and to solve phase 1 and phase 2 for every instance.
The output of the secondary model for phase 1 is stored in the decision variable array.
The decision variable array is then passed to the secondary model for phase 2 by the primary
model. Phase 2 model defines a new decision variable array and uses the solution of phase
1. The output of the secondary model for phase 2 is the value which denotes how much
interference can be improved after phase 1. Secondary model for phase 2 also produces
a decision variable array as output. Based on the two decision variable array generated in
both phase, the primary model determines the edges which are the part of the improved
matching. The primary model outputs the final matching. We define a decision matrix to
determine whether an edge is a part of the final matching. The decision matrix for an edge
describes the overall states of the edge at multiple stages. The decision matrix for an edge
is shown in Figure 6.1.
According to the decision matrix shown in Figure 6.1, if an edge has value (0,0), then
the edge has never been selected. The edge e is in the final matching if the value are either
(1,0) or (0,1). The value (1,1) for an edge means that the edge has been selected in phase
1 and later removed in phase 2.
We tested the ILP for the second phase with many instances. We show some results
48
6.3. EVALUATION AND DISCUSSION
          ye 
0 1
ce
0 (0,0) 
(0,1) 
never selected added in the  second phase
(1,0) (1,1) 
1 added in the  removed in the 
first phase second phase
Figure 6.1: Decision matrix for an edge
in the next section. The results shown are significant in terms of determining the level of
improvement possible after the first phase in less time.
6.3 Evaluation and discussion
We conduct simulations to verify the ILP for the second phase. We generate the graph
instance for a number of vertices within the range [4−500]. We assign the zero-sum rate to
some of the edges at random. We perform the simulation three times for each possible case
in the vertex range. The integer linear program always returns the matching which has the
minimum number of edges.
Table 6.1: Time needed for the second phase integer linear program
| | | | I in Time: phase 2 Time:V E Total time
phase 1 phase 1 result phase 2
16 64 8 0.0119 1 0.0212 0.0337
50 625 25 0.0344 2 0.1263 0.1637
84 1764 42 0.1807 4 0.8359 1.0282
200 10000 100 1.4511 4 6.5637 8.0591
324 26244 160 4.3432 5 43.8423 48.1855
410 42025 203 12.0597 6 61.8751 74.1621
500 62500 247 24.3488 6 123.6714 148.2288
Notes: All times are given in second. |V |= number of vertex, |E|= number of edge.
I = number of edge in matching.
Phase 2 result indicates the possible improvement after phase 1.
We noted the time taken by the ILP to improve the solution. Results are shown in Table
6.1.
49
6.4. CORS CONFERENCE PRESENTATION
We generated graph instances such that a solution exists. If an instance, the matching
in the first phase is not optimal regarding the number of edges, then it can be improved in
the second phase. We observe that in all the cases the ILP returns the optimal solution as
shown in Table 6.1.
6.4 CORS conference presentation
We presented the model formulation for the second phase at the 60th Annual Canadian
Operational Research Society Conference, June 4−6, 2018, Halifax, Nova Scotia, Canada.
The purpose of the optimization cluster of CORS conference is to share the ongoing re-
search results.
We observed that the overall time needed for the model to improve the solution in the
second phase is less than the time needed for the ILP of the interference minimization
problem that we show in section 3.3. Algorithm 1 gives a combinatorial way to solve the
formulation in this chapter.
50
Chapter 7
Conclusion and Future Work
We discuss the main contribution of the thesis and outline some directions for further re-
search.
7.1 Summary
D2D communication underlaying cellular networks enables a large number of location-
based modern services, which is important for future generation networks. We study the
problem of allocation of resources to minimize interference while maintaining a target sum
rate in D2D communication underlaying cellular networks. We give a two-phase combina-
torial algorithm for the general problem. For the case when the interferences are uniform,
we prove that the algorithm generates an optimal solution in polynomial time. We evaluate
the performance of our algorithm empirically both on synthetic and random data.
7.2 Future work
The interference minimization problem that arises while allocating resources in cellular
networks to communicate D2D pairs has some interesting open questions which are not
answered yet. There are some possibilities to extend the work in the thesis. Our two-phase
combinatorial algorithm (Algorithm 1) generates an optimal solution for the uniform inter-
ference scenario. However, the time needed for the algorithm is significant than the time
needed to solve the integer linear program (ILP) of the interference minimization prob-
lem. A faster algorithm to detect special triples which are the part of some decrementing
51
7.2. FUTURE WORK
structures can be an important extension of the thesis. An empirical study can be done
for the general case when the interferences are nonuniform. An approximation algorithm
can be developed for this case based on scaling and our two-phase combinatorial algorithm
(Algorithm 1).
52
Bibliography
[1] Marko Hoyhtya, Olli Apilo, and Mika Lasanen. Review of Latest Advances in 3GPP
Standardization: D2D Communication in 5G Systems and Its Energy Consumption
Models. Future Internet, 10(1)(3), 2018.
[2] Gabor Fodor. D2D Communications What Part Will It Play in 5G? https://www.
ericsson.com/research-blog/device-device-communications/.
[3] M. T. Islam, A. M. Taha, S. Akl, and S. Choudhury. A Local Search Algorithm for
Resource Allocation for Underlaying Device-to-Device Communications. In IEEE
Global Communications Conference (GLOBECOM), pages 1–6, Dec 2015.
[4] Y. Hassan, F. Hussain, S. Hossen, S. Choudhury, and M. M. Alam. Interference
Minimization in D2D Communication Underlaying Cellular Networks. IEEE Access,
5:22471–22484, 2017.
[5] Md Yeakub Hassan, Faisal Hussain, Md Sakhawat Hossen, Salimur Choudhury, and
Muhammad Mahbub Alam. A near optimal interference minimization resource al-
location algorithm for D2D communication. In IEEE International Conference on
Communications, ICC 2017, Paris, France, May 21-25, 2017, pages 1–6.
[6] M. T. Islam, A. M. Taha, S. Akl, and S. Choudhury. A two-phase auction-based
fair resource allocation for underlaying D2D communications. In IEEE International
Conference on Communications (ICC), pages 1–6, May 2016.
[7] C. Xu, L. Song, Z. Han, Q. Zhao, X. Wang, X. Cheng, and B. Jiao. Efficiency Re-
source Allocation for Device-to-Device Underlay Communication Systems: A Re-
verse Iterative Combinatorial Auction Based Approach. IEEE Journal on Selected
Areas in Communications, 31(9):348–358, September 2013.
[8] C. Xu, L. Song, Z. Han, Q. Zhao, X. Wang, and B. Jiao. Interference-aware resource
allocation for device-to-device communications as an underlay using sequential sec-
ond price auction. In IEEE International Conference on Communications (ICC),
pages 445–449, June 2012.
[9] Tae-Sub Kim, Kwang-Hyung Lee, Seungwan Ryu, and Choong-Ho Cho. Re-
source Allocation and Power Control Scheme for Interference Avoidance in an LTE-
Advanced Cellular Networks with Device-to-Device Communication. SERSC : Inter-
national journal of control and automation, 6(1):181–190, February 2013.
53
BIBLIOGRAPHY
[10] Rongqing Zhang, Xiang Cheng, Liuqing Yang, and Bingli Jiao. Interference-aware
graph based resource sharing for device-to-device communications underlaying cel-
lular networks. In IEEE Wireless Communications and Networking Conference
(WCNC), pages 140–145, April 2013.
[11] Pimmy Gandotra and Rakesh Kumar Jha. Device-to-Device Communication in Cellu-
lar Networks: A Survey. Journal of Network and Computer Applications, 71:99–117,
August 2016.
[12] CORS. Canadian Operational Research Society. https://www.cors.ca/?q=
content/cors-annual-conferences.
[13] H.W. Kuhn. The Hungarian method for the assignment problem. Naval Research
Logistic, 2(1-2):83–97, 1955.
[14] John E. Hopcroft and Richard M. Karp. An n5/2 Algorithm for Maximum Matchings
in Bipartite Graphs. SIAM Journal on Computing, 2(4):225231, July 2006.
[15] Robert J. Vanderbei. Linear Programming: Foundations and Extensions. International
Series in Operations Research and Management Science, Springer, 1996.
[16] Michael R. Garey and David S. Johnson. Computers and Intractability: A Guide to
the Theory of NP-Completeness. W. H. Freeman and Company, 1979.
[17] Holger Boche and Martin Schubert. A unifying approach to interference modeling for
wireless networks. IEEE Transactions on Signal Processing, 58(6):3282–3297, June
2010.
[18] Yue Li and Lin Cai. Cooperative Device-to-Device Communication for Uplink Trans-
mission in Cellular System. IEEE Transactions on Wireless Communications, 2018.
[19] Aravind Iyer, Catherine Rosenberg, and Aditya Karnik. What is the right model for
wireless channel interference? IEEE Transactions on Wireless Communications, 8(5),
2009.
[20] Massimo Franceschetti and Ronald Meester. Random Networks for Communication,
From Statistical Physics to Information Systems. Cambridge Series in Statistical and
Probabilistic Mathematics, 24, 2008.
[21] F. Hussain, M. Y. Hassan, M. S. Hossen, and S. Choudhury. An optimal resource allo-
cation algorithm for D2D communication underlaying cellular networks. In 14th IEEE
Annual Consumer Communications Networking Conference (CCNC), pages 867–872,
Jan 2017.
[22] M. Zulhasnine, C. Huang, and A. Srinivasan. Efficient resource allocation for device-
to-device communication underlaying LTE network. In IEEE 6th International Con-
ference on Wireless and Mobile Computing, Networking and Communications, pages
368–375, Oct 2010.
54
BIBLIOGRAPHY
[23] Christos H. Papadimitriou and Kenneth Steiglitz. Combinatorial Optimization: Algo-
rithm and Complexity. Prentice-Hall, Englewood Cliffs, NJ, 1982.
[24] F. Hussain, M. Y. Hassan, M. S. Hossen, and S. Choudhury. System Capacity Max-
imization With Efficient Resource Allocation Algorithms in D2D Communication.
IEEE Access, 6:32409–32424, 2018.
[25] D. Feng, L. Lu, Y. Yuan-Wu, G. Y. Li, G. Feng, and S. Li. Device-to-Device Commu-
nications Underlaying Cellular Networks. IEEE Transactions on Communications,
61(8):3541–3551, August 2013.
[26] P. Janis, V. Koivunen, C. Ribeiro, J. Korhonen, K. Doppler, and K. Hugl. Interference-
Aware Resource Allocation for Device-to-Device Radio Underlaying Cellular Net-
works. In VTC Spring 2009 - IEEE 69th Vehicular Technology Conference, pages
1–5, April 2009.
[27] K. Doppler, M. Rinne, C. Wijting, C. B. Ribeiro, and K. Hugl. Device-to-device
communication as an underlay to LTE-advanced networks. IEEE Communications
Magazine, 47(12):42–49, Dec 2009.
[28] Xu C. Song L., Han Z. Radio resource management. In: Resource Management for
Device-to-Device Underlay Communication, Springer, New York, NY, 2013.
[29] Pekka Janis, Chia-Hao Yu, Klaus Doppler, Cssio Ribeiro, Carl Wijting, Klaus Hugl,
Olav Tirkkonen, and Visa Koivunen. Device-to-Device Communication under Laying
Cellular Communications Systems. International journal of communications, Net-
work and System Sciencese, 2(3):169–178, June 2009.
[30] C. Yu, O. Tirkkonen, K. Doppler, and C. Ribeiro. On the Performance of Device-to-
Device Underlay Communication with Simple Power Control. In VTC Spring 2009 -
IEEE 69th Vehicular Technology Conference, pages 1–5, April 2009.
[31] M. Tauhidul Islam, A. M. Taha, and S. Akl. Reducing the complexity of Resource
Allocation for underlaying Device-to-Device communications. In International Wire-
less Communications and Mobile Computing Conference (IWCMC), pages 61–66,
Aug 2015.
[32] Saha Peash, Choudhury Salimur, and Gaur Daya. Interference Minimization for
Device-to-Device Communications: A Combinatorial Approach. In preparation.
[33] Aleksejs Voroncovs Mark J. Becker. Hungarian Method. https://www-m9.ma.tum.
de/graph-algorithms/matchings-hungarian-method/index_en.html.
[34] Zvi Galil. Efficient algorithms for finding maximum matching in graphs. ACM Com-
puting Surveys (CSUR), 18(1):23–38, 1986.
[35] P. Sun, K. G. Shin, H. Zhang, and L. He. Transmit Power Control for D2D-Underlaid
Cellular Networks Based on Statistical Features. IEEE Transactions on Vehicular
Technology, 66(5):4110–4119, May 2017.
55
7.2. FUTURE WORK
[36] IBM analytics. CPLEX Optimizer. https://www.ibm.com/analytics/
cplex-optimizer.
[37] George Dantzig, Randolph Ford Jr, and Delbert Fulkerson. A primal-dual algorithm.
Technical report, Rand Corp Santa Monica CA, 1956.
56
Appendix A
Table for different δ values
This appendix includes the tables on the performance of algorithm 1 on random data for
different δ values. We mention some instances where the improvement in phase 2 is signif-
icant for some values of δ in the range of [0−1].
Notes in Table (5.1) are also applicable to Tables (A.1) - (A.6).
Table A.1: Performance on random data for δ = 0.2
I I Time: Time: Time:
I Time:
| | | | Time: after after phase 2 phase 2 Algo-V E in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
20 100 9 0.067 10 9 0.0326 2.1059 0.1726 3.4524
44 484 21 0.888 22 21 0.24 4.1955 0.4060 6.767
120 3600 58 3.6749 60 58 0.883 16.2066 1.0289 25.7249
268 17956 130 51.7430 134 130 2.0697 228.574 9.507 388.072
Table A.2: Performance on random data for δ = 0.3
I I Time: Time: Time:
I Time:
| Time: after after phase 2 phase 2 Algo-V | |E| in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
34 289 16 0.616 17 16 0.0924 2.214 0.160 3.5604
68 1156 32 0.969 34 32 0.5426 2.9066 0.2966 4.7471
200 10000 97 6.931 100 97 1.4050 27.237 4.5953 41.776
312 24336 152 57.3633 156 152 2.2945 252.6795 8.0624 336.01
57
A. TABLE FOR DIFFERENT δ VALUES
Table A.3: Performance on random data for δ = 0.5
I I Time: Time: Time:
I Time:
| | | | Time: after after phase 2 phase 2 Algo-V E in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
24 144 11 0.0579 12 11 0.0212 0.5836 0.0486 0.9728
72 1296 34 0.7825 36 34 0.4226 2.029 0.2276 3.2520
180 8100 235 3.9780 240 235 1.139 13.829 1.4804 20.336
240 14400 116 9.326 120 116 2.0305 81.7943 3.4080 136.3238
356 31864 172 64.1358 178 172 3.2438 292.0130 8.920 400.0179
Table A.4: Performance on random data for δ = 0.6
I I Time: Time: Time:
I Time:
| | | | Time: after after phase 2 phase 2 Algo-V E in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
20 100 9 0.0496 10 9 0.0068 0.5336 0.0444 0.8893
64 1024 30 0.8131 32 30 0.3507 2.9595 0.4097 4.5532
150 5625 71 3.722 75 71 0.7711 15.1595 1.9707 21.6565
274 18769 131 22.6669 136 131 2.3104 116.8854 3.7704 188.5248
310 24025 150 57.230 155 150 3.2912 269.944 9.6557 355.190
Table A.5: Performance on random data for δ = 0.7
I I Time: Time: Time:
I Time:
| | | | Time: after after phase 2 phase 2 Algo-V E in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
16 64 6 0.0183 7 6 0.0037 0.0955 0.0079 0.1592
70 1225 33 0.8961 35 33 0.359 2.8827 0.4804 4.8045
120 3600 54 2.1616 58 54 0.5203 11.193 1.691 17.91
320 25600 153 62.423 158 153 3.0064 241.4716 10.6051 326.313
58
A. TABLE FOR DIFFERENT δ VALUES
Table A.6: Performance on random data for δ = 0.8
I I Time: Time: Time:
I Time:
| | | | Time: after after phase 2 phase 2 Algo-V E in phase
ILP phase phase with no with rithm
ILP 1
1 2 update update 1
18 81 8 0.0445 9 8 0.0050 0.3339 0.0279 0.5566
74 1369 34 0.7522 36 34 0.2535 2.4219 0.3440 3.8876
210 11025 100 9.0854 103 100 2.1814 42.4252 4.8486 60.6075
392 38416 186 105.559 192 186 5.6918 425.465 20.422 567.2867
59