AA_research_paper

Journal of Computer and System Sciences 130 (2022) 43–70
Contents lists available at ScienceDirect
Journal of Computer and System Sciences

www.elsevier.com/locate/jcss
Approximation algorithms for clustering with dynamic

points ✩
Shichuan Deng a,∗ , Jian Li a , Yuval Rabani b
a
Institute for Interdisciplinary Information Sciences, Tsinghua University, Haidian District, Beijing, 100084, China
b
The Rachel and Selim Benin School of Computer Science and Engineering, The Hebrew University of Jerusalem, Jerusalem, Israel
a r t i c l e i n f o a b s t r a c t
Article history: We study two generalizations of classic clustering problems called dynamic ordered k-
Received 31 May 2021 median and dynamic k-supplier, where the points that need clustering evolve over time,
Received in revised form 16 April 2022 and we are allowed to move the cluster centers between consecutive time steps. In these
Accepted 7 July 2022
dynamic clustering problems, the general goal is to minimize certain combinations of the
Available online 22 July 2022
service cost of points and the movement cost of centers, or to minimize one subject to some
Keywords: constraints on the other. We obtain a constant-factor approximation algorithm for dynamic
Clustering ordered k-median under mild assumptions on the input. We give a 3-approximation for
Facility location dynamic k-supplier and a multi-criteria approximation for its outlier version where some
Dynamic points points can be discarded, when the number of time steps is two. We complement the
Multi-objective optimization algorithms with almost matching hardness results.
Approximation algorithms © 2022 Elsevier Inc. All rights reserved.
1. Introduction
Clustering a data set of points in a metric space is a fundamental abstraction of many practical problems of interest and
has been subject to extensive study as a fundamental problem of both machine learning and combinatorial optimization. In
particular, cluster analysis is one of the main methods of unsupervised learning, and clustering often models facility location
problems. More specifically, some of the most well-studied clustering problems involve the following generic setting. We
are given a set C of points in a metric space, and our goal is to compute a set of k centers that optimizes a certain objective
function which involves the distances between the points in C and the computed centers. Two prominent examples are
the k-median problem and the k-center problem. They are formally defined as follows. Let S denote the computed set of k
cluster centers, d( j , S ) = mini ∈ S d(i , j ) be the minimum distance from a pointj ∈ C to S, and D = (d( j , S )) j ∈C be the service
cost vector. The k-median problem aims to minimize the L 1 objective D 1 = j ∈C d( j , S ) over the choices of S, and k-center
aims to minimize the L ∞ objective D ∞ = max j ∈C d( j , S ). In general metric spaces and when k is not a fixed constant, both
problems are APX-hard [2,3] and exhibit constant factor approximation algorithms [4–9]. An important generalization is the
|C |
ordered k-median problem. Here, in addition to C and k, we are given also a non-increasing weight vector w ∈ R≥0 . Letting
D ↓ denote the sorted version of D in non-increasing order, the objective of ordered k-median is to minimize w D ↓ . This
problem generalizes both k-center and k-median and has attracted significant attention recently. Several constant factor
✩
A preliminary version of this paper [1] was presented at the 28th Annual European Symposium on Algorithms (ESA 2020). This version contains full
proof details and several new hardness results on the supplier problem.
*
Corresponding author.
E-mail addresses: moleculecloud@gmail.com (S. Deng), lijian83@mail.tsinghua.edu.cn (J. Li), yrabani@cs.huji.ac.il (Y. Rabani).
https://doi.org/10.1016/j.jcss.2022.07.001
0022-0000/© 2022 Elsevier Inc. All rights reserved.
S. Deng, J. Li and Y. Rabani Journal of Computer and System Sciences 130 (2022) 43–70
approximation algorithms have been developed [10–13]. We note that in the facility location literature, points are called
clients and centers are called facilities, and we will use these terms interchangeably.
In this paper, we study several dynamic versions of the classic clustering problems, in which the points that need clus-
tering may change for each time step, and we are allowed to move the cluster centers in each time step, either subject to
a constraint on the distance moved, or by incurring a cost proportional to that distance. These versions are motivated in
general by practical applications of clustering, where the data set evolves over time, reflecting an evolution of the under-
lying clustering model. Consider, for instance, a data set representing the active users of a web service, and a clustering
representing some meaningful segmentation of the user base. The segmentation should be allowed to change over time, but
if it is changed drastically between time steps, then it is probably meaningless. For a more concrete example, consider the
following application scenario. A giant construction company with several construction teams is working in a city. The com-
pany has k movable wireless base stations for their private radio communication, and each team also has a terminal device.
The teams need to put their devices at a certain energy level, in order to maintain the communication channel between the
device and the nearest base station. A team may finish their current project and move to another place to carry out new
tasks at some time. Note that the wireless base stations are also movable at a certain expense. Our high level objective is to
have all teams covered by the base stations at all times, meanwhile minimizing the energy cost of all teams plus the cost
of moving these base stations.
We study two problems of this flavor. The first problem, called dynamic ordered k-median (DOkMed), is a very general
model that captures a wide range of dynamic clustering problems where the objective is to minimize the sum of service
cost and movement cost. In particular, it generalizes dynamic versions of k-center and k-median. The problem is defined as
follows. We are given a metric space and T time steps. In each time step t, there is a set C t of clients that we need to serve,
and we can choose the locations for k mobile facilities to serve the clients (each client is served by its closest facility). Our
goal is to minimize the total ordered service cost (i.e., the ordered k-median objective), summed over all time steps, plus
the total cost of moving the k facilities. We define the problem formally as follows.
Definition 1. (DOkMed). Given a metric space ( X , d), an instance of T -DOkMed, T ∈ Z+ is specified by {C t }tT=1 , { F t }tT=1 ,
|C |
non-increasing weight vectors { w t ∈ R≥0t }tT=1 , and γ > 0, k ∈ Z+ . T ≥ 2 and γ are constants. C t ⊆ X is the set of clients at
time t, and F t ⊆ X is the set of candidate locations where we can place facilities. We are required to compute a sequence of
multi-sets of facilities { A t }tT=1 with A t ⊆ F t , | A t | = k, minimizing the following sum of ordered service cost and movement
cost,

T
T −1
↓
w t (d( j , A t )) j ∈C + γ · m( A t , A t +1 ),
t
t =1 t =1

where for a vector v, v ↓ is the non-increasingly sorted version of v , and m(Y , Z ) = min M 0 ∈ M (Y , Z ) (i ,i )∈ M 0 d(i , i ) is the
minimum total distance among perfect matchings between two equal-sized multi-sets Y and Z .
In DOkMed, the second term in the objective (i.e., the total distance traveled by all facilities) can also be motivated by
the online optimization problems of k-server [14,15] and k-taxi [16,17]. In these problems, k servers are present in a metric
space; at each time step, one client is revealed with a request, and some server needs to travel in the metric space to fulfill
the client’s request, incurring a movement cost. The unweighted version of DOkMed (i.e., each w t is an all-one vector) can
be roughly regarded as an offline version of k-server, except that many clients show up at each time step, and the servers
need to fix their locations and serve all clients simultaneously.
DOkMed is also related to the stochastic k-server problem, first studied by Dehghani et al. [18]. In this problem, we have
T time steps and T distributions { P t }t ∈[ T ] given in advance. The t-th client is drawn from P t , and we can use k movable
servers. One variant they consider is that, after a client shows up, its closest server goes to the client’s location and comes
back, and the optimization objective is the total distance traveled by all servers; in expectation, this objective is the same as
in DOkMed, if we consider non-ordered and weighted clients and all weights sum up to 1 for each time slot. We will further
discuss the relation between DOkMed and stochastic k-server in Section 1.2.
It is also natural to formulate dynamic clustering problems where the objective is to minimize just the service cost,
subject to some constraints on the movement cost. This turns out to be technically very different from DOkMed. Our second
problem, which we call dynamic k-supplier (DkSup), is such a concrete problem, where the service cost is the k-supplier
objective, i.e., the maximum client-facility connection distance over all time steps, and the constraint is that any facility
cannot be moved further than a fixed distance B ≥ 0 between any two consecutive time steps. More formally:
Definition 2. (DkSup). Given a metric space ( X , d), an instance of T -DkSup, T ∈ Z+ is specified by {C t }tT=1 , { F t }tT=1 , and
B ≥ 0, k ∈ Z+ . T ≥ 2 is the number of time steps, C t ⊆ X is the set of clients for time t, and F t ⊆ X is the set of candidate
locations where we can place facilities. We are required to compute a sequence of multi-sets of facilities { A t }tT=1 , with
A t ⊆ F t , | A t | = k, minimizing the maximum service cost of any client, maxt max j ∈C t d( j , A t ), subject to the constraint that
there must exist a perfect matching between A t and A t +1 for each t ∈ [ T − 1], and the distance between each matched pair
is at most B.
44
In the outlier version (DkSupOut), we are additionally given the outlier constraints {lt ∈ Z+ }tT=1 . We are required to
identify a sequence of multi-sets of facilities { A t }tT=1 and a sequence of subsets of served clients { S t ⊆ C t }tT=1 , with A t ⊆
F t , | A t | = k, | S t | ≥ lt . The goal is to minimize the maximum service cost of any served client, maxt max j ∈ S t d( j , A t ), with
the constraint that there must exist a perfect matching between A t and A t +1 for each t ∈ [ T − 1], and the distance between
each matched pair is at most B.
For the outlier problem T -DkSupOut, a multi-criteria (α0 , α1 , . . . , α T )-approximation is a polynomial-time algorithm that
always outputs a solution with objective value at most α0 times the optimum, while the number of served clients is at least
αt lt at time t, t ∈ [ T ]. Because standard k-supplier is APX-hard and prohibits any polynomial time approximation schemes
(PTAS) unless P = NP, we are interested in obtaining multi-criteria (α0 , 1 − , . . . , 1 − )-approximation algorithms for some
non-trivial α0 and any constant > 0, or pure approximations where α0 is non-trivial and all other αt ’s are 1.
Remark 1. The solutions to both DOkMed and DkSup are allowed to be multi-sets, since we regard the facilities as mobile
ones and it is natural for them to be co-located. We note that all of our hardness results also apply if we only allow the
solution to consist of subsets instead of multi-sets (See Section 3).
1.1. Our results
1.1.1. Dynamic ordered k-median

We first study DOkMed. When all entries of the weight vectors are larger than some small constant > 0, we present a
constant-approximation on general metrics. Moreover, for 2-DOkMed, i.e., T = 2, we present a constant-factor approximation
algorithm with no assumptions on the weight vectors.
Theorem 1. (Informal; see Theorem 8). When T = 2, there exists a polynomial-time constant-factor approximation algorithm for
2-DOkMed. When T ≥ 3 is a constant, and the smallest entry in { w t }tT=1 is at least some constant > 0, there exists a polynomial-
time O (γ / )-approximation algorithm for T -DOkMed.
Our techniques The key idea in our algorithm is to design a surrogate LP relaxation to approximate the ordered objective,
and embed its fractional solution in a network flow instance. We proceed to round the fractional flow to an integral flow,
thus obtaining the induced integral solution to the original problem. The network is constructed based on a filtering process
introduced by Charikar and Li [7]. When estimating the service cost of each client, we also adapt the oblivious clustering
arguments by Byrka et al. [11], with a slight increase in the approximation factor due to the structure of the network.
One notable difficulty we manage to overcome when translating an integral flow into an integral solution, is that the flow
oftentimes indicates the opening of more than k facilities, and we need to remove some of them without incurring a cost
that is unbounded compared to the LP objective.
1.1.2. Dynamic k-supplier

We first obtain inapproximability results for DkSup and its outlier version DkSupOut. These two problems turn out to be
much harder when T ≥ 3 than the case of T = 2.
Theorem 2. One has the following hardness results.
• Let T ≥ 3. There are no polynomial-time algorithms for T -DkSup with non-trivial approximation factors, unless P = NP (Theo-
rem 9). There exists a constant 0 ∈ (0, 1), such that T -DkSupOut admits no multi-criteria (α , 1 − 0 , . . . , 1 − 0 )-approximations
for any non-trivial factor α , unless P = NP (Theorem 10).
• There are no polynomial-time multi-criteria (α , 1, 1)-approximation algorithms for 2-DkSupOut for any non-trivial factor α ,
unless P = NP (Theorem 13).
On the positive side, we present a flow-based 3-approximation for 2-DkSup and a matching-based multi-criteria ap-
proximation for 2-DkSupOut. The approximation guarantee for 2-DkSup is optimal since vanilla k-supplier is NP-hard to
approximate within a factor of (3 − ) for any > 0 [19]. The multi-criteria approximation guarantee for 2-DkSupOut is also
nicely complemented by the aforementioned hardness result in Theorem 13.
Theorem 3. There exists a 3-approximation for 2-DkSup (Theorem 12). For every constant > 0, there exists a multi-criteria (3, 1 −
, 1 − )-approximation for 2-DkSupOut (Theorem 21).
Our techniques Our algorithm for 2-DkSup first guesses the optimum R , and uses a standard greedy algorithm (see, e.g.,
[4]) to create clusters for both time steps t = 1, 2. Since the maximum distance any facility can travel is B, we represent all
such possible movements using edges of a bipartite graph. The bipartite graph is then embedded into a network by append-
ing the aforementioned clusters on both sides, forming a network flow instance. The link capacities in the network are all
45
integers, thus we can directly obtain an integral flow (because the corresponding coefficient matrix is totally unimodular),
which in turn induces a 3-approximation to the original problem.
For 2-DkSupOut, we first guess a constant-size portion of “heavy” facilities in the optimal solution, properly modify the
instance and solve an LP relaxation on the remaining problem. This guessing step is standard in multi-objective optimization
[20]. According to the LP solution, we form clusters using the filtering algorithm by Harris et al. [21], and create a bipartite
matching instance where each vertex either has no contribution to the coverage of clients, or represents a cluster and
provides a certain number of nearby distinct clients to cover. By assigning vertex weights, the problem of covering some
specified numbers of clients for t = 1, 2 becomes finding a matching that satisfies a lower bound of total weights on both
sides of the bipartite graph. We round the LP-induced fractional matching to an integral one using the iterative rounding
methods by Grandoni et al. [20].
1.2. Related work
The ordered k-median problem generalizes a number of classic clustering problems like k-center, k-median, k-facility
l-centrum, and has been studied extensively in the literature. There are numerous approximation algorithms known for its
special cases. We survey here only the results most relevant to our work (ignoring, for instance, results regarding restricted
metric spaces or fixed k). Constant approximations for k-median can be obtained via local search, Lagrangian relaxations
and the primal-dual schema, and LP-rounding [6,9,22,23]. Constant approximations for k-center are obtained via greedy
algorithms [4,5]. Aouad and Segev [10] employ the idea of surrogate models and give the first O (log n)-approximation
for ordered k-median. Later, Byrka et al. [11] and Chakrabarty and Swamy [12] both successfully devise constant-factor
approximations for k-facility l-centrum and ordered k-median. Chakrabarty and Swamy [13] subsequently improve the ap-
proximation factor for ordered k-median to (5 + ), using deterministic rounding in a unified framework.
In the closely-related k-center with outliers problem, a.k.a. robust k-center, we are required to select k open facilities
F ⊆ X in the finite metric space ( X , d), m served clients S ⊆ X , and the goal is to minimize max j ∈ S d( j , F ). This problem
is introduced by Charikar et al. in [24], where they give a greedy algorithm that achieves an approximation factor of 3. A
best-possible 2-approximation is given independently by Chakrabarty et al. [25] and Harris et al. [21]. Many of its variants
are also studied. In matroid center with outliers, the input is the same as k-center with outliers, except that the cardinality
constraint is replaced with a given matroid, and the set of open facilities is required to be an independent set of the
matroid. Chen et al. [26] give the first 7-approximation for this problem, and a tight 3-approximation is later obtained by
Harris et al. [21]. In knapsack center with outliers, the input is the same as k-center with outliers with the cardinality
constraint removed, every facility has a non-negative weight, and the total weight of open facilities is required to be no
more than a given threshold. Chen et al. [26] give a 3-approximation that violates the knapsack constraint by a factor of
(1 + ) for this problem, and Chakrabarty and Negahbani [27] give the first pure 3-approximation.
Recently, Bandyapadhyay et al. [28] introduce the T -colorful k-center problem (T CkC), as a generalization of k-center with
outliers. In this problem, given a finite metric space ( X , d), there are T = O (1) subsets X t ⊆ X and T thresholds mt ∈ Z+
for t ∈ [ T ]. We are then asked to select k open facilities F ⊆ X , T served client subsets S t ⊆ X t satisfying | S t | ≥ mt for t ∈ [ T ],
and the goal is to minimize maxt ∈[ T ] max j ∈ S t d( j , F ). Evidently, T CkC recovers k-center with outliers by setting T = 1.
Bandyapadhyay et al. [28] give a pseudo 2-approximation for T CkC by opening at most k + T − 1 facilities. Anegg et al.
[29] and Jia et al. [30] independently obtain pure constant factor approximations for T CkC. It is fairly easy to see that
T -DkSupOut exactly recovers T CkC by setting the movement constraint B = 0, regarding the T client sets at different time
steps as having T colors, and removing redundant co-located open facilities from the solution. Unfortunately, unlike T CkC
which has pure constant approximations for T = O (1), our formulation of T -DkSupOut is seemingly much harder. In many
cases, pure approximations with any non-trivial factors, or even multi-criteria approximations are impossible unless P = NP.
Our problems are closely related to the mobile facility location problems (MFL), introduced by Demaine et al. [31]. In
these problems, a static set of clients has to be served by a set of facilities that are given initial locations and can be moved
to improve the service cost at the expense of incurring a facility movement cost. For the minimum total movement MFL
problem (TM-MFL), Friggstad and Salavatipour [32] give an 8-approximation using LP-rounding, where all facilities have unit
weights. Ahmadian et al. [33] give a local search algorithmfor TM-MFL with weighted facilities and proportional movement
costs via p-swaps with an approximation factor of 3 + O ( log log p / log p ), and specifically show that the factor is at most
499 for p = 1. Swamy [34] obtains an 8-approximation for the case of arbitrary movement costs using the reduction to the
matroid median problem. Krishnaswamy et al. [35] later improve the approximation factor of matroid median to 7.081.
The dynamic formulations of our problems are closely related to the facility location problem with evolving metrics,
proposed by Eisenstat et al. [36]. In this problem, there are also T time steps. While the facilities and clients are fixed,
the underlying metric is changing. The total cost is the sum of facility-opening cost, client-serving cost and additional
switching costs for each client. The switching cost is paid whenever a client switches facility between adjacent time steps.
In comparison, our problem DkSup considers the cost of moving facilities instead of opening costs, and allows the number
of clients to change over time. Eisenstat et al. [36] consider the problem when the open facility set A is fixed, and give a
O (log(nT ))-approximation, where n is the number of clients. They also show a hardness result on o(log T )-approximations.
An et al. [37] consider the case when the open facilities are allowed to evolve as well, and give a 14-approximation.
Our problems are also related to stochastic k-server [18] and the page migration problem [38,39]. As mentioned before,
the expected objective of stochastic k-server is the same as in DOkMed, if we consider non-ordered and weighted clients
46
and all weights sum up to 1 for each time slot. Dehghani et al. [18] provide an O (log n)-approximation for stochastic
k-server in general metrics, where n is the size of the distribution support. Our result in Theorem 1 does not imply a
constant approximation for stochastic k-server. The difficulty is that if one maps the stochastic k-server problem to ours,
the corresponding weight coefficient γ is not necessarily a constant and our approximation factor is linear in γ . Obtaining
a constant-factor approximation algorithm for stochastic k-server is still an interesting open problem.
Another related dynamic formulation of clustering problems is the fully dynamic model [40]. In these dynamic clustering
problems, there is an adversary with a hidden sequence of operations. At each time step t, the adversary inserts a point xt
or deletes a point xt , based on the t-th operation in its sequence. We are required to maintain a good approximate solution
after each adversarial operation, such that small running time and space are achieved in the long run. For classic k-clustering
objectives such as k-center, k-median and k-means, constant-factor approximation algorithms under fully dynamic models
are developed in [40,41].
1.3. Organization
The remainder of this paper is organized as follows. In Section 2, we give a polynomial-time approximation algorithm for
T -DOkMed based on LP rounding and a network flow instance. In Section 3, we first prove the hardness of approximation
for T -DkSup and T -DkSupOut when T ≥ 3, then provide a 3-approximation for 2-DkSup; we show another hardness result
on pure approximations for 2-DkSupOut, and complete the results with a multi-criteria (3, 1 − , 1 − )-approximation for
2-DkSupOut. Finally, we list some future directions and open problems in Section 4.
2. A constant approximation for DOkMed
In this section, we devise an LP-based algorithm for DOkMed. The ordered objective is estimated using reduced cost
functions as introduced in [11] (see Section 2.4 for the formal definitions). At the center of our algorithm, we construct a
network flow instance by applying a modified version of the filtering algorithm by Charikar and Li [7], and use an integral
flow to induce the output solution. When T ≥ 3, the integral flow may open more than k facilities in each time step, thus
we apply a crucial subroutine in Section 2.2 to reroute some part of the flow and prune the extra open facilities. We analyze
the approximation factor by adapting the oblivious rounding analysis by Byrka et al. [11] in Section 2.3, and provide the
missing proofs in Section 2.4.
2.1. LP relaxation
(t )
We first give the LP relaxation. By adding a superscript to every variable to indicate the time step, denote xi j ∈ [0, 1]
(t )
the extent of connection between client j and facility i, and y i ≥ 0 the extent of opening facility location i at time step t.
(t )
Moreover, denote zii the fractional movement from facility i to facility i between neighboring time steps t and t + 1. We
use the cost reduction trick by Byrka et al. [11]. Call d : X × X → R ≥0 a reduced cost function (not necessarily a metric) of
metric d, if for any x, y ∈ X , one has d (x, y ) ≥ 0, d (x, y ) = d ( y , x), and d(x1 , y 1 ) ≤ d(x2 , y 2 ) ⇒ d (x1 , y 1 ) ≤ d (x2 , y 2 ). For a
sequence of reduced cost functions D = {d(t ) }tT=1 of d, the relaxation is defined as follows.

T
T −1
(t ) (t )
min d(t ) (i , j )xi j + γ d(i , i ) zii (LP(D))
t =1 j ∈ C t i ∈ F t t =1 i ∈ F t i ∈ F t +1
(t )
s.t. xi j = 1 ∀ j ∈ Ct , t ∈ [T ]
i∈ F t
(t )
yi = k ∀t ∈ [ T ]
i∈ F t
(t ) (t )
0 ≤ xi j ≤ y i ∀i ∈ F t , j ∈ C t , t ∈ [ T ]
(t ) (t )
zii = y i ∀ i ∈ F t , t ∈ [ T − 1]
i ∈ F t +1
(t ) (t +1)
zii = y i . ∀ i ∈ F t +1 , t ∈ [ T − 1]
i∈ F t
(t ) (t ) (t )
We solve (LP(D)) and obtain an optimal solution (x, y , z), assuming that whenever xi j > 0, we have xi j = y i , via the
standard duplication technique on facility locations (for example, see [7]). Strictly speaking, whenever we split i ∈ F t into
co-located copies, we also need to split the corresponding variables in y (t ) and z(t −1) , z(t ) in order for it to remain feasible
(t ) (t )
to the LP relaxation. Here, we first split y (t ) s.t. xi j ∈ {0, y i }, then arbitrarily split the related z variables such that the
last two constraints in (LP(D)) are still satisfied. Denote B ( j , R ) = {x ∈ X : d(x, j ) < R } the open ball centered at j with
o
47
Algorithm 1: FILTER&MATCH.
Input : ( X , d), (x, y , z), {C t }tT=1 , { F t }tT=1
Output: filtered client subsets with a bundle for each filtered client and a partition on each subset
1 for t ∈ [ T ] do
2 C t ← ∅, C t ← C t
3 while C t is nonempty do // filtering phase
(t )
4 choose j ∈ C t s.t. dav ( j ) is minimized
(t )
5 C t ← C t ∪ { j }, C t ← C t \ { j }, delete each j ∈ C t s.t. d( j , j ) ≤ 4dav ( j )
6 for j ∈ C t do
(t ) (t ) 1 (t ) (t ) (t ) (t )
7 n j ← arg min j ∈C , j = j d( j , j ), R j ← 2
d( j , n j ), U j ← E j ∩ B o ( j, R j )
t
8 P t ← ∅, C t ← C t
(t )
9 while ∃ j ∈ C t s.t. n j ∈ C t do // matching phase
(t )
10 choose such j ∈ C t s.t. d( j , n j ) is minimized
(t ) (t )
11 P t ← P t ∪ {( j , n j )}, C t ← C t \ { j , n j }
12 for j ∈ C t do
13 P t ← P t ∪ {( j )}, C t ← C t \ { j }
(t )
14 return {C t }tT=1 , {U j : j ∈ C t }tT=1 , { P t }tT=1
(t ) (t ) (t ) (t )
radius R, and E j = {i ∈ F t : xi j > 0} the relevant facilities for client j. Denote dav ( j ) = i∈ F t d(i , j )xi j the average service
(t )
cost of client j and y (t ) ( S ) = i∈ S yi the “volume” of fractional facilities in S ⊆ F t , with respect to y (t ) . We perform a
(t )
filter-and-match algorithm (see Algorithm 1) to obtain a subset C t ⊆ C t for each t, a so-called “bundle” U j ⊆ F t for each
j ∈ C t and a partition P t of C t , where
• C t is a subset of “well-separated” clients of C t , and for any client j ∈ C t \ C t , there exists a relatively close client in C t ,
acting like a proxy for j ;
(t )
• Uj is a subset of fractionally open facility locations that are relatively close to client j ∈ C t , and the bundles in
(t )
{U j : j ∈ C t } are pair-wise disjoint;
• P t is a judiciously created partition of C t , where every subset is either a pair of clients (a.k.a. a normal pair), or a single
client (a.k.a. a singleton pair). Each normal pair { j , j } in P t is chosen such that either j or j is the nearest neighbor of
(t ) (t )
the other in C t , and at least one facility is opened in U j ∪ U j .
(t )
We first provide some basic properties of C t and U j (also see the results in [7]).
Lemma 4. (Charikar and Li [7]). Fix t ∈ [ T ]. The following statements hold.
(t ) (t )
(1) For any j , j ∈ C t , j = j , one has d( j , j ) > 4 max{dav ( j ), dav ( j )}.
(t ) (t ) (t )
(2) For any j ∈ C t \ C t , there exists j ∈ C t s.t. dav ( j ) ≤ dav ( j ), d( j , j ) ≤ 4dav ( j ).
(t )
(3) For any j ∈ C t , 1/2 ≤ y (t ) (U j ) ≤ 1.
(t ) (t )
(4) For any j , j ∈ C t , j = j , U j ∩ U j = ∅.
Proof. The first two assertions are obvious given the order we filter the clients in Algorithm 1. For the third one, we have
(t ) (t ) (t ) (t ) (t )
y (t ) (U j ) ≤ y (t ) ( E j ) = 1, since we assume xi j = y i whenever i ∈ E j . For the other inequality, assume otherwise and
(t ) (t ) (t ) (t ) (t ) (t ) (t ) (t ) (t )
we obtain y (t ) ( E j \ B o ( j , R j )) = 1 − y (t ) (U j ) > 1/2. One also has R j = 0.5d( j , n j ) > 2 max{dav ( j ), dav (n j )} ≥ 2dav ( j )
(t )
because j and n j are both in C t . This puts the average service cost of j at least
(t ) (t ) (t ) (t )
d(avt ) ( j ) ≥ y i d(i , j ) ≥ R j · y (t ) ( E j \ B o ( j , R j )) > d(avt ) ( j ),
(t ) (t )
i∈ E j \ B o ( j, R j )
(t ) (t ) (t ) (t )
which is a contradiction. For the last proposition, we simply notice that U j ⊆ B o ( j , R j ) and U j ⊆ B o ( j , R j ). Since the
(t ) (t ) (t ) (t )
sum of two radii is at most R j + R j = 0.5(d( j , n j ) + d( j , n j )) ≤ d( j , j ), the two open balls must be disjoint.
On the partition P t , we discuss the differences between our Algorithm 1 and the algorithm used by Charikar and Li [7]. In
each time step t, Charikar and Li use a simple greedy algorithm on the filtered client set C t , matching the closest unmatched
pair in C t whenever possible, thus leaving at most one client in C t unmatched. In our algorithm, we first compute for every
48
Fig. 1. Several intermediate layers of N representing a single time step t. Some nodes and links are left out for simplicity.
(t ) (t )
j ∈ C t its nearest neighbor n j ∈ C t ; whenever there exists an unmatched j such that n j is also unmatched, we choose j
(t )
which minimizes d( j , n j ) among such choices and match the two. Notice that this process may leave an arbitrary number
of clients in C t unmatched and they all end up in singleton pairs.
The motivation in our algorithm to leave many clients unmatched, is that by restricting them to singleton pairs, the
additional rerouting cost incurred in the post-processing phase (see the next section) is much easier to bound. In fact,
we are not certain whether the original matching in [7] can result in a good approximate solution under our rounding
framework. It is also worth noting that, while we define the objective of (LP(D)) using reduced cost functions in D that
simulate the ordered service cost, Algorithm 1 is completely oblivious of them and only uses the underlying metric d.
2.2. Flow-based LP rounding
We construct an instance of network flow N , and embed the LP solution as a fractional flow f̃ . The network N consists
of a source s, a sink t and 6T intermediate layers L 1 , L 2 , . . . , L 6T arranged in a linear fashion. For each time step t ∈ [ T ],
(t )
we create two nodes for every pair p ∈ P t , every bundle U j and every candidate facility location i ∈ F t . All these nodes
are contained in the layers L 6t −5 , . . . , L 6t . To distinguish between the two mirror nodes, we use L (·) and R (·) to represent
the nodes in { L 6t −5 , L 6t −4 , L 6t −3 } (on the left) and the nodes in { L 6t −2 , L 6t −1 , L 6t } (on the right), respectively. The network
is constructed as follows, and an illustration is given in Fig. 1.
Step 1. For each i ∈ F t , t ∈ [ T ], add nodes L (i ) to L 6t −5 and R (i ) to L 6t .

(t ) (t ) (t )
Step 2. For each U j , t ∈ [ T ], add nodes L (U j ) to L 6t −4 and R (U j ) to L 6t −1 .
Step 3. For each p ∈ P t , t ∈ [ T ], add nodes L ( p ) to L 6t −3 and R ( p ) to L 6t −2 .
(t ) (t )
Step 4. For each j ∈ C t , p ∈ P t s.t. j ∈ p , t ∈ [ T ], add links (L (U j ), L ( p )), (R ( p ), R (U j )) both with capacity range

(t ) (t ) (t ) (t ) (t )
y (t ) (U j ), y (t ) (U j ) and initial fractional flow values f̃ (L (U j ), L ( p )) = f̃ (R ( p ), R (U j )) = y (t ) (U j ).
Note that the capacity is either [0, 1] or {1}.
Step 5. For each p ∈ P t , t ∈ [ T ], add the link (L ( p ), R ( p )) with capacity [ y (t ) ( p ), y (t ) ( p )] and f̃ (L ( p ), R ( p )) =
(t )
y (t ) ( p ) = j ∈ p y (U j ). According to Lemma 4, if p is a normal pair, the capacity is either [1, 2] or {1} or {2},
(t )
since y ( p ) ≥ 2 × (1/2) = 1; if p is a singleton pair, the capacity is either [0, 1] or {1}.

(t)
(t ) (t ) (t )
Step 6. For each j ∈ C t and i ∈ U j , t ∈ [ T ], add links (L (i ), L (U j )), (R (U j ), R (i )) with unit capacity [0, 1]. Let the
(t ) (t ) (t )
initial fractional flows be f̃ (L (i ), L (U j )) = f̃ (R (U j ), R (i )) = y i .
(t ) (t ) (t )
Step 7. For each i ∈ F t \ ( j ∈C U j ), t ∈ [ T ], add the link (L (i ), R (i )) with capacity [ y i , y i ] (across intermediate
t
(t )
layers L 6t −4 , . . . , L 6t −1 ). Let its initial fractional flow be f̃ (L (i ), R (i )) = y i . Note that the flow may be larger than
(t )
one, since we allow multi-sets and do not impose the constraint y i ≤ 1 in (LP(D)).
(t ) (t )
Step 8. For each i ∈ F t , i ∈ F t +1 , t ∈ [ T − 1], add the link (R (i ), L (i )) with capacity [ zii , zii ]. Let its initial fractional
(t )
flow be f̃ (R (i ), L (i )) = zii .
(1) (1) (1)
Step 9. For each i ∈ F 1 , add the link (s, L (i )) with capacity [ y i , y i ] and initial flow f̃ (s, L (i )) = y i . For each
(T ) (T ) (T )
i ∈ F T , add the link (R (i ), t) with capacity [ y i , y i ] and initial flow f̃ (R (i ), t) = y i .
Notice f̃ is naturally a flow with value k, since the flow conservation constraints are directly satisfied by the last two
constraints of (LP(D)). Because the flow polytope is defined by a totally unimodular matrix, and our capacity constraints
49
Algorithm 2: REROUTE.
(t )
Input : N, f̄ , {C t }tT=1 , {U j : j ∈ C t }tT=1 , { P t }tT=1
Output: a feasible solution to the original T -DOkMed instance
1 for t ∈ [ T ] do
2 At ← ∅
3 for p ∈ P t do
4 if f̄ (L ( p ), R( p )) = 2 then
5 pick the left-activated i 1 , i 2 , A t ← A t ∪ {i 1 , i 2 }
6 else if f̄ (L ( p ), R( p )) = 1 then
(t ) (t )
7 if p = ( j 1 , j 2 ), n j = j 2 , n j = j 1 then // j 1 and j 2 are closest to each other
1 2
8 pick the left-activated i, A t ← A t ∪ {i }
(t ) (t )
9 else if p = ( j 1 , j 2 ), n j = j 2 , n j = j 1 then
1 2
10 if f̄ passes through the same bundle in L 6t −5 and L 6t then
12 else
(t )
13 pick the (left or right) activated i ∈ U j , A t ← A t ∪ {i }
2
14 else if p = ( j ) then
(t )
16 for i ∈ F t \ j ∈C t U j do
17 if f̄ (L (i ), R(i )) = m ≥ 1 then
18 add m copies of i to A t
19 return { A t }tT=1
are all integers, it is well-known (see, e.g., [42]) that the flow polytope has integral extreme point solutions. Using the
dependent rounding algorithm by Kumar et al. [43], one can efficiently sample an integral flow f̄ corresponding to an
integral extreme point solution, such that f̄ is guaranteed to have value k, and E[ f̄ ] = f̃ holds for all links (cf. [44]). Next,
given the integral flow f̄ , we deterministically construct the solution { A t }tT=1 as follows,
• If T = 2, there are 12 layers L 1 , L 2 , . . . , L 12 in the network. For each link e = (R (i 1 ), L (i 2 )) between L 6 and L 7 such
that f̄ (e ) = m ≥ 1, add m copies of i 1 to A 1 and m copies of i 2 to A 2 .
• If T ≥ 3, f̄ does not immediately reveal a feasible solution. To see this, focus on any unit flow in f̄ . It may enter L 7
and exit from L 12 (both for t = 2) through nodes that represent different facility locations. We design the following
Algorithm 2 to resolve this issue.
For a facility i ∈ F t , if there is (at least) one unit of flow through L (i ) or R (i ), we call the facility i left-activated or right-
activated correspondingly. The algorithm looks at each pair ( j 1 , j 2 ) = p ∈ P t independently, and considers the 1 or 2
units of flow f̄ on the link (L ( p ), R ( p )). There are two cases: The first one is when f̄ (L ( p ), R ( p )) = 2. Since the links
(t ) (t ) (t )
(L (U j ), L ( p )) and (L (U j ), L ( p )) both have capacities at most one, one must also have f̄ (L (U j ), L ( p )) =
1 2 1
(t )
f̄ (L (U j ), L ( p )) = 1. This happens for the nodes on the right (i.e., in layers L 6t −2 , L 6t −1 , L 6t ) as well. In this case,
2
we simply ignore the activated facilities in L 6t , and add the activated facilities in L 6t −5 to A t ; the other case is when
f̄ (L ( p ), R ( p )) = 1. Now, the two activated facilities in L 6t −5 and L 6t may not even be in the same bundle. Here, we
exploit the crucial property of p that either j 1 or j 2 is the nearest neighbor of the other in C t , and are able to obtain
a simple criterion of choosing whether the left-activated or right-activated facility location. The cases of singleton pairs
(t )
and facility locations in F t \ j ∈C U j are similar and much easier to handle, since the integral flow is less ambiguous
t
for them.
One may also regard Algorithm 2 as implicitly rerouting some flows of f̄ to make it symmetric on the six layers
L 6t −5 , L 6t −4 , . . . , L 6t , for every t ∈ [ T ] (see Fig. 2). This results in another integral flow that has marginal distributions
different from f̄ , in particular the marginal distributions of A t , t ∈ [ T ]. We first estimate the movement cost of the solution
{ A t }tT=1 . To this end, we have the following lemma.
Lemma 5. Recall that m( A , A ) is the minimum total distance of perfect matching between A and A . If T = 2, the expected movement
cost of solution { A 1 , A 2 } satisfies
(1 )
E[m( A 1 , A 2 )] ≤ d(i , i ) zii .
i∈ F 1 i ∈ F 2
50
(t ) (t )
Fig. 2. An illustration of Algorithm 2, where p = ( j 1 , j 2 ), n j = j 2 and n j = j 1 . We highlight the unit flow on (L ( p ), R( p )) before and after rerouting
1 2
using solid links, and use dashed lines to represent links with no flow on them. The algorithm resolves the ambiguity in f̄ by choosing the right-activated
(t )
facility location in U j , which is equivalent to rerouting the flow through its mirror node on the left in L 6t −5 .
2
If T ≥ 3, the expected movement cost of solution { A t }tT=1 satisfies

⎡ ⎤

E⎣ m( A t , A t +1 )⎦ ≤
(t )
d(i , i ) zii + 6 d(avt ) ( j ).
t ∈[ T −1] t ∈[ T −1] i ∈ F t i ∈ F t +1 t ∈[ T ] j ∈C t
Proof. If T = 2, recall that we choose A 1 , A 2 solely based on f̄ and the links between L 6 and L 7 , thus the matching cost of

A 1 , A 2 is obviously at most i ∈ F 1 ,i ∈ F 2 f̄ (R (i ), L (i ))d(i , i ). Because E[ f̄ ] = f̃ , each flow on link e = (R (i ), L (i )) between
(1)
L 6 and L 7 has expectation zii . Hence the expectation of total movement cost is at most,
(1 )
E[m( A 1 , A 2 )] ≤ E[ f̄ (R (i ), L (i ))] · d(i , i ) = d(i , i ) zii .
i ∈ F 1 ,i ∈ F 2 i∈ F 1 i ∈ F 2
Now we consider the case where T ≥ 3. Denote f̄ ( L ) the multi-set of facilities that f̄ activates in layer L (where
L has an index of 6t − 5 or 6t, t ∈ [ T ]). According to Algorithm 2, | f̄ ( L )| = k for each L ∈ { L 6t −5 , L 6t }, and we define
m P t ( f̄ ( L 6t −5 ), f̄ ( L 6t )) the total distance of a perfect matching between f̄ ( L 6t −5 ) and f̄ ( L 6t ), restricted to the matchings
(t )
within each pair p ∈ P t (note that this does not affect F t \ j ∈C t U j , because f̄ is naturally consistent for them). Using the
triangle inequality, one has
⎡ ⎤ ⎡ ⎤

E⎣ m( A t , A t +1 )⎦ ≤ E ⎣ m( A t , f̄ ( L 6t )) + m( f̄ ( L 6t ), f̄ ( L 6t +1 )) + m( f̄ ( L 6t +1 ), A t +1 )⎦
t ∈[ T −1] t ∈[ T −1]

≤ E m( A t , f̄ ( L 6t )) + m( f̄ ( L 6t −5 ), A t ) + E m( f̄ ( L 6t ), f̄ ( L 6t +1 ))
t ∈[ T ] t ∈[ T −1]
(t )
≤ E m( A t , f̄ ( L 6t )) + m( f̄ ( L 6t −5 ), A t ) + d(i , i ) zii
t ∈[ T ] t ∈[ T −1] i ∈ F t i ∈ F t +1
(t )
≤ E m P t ( f̄ ( L 6t −5 ), f̄ ( L 6t )) + d(i , i ) zii , (1)
t ∈[ T ] t ∈[ T −1] i ∈ F t i ∈ F t +1
where the second to last inequality follows from the same analysis as the case when T = 2. To see the last inequality, notice
that by definition of Algorithm 2, we have A t ⊆ f̄ ( L 6t −5 ) ∪ f̄ ( L 6t ) and
m( A t , f̄ ( L 6t )) + m( f̄ ( L 6t −5 ), A t ) ≤ m P t ( A t , f̄ ( L 6t )) + m P t ( f̄ ( L 6t −5 ), A t ),
when the matchings are also restricted within pairs of P t . Then consider each pair p ∈ P t separately. When f̄ (L ( p ),
R ( p )) = 2, A t always agrees with the left-activated facilities in f̄ ( L 6t −5 ); when f̄ (L ( p ), R ( p )) = 1, A t must agree with an
activated facility either in f̄ ( L 6t −5 ) or in f̄ ( L 6t ). Summing over all pairs, one has the following equality, which verifies (1),
m P t ( A t , f̄ ( L 6t )) + m P t ( f̄ ( L 6t −5 ), A t ) = m P t ( f̄ ( L 6t −5 ), f̄ ( L 6t )).
51
In what follows, we fix t and leave out the superscript for convenience. Denote the random variable t =

m P t ( f̄ ( L 6t −5 ), f̄ ( L 6t )) and t , p the partial matching cost within pair p ∈ P t , thus t = p ∈ P t t , p . To further obtain
an upper bound on t , p , we prioritize the facilities in the same bundle, and only make a cross-bundle matching if we
have to. For example for p = ( j 1 , j 2 ), if f̄ (L ( p ), R ( p )) = 2, with left-activated facilities i 1 ∈ U j 1 , i 2 ∈ U j 2 and right-
activated ones i 1 ∈ U j 1 , i 2 ∈ U j 2 , we (possibly suboptimally) match (i 1 , i 1 ) and (i 2 , i 2 ). We only consider cross-bundle pairs
when the unit flow on the link (L ( p ), R ( p )) passes through different bundles in L 6t −5 and L 6t , in which case we have
d(i 1 , i 2 ) ≤ d(i 1 , j 1 ) + d( j 1 , j 2 ) + d( j 2 , i 2 ) and pay the cost d( j 1 , j 2 ). It is then easy to obtain the following using triangle
inequality (we omit the singleton case here because it is obviously easier),

E[t , p ] ≤ Pr[ f̄ left-activates i 1 , right-activates i 1 ]d(i 1 , i 1 )
i 1 ,i 1 ∈U j 1

+ Pr[ f̄ left-activates i 2 , right-activates i 2 ]d(i 2 , i 2 )
i 2 ,i 2 ∈U j 2

+ Pr[ f̄ (L ( p ), R ( p )) = 1, f̄ activates i 1 , i 2 ]d(i 1 , i 2 )
i 1 ∈U j 1 ,i 2 ∈U j 2

≤ d(i , j 1 ) · Pr[ f̄ (L (i ), L (U j 1 )) = 1] + d(i , j 1 ) · Pr[ f̄ (R (U j 1 ), R (i )) = 1]
i ∈U j 1

+ d(i , j 2 ) · Pr[ f̄ (L (i ), L (U j 2 )) = 1] + d(i , j 2 ) · Pr[ f̄ (R (U j 2 ), R (i )) = 1]
i ∈U j 2
+ Pr[ f̄ (L ( p ), R ( p )) = 1, L 6t −5 , L 6t disagree on U j 1 , U j 2 ] · d( j 1 , j 2 )
≤ 2dav ( j 1 ) + 2dav ( j 2 ) + Pr[ f̄ (L ( p ), R ( p )) = 1, L 6t −5 , L 6t disagree on U j 1 , U j 2 ] · d( j 1 , j 2 ). (2)
Since ( j 1 , j 2 ) is a pair, w.l.o.g. let j 2 be the nearest neighbor of j 1 among C t , R j 1 = 0.5d( j 1 , j 2 ) and dav ( j 1 ) ≥ (1 −
y (U j 1 )) · R j 1 = 0.5(1 − y (U j 1 )) · d( j 1 , j 2 ), hence (1 − y (U j 1 ))d( j 1 , j 2 ) ≤ 2dav ( j 1 ). Furthermore, we can obtain a simple
bound on the probability in (2),
Pr[ f̄ (L ( p ), R ( p )) = 1, L 6t −5 , L 6t disagree on U j 1 , U j 2 ]
= Pr[ f̄ (L (U j 1 ), L ( p )) = 1, f̄ (L ( p ), R ( p )) = 1, f̄ (R ( p ), R (U j 2 )) = 1]
+ Pr[ f̄ (L (U j 2 ), L ( p )) = 1, f̄ (L ( p ), R ( p )) = 1, f̄ (R ( p ), R (U j 1 )) = 1]
≤ Pr[ f̄ (R ( p ), R (U j 1 )) = 0] + Pr[ f̄ (L (U j 1 ), L ( p )) = 0]
≤ 2(1 − y (U j 1 )).
Therefore, the expectation above can be further bounded as
E[t , p ] ≤ 2dav ( j 1 ) + 2dav ( j 2 ) + 2d( j 1 , j 2 ) · (1 − y (U j 1 )) ≤ 6dav ( j 1 ) + 2dav ( j 2 ).

Notice that for each t, P t is a partition of C t ⊆ C t . By summing over all pairs and all time steps, we have

E[t ] ≤ 6 d(avt ) ( j ) ≤ 6 d(avt ) ( j ),
t ∈[ T ] t ∈[ T ] j ∈C t t ∈[ T ] j ∈C t
which yields the lemma when combined with (1).
2.3. Analysis
For the solution { A t : t ∈ [ T ]} given by Algorithm 2, recall that the service cost vector at time t is defined as d(C t , A t ) =
(d( j , A t )) j ∈C t . We first provide a lemma that bounds the expectation of Top norm of the service cost vector at time t,
for any ∈ [|C t |], where the Top norm of a non-negative vector is the sum of its largest entries. Consequently, the
more general ordered cost can be written as a conic combination of Top norms, and easily bounded. A simpler version is
presented as Lemma 4.7 in [11], which is the result of dependent rounding by Charikar and Li [7]. To obtain our following
lemma, we need to further examine the rerouting procedures in Algorithm 2, and consider the new marginal distributions
of A t , t ∈ [ T ].
52
Lemma 6. (Adapted from [11]). Fix t ∈ [ T ] and let ∈ [|C t |], h > 0 be arbitrary. Define d−h ( j , j ) = 0 if d( j , j ) < h and d−h ( j , j ) =
d( j , j ) otherwise. One has
(t )
E[Top (d(C t , A t ))] ≤ 41.33 · h + 41.33 d−h (i , j )xi j .
j ∈C t i ∈ F t
We prove Lemma 6 in Section 2.4. For now, we turn to the case with general weights. We use an argument by Byrka et al.
[11]. Suppose all pair-wise distances are distinct (via small perturbation) and the weight vector w t has N t distinct entries
{ w̄ tr : r = 1, . . . , N t } in decreasing order. For each distinct entry, we guess the exact value T r(t ) which is the smallest distance
(t ) (t )
that is multiplied with w̄ tr in some fixed optimum. Also let T 0 = ∞ and T Nt +1 = 0. Now, the t-th reduced cost function
(t ) (t )
in (LP(D)) is defined as d(t ) (i , j ) = d(i , j ) w̄ ts , where T s ≤ d(i , j ) < T s−1 . After solving the corresponding relaxation (LP(D))
and using Algorithm 2 to obtain the solution { A t : t ∈ [ T ]}, we have the following lemma. The proof is in spirit similar to
Lemma 5.1 in [11], which we provide in Section 2.4.
Lemma 7. When T = 2, { A t : t ∈ [ T ]} is an 82.66-approximation for 2-DOkMed. If T ≥ 3 is a constant and the smallest entry in
{ w t }tT=1 is at least some constant > 0, { A t : t ∈ [ T ]} is an (82.66 + 6γ / )-approximation for T -DOkMed. In both cases, the algorithm
T
t =1 (| F t | · |C t |) time, where N t is the number of distinct entries in the weight vector w t , t ∈ [ T ].
O (Nt )
runs in
Remark 2. We remark that when T ≥ 3, our algorithm only works if the smallest entry in the weight vectors is at least
> 0. This means that our algorithm cannot handle 0-1 vectors. This technical difficulty arises because for 0-1 vectors, the
contribution of a filtered client in C t may be zero (in the inner products) and it is unclear how to bound the extra rerouting
cost in Lemma 5 in terms of the actual contributions of clients. We leave the more general case as an interesting open
question.
One obvious limitation of the algorithm above is the prohibitive running time when N t = ω(1). Here, we use the follow-
ing logarithmic bucketing trick by Aouad and Segev [10] and Byrka et al. [11]: To obtain a constant-factor approximation in
polynomial time, one does not have to strictly adhere to the original weight vectors, and the guessed thresholds only need
to be approximate. Since T is a constant, it is possible to guess the largest service distance for each time step by losing a
polynomial factor in the running time. The connection distances are grouped into logarithmically many buckets thus losing
only a factor of 1 + δ
. For each bucket,
its average weight is also guessed up to a small multiplicative error of δ . There are
at most O log1+δ nδ = O 1δ log nδ buckets for each time step, where n = | F t | + |C t |, therefore guessing a non-increasing

sequence of the average weights only causes another polynomial factor exp O 1δ log nδ = (n/δ) O (1/δ) in the running time.
Finally, because T is a constant, the overall number of guesses is still bounded by a polynomial (n/δ) O ( T /δ) . Formally, we
have the following main theorem adapted from Theorem 6.1 in [11], and provide a simplified proof in Section 2.4.
Theorem 8. When T = 2, for any δ > 0 there exists an 82.66(1 + δ)-approximation algorithm for 2-DOkMed, with running
time (| F 1 | + |C 1 |) O (1/δ) · (| F 2 | + |C 2 |) O (1/δ) . When T ≥ 3 is a constant, and the smallest entry in { w t }tT=1 is at least some con-
stant > 0, for any δ > 0 there exists an (82.66 + 6γ / ) (1 + δ)-approximation algorithm for T -DOkMed, with running time
T
t =1 (| F t | + |C t |)
O (1/δ)
.
2.4. Missing proofs
Proof of Lemma 6. We fix and leave out the superscript t ∈ [ T ] in the following. The proof follows the one of Lemma 4.7
by Byrka et al. in [11], by splitting the service cost of j ∈ C t into a deterministic part D j and another stochastic part
X j , s.t. Pr[d( j , A t ) ≤ D j + X j ] = 1. One major difficulty we need to surmount is the fact that Algorithm 2 displaces some
facilities, and A t no longer follows the original marginal distributions described by Charikar and Li [7]. We emphasize that
our modified matching phase in Algorithm 1 and the structure of the network enable our analysis.
Like the definition of average service cost d av ( j ) = i ∈ F t d(i , j )xi j , we define the average reduced service cost d−h,av ( j ) =
i ∈ F t d−h (i , j )xi j , where d−h is defined as in the lemma. In the following, for any fixed client j, we progressively charge
parts of d( j , A t ) to either D j or X j s.t. d( j , A t ) ≤ D j + X j always holds, where we initially charge 0 to D j and X j . By
Lemma 4, there must exist j ∈ C t (probably j = j when j ∈ C t ) such that d( j , j ) ≤ 4dav ( j ) and dav ( j ) ≤ dav ( j ). Notice by
definition of d−h , it is obvious that dav ( j ) ≤ d−h,av ( j ) + h. Using triangle inequality, one has d( j , A t ) ≤ d( j , A t ) + d( j , j ) ≤
d( j , A t ) + 4dav ( j ) ≤ d( j , A t ) + 4d−h,av ( j ) + 4h. We charge 4h to D j and 4d−h,av ( j ) (a fixed value) with probability 1 to X j ,
so D j = 4h and E[ X j ] = 4d−h,av ( j ) at the moment.
Next, we charge the stochastic service cost d( j , A t ) to D j or X j , and consider the case where j and its nearest neighbor
j ∈ C t are not matched in Algorithm 1 (the case where they are matched is simpler but somewhat different, which will be
explained later). Fix a constant β > 5 which is determined afterwards. We have the following cases.
Case 1: A t ∩ B o ( j , β h) = ∅. We have d( j , A t ) ≤ β h and charge β h to D j .
53

Case 2: A t ∩ U j \ B o ( j , β h) = ∅. We charge this (stochastic) cost to X j . According to Algorithm 2, if j is in a singleton
pair, the marginal distribution over bundle U j is not changed, hence the expectation of this stochastic cost is at
most

xi j d ( i , j ) = xi j d−h (i , j )
i ∈U j \ B o ( j ,β h) i ∈U j \ B o ( j ,β h)

≤ xi j d−h (i , j ) d−h,far ( j ) ≤ d−h,av ( j ), (3)
i ∈ E j \ B o ( j ,β h)
where the last inequality is due to the definition of d−h,av ( j ).

When j is not in a singleton pair, since j is not matched with its nearest neighbor j , it is matched with j̃
in p = ( j , j̃ ), and j is the nearest neighbor of j̃ in C t . The marginal distribution over U j is modified in this case,
because whenever the rounded f̄ disagrees on the bundles U j versus U j̃ in layers L 6t −5 and L 6t (note that this
only happens when f̄ (L ( p ), R ( p )) = 1), Algorithm 2 always chooses the activated facility in U j . Comparing to
the singleton case in (3), there is some extra cost corresponding to the event that f̄ right-activates some facility in
U j \ B o ( j , β h) in L 6t and left-activates some facility in U j̃ in L 6t −5 , and the expectation is obviously no more than
the amount in (3). In summary, the increase of E[ X j ] in this case is at most 2d−h,far ( j ).
Case 3: A t ∩ U j = ∅. Because j is not matched with its nearest neighbor j , j has to be matched with another j in
another pair p = ( j , j ) s.t. d( j , j ) ≤ d( j , j ). Denote R = R j = 12 d( j , j ), since our algorithm guarantees that
there exists an open facility in U j ∪ U j (since we have initial flow f̃ (L ( p ), R ( p )) ≥ 1), j can always be served
by such a facility at a distance
d( j , j ) + d( j , j ) + max{ R j , R j } ≤ 4R + max{ R j , R j } ≤ 5R ,
where the last inequality holds when j is the nearest neighbor of j , in which case R j ≤ R j = 12 d( j , j ) ≤ R, or
1
vice versa and R j ≤ R j = , j ) ≤ R. If R ≤ β h, we simply charge 5R ≤ 5β h to D j , hence we assume R > β h
2
d( j
in the following.
When j is in a singleton pair, we have E j ∩ B o ( j , β h) ⊆ U j from R > β h, and

d−h,far ( j ) = xi j d−h (i , j ) ≥ xi j d−h (i , j ) ≥ R xi j = R (1 − y (U j )); (4)
i ∈ E j \ B o ( j ,β h) i ∈ E j \U j i ∈ E j \U j
according to our rounding algorithm, there is a probability exactly 1 − y (U j ) that none of the facilities in U j is
chosen, in which case we charge 5R to X j , resulting in an increase of expectation of at most 5R (1 − y (U j )) ≤
5d−h,far ( j ) to E[ X j ]. When j is in a normal pair with j̃, j must be the nearest neighbor of j̃, and the marginal
probability of not selecting any facility in U j is even smaller, hence the argument above still holds.
Up till now, we have charged costs to D j and X j in three different cases, such that D j ≤ (4 + 5β)h (the deterministic
part from d( j , A t ) is always at most 5β h) and E[ X j ] ≤ 4d−h,av ( j ) + 7d−h,far ( j ) (by taking the sum over all cases), only for
the case when j is not matched with its nearest neighbor j .
Next, we briefly discuss the case when j is actually matched with j in p = ( j , j ). First, if j is not the nearest neighbor
of j , Case 1 has the same analysis; for Case 2, since the probability that none of the facilities in U j is chosen is in fact
increased, the analysis is similar to (3) with an increase of d−h,far ( j ) to E[ X j ] instead of 2d−h,far ( j ); for Case 3, we must
have f̄ (L ( p ), R ( p )) = 1, and the event A t ∩ U j = ∅ happens with probability at most 3(1 − y (U j )) (considering the event
that f̄ naturally activates only U j in both L 6t −5 and L 6t with prob. at most Pr[ f̄ (L (U j ), L ( p )) = 0] = 1 − y (U j ), and the
event that the disagreement happens with prob. at most 2(1 − y (U j )) using the same argument in the proof of Lemma 5).
Since at least one facility is opened in the bundle pair (U j , U j ) at most d( j , j ) + R j ≤ 3R away from j , we (i) either
charge 3β h to D j if R ≤ β h, (ii) or charge 3R more to X j with probability at most 3(1 − y (U j )) if R > β h, hence an increase
of 9d−h,far ( j ) in expectation by (4). In total, D j ≤ (4 + 3β)h and X j has expectation at most 4d−h,av ( j ) + 10d−h,far ( j ) in this
case. On the other hand, when j and j are matched and they are nearest neighbors of each other, the analysis and upper
bounds are the same.
With all cases taken into consideration, we take the maximum on D j and E[ X j ] and obtain
D j = (4 + 5β)h, E[ X j ] ≤ 4d−h,av ( j ) + 10d−h,far ( j ). (5)
We first assume d( j , j ) > α h for another parameter α ∈ (4, β − 1] to be determined later (recall that β > 5). Then from
αh < d( j , j ) ≤ 4dav ( j ) and d−h,far ( j ) ≤ d−h,av ( j ) ≤ dav ( j ) ≤ dav ( j ) ≤ d−h,av ( j ) + h, it is easy to see that
4 α
h< · d−h,av ( j ) ⇒ d−h,far ( j ) ≤ ·d ( j ),
α−4 α − 4 −h,av
54
so in this case, E[ X j ] is bounded by

10α 14α − 16
E[ X j ] ≤ 4 + d−h,av ( j ) = · d−h,av ( j ). (6)
α−4 α−4
Second, if d( j , j ) ≤ α h, we claim that in the fractional assignment x, j is served to an extent by facilities in F t \ B o ( j , β h)
no less than j does, otherwise we could shift a positive amount of assignment from F t \ B o ( j , β h) into B o ( j , β h) for j ,
resulting in a better assignment of j , which is a contradiction. Therefore, we know for sure that

xi j ≤ xi j .
i ∈ F t \ B o ( j ,β h) i ∈ F t \ B o ( j ,β h)
We then try to modify the assignment of j and get a sub-optimal solution, which helps us relate d−h,far ( j ) to d−h,av ( j ).
In the altered assignment x , only assignments of j are changed, where xi j = xi j for every i ∈ F t \ B o ( j , β h). It is easy to
see that this is possible, simply by shifting some fractional assignment of j in B o ( j , β h) to the outside. Using triangle
inequality, for any i ∈ F t \ B o ( j , β h), d(i , j ) ≥ d(i , j ) − d( j , j ) ≥ (β − α )h ≥ h, so d−h (i , j ) = d(i , j ) and
d−h (i , j ) d (i , j ) d (i , j ) βh β
= ≤ ≤ = .
d−h (i , j ) d (i , j ) d (i , j ) − α h β h − αh β − α
Because x may not be the optimal assignment for j , one has
β β
d−h,far ( j ) ≤ xi j d−h (i , j ) ≤ xi j d−h (i , j ) ≤ · d−h,av ( j ). (7)
β −α β −α
i ∈ F t \ B o ( j ,β h) i ∈ F t \ B o ( j ,β h)
Combining (5)(6)(7), we have

14α − 16 14β − 4α
D j ≤ (4 + 5β)h, E[ X j ] ≤ max , d−h,av ( j ), where β > 5, 4 < α ≤ β − 1.
α−4 β −α
√ √
Plugging in α = 2 + 2 3, β = 4 + 2 3, one has
√ √
D j ≤ (24 + 10 3)h, E[ X j ] ≤ (24 + 10 3)d−h,av ( j ), Pr[d( j , A t ) ≤ D j + X j ] = 1.
Notice that this holds for all clients j ∈ C t . In this case, for the expectation of Top (d(C t , A t )), since we are only paying
for at most of the deterministic values (i.e., D j ’s), it follows that,

E [Top (d(C t , A t ))] ≤ · max D j + E[ X j ] ≤ 41.33 · h + 41.33 d−h,av ( j ).
j
j ∈C t j ∈C t
Proof of Lemma 7. There are | F t | · |C t | possible distinct distances at time t, hence the number of guesses is obvious, and
each guess is associated with a family of reduced cost functions D, which fully dictates (LP(D)) and the running (the
non-stochastic part) of our algorithm.
Assume our guessed thresholds are exactly those in the optimal solution in the following. Let LP be the optimum of
our corresponding (LP(D)), { A t }tT=1 be the stochastic output solution, and { O t }tT=1 be the optimal solution with a total cost
OPT = OPTservice + OPTmove , where OPTservice is the total client service cost and OPTmove is the total facility movement
cost (γ -scaled). We can easily see OPT ≥ LP by considering the LP solution induced by the optimum. Also define SOL =
SOLservice + SOLmove , where the two parts represent the total client service cost and total facility movement cost (γ -scaled)
in our solution.
(t ) (t )
For each t ∈ [ T ], r ∈ [ N t ], let I r be the largest index in w t that has value w t I (t ) = w̄ tr . Using Lemma 6 with h = T r and
r
(t )
= I r , we have

(t ) (t ) (t )
E Top I (t ) d(C t , A t )) ≤ 41.33I r T r + 41.33 d− T (t ) (i , j )xi j , (8)
r r
j ∈C t i ∈ F t
and we decompose the ordered cost w t d(C t , A t )↓ as a conic combination of Top norms (see, e.g., [11]),

Nt
E[ w t d(C t , A t )↓ ] = ( w̄ tr − w̄ t (r +1) ) · E Top I (t ) d(C t , A t )
r
r =1
⎛ ⎞

Nt
≤ ( w̄ tr − w̄ t (r +1) ) ⎝41.33I r T r + 41.33
(t ) (t )
d− T (t ) (i , j )xi j ⎠
( t )
r
r =1 j ∈C t i ∈ F t
55

Nt
Nt
= 41.33 ( w̄ tr − w̄ t (r +1) ) I r(t ) T r(t ) + 41.33 ( w̄ tr − w̄ t (r +1) ) (t )
d− T (t ) (i , j )xi j . (9)
r
r =1 r =1 j ∈C t i ∈ F t
Since our guessed thresholds are correct, in the inner product OPTt = w t d(C t , O t )↓ of the optimum, every weight that
(t )
is equal to w̄ tr should be multiplied by a distance at least T r , thus one has

Nt
(t ) (t ) (t )

Nt
(t ) (t )

Nt
(t ) (t )
OPTt ≥ w̄ tr ( I r − I r −1 ) T r ≥ w̄ tr I r T r − w̄ tr I r −1 T r −1
r =1 r =1 r =1
Nt

Nt
(t ) (t ) (t ) (t ) (t ) (t )
≥ w̄ tr I r T r − w̄ t (r +1) I r T r = ( w̄ tr − w̄ t (r +1) ) I r T r , (10)
r =1 r =1
(t )
where I 0 = 0, OPTt denotes the service cost of all clients in C t , and thus t ∈[ T ] OPTt = OPTservice .
Further, for the second term in (9), we have

Nt (t )
(t )

Nt
( w̄ tr − w̄ t (r +1) ) d− T (t ) (i , j )xi j = xi j · ( w̄ tr − w̄ t (r +1) )d−T (t ) (i , j )
r r
r =1 j ∈C t i ∈ F t j ∈C t i ∈ F t r =1
(t )

Nt
= xi j · d ( i , j ) ( w̄ tr − w̄ t (r +1) )
j ∈C t i ∈ F t (t )
r : T r ≤d(i , j )
(t ) (t )
= xi j · d (i , j ). (11)
j ∈C t i ∈ F t
(t )
To see the last equality, notice that the sum is taken over all r s.t. T r ≤ d(i , j ). Since T r(t ) is decreasing and T N(tt)+1 = 0, this
(t )
sum is equal to w̄ tr0 s.t. r0 = min{r ∈ [ N t ] : T r ≤ d(i , j )}, and w̄ tr0 d(i , j ) = d(t ) (i , j ) by definition.
Combining (9)(10)(11), our algorithm outputs a stochastic solution { A t }t ∈[ T ] such that the expected service cost of all
clients is at most
T

T
T
↓ (t )
E[SOLservice ] = E w t d(C t , A t ) ≤ 41.33 OPTt + 41.33 xi j · d(t ) (i , j ). (12)
t =1 t =1 t =1 j ∈ C t i ∈ F t
Meanwhile by Lemma 5, the movement cost of all facilities in { A t }t ∈[ T ] has an expectation of at most

T −1
(t )
E[SOLmove ] ≤ γ d(i , i ) zii + 6γ · 1[ T ≥ 3] d(avt ) ( j ). (13)
t =1 i ∈ F t i ∈ F t +1 t ∈[ T ] j ∈C t
By combining (12), (13), if T = 2, our overall cost is

T
T −1
(t ) (t )
E[SOL] ≤ 41.33OPTservice + 41.33 xi j · d(t ) (i , j ) + γ d(i , i ) zii
t =1 j ∈ C t i ∈ F t t =1 i ∈ F t i ∈ F t +1

T
T −1
(t ) (t )
≤ 41.33OPT + 41.33 xi j · d(t ) (i , j ) + 41.33γ d(i , i ) zii
t =1 j ∈ C t i ∈ F t t =1 i ∈ F t i ∈ F t +1
≤ 41.33OPT + 41.33LP ≤ 82.66OPT.

For the case where T ≥ 3, recall that we assume the smallest entry in { w t }t ∈[ T ] is at least > 0. Using Lemma 5, we
have the following,

T
T −1
(t ) (t )
E[SOL] ≤ 41.33OPTservice + 41.33 xi j · d(t ) (i , j ) + γ d(i , i ) zii
t =1 j ∈ C t i ∈ F t t =1 i ∈ F t i ∈ F t +1
(t )
+ 6γ xi j · d ( i , j )
t ∈[ T ] j ∈C t i ∈ F t

T
T −1
(t ) (t )
≤ 41.33OPT + 41.33 xi j · d(t ) (i , j ) + γ d(i , i ) zii
t =1 j ∈ C t i ∈ F t t =1 i ∈ F t i ∈ F t +1
56
6γ (t ) (t )
+ xi j · d ( i , j )
t ∈[ T ] j ∈C t i ∈ F t

6γ 6γ
≤ 41.33OPT + 41.33 + LP ≤ 82.66 + OPT,

where the second inequality is because the smallest entry in { w t }t ∈[ T ] is at least , thus our reduced cost function d(t )
always satisfies d(t ) ≥ d for any t, since d(t ) (i , j ) = d(i , j ) w̄ ts ≥ d(i , j ) by definition.
Proof of Theorem 8. To start with, we fix a small constant δ > 0 and consider a slightly modified instance, where instead of
using the weight vector w t , we define w t where w ts = max{ w ts , δ w t1 /|C t |}. It is easy to see that for any non-negative
|C |
w t v ≥ w t v, and
non-increasing vector v ∈ R≥0t , one has w t v − w t v ≤ δ w t1 v 1 ≤ δ w t v, which can be rewritten as
w t v ∈ [ w t v , (1 + δ) w t v ]. Therefore, by solving T -DOkMed on these new weight vectors, we lose a factor of at most

1 + δ in the approximation ratio. Suppose we do this in the sequel.
(t )
Since T = O (1), we also assume that the maximum connection distance T max ≥ 0 is known to us for any t ∈ [ T ] in a
fixed optimal solution. This is done via simple exhaustive search, and we only lose a polynomial factor in the running time.
Fix t and let M t = log1+δ (|C t |/δ). For r = 1, 2, . . . , M t , define the real intervals

J r = (1 + δ)−( Mt −r +1) T max , (1 + δ)−( Mt −r ) T max ,
(t ) (t )
and let J 0 = [0, (1 + δ)− Mt T max ]. We have that J = { J 0 , . . . , J Mt } is a partition of the interval [0, T max ], i.e., the interval that
(t ) (t )
all connection distances at t fall in, with respect to the optimum.

Denote o the sorted service cost vector at time t in the optimum, and I J = {s ∈ [|C t |] : o s ∈ J } the indices of entries of o
(t ) (t )
in J . Obviously, I J Mt is non-empty, because our guess T max is correct and T max ∈ J Mt . For r = M t , M t − 1, . . . , 1, 0, iteratively
define the average weight in the following way (we omit the superscript (t ) here for simplicity),

1
|I Jr | r∈I Jr
w tr I J r = ∅
wr = (14)
w r +1 I J r = ∅,
and this is well-defined since I J Mt is always non-empty. Additionally, note that w = { w 0 , . . . , w Mt } is non-decreasing and
they are all in the interval [δ w t1 /|C t |,
w t1 ]. Although we have no knowledge of the optimum or the exact values of these
average weights, it is possible for us to approximately guess their values. To achieve this, recall that M t = log1+δ (|C t |/δ)
and all entries of w are in the interval [δ w t1 /|C t |,
w t1 ]. Therefore, if we guess each w r to its closest and no smaller power
of 1 + δ (hence w r ∈ [ w r , (1 + δ) w r ]), there are at most 2 + log1+δ (|C t |/δ) possible candidates, and they must also be
non-decreasing as w does. Using a basic counting method, the number of such non-decreasing sequences is at most
exp( O (log1+δ (|C t |/δ))) = (|C t |/δ) O (1/δ) , thus bounded by a polynomial.
(t )
We define the reduced cost functions now. For each t ∈ [ T ], we guess the maximum connection distance T max and the
(t ) (t )
approximate values of average weights w0 , . . . ,
w (t ) = { w Mt }. Suppose the guesses are correct in what follows. The reduced
(t )
cost function d(t ) is defined as d(t ) (i , j ) =
w r · d(i , j ) where d(i , j ) ∈ J r , r = 0, . . . , M t . We notice that d(t ) is undefined
(t ) (t )
for distances larger than T max , and this is easily handled by explicitly adding to (LP(D)) the constraints xi j = 0 for all
(t )
i ∈ F t , j ∈ C t s.t. d(i , j ) > T max . Fix t in the following, and for each s ∈ [|C t |], define T s = sup( J ) where o s ∈ J and J ∈ J . We
(t )
obviously have (1 + δ)o s ≥ T s if s ≥ 1, and otherwise we have o s ∈ J 0 and thus o s ≤ T s ≤ δ T max /|C t |. Using Lemma 6 with
h = T s and = s, s = 1, 2, . . . , |C t |,
|C t |

w t d(C t , A t )↓ ] =
E[ w ts −
( w t (s+1) )E[Tops (d(C t , A t ))]
s =1
⎛ ⎞
|C t |

w t (s+1) ) ⎝41.33s · T s + 41.33 d − T s ( i , j ) xi j ⎠
(t )
≤ w ts −
(
s =1 j ∈C t i ∈ F t
(t )
δ T max
≤ 41.33(1 + δ) w ts −
( w t (s+1) )s · o s + 41.33 w ts −
( w t ( s +1 ) ) s ·
|C t |
s:o s ∈ J ≥1 s:o s ∈ J 0
|C t |
(t )
+ 41.33 w ts −
( w t ( s +1 ) ) d − T s ( i , j ) xi j
s =1 j ∈C t i ∈ F t
|C t |
(t )
≤ 41.33(1 + δ) w ts −
( w t (s+1) )Tops (o) + 41.33δ
w t1 T max
s =1
57
|C t |
(t )
+ 41.33 w ts −
( w t ( s +1 ) ) d − T s ( i , j ) xi j
s =1 j ∈C t i ∈ F t
|C t |

w t o + 41.33
(t )
≤ 41.33(1 + 2δ) w ts −
( w t ( s +1 ) ) d − T s ( i , j ) xi j . (15)
s =1 j ∈C t i ∈ F t
To bound (15), we notice that its first part is exactly 41.33(1 + 2δ)OPTt , where OPTt is the total ordered service cost
paid in the optimum at time t. For the second part, we obtain
|C t |
|C t |

(t ) (t )
w ts −
( w t ( s +1 ) ) d − T s ( i , j ) xi j = w ts −
( w t ( s +1 ) ) x i j d − T s ( i , j )
s =1 j ∈C t i ∈ F t j ∈ C t i ∈ F t s =1
(t )

= xi j d ( i , j ) w ts −
( w t (s+1) ). (16)
j ∈C t i ∈ F t s∈[|C t |], T s ≤d(i , j )
(t ) (t )
We group all d(i , j ) s.t. d(i , j ) ∈ J 0 together, and notice that d(i , j ) ≤ δ T max /|C t |. Since we have i∈ F t xi j = 1 for all
(t )
j ∈ C t , the total contribution of all such pairs (i , j ) in (16) is at most δ w t1 T max ≤ δ OPTt . Now suppose d(i , j ) ∈ J r0 , r0 ≥ 1,
and for a technical reason, we assume that none of the input distances coincides with any boundary of the intervals. This is
achieved by slightly increasing the boundaries of all intervals (also see [11] about this trick). If there is no such T s ≤ d(i , j ),
the corresponding contribution in (16) is zero, while d(t ) (i , j ) = w r0 · d(i , j ) is non-negative; otherwise, the sum of weights
is the maximum w ts s.t. o s ∈ J <r0 , which is at most the average weight on J r0 and hence at most w r0 ≤ w r0 . This implies
that (16) is at most
|C t |

(t ) (t )
w ts −
( w t ( s +1 ) ) d − T s ( i , j ) xi j ≤ xi j d(t ) (i , j ) + δ OPTt , (17)
s =1 j ∈C t i ∈ F t j ∈C t i ∈ F t
and the first part is exactly the service cost of time step t in the objective of (LP(D)).
Finally, via considering the solution (x , y , z ) induced by the optimum of the original problem, the optimum of the
relaxation (LP(D)) under the aforementioned reduced cost functions is at most,

T
T −1
(t ) (t )
LP ≤ d(t ) (i , j )xi j +γ d(i , i ) zii
t =1 j ∈ C t i ∈ F t t =1 i ∈ F t i ∈ F t +1
|C t |
T
(t ) (t )
≤
w r s · o s + OPTmove
t =1 s =1
| Mt |
T (t )
≤ (1 + δ) w r(t ) o s + OPTmove
t =1 r =0 (t )
s∈ I J
r
⎞⎛ ⎛ ⎞
T ⎜ ⎟⎜ 1 (t ) ⎟
≤ (1 + δ) ⎝
w ts ⎠ ⎝ (t ) o s ⎠ + OPTmove
t =1 (t ) (t ) | I J r | (t )
r :| I J |=∅ s∈ I J s∈ I J
r r r

T (t )

T
(t )
≤ (1 + δ)2
w ts o s + (1 + δ) w t1 T max + OPTmove
δ
t =1 r :| I (t ) |=∅ s∈ I (t ) t =1
J r J r

T
≤ (1 + 3δ + 2δ 2 ) w t o(t ) + OPTmove ≤ (1 + 3δ + 2δ 2 )OPT,
(18)
t =1
(t ) (t )
where
w r s denotes the guessed average weight where o s ∈ J r s . We note that we again consider J 0 and J >0 separately.
(t )
For the former, the total contribution is easily bounded using sup( J 0 ) ≤ δ T max /|C t |; for the latter, since the distances falling
into the same interval differ by a multiplicative factor of at most δ , we have the above inequality.
By combining (15), (17), (18) and taking the sum over t ∈ [ T ], the rest follows from the same analysis in Lemma 7, thus
omitted here.
58
3. Approximating DkSup
In this section, we present our various results on DkSup. We start by showing two hardness results for T -DkSup and
T -DkSupOut when T ≥ 3, then present a simple 3-approximation for 2-DkSup and another more involved multi-criteria
approximation algorithm for 2-DkSupOut in the more nuanced outlier setting. We finish this section with a hardness result
on pure approximations for 2-DkSupOut.
3.1. The hardness of approximating 3-DkSup and 3-DkSupOut
As a warm-up, we use a simple argument to show that when the number of time steps is at least 3, T -DkSup admits
no non-trivial approximations unless P = NP. The proof uses the reduction from the perfect 3D-matching problem, which is
known to be NP-complete [45].
Theorem 9. If T ≥ 3, there are no polynomial-time algorithms for T -DkSup with non-trivial approximation factors, unless P = NP.
Proof. Notice that we only need to prove the hardness of 3-DkSup, since this is a special case of T > 3 by setting C t = ∅, t ≥
4. We reduce an arbitrary instance of perfect 3D-matching to 3-DkSup. Recall that in an instance of perfect 3D-matching,
we are given three finite ground sets A , B , C with | A | = | B | = |C |, and a triplet set T ⊆ A × B × C . Suppose | A | = n and
|T | = m, and one needs to decide whether there exists a subset S ⊆ T , such that |S| = n, and each element in A ∪ B ∪ C
appears exactly once in some triplet in S. W.l.o.g., we assume that each element in A ∪ B ∪ C appears in at least one triplet,
otherwise the answer is trivially negative.
Assume there is an α -approximation for 3-DkSup with factor α > 1. We construct the following graph G = ( V A ∪ V B ∪
V C , E ), where the vertex set and edge set are initially empty.
Step 1. For each triplet g = (a, b, c ) ∈ T where a ∈ A , b ∈ B , c ∈ C , add three new vertices V A ← V A ∪ {a g }, V B ← V B ∪
{b g }, V C ← V C ∪ {c g }. Connect (a g , b g ) with an edge of length α and add it to E. Connect (b g , c g ) with an edge of
length α and add it to E.
Step 2. For any two vertices a g ∈ V A , a g ∈ V A corresponding to the same element a ∈ A, connect them with an edge of
length 1 and add the edge to E, thus forming a clique K a with unit-length edges. Repeat the same procedure for
V B, VC .
An illustration is shown in Fig. 3. We approximately solve 3-DkSup on G with its graph metric d G (set d G (s, t ) = ∞ if
they are in different connected components), with k = n and the movement constraint B = α , where the client sets are
C 1 = V A , C 2 = V B , C 3 = V C and facility sets are F 1 = V A , F 2 = V B , F 3 = V C . W.l.o.g., this 3-DkSup instance has optimum at
least 1, otherwise one has m = n and the original instance is again trivial.
Now, if the original perfect 3D-matching instance is feasible, by letting each of the n mobile facilities move along the n
trajectories induced by the feasible solution, we have a solution to the 3-DkSup instance on G with objective exactly 1, since
each clique K u corresponding to some u ∈ A ∪ B ∪ C is formed using unit-length edges, and there is exactly one open facility
in every such clique. Conversely, if the 3-DkSup instance has optimum 1, it is easy to see that there must be exactly one
vertex chosen as an open facility in each clique K u , u ∈ A ∪ B ∪ C , and since the movement constraint is B = α , each open
facility must follow the edges (a g , b g ) and (b g , c g ) for some g ∈ T when it moves. Thus, this optimum directly induces a
feasible solution to the original perfect 3D-matching instance.
By noticing that for any s = t ∈ V A , either d G (s, t ) = 1 or d G (s, t ) ≥ 2α , one has that the optimum of the 3-DkSup
instance is either equal to 1 or no less than 2α . Hence any α -approximation for 3-DkSup implies an efficient algorithm that
determines the feasibility of the perfect 3D-matching instance, which is prohibited by its NP-completeness.
From the result above, the hardness of pure approximations for T -DkSupOut with T ≥ 3 directly follows, since T -DkSup
is only its special case. We strengthen this result by showing that, even by allowing multi-criteria approximations, i.e.,
violating the outlier constraints by some small -fraction, one is still unable to obtain any non-trivial approximation factors
for T -DkSupOut, T ≥ 3. The result is formally given in the following theorem, where we use the reduction from maximum
3D-matching, which is known to be APX-complete [46,47].
Theorem 10. If T ≥ 3, there exists a constant 0 ∈ (0, 1), such that T -DkSupOut admits no multi-criteria (α , 1 − 0 , . . . , 1 − 0 )-
approximations for any non-trivial factor α > 1, unless P = NP.
Proof. Again, we only need to prove the case of 3-DkSupOut, since this is a special case of larger T ’s. For a maximum
3D-matching instance on ground sets A , B , C (not necessarily with identical cardinalities) and a triplet set T ⊆ A × B × C ,
we are asked to find S ⊆ T such that the triplets in S are pairwise disjoint and |S| is maximized.
By way of contradiction, assume we have for some α > 1 and any constant > 0, an efficient multi-criteria (α , 1 −
, 1 − , 1 − )-approximation for 3-DkSupOut. Given any maximum 3D-matching instance, we first create the same graph
59
Fig. 3. An illustration of the reduction from perfect 3D-matching to 3-DkSup, with the instance A = {a1 , a2 , a3 }, B = {b1 , b2 , b3 }, C = {c 1 , c 2 , c 3 }, and
T = {(a1 , b1 , c 1 ), (a2 , b3 , c 2 ), (a3 , b2 , c 3 ), (a1 , b1 , c 2 ), (a3 , b2 , c 1 ), (a3 , b3 , c 2 )}. Every intra-clique edge has length 1, and we highlight a feasible solution with
objective value 1.
G = ( V A ∪ V B ∪ V C , E ) as in the proof for Theorem 9. Let m = maxu ∈ A ∪ B ∪C | K u | be the maximum size of the cliques in G,
i.e., the maximum number of occurrences of elements in the triplets, and add additional dummy vertices to each clique
so that each clique has size exactly m. Evidently, any intra-clique distance is 1 and any inter-clique distance is at least 2α
(within the same time step).
Now, we run the approximation algorithm for 3-DkSupOut on G with C 1 = F 1 = V A , C 2 = F 2 = V B and C 3 =
F 3 = V C and movement constraint B = α as the input. The algorithm is run for multiple times, for the parameter
k = 1, 2, . . . , min{| A |, | B |, |C |} and outlier constraints l1 = l2 = l3 = mk each time.
Suppose the optimum for the maximum 3D-matching instance is k0 ∈ Z+ , and it is easy to see that this optimum induces
a solution to the 3-DkSupOut instance on G that utilizes k0 open facilities and covers exactly mk0 clients within distance 1
at each time step, using intra-clique unit-length edges. Because our algorithm is an (α , 1 − , 1 − , 1 − )-approximation,
when we run it with parameter k = k0 , it must be able to output a solution that uses k0 open facilities and covers at least
(1 − )mk0 clients within distance α at each time step, and this solution naturally induces a subset S ⊆ T with |S| = k0 ,
according to the definition of G. Recall that inter-clique distances are at least 2α , hence this solution must use only intra-
clique edges. We say that u ∈ A ∪ B ∪ C is “matched” in this solution, if and only if some vertex in the clique K u is selected
as an open facility location. Thus the numbers of elements in A , B and C that are matched in this 3-DkSupOut solution are
all at least (1 − )k0 m/m = (1 − )k0 , respectively.
This means that S with |S| = k0 , which may have intersecting triplets, matches at least (1 − )k0 elements in A , B and
C . To further distill the induced subset S ⊆ T and remove intersecting triplets, we iteratively remove some arbitrary g ∈ S
and assign S ← S \ { g }, whenever g intersects some other triplet in S. It is easy to see that, at most 3 k0 triplets are
removed in this process, and we obtain |S| ≥ (1 − 3 )k0 as a feasible 3D-matching solution. This happens for any constant
> 0, hence a de facto PTAS for maximum 3D-matching. But unless P = NP, this is impossible since maximum 3D-matching
is APX-complete. Therefore, our initial assumption is incorrect, which yields the theorem as its contrapositive assertion.
Remark 3. One can easily see that the two hardness results above are also valid if we only allow the DkSup and DkSupOut
solutions to consist of subsets of open facilities at each time step, instead of multi-sets.
3.2. A 3-approximation for 2-DkSup
In contrast to the hardness of approximating T -DkSup for T ≥ 3, we consider 2-DkSup on general metrics and present a
simple flow-based 3-approximation. Suppose that we have successfully guessed the optimum R (using binary search). We
construct the following network flow instance G(V, E). V consists of 4 layers of vertices (two layers L11 , L12 for t = 1,
two layers L21 , L22 for t = 2), a source s and sink t. We define G as follows:
Step 1. For each i ∈ F 1 , add a vertex in L12 . For i ∈ F 2 , add a vertex in L21 .
Step 2. Repeatedly pick an arbitrary client j ∈ C 1 and remove from C 1 every client within distance 2R from j. Denote these
clients a new cluster corresponding to j and call j a cluster center. Since R is optimal, it is easy to see we obtain
at most k such clusters (otherwise there exist k + 1 clients with pair-wise distance > 2R , and every open facility
can obviously cover at most one of them, which is a contradiction). And if there are less than k clusters, we create
some extra dummy clusters to obtain exactly k clusters, while dummy clusters do not represent any clients. For each
cluster, add a vertex to L11 . Repeat this for C 2 and form L22 .
Step 3. The four layers are arranged in order as L11 , L12 , L21 , L22 . With a slight abuse of notation, for a non-dummy
cluster center u ∈ L11 and facility location v ∈ L12 , connect them using a link with unit capacity if d(u , v ) ≤ R ; for
facility location w ∈ L21 and a non-dummy cluster center z ∈ L22 , connect them using a link with unit capacity if
d( w , z) ≤ R . For v ∈ L12 , w ∈ L21 (both are facility locations), connect them using a link with unbounded capacity
if d( v , w ) ≤ B.
60
Step 4. Connect every dummy cluster in L11 with every facility location vertex in L12 . Connect every dummy cluster in
L22 with every facility location vertex in L21 . Every such link has unit capacity.
Step 5. Finally, the source s is connected to every vertex in L11 and the sink t is connected to every vertex in L22 , with
every edge having unit capacity.
Lemma 11. G(V, E) admits a flow of value k. Moreover, we can obtain a feasible solution of cost at most 3R from an integral flow of
value k in G(V, E).
Proof. Consider an optimal solution A 1 ⊆ F 1 , A 2 ⊆ F 2 with objective R , and there exists a perfect matching between the
two. For any i ∈ F 1 , i ∈ F 2 , if the pair (i , i ) appears m times in the perfect matching, define a flow value f (i , i ) = m over
link (i , i ). Notice that this is well-defined, since one must have d(i , i ) ≤ B and the corresponding vertices are linked in the
network.
Consider the first time step. For any facility location i and two non-dummy cluster centers j , j , either d(i , j ) or d(i , j )
is larger than R , otherwise d( j , j ) ≤ 2R using triangle inequality, contradicting with our construction. This means that the
subsets of facility locations linked with the non-dummy clusters are disjoint. Next, since A 1 covers all j ∈ C 1 with radius
R , for every cluster center j ∈ L11 , we can always find a distinct open facility i ∈ A 1 such that d(i , j ) ≤ R , and add a unit
flow as f ( j , i ) = 1. The same process is repeated for L22 and A 2 .
The total flow between L12 and L21 is now obviously k, since the perfect matching between A 1 and A 2 has size k. After
the construction of unit flows for non-dummy clusters, we need to satisfy the flow conservation constraint at all facility
vertices. Arbitrarily direct the remaining flows from facility vertices to dummy clusters, one unit each time, and this always
satisfies the flow conservation constraint at facility vertices since the total numbers of clusters are both k in L11 and L22 .
Finally, for every cluster with unit flow, define the flow value between it and the source/sink as 1. This completes an integral
flow of value k on G.
For the second assertion in the lemma, suppose we have an integral flow f̄ of value k on G. For any facility location
i ∈ F 1 , denote g (i ) the total flow through i. We place g (i ) facilities at location i, and repeat the same procedures for i ∈ F 2 .
If f̄ (i , i ) = m for i ∈ F 1 , i ∈ F 2 , move m facilities from i to i in the transition between the two time steps.
For every j ∈ C 1 , if j is the cluster center it belongs to, using triangle inequality and the fact that j has unit flow on
its vertex, the nearest open facility for j is at a distance at most R , and there exists an open facility at most d( j , i ) ≤
d( j , j ) + d( j , i ) ≤ 3R away from j .
Theorem 12. There exists a 3-approximation for 2-DkSup.
Proof. Consider the aforementioned network flow instance. It only has integer constraints and the coefficient matrix is
totally unimodular. Moreover, there exists a flow of value k due to Lemma 11, hence we can efficiently compute an integral
flow f̄ of value k, thus obtaining a 3-approximate solution.
3.3. The hardness of pure approximations for 2-DkSupOut
In this section, we show that 2-DkSupOut admits no pure approximations for any non-trivial factor α > 1, i.e., multi-
criteria (α , 1, 1)-approximation algorithms, complementing our multi-criteria (3, 1 − , 1 − )-approximation in the next
section for every > 0.
Theorem 13. There are no polynomial-time multi-criteria (α , 1, 1)-approximation algorithms for 2-DkSupOut for any non-trivial
factor α , unless P = NP.
Proof. We use the reduction from maximum satisfiability (MAX-SAT), which is known to be NP-hard. By way of contra-
diction, assume there exists a pure α -approximation for 2-DkSupOut. For an arbitrary SAT instance with boolean variables
x1 , x2 , . . . , xn and m clauses s1 ∧ s2 ∧ · · · ∧ sm in conjunctive normal form (CNF), we create a quadripartite graph G and
subsequently define the 2-DkSupOut instance on the graph using its graph metric d G .
G is initially empty. First create the client sets by adding n + m vertices to G, where C 1 = {x1 , . . . , xn } and
C 2 = {s1 , . . . , sm }, with a slight abuse of notation. For the facility sets, add 4n new vertices to G, where F 1 =
{ y 1,0 , y 1,1 , . . . , yn,0 , yn,1 } and F 2 = { z1,0 , z1,1 , . . . , zn,0 , zn,1 }. Here, y i ,0 and zi ,0 represent the choice of setting the variable
xi = 0, and y i ,1 , zi ,1 are similarly defined.
With the 5n + m vertices added to G, we proceed to define the edges. For each variable xi , link (xi , y i ,0 ) and (xi , y i ,1 )
using an edge with length 1. For each variable xi , also link ( y i ,0 , zi ,0 ) and ( y i ,1 , zi ,1 ) using an edge with length α . Finally,
for each clause s j and each literal in s j , link s j to the corresponding “literal vertex” in F 2 that satisfies s j , using an edge
with length 1. For example, for clause s j = x1 ∨ x2 ∨ ¬x3 , we link (s j , z1,1 ), (s j , z2,1 ) and (s j , z3,0 ). Also see an illustration in
Fig. 4. Evidently, any facility-client connection distance is either 1 or ≥ 2α + 1 on this graph metric, within the same time
step.
61
Fig. 4. An illustration of the reduction from MAX-SAT to 2-DkSupOut, where the boolean formula in CNF is (¬x1 ∨ ¬x2 ) ∧ (¬x2 ∨ x3 ) ∧ (x1 ∨ x2 ∨ ¬x3 ). We
use solid lines to represent edges with length 1, and dashed lines for edges with length α .
We run our pure α -approximation on G, with the client sets and facility sets defined as above, and k = n, B = α . For the
outlier constraints, we fix l1 = n and try l2 = 1, 2, . . . , m. Suppose the optimum for the MAX-SAT instance is m0 , and it is
easy to see that, this optimum induces a solution to the 2-DkSupOut instance with objective 1 when we have l2 = m0 . Since
any facility-client connection is either 1 or ≥ 2α + 1, our pure α -approximation in fact solves MAX-SAT exactly by trying
all possible values of l2 , in particular when l2 = m0 . But this is impossible unless P = NP, hence our initial assumption is
incorrect, yielding the theorem.
As a simple corollary, one can also show the hardness of guaranteed multi-criteria (α , 1, 1 − ) or (α , 1 − , 1)-
approximations for 2-DkSupOut with some non-trivial α and any constant > 0, using the identical reduction from the
above and the fact that MAX-SAT is APX-hard [48]. Here, a guaranteed multi-criteria approximation means that the factors
are strictly aligned with the time steps, regardless of the given instance.
Corollary 14. There exists a constant 0 ∈ (0, 1), such that 2-DkSupOut admits no guaranteed multi-criteria (α , 1, 1 − 0 ) or (α , 1 −
0 , 1)-approximations for any non-trivial factor α > 1, unless P = NP.
3.4. A multi-criteria approximation for 2-DkSupOut
In this section, we present our more involved multi-criteria (3, 1 − , 1 − )-approximation for 2-DkSupOut. At a high
level, the algorithm first guesses the optimal objective via binary search. It then guesses a small portion of the unknown
optimal solution that “covers” as many clients as possible (Section 3.4.1), and attempts to formulate the remaining problem
as finding a bipartite matching that attains good coverage of clients on both sides (Section 3.4.2).
Before we start, we formulate the problem in an alternative and more convenient way. In the original instance, the client
sets C 1 , C 2 and candidate facility locations F 1 , F 2 all live in the same finite metric space ( X , d), while a facility can move
to any other location as long as the distance is no more than some given threshold B. In this section, we create two copies
( X (1) , d) and ( X (2) , d) of the original space ( X , d), and consider C 1 , F 1 ⊆ X (1) , C 2 , F 2 ⊆ X (2) instead. Now, a facility i ∈ F 1 can
move to i ∈ F 2 if and only if their corresponding locations in X satisfy d(i , i ) ≤ B. Obviously, this formulation is equivalent
to the original one. In the sequel, we also assume that we have correctly guessed the optimum R via binary search, since
it only has a polynomial number of possibilities.
3.4.1. Guessing and modifying the instance

We start with some pre-processing, including some must-have facility locations and clients and excluding some others.
The intuition is to guess a constant-size “heavy” subset of the optimal solution, and consider the problem on the reduced
instance. Denote B 1 (i , R ) = { j ∈ C 1 : d(i , j ) ≤ R }, B 2 (i , R ) = { j ∈ C 2 : d(i , j ) ≤ R }. Fix some constant γ > 0 (its value will be
determined later), we do the following to modify the original instance.
Step 1. Enumerate all possible choices of γ / distinct facilities in F 1 and γ / distinct facilities in F 2 . Denote the two chosen
sets of facilities T 1 , T 2 . Additionally, enumerate all possible movements associated with T 1 and T 2 , i.e., g : T 1 → F 2
and h : T 2 → F 1 such that ∀i ∈ T 1 , d(i , g (i )) ≤ B , ∀i ∈ T 2 , d(i , h(i )) ≤ B.
Step 2. Recursively sort T 1 , each time by choosing the unchosen i ∈ T 1 such that B 1 (i , 3R ) covers the largest number
(1) (1)
of uncovered clients in C 1 , and letting this number be u i . Denote u 0 the number of remaining clients that
(2) (2)
are covered by i ∈h( T 2) B 1 (i , 3R ). Same for T 2 and C 2 , and we have u 0 and u i for i ∈ T 2 . Denote C 1 = C 1 \

i ∈ T 1 ∪h( T 2 ) B 1 (i , 3R ) and C 2 = C 2 \
i ∈ T 2 ∪ g ( T 1 ) B 2 (i , 3R ) , i.e., the clients that are not covered by any of the
closed balls above.
62
(1) (2)
Step 3. Remove any i ∈ F 1 s.t. | B 1 (i , 3R ) ∩ C 1 | > mini ∈ T 1 {u i } and any i ∈ F 2 s.t. | B 2 (i , 3R ) ∩ C 2 | > mini ∈ T 2 {u i }. Denote
the remaining two facility sets F 1 and F 2 . Notice that each removed i ∈ F 1 is not co-located with T 1 ∪ h( T 2 ), and
each removed i ∈ F 2 is not co-located with T 2 ∪ g ( T 1 );
Step 4. Denote the reduced problem P , where the client sets are C 1 , C 2 , the facility sets are F 1 , F 2 and the new outlier
! !
(1) (1) (2) (2)
constraints are l1 = max l1 − u 0 − i∈T 1 u i , 0 and l2 = max l2 − u 0 − i∈T 2 u i , 0 , respectively.

(1)
Lemma 15. In the reduced problem P , every facility location i ∈ F 1 covers at most γ |C 1 \ C 1 | − u 0 clients with radius 3R , and

(2)
every facility location i ∈ F 2 covers at most γ |C 2 \ C 2 | − u 0 clients with radius 3R .
Proof. First, if i ∈ F 1 is co-located with some location in T 1 ∪ h( T 2 ), the number of clients it covers in C 1 is obviously zero,
since they are all removed from C 1 . Otherwise, according to our construction of F 1 , for any i ∈ F 1 that are not co-located
(1)
with T 1 ∪ h( T 2 ), it must cover at most mini ∈ T 1 {u i } remaining clients in C 1 with radius 3R , otherwise it would have been
(1)
removed during the construction of P . Of course, this number is also smaller than the average of all u i ’s, which is exactly
(1 )
ui = |C 1 \ C 1 | − u (01) ,
γ i∈T γ
1
where we recall C 1 \ C 1 contains all clients covered by T 1 ∪ h( T 2 ) and | T 1 | = γ / . The case with F 2 is identical.
Lemma 16. There exists a guess ( T 1 , T 2 , g , h) such that the original instance reveals a solution ( A 1 ⊆ F 1 , A 2 ⊆ F 2 ) with objective at
most R , satisfying T 1 ∪ h( T 2 ) ⊆ A 1 , T 2 ∪ g ( T 1 ) ⊆ A 2 in terms of multi-sets, facility
at i ∈ T 1 moves to g (i ) ∈ g ( T 1 ) and facility at
i ∈ T 2 comes from h(i ) ∈ h( T 2 ). Moreover, there are at most O (| F 1 | · | F 2 |)2γ / such different guesses.
Proof. The number of possible guesses is easy to see from the definition. Consider the optimal solution (U 1 ⊆ F 1 , U 2 ⊆ F 2 )
to the original problem. Choose γ / facilities U 1 ⊆ U 1 , such that the total number of clients they cover in radius 3R is
maximized. Define U 2 ⊆ U 2 similarly. Assume that we have made the correct guesses T 1 = U 1 , T 2 = U 2 (see that their
cardinalities are the same), and our guesses g and h correctly depict their movement in the optimal solution.
Now we only need to prove that U 1 ⊆ F 1 and U 2 ⊆ F 2 . Since U 1 = T 1 , we can sort U 1 in the same way as T 1 , and it is
(1)
easy to see that any i ∈ U 1 \ U 1 can only cover at most mini ∈ T 1 {u i } new clients in C 1 \ i ∈ T 1 B 1 (i , 3R ), otherwise
we may
always replace the last one in U 1 with i and cover more clients, a contradiction. Since C 1 is a subset of C 1 \ i ∈ T 1 B 1 (i , 3R ),
each i ∈ U 1 \ U 1 is not removed from F 1 when F 1 is created, thus one has U 1 \ U 1 ⊆ F 1 and U 1 ⊆ F 1 as a result. The same
argument also shows U 2 \ U 2 ⊆ F 2 and U 2 ⊆ F 2 . We recover the lemma by setting A 1 = U 1 and A 2 = U 2 .
Remark 4. We remark that, while we allow g ( T 1 ) and h( T 2 ) to induce multi-sets, we only consider the cases where T 1 ⊆ F 1
and T 2 ⊆ F 2 are subsets. We briefly explain the reason here. From the proof of the above Lemma 16, one can easily see that
if U 1 in the optimum (as a multi-set itself) has less than γ / distinct members of F 1 , we can simply guess these members,
delete every other facility location that is not co-located with them, and still have an optimal solution with objective R
to the remaining instance. This process takes at most | F 1 | + | F 1 |2 + · · · + | F 1 |γ / = O (| F 1 |1+γ / ) guesses in total, and the
remaining instance is much easier to solve, since the guessed distinct facilities in U 1 should already cover at least l1 clients
with radius R . Therefore, we omit these scenarios here for simplicity.
From now on, further assume that we have made the correct guess ( T 1 , T 2 , g , h) as shown in Lemma 16 and reached
the reduced instance P . We define the following natural LP relaxation. By adding a superscript to every variable to indicate
(t ) (t )
the time step, denote xi j ∈ [0, 1] the partial assignment of client j to facility i and y i ≥ 0 the extent of opening facility
location i at time step t. Moreover, denote zii the extent of movement from facility i to facility i , between neighboring
time steps t = 1 and t = 2.
(t )
xi j ≥ lt ∀t = 1, 2 (LP(P ))
j ∈C t i ∈ F t ,d(i , j )≤ R
(t )
xi j ≤ 1 ∀ j ∈ C t , t = 1, 2
i ∈ F t ,d(i , j )≤ R
(t )
yi = k ∀t = 1, 2
i∈ F t
(t ) (t )
0 ≤ xi j ≤ y i ∀ i ∈ F t , j ∈ C t , t = 1, 2
63
Algorithm 3: GREEDYFILTER.
Input : P , R , (x, y , z)
Output: two subsets of filtered clients for t = 1, 2, with each client having a certain profit value
1 for t = 1, 2 do
2 Ct ← ∅
(t )
3 for unmarked j ∈ C t in non-increasing order of s j do
4 Ct ← Ct ∪ { j}
5 set each unmarked j ∈ C t s.t. d( j , j ) ≤ 2R as marked
(t )
6 let c j be the number of clients marked in this iteration // the “profit” of j ∈ C t

(t )
7 c (t ) ← c j : j ∈ C t
8 return (C t , c (t ) ), t = 1, 2
(t )
xi j = 0 ∀ i ∈ F t , j ∈ C t , d ( i , j ) > R , t = 1, 2
(2 )
zii = y i ∀i ∈ F 2
i ∈ F 1 ,d(i ,i )≤ B
(1 )
zii = y i ∀i ∈ F 1
i ∈ F 2 ,d(i ,i )≤ B
zii ≥ 0 ∀i ∈ F 1 , i ∈ F 2
zii = 0 ∀i ∈ F 1 , i ∈ F 2 , d (i , i ) > B
zig (i ) ≥ 1 ∀i ∈ T 1
z h ( i ) i ≥ 1. ∀i ∈ T 2
Lemma 17. (LP(P )) is feasible.
(t ) (t )
Proof. Consider the optimal solution U 1 ⊆ F 1 , U 2 ⊆ F 2 , and define the variables xi j , y i and zii which are restricted to
the reduced instance P accordingly. Using Lemma 16, if F 1 , F 2 is computed according to T 1 and T 2 , we have U 1 ⊆ F 1 , U 2 ⊆
F 2 , and it is easy to check that all but the first constraint of (LP(P )) are satisfied by (x , y , z ).
Now consider the first constraint. In the optimal solution, U 1 covers at least l1 clients in C 1 with radius R and T 1 ⊆
(1) (1)
U 1 , h( T 2 ) ⊆ U 1 . The facilities in T 1 ∪ h( T 2 ) cover exactly u 0 + i ∈ T 1 u i clients with radius 3R > R , which are all
included in C 1 \ C 1 . Evidently, for those in C 1 that are not covered by the previous larger balls, the remaining facilities in
(1) (1)
U 1 \ ( T 1 ∪ h( T 2 )) have to cover at least max{l1 − u 0 − i ∈ T 1 u i , 0} = l1 of them with radius R . The same holds for C 2
and l2 , thus the first constraint is also satisfied by (x , y , z ), and (LP(P )) is feasible.
3.4.2. Matching-based LP rounding

(t ) (t )
Given a fractional solution (x, y , z) to (LP(P )), define s j = i ∈ F ,d(i , j )≤ R xi j the extent of connection of client j ∈ C t .
t
We conduct the standard filtering algorithm (see, e.g., [21,27]) to filter the clients, defined as in Algorithm 3. According to
the definition of the algorithm, for any j ∈ C t , there exists some j ∈ C t such that d( j , j ) ≤ 2R . Therefore, if we can cover
all the clients in C t with a radius of at most κ R , using triangle inequality, every client in C t can be covered within radius
(κ + 2) R .
The following two lemmas illustrate the relative sparsity of C 1 around (almost) each client j ∈ C 1 and provide an inequal-
ity which will be useful for achieving our final approximate solution. The lemmas for t = 2 are the same, hence omitted
here.

(1)
Lemma 18. For any j ∈ C 1 with s j > 0, c (j1) ≤ γ |C 1 \ C 1 | − u (01) .
(1) (1)
Proof. Evidently, c j is at most the number of clients in C 1 that are ≤ 2R away from j. Since s j > 0, there exists i ∈ F 1
(1)
with d(i , j ) ≤ R , and these c j clients are at most 3R away from i
using triangle inequality. One also has that i is
not co-located with T 1 ∪ h( T 2 ), since otherwise j would have been removed from C 1 during pre-processing, which is a
(1)
contradiction. Using Lemma 15, c j is at most the number of clients in C 1 covered by i with radius 3R , hence at most

(1)
|C 1 \ C 1 | − u 0 /γ .
64
Algorithm 4: SPLIT&MERGE.
Input : P , R , (x, y , z), C 1 , C 2 , c (1) , c (2) , T 1 , T 2 , g , h
Output: a modified instance of 2-DkSupOut and a feasible fractional solution
1 F 1 ← F 1 , F 2 ← F 2 , ŷ ← y , ẑ ← z
2 for each distinct pair (i , i ) ∈ {(i , g (i )) : i ∈ T 1 } ∪ {(h(i ), i ) : i ∈ T 2 } do
3 create co-located copies i 1 , i 1 of i , i , respectively, F 1 ← F 1 ∪ {i 1 }, F 2 ← F 2 ∪ {i 1 }
(1) (2) (1) (1) (2) (2)
4 ŷ i
1
← 1, ŷ i ← 1, ẑi 1 i 1 ← 1, ŷ i ← ŷ i − 1, ŷ i ← ŷ i − 1, ẑii ← ẑii − 1
1
5 for t = 1, 2 and j ∈ C t do
(t )
6 E j ← {i ∈ F t : d (i , j ) ≤ R }
(t ) (t ) (t )
7 while there exists i ∈ E j s.t. ŷ i > xi j do
(t ) (t ) (t ) (t ) (t )
8 split i into co-located copies i 1 , i 2 , ŷ i ← xi j , ŷ i 2 ← ŷ i − ŷ i 1 , split corresponding movement variables ẑ accordingly
1
(t ) (t )
9 F t ← F t \ {i } ∪ {i 1 , i 2 } , E j ← E j \ {i } ∪ {i 1 }
(t ) (t ) (t ) (t )
10 merge all locations in E j into a single one f j , ŷ (t ) ← s j , merge corresponding ẑ values accordingly
fj

(t ) (t )
11 Ft ← Ft \ E j ∪{fj }
(t )
12 while there exists i ∈ F t s.t. ŷ i > 1 do
(t ) (t ) (t ) (t )
13 split i into ŷ i co-located copies, such that the first ŷ i of them have ŷ (t ) values 1, and the last one (if any) has ŷ (t ) value ŷ i − ŷ i ,
split corresponding movement variables ẑ accordingly
14 return ( ŷ , ẑ), F 1 , F 2
(1)
We note that if s j = 0, the lemma above does not hold. Nevertheless, our algorithm handles such clients by simply
ignoring them in the rounding process.
(1) (1)
Lemma 19. j ∈C 1 cj sj ≥ l1 .
(1) (1)
Proof. In Algorithm 3, the subsets are processed in non-increasing order of s j , and c j denotes the number of clients that
are marked in each step, hence one has
(1 ) (1 )
(1 )
cj sj ≥ sj ≥ l1 ,
j ∈C 1 j ∈C 1
(1) (1) (1)

where we recall that s j = i ∈ F 1 ,d(i , j )≤ R xi j , and the last inequality is due to the feasibility of xi j for (LP(P )).
We first use the following Algorithm 4 to modify the LP solution. The main goal of this algorithm is to duplicate, and in
some cases, merge facilities, such that the modified solution ( ŷ , ẑ) is easier to round on the new instance.
In Algorithm 4, the first loop is used to reserve the variables for our pre-selected facilities T 1 , T 2 and g ( T 1 ), h( T 2 ). We
note that this is always possible, since in (LP(P )), we explicitly set zig (i ) ≥ 1 and zh(i )i ≥ 1 for i ∈ T 1 and i ∈ T 2 , which
(t )
also implicitly sets up lower bounds for relevant y variables. Secondly, we define a subset E j of facilities that are close to
(t ) (t ) (t )
j ∈ C t , t = 1, 2. Using the standard facility duplication technique, we make sure that ∀i ∈ E j , one has ŷ i = xi j , which in
(t ) (t ) (t )
turn implies that ŷ (t ) ( E j ) = s j . This even holds for j ∈ C t with s j = 0. Using this property, we simply merge the facilities
(t ) (t ) (t ) (t )
in E j into a single one, called f j , and assign it ŷ (t ) ← s j (recall that F 1 and F 2 live in different metric spaces now, so
fj
this merging process has no effect on the other time step). We also clarify that, whenever a facility location i is split into
multiple copies, we always need to respect the fractional movement variables ẑ and make sure that they are also properly
split in order to satisfy the “conservation” of fractional facilities. Here, we do not need any particular splitting criteria, and
simply use an arbitrary one.
From now on, we view the facilities in F 1 and F 2 as vertices in a bipartite graph, and the variables in ẑ as a fractional
matching on it. Roughly speaking, if some ẑii is rounded to one, it means that a facility is moved from i ∈ F 1 to i ∈ F 2 .
(1) (1)
Further, if i is equal to some merged facility f j with j ∈ C 1 , it means that some facility in E j (before merging) can be
(1)
opened, and we can use this facility to cover at least c j distinct clients in C 1 with radius 3R , using triangle inequality
and the fact that any client j marked by j satisfies d( j , j ) ≤ 2R . This observation motivates us to assign a “profit” value
to each merged facility, and define a multi-objective optimization problem in the following lemma.
Lemma 20. For the modified LP solution ( ŷ , ẑ) by Algorithm 4, create a bipartite graph G = ( V 1 ∪ V 2 , E ), where V 1 , V 2 represent
the entries in F 1 , F 2 , respectively, and edge (i , i ) ∈ E is defined for every non-zero ẑii > 0. ẑ is a fractional k-cardinality bipartite
matching over G.
65
(1) (1) (2)

Moreover, assign profit p 1 (e ) = c j to every edge e incident on merged facility f j , j ∈ C 1 , profit p 2 (e ) = c j to every edge e on
(2)
merged facility f j , j ∈ C 2 , and set all other profit values to zero. ẑ satisfies the following “profitability constraints”,

P 1 (ẑ) = ẑe p 1 (e ) ≥ l1 , P 2 (ẑ) = ẑe p 2 (e ) ≥ l2 .
e∈ E e∈ E
(t ) (t ) (t )
Proof. In Algorithm 4, for every merged facility f j , it satisfies ŷ (t ) = s j ≤ 1. For other facilities, we split those with
fj
ŷ i > 1 into co-located copies, and every copy is matched up to an extent of 1. The total extent of matched edges is directly
(1)
from the constraint i ∈ F y i = k in (LP(P )). Hence ẑ is a fractional k-cardinality matching on G.
1
To see that the profitability constraints are satisfied, we focus on t = 1. Since each edge that is not incident on any
(1)
f j , j ∈ C 1 has p 1 (e ) = 0, we may rewrite the total profit as,
⎛ ⎞
⎜ ⎟
cj ⎜ ẑe ⎟
(1 ) (1 ) (1 ) (1 ) (1 )
P 1 (ẑ) = ẑe p 1 (e ) = ⎝ ⎠= c j ŷ (1) = cj sj ≥ l1 ,
fj
j ∈ C 1 e = f ,i(1) j ∈C 1 (1)
e = f j ,i j ∈C 1 j ∈C 1
j
where the last inequality is due to Lemma 19.
(1) (1)
We also notice that if s j = 0, j ∈ C 1 , the merged facility satisfies ŷ (1 ) = 0 and every associated movement variable also
fj
(1)
satisfies ẑ f (1) i = 0. This means that no edges are incident on f j , and we can simply remove such vertices from the graph.
j
With the feasibility lemma above, we present the following main theorem, in which we round the aforementioned fractional
bipartite matching on G to an integral matching, and directly obtain a solution to the original 2-DkSupOut instance.
Theorem 21. For any constant > 0, there exists a multi-criteria (3, 1 − , 1 − )-approximation for 2-DkSupOut.
Proof. We consider the fractional matching induced by ẑ in Lemma 20 and adapt the rounding procedures in Section 4
of [20], with some necessary modifications. First, since Algorithm 4 explicitly sets ẑig (i ) = 1 and ẑh(i )i = 1 for each i ∈
T 1 , i ∈ T 2 , these edges are always matched, and we can remove them in advance. Suppose there are κ ≥ 0 such edges, and
we remove these edges and their endpoints.
Notice that p 1 (e ) = p 2 (e ) = 0 if e is removed in this process. To see this, both endpoints of e are in T 1 ∪ h( T 2 ) and
T 2 ∪ g ( T 1 ), and every client within 3R from them is already removed from C 1 and C 2 . W.l.o.g., assume p 1 (e ) > 0 for some
(1)
e removed, which means that there exists j ∈ C 1 such that E j contained some i ∈ T 1 ∪ h( T 2 ) before being merged. This
puts d(i , j ) ≤ R , and we would have removed j from C 1 during pre-processing, a contradiction.
Let P M be the matching polytope of the remaining graph, and consider the following LP,
"

min 1 z : z ∈ P M , P 1 ( z) = ze p 1 (e ) ≥ l1 , P 2 ( z) = ze p 2 (e ) ≥ l2 . (19)
e∈ E e∈ E
According to Lemma 20, the solution ẑ given by Algorithm 4 is feasible to the above LP with objective k − κ (after the
removal of those must-have edges, with each of them having zero profits), hence the optimum is at most k − κ . Let z0 be
such an optimal basic solution with objective 1 z0 ≤ k − κ .
Now that z0 is a basic solution, it lies on a face of P M of dimension at most 2, thus using Carathéodory theorem,
it can be written as the convex combination of 3 (integral) basic solutions of P M , say z0 = α1 z1 + α2 z2 + α3 z3 , where
αi ∈ [0, 1], α1 + α2 + α3 = 1, and z1 , z2 , z3 are three basic solutions to P M . W.l.o.g., we assume that α1 , α2 , α3 are all
positive numbers.
Construction of an intermediate matching ẑ2 . We create an almost-matching that fractionally combines z1 , z2 . To be more
precise, let z1,2 = α α+1α z1 + α α+2α z2 be the convex combination of z1 and z2 thus a fractional matching, we want
1 2 1 2
to find another almost-integral matching z2 ∈ [0, 1] E s.t.
P 1 ( z2 ) = P 1 ( z1,2 ), P 2 ( z2 ) = P 2 ( z1,2 ), 1 z2 = 1 z1,2 , (20)

and it is possible to set at most 4 variables in z2 to 0 and obtain an integral matching. Using Corollary 4.10 in [20],
such an almost-matching z2 exists and can be efficiently computed. We then set at most 4 variables in z2 to 0 and
obtain a matching ẑ2 . It is obvious that
P 1 (ẑ2 ) ≥ P 1 ( z2 ) − 4 max p 1 (e ), P 2 (ẑ2 ) ≥ P 2 ( z2 ) − 4 max p 2 (e ), 1 ẑ2 ≤ 1 z2 . (21)

e∈ E e∈ E
66
Construction of the final matching ẑ3 . Using Corollary 4.10 in [20] again, let z2,3 = (α1 + α2 )ẑ2 + α3 z3 , we can efficiently
find z3 ∈ [0, 1] E s.t.
P 1 ( z3 ) = P 1 ( z2,3 ), P 2 ( z3 ) = P 2 ( z2,3 ), 1 z3 = 1 z2,3 , (22)
and it is possible to set at most 4 variables in z3 to 0 and obtain a matching ẑ3 . It is obvious that
P 1 (ẑ3 ) ≥ P 1 ( z3 ) − 4 max p 1 (e ), P 2 (ẑ3 ) ≥ P 2 ( z3 ) − 4 max p 2 (e ), 1 ẑ3 ≤ 1 z3 . (23)

e∈ E e∈ E
Combining (20)(21)(22)(23), one has
P 1 (ẑ3 ) ≥ P 1 ( z3 ) − 4 max p 1 (e ) = P 1 ( z2,3 ) − 4 max p 1 (e )

e∈ E e∈ E
= (α1 + α2 ) P 1 (ẑ2 ) + α3 P 1 ( z3 ) − 4 max p 1 (e )

e∈ E

≥ (α1 + α2 ) P 1 ( z2 ) − 4 max p 1 (e ) + α3 P 1 ( z3 ) − 4 max p 1 (e )
e∈ E e∈ E
≥ (α1 + α2 ) P 1 ( z1,2 ) + α3 P 1 ( z3 ) − 8 max p 1 (e )

e∈ E
8 (1 )

= P 1 (α1 z1 + α2 z2 + α3 z3 ) − 8 max p 1 (e ) ≥ l1 − |C 1 \ C 1 | − u 0 ,
e∈ E γ
where the last inequality is due to Lemma 18 and the fact that any non-zero profit
must be defined
for an edge
(1) (1) (1) (2)
incident on f j , j ∈ C 1 with 0 < ŷ (1 ) = s j . Similarly, one has P 2 (ẑ3 ) ≥ l2 − 8γ |C 2 \ C 2 | − u 0 and easily sees
fj
1 ẑ
3 ≤ 1 z
3 = 1 z
2,3 = (α1 + α2 2 + α3 )1 ẑ
3 ≤ (α1 + α2 1 z
1, 2 + α 3 3 )1 z
0 ≤ k − κ . Hence the cardinality 1 z = 1 z
of ẑ3 is at most k − κ .
Construction of the output solution. Let γ = 8, M be the set of edges matched in ẑ3 plus the κ edges removed in the
beginning, and A 1 , A 2 be two multi-sets that are initially empty. We have | M | ≤ k. For each e ∈ M, consider the
following cases.
e is among the κ edges removed in the beginning. Add the two endpoints to A 1 and A 2 , respectively.
e is not incident on any merged facility. Add the two endpoints to A 1 and A 2 , respectively.
(1)
e is incident on exactly one merged facility. W.l.o.g., let this endpoint be f j , j ∈ C 1 and the other be i ∈ F 2 . No-
tice that according to Lemma 20 and the merging process in Algorithm 4, e ∈ E indicates that ẑ f (1) i > 0
j
(1)
when we construct G in the first place, thus there must exist i ∈ E j s.t. ẑii > 0 before merging and hence
d(i , i ) ≤ B. Add i to A 1 and i to A 2 .
(1) (2)
e is incident on two merged facilities. Let them be f j , j ∈ C 1 and f j , j ∈ C 2 . Again, e ∈ E implies that ẑ f (1) f (2) >
j j
(1) (2)
0, thus there must exist i ∈ E j and i ∈ E j s.t. ẑii > 0 and hence d(i , i ) ≤ B. Add i to A 1 and i to A 2 .
The construction of ( A 1 , A 2 ) naturally induces a feasible perfect matching between them, such that any matched pair is
no more than B away from each other. We already know ( T 1 ∪ h( T 2 )) ⊆ A 1 , and it covers all clients in C 1 \ C 1 with radius

(1)
3R . From the value of P 1 (ẑ3 ), at least l1 − |C 1 \ C 1 | − u 0 distinct clients in C 1 are covered additionally with radius
3R , so the total number of clients covered by A 1 (using radius 3R ) is at least

(1 )
P 1 (ẑ3 ) + |C 1 \ C 1 | ≥ l1 − |C 1 \ C 1 | − u 0 + u (01) + (1 )
ui
i∈T 1
⎧ ⎫
⎨ ⎬
≥ max l1 − u (01) − (1 )
ui , 0 + u (01) + ui
(1 )
− |C 1 \ C 1 | − u (01)
⎩ ⎭
i∈T 1 i∈T 1
≥ l1 − l1 = (1 − )l1 ,
(1)
where we assume that |C 1 \ C 1 | − u 0 < l1 , otherwise we would have already covered ≥ l1 clients using only T 1 . The proof
is the same for A 2 and the second time step. Arbitrarily add more facilities to A 1 and A 2 until | A 1 | = | A 2 | = k, and this
gives our final multi-criteria (3, 1 − , 1 − )-approximate solution.
67
4. Future work
It would be very interesting to remove the dependency of γ (the coefficient of movement cost) and (the lower bound
of the weight) from the approximation factor for DOkMed in Theorem 8, or to show that such dependency is inevitable. We
leave it as an important open problem. We note that a constant approximation factor for DOkMed without depending on γ
would imply a constant approximation for stochastic k-server, for which only a logarithmic-factor approximation algorithm
is previously known [18].
Our approximation algorithm for DOkMed is based on the technique developed by Aouad and Segev [10] and Byrka et al.
[11]. The original ordered k-median problem has subsequently seen improved approximation results in [12,13]. We did not
try hard to optimize the constant factors. Nevertheless, it is an interesting future direction to improve the constant factors
by leveraging new techniques and ideas. On the other hand, it would also be very interesting to obtain better lower bounds
than the trivial lower bounds of k-supplier [19] and k-median [3].
It is possible to formulate other problems that naturally fit into the dynamic clustering theme and are well motivated
by realistic applications. Besides the obvious variants of DOkMed and DkSup via changing the clustering and movement
objectives, we list some examples that are of particular interest.
• In each time step, the constraint of opening at most k facilities can be replaced by a matroid or knapsack constraint;
that is, we require the open facilities to form an independent set of a given matroid, or to have a total weight no
more than a given threshold. These formulations generalize the k-clustering setting we adopt in this paper. Many such
clustering problems are studied in the literature, e.g., matroid center [26], matroid median [34,35,49], knapsack center
[19] and knapsack median [34,35,50,51], and we can study the dynamic versions of these problems.
• We can consider fair dynamic clustering problems. For example, at each time step t, the clients have different colors
representing the demographics they belong to, and the goal is to cluster the clients such that for each t, the proportions
of different colors in each cluster are roughly the same as their global proportions at time t (see, e.g., [52–54]).
• One may also study other clustering criteria with additional constraints under the dynamic setting, e.g., the fault-
tolerant versions [44,55,56] and the capacitated versions [57–59].
CRediT authorship contribution statement
Shichuan Deng: Methodology, Formal analysis, Writing. Jian Li: Conceptualization, Methodology, Writing. Yuval Rabani:
Conceptualization, Methodology, Writing.
Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have
appeared to influence the work reported in this paper.
Acknowledgments
Shichuan Deng and Jian Li were supported by the National Natural Science Foundation of China Grant 61822203,
61772297, 61632016, 61761146003, the Zhongguancun Haihua Institute for Frontier Information Technology, Turing AI In-
stitute of Nanjing, and Xi’an Institute for Interdisciplinary Information Core Technology. Yuval Rabani was supported by ISF
grant number 2553-17.
We thank Chaitanya Swamy for kindly pointing out studies relevant to our results. We also thank the anonymous referees
for their insightful and constructive comments.
References
[1] S. Deng, J. Li, Y. Rabani, Approximation algorithms for clustering with dynamic points, in: 28th Annual European Symposium on Algorithms, in: LIPIcs,
vol. 173, 2020, 37.
[2] W. Hsu, G.L. Nemhauser, Easy and hard bottleneck location problems, Discrete Appl. Math. 1 (3) (1979) 209–215, https://doi.org/10.1016/0166-218X(79)
90044-1.
[3] K. Jain, M. Mahdian, A. Saberi, A new greedy approach for facility location problems, in: Proceedings on 34th Annual ACM Symposium on Theory of
Computing, 2002, pp. 731–740.
[4] T.F. Gonzalez, Clustering to minimize the maximum intercluster distance, Theor. Comput. Sci. 38 (1985) 293–306, https://doi.org/10.1016/0304-3975(85)
90224-5.
[5] D.S. Hochbaum, D.B. Shmoys, A best possible heuristic for the k-center problem, Math. Oper. Res. 10 (2) (1985) 180–184, https://doi.org/10.1287/moor.
10.2.180.
[6] V. Arya, N. Garg, R. Khandekar, A. Meyerson, K. Munagala, V. Pandit, Local search heuristics for k-median and facility location problems, SIAM J.
Comput. 33 (3) (2004) 544–562, https://doi.org/10.1137/S0097539702416402.
[7] M. Charikar, S. Li, A dependent LP-rounding approach for the k-median problem, in: Automata, Languages, and Programming - 39th International
Colloquium, vol. 7391, 2012, pp. 194–205.
[8] S. Li, O. Svensson, Approximating k-median via pseudo-approximation, SIAM J. Comput. 45 (2) (2016) 530–547, https://doi.org/10.1137/130938645.
68
[9] J. Byrka, T.W. Pensyl, B. Rybicki, A. Srinivasan, K. Trinh, An improved approximation for k-median and positive correlation in budgeted optimization,
ACM Trans. Algorithms 13 (2) (2017) 23, https://doi.org/10.1145/2981561.
[10] A. Aouad, D. Segev, The ordered k-median problem: surrogate models and approximation algorithms, Math. Program. 177 (1–2) (2019) 55–83, https://
doi.org/10.1007/s10107-018-1259-3.
[11] J. Byrka, K. Sornat, J. Spoerhase, Constant-factor approximation for ordered k-median, in: Proceedings of the 50th Annual ACM SIGACT Symposium on
Theory of Computing, 2018, pp. 620–631.
[12] D. Chakrabarty, C. Swamy, Interpolating between k-median and k-center: approximation algorithms for ordered k-median, in: 45th International Col-
loquium on Automata, Languages, and Programming, in: LIPIcs, vol. 107, 2018, 29.
[13] D. Chakrabarty, C. Swamy, Approximation algorithms for minimum norm and ordered optimization problems, in: Proceedings of the 51st Annual ACM
SIGACT Symposium on Theory of Computing, 2019, pp. 126–137.
[14] M.S. Manasse, L.A. McGeoch, D.D. Sleator, Competitive algorithms for server problems, J. Algorithms 11 (2) (1990) 208–230, https://doi.org/10.1016/
0196-6774(90)90003-W.
[15] E. Koutsoupias, C.H. Papadimitriou, On the k-server conjecture, J. ACM 42 (5) (1995) 971–983, https://doi.org/10.1145/210118.210128.
[16] C. Coester, E. Koutsoupias, The online k-taxi problem, in: Proceedings of the 51st Annual ACM SIGACT Symposium on Theory of Computing, 2019,
pp. 1136–1147.
[17] N. Buchbinder, C. Coester, J.S. Naor, Online k-taxi via double coverage and time-reverse primal-dual, in: Integer Programming and Combinatorial
Optimization - 22nd International Conference, in: Lecture Notes in Computer Science, vol. 12707, 2021, pp. 15–29.
[18] S. Dehghani, S. Ehsani, M. Hajiaghayi, V. Liaghat, S. Seddighin, Stochastic k-server: how should Uber work?, in: 44th International Colloquium on
Automata, Languages, and Programming, in: LIPIcs, vol. 80, 2017, 126.
[19] D.S. Hochbaum, D.B. Shmoys, A unified approach to approximation algorithms for bottleneck problems, J. ACM 33 (3) (1986) 533–550, https://doi.org/
10.1145/5925.5933.
[20] F. Grandoni, R. Ravi, M. Singh, R. Zenklusen, New approaches to multi-objective optimization, Math. Program. 146 (1–2) (2014) 525–554, https://
doi.org/10.1007/s10107-013-0703-7.
[21] D.G. Harris, T.W. Pensyl, A. Srinivasan, K. Trinh, A lottery model for center-type problems with outliers, ACM Trans. Algorithms 15 (3) (2019) 36,
https://doi.org/10.1145/3311953.
[22] K. Jain, V.V. Vazirani, Approximation algorithms for metric facility location and k-median problems using the primal-dual schema and Lagrangian
relaxation, J. ACM 48 (2) (2001) 274–296, https://doi.org/10.1145/375827.375845.
[23] V. Cohen-Addad, A. Gupta, L. Hu, H. Oh, D. Saulpic, An improved local search algorithm for k-median, in: Proceedings of the 2022 Annual ACM-SIAM
Symposium on Discrete Algorithms, 2022, pp. 1556–1612.
[24] M. Charikar, S. Khuller, D.M. Mount, G. Narasimhan, Algorithms for facility location problems with outliers, in: Proceedings of the Twelfth Annual
Symposium on Discrete Algorithms, 2001, pp. 642–651, http://dl.acm.org/citation.cfm?id=365411.365555.
[25] D. Chakrabarty, P. Goyal, R. Krishnaswamy, The non-uniform k-center problem, in: 43rd International Colloquium on Automata, Languages, and Pro-
gramming, in: LIPIcs, vol. 55, 2016, 67.
[26] D.Z. Chen, J. Li, H. Liang, H. Wang, Matroid and knapsack center problems, Algorithmica 75 (1) (2016) 27–52, https://doi.org/10.1007/s00453-015-0010-
1.
[27] D. Chakrabarty, M. Negahbani, Generalized center problems with outliers, ACM Trans. Algorithms 15 (3) (2019) 41, https://doi.org/10.1145/3338513.
[28] S. Bandyapadhyay, T. Inamdar, S. Pai, K.R. Varadarajan, A constant approximation for colorful k-center, in: 27th Annual European Symposium on
Algorithms, in: LIPIcs, vol. 144, 2019, 12.
[29] G. Anegg, H. Angelidakis, A. Kurpisz, R. Zenklusen, A technique for obtaining true approximations for k-center with covering constraints, Math. Program.
192 (1) (2022) 3–27, https://doi.org/10.1007/s10107-021-01645-y.
[30] X. Jia, K. Sheth, O. Svensson, Fair colorful k-center clustering, in: Integer Programming and Combinatorial Optimization - 21st International Conference,
in: Lecture Notes in Computer Science, vol. 12125, 2020, pp. 209–222.
[31] E.D. Demaine, M.T. Hajiaghayi, H. Mahini, A.S. Sayedi-Roshkhar, S.O. Gharan, M. Zadimoghaddam, Minimizing movement, ACM Trans. Algorithms 5 (3)
(2009) 30, https://doi.org/10.1145/1541885.1541891.
[32] Z. Friggstad, M.R. Salavatipour, Minimizing movement in mobile facility location problems, ACM Trans. Algorithms 7 (3) (2011) 28, https://doi.org/10.
1145/1978782.1978783.
[33] S. Ahmadian, Z. Friggstad, C. Swamy, Local-search based approximation algorithms for mobile facility location problems, in: Proceedings of the Twenty-
Fourth Annual ACM-SIAM Symposium on Discrete Algorithms, 2013, pp. 1607–1621.
[34] C. Swamy, Improved approximation algorithms for matroid and knapsack median problems and applications, ACM Trans. Algorithms 12 (4) (2016) 49,
https://doi.org/10.1145/2963170.
[35] R. Krishnaswamy, S. Li, S. Sandeep, Constant approximation for k-median and k-means with outliers via iterative rounding, in: Proceedings of the 50th
Annual ACM SIGACT Symposium on Theory of Computing, 2018, pp. 646–659.
[36] D. Eisenstat, C. Mathieu, N. Schabanel, Facility location in evolving metrics, in: Automata, Languages, and Programming - 41st International Colloquium,
vol. 8573, 2014, pp. 459–470.
[37] H. An, A. Norouzi-Fard, O. Svensson, Dynamic facility location via exponential clocks, ACM Trans. Algorithms 13 (2) (2017) 21, https://doi.org/10.1145/
2928272.
[38] D.L. Black, A. Gupta, W. Weber, Competitive management of distributed shared memory, in: Thirty-Fourth IEEE Computer Society International Confer-
ence, 1989, pp. 184–190.
[39] J.R. Westbrook, Randomized algorithms for multiprocessor page migration, SIAM J. Comput. 23 (5) (1994) 951–965, https://doi.org/10.1137/
S0097539791199796.
[40] T.H. Chan, A. Guerqin, M. Sozio, Fully dynamic k-center clustering, in: Proceedings of the 2018 World Wide Web Conference on World Wide Web,
2018, pp. 579–587.
[41] V. Cohen-Addad, N. Hjuler, N. Parotsidis, D. Saulpic, C. Schwiegelshohn, Fully dynamic consistent facility location, in: Advances in Neural Information
Processing Systems 32: Annual Conference on Neural Information Processing Systems 2019, 2019, pp. 3250–3260, https://proceedings.neurips.cc/paper/
2019/hash/fface8385abbf94b4593a0ed53a0c70f-Abstract.html.
[42] A. Schrijver, Combinatorial Optimization: Polyhedra and Efficiency, vol. 24, Springer, 2003.
[43] V. Kumar, M.V. Marathe, S. Parthasarathy, A. Srinivasan, A unified approach to scheduling on unrelated parallel machines, J. ACM 56 (5) (2009) 28,
https://doi.org/10.1145/1552285.1552289.
[44] M.T. Hajiaghayi, W. Hu, J. Li, S. Li, B. Saha, A constant factor approximation algorithm for fault-tolerant k-median, ACM Trans. Algorithms 12 (3) (2016)
36, https://doi.org/10.1145/2854153.
[45] R.M. Karp, Reducibility among combinatorial problems, in: 50 Years of Integer Programming 1958-2008 - from the Early Years to the State-of-the-Art,
Springer, 2010, pp. 219–241.
[46] M. Chlebík, J. Chlebíková, Complexity of approximating bounded variants of optimization problems, Theor. Comput. Sci. 354 (3) (2006) 320–338,
https://doi.org/10.1016/j.tcs.2005.11.029.
69
[47] V. Kann, Maximum bounded 3-dimensional matching is MAX SNP-complete, Inf. Process. Lett. 37 (1) (1991) 27–35, https://doi.org/10.1016/0020-
0190(91)90246-E.
[48] S. Arora, C. Lund, R. Motwani, M. Sudan, M. Szegedy, Proof verification and the hardness of approximation problems, J. ACM 45 (3) (1998) 501–555,
https://doi.org/10.1145/278298.278306.
[49] R. Krishnaswamy, A. Kumar, V. Nagarajan, Y. Sabharwal, B. Saha, The matroid median problem, in: Proceedings of the Twenty-Second Annual ACM-SIAM
Symposium on Discrete Algorithms, 2011, pp. 1117–1130.
[50] A. Kumar, Constant factor approximation algorithm for the knapsack median problem, in: Proceedings of the Twenty-Third Annual ACM-SIAM Sympo-
sium on Discrete Algorithms, 2012, pp. 824–832.
[51] A. Gupta, B. Moseley, R. Zhou, Structural iterative rounding for generalized k-median problems, in: 48th International Colloquium on Automata, Lan-
guages, and Programming, vol. 198, 2021, 77.
[52] F. Chierichetti, R. Kumar, S. Lattanzi, S. Vassilvitskii, Fair clustering through fairlets, in: Advances in Neural Information Processing Systems
30: Annual Conference on Neural Information Processing Systems 2017, 2017, pp. 5029–5037, https://proceedings.neurips.cc/paper/2017/hash/
978fce5bcc4eccc88ad48ce3914124a2-Abstract.html.
[53] S.K. Bera, D. Chakrabarty, N. Flores, M. Negahbani, Fair algorithms for clustering, in: Advances in Neural Information Processing Systems
fc192b0c0d270dbf41870a63a8c76c2f-Abstract.html.
[54] L. Huang, S.H. Jiang, N.K. Vishnoi, Coresets for clustering with fairness constraints, in: Advances in Neural Information Processing Systems
810dfbbebb17302018ae903e9cb7a483-Abstract.html.
[55] S. Khuller, R. Pless, Y.J. Sussmann, Fault tolerant k-center problems, Theor. Comput. Sci. 242 (1–2) (2000) 237–245, https://doi.org/10.1016/S0304-
3975(98)00222-9.
[56] C. Swamy, D.B. Shmoys, Fault-tolerant facility location, ACM Trans. Algorithms 4 (4) (2008) 51, https://doi.org/10.1145/1383369.1383382.
[57] S. Khuller, Y.J. Sussmann, The capacitated k-center problem, SIAM J. Discrete Math. 13 (3) (2000) 403–418, https://doi.org/10.1137/S0895480197329776.
[58] J. Chuzhoy, Y. Rabani, Approximating k-median with non-uniform capacities, in: Proceedings of the Sixteenth Annual ACM-SIAM Symposium on Discrete
Algorithms, 2005, pp. 952–958, http://dl.acm.org/citation.cfm?id=1070432.1070569.
[59] S. Li, Approximating capacitated k-median with (1 + )k open facilities, in: Proceedings of the Twenty-Seventh Annual ACM-SIAM Symposium on
Discrete Algorithms, 2016, pp. 786–796.
70

AA_research_paper

Uploaded by

Document Informationclick to expand document information

Copyright:

Available Formats

AA_research_paper

Uploaded by

Document Information

Original Description:

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

AA_research_paper

Uploaded by

Copyright:

Available Formats

Journal of Computer and System Sciences 130 (2022) 43–70

Contents lists available at ScienceDirect

Journal of Computer and System Sciences

Approximation algorithms for clustering with dynamic

1.1. Our results

1.1.1. Dynamic ordered k-median

1.1.2. Dynamic k-supplier

Theorem 2. One has the following hardness results.

1.2. Related work

2. A constant approximation for DOkMed

Lemma 4. (Charikar and Li [7]). Fix t ∈ [ T ]. The following statements hold.

2.2. Flow-based LP rounding

Step 1. For each i ∈ F t , t ∈ [ T ], add nodes L (i ) to L 6t −5 and R (i ) to L 6t .

since y ( p ) ≥ 2 × (1/2) = 1; if p is a singleton pair, the capacity is either [0, 1] or {1}.

If T ≥ 3, the expected movement cost of solution { A t }tT=1 satisﬁes

E[t , p ] ≤ 2dav ( j 1 ) + 2dav ( j 2 ) + 2d( j 1 , j 2 ) · (1 − y (U j 1 )) ≤ 6dav ( j 1 ) + 2dav ( j 2 ).

which yields the lemma when combined with (1).

2.4. Missing proofs

Case 1: A t ∩ B o ( j , β h) = ∅. We have d( j , A t ) ≤ β h and charge β h to D j .

where the last inequality is due to the deﬁnition of d−h,av ( j ).

D j = (4 + 5β)h, E[ X j ] ≤ 4d−h,av ( j ) + 10d−h,far ( j ). (5)

so in this case, E[ X j ] is bounded by

Combining (5)(6)(7), we have

By combining (12), (13), if T = 2, our overall cost is

≤ 41.33OPT + 41.33LP ≤ 82.66OPT.

all connection distances at t fall in, with respect to the optimum.

3.1. The hardness of approximating 3-DkSup and 3-DkSupOut

3.2. A 3-approximation for 2-DkSup

Theorem 12. There exists a 3-approximation for 2-DkSup.

3.3. The hardness of pure approximations for 2-DkSupOut

3.4. A multi-criteria approximation for 2-DkSupOut

3.4.1. Guessing and modifying the instance

Lemma 17. (LP(P )) is feasible.

3.4.2. Matching-based LP rounding

(1) (1) (1)

(1) (1) (2)

where the last inequality is due to Lemma 19.

P 1 ( z2 ) = P 1 ( z1,2 ), P 2 ( z2 ) = P 2 ( z1,2 ), 1 z2 = 1 z1,2 , (20)

P 1 (ẑ2 ) ≥ P 1 ( z2 ) − 4 max p 1 (e ), P 2 (ẑ2 ) ≥ P 2 ( z2 ) − 4 max p 2 (e ), 1 ẑ2 ≤ 1 z2 . (21)

P 1 ( z3 ) = P 1 ( z2,3 ), P 2 ( z3 ) = P 2 ( z2,3 ), 1 z3 = 1 z2,3 , (22)

P 1 (ẑ3 ) ≥ P 1 ( z3 ) − 4 max p 1 (e ), P 2 (ẑ3 ) ≥ P 2 ( z3 ) − 4 max p 2 (e ), 1 ẑ3 ≤ 1 z3 . (23)

Combining (20)(21)(22)(23), one has

P 1 (ẑ3 ) ≥ P 1 ( z3 ) − 4 max p 1 (e ) = P 1 ( z2,3 ) − 4 max p 1 (e )

= (α1 + α2 ) P 1 (ẑ2 ) + α3 P 1 ( z3 ) − 4 max p 1 (e )

≥ (α1 + α2 ) P 1 ( z1,2 ) + α3 P 1 ( z3 ) − 8 max p 1 (e )

CRediT authorship contribution statement

Declaration of competing interest

You might also like

E[t , p ] ≤ 2dav ( j 1 ) + 2dav ( j 2 ) + 2d( j 1 , j 2 ) · (1 − y (U j 1 )) ≤ 6dav ( j 1 ) + 2dav ( j 2 ).

Case 1: A t ∩ B o ( j , β h) = ∅. We have d( j , A t ) ≤ β h and charge β h to D j .

P 1 ( z2 ) = P 1 ( z1,2 ), P 2 ( z2 ) = P 2 ( z1,2 ), 1 z2 = 1 z1,2 , (20)

P 1 (ẑ2 ) ≥ P 1 ( z2 ) − 4 max p 1 (e ), P 2 (ẑ2 ) ≥ P 2 ( z2 ) − 4 max p 2 (e ), 1 ẑ2 ≤ 1 z2 . (21)

P 1 ( z3 ) = P 1 ( z2,3 ), P 2 ( z3 ) = P 2 ( z2,3 ), 1 z3 = 1 z2,3 , (22)

P 1 (ẑ3 ) ≥ P 1 ( z3 ) − 4 max p 1 (e ), P 2 (ẑ3 ) ≥ P 2 ( z3 ) − 4 max p 2 (e ), 1 ẑ3 ≤ 1 z3 . (23)