Spherical coverage verification
We consider the problem of covering hypersphere by a set of spherical hypercaps. This sort of problem has numerous practical applications such as error correcting codes and reverse k-nearest neighbor problem. Using the reduction of non degenerated co…
Authors: Marko D. Petkovic, Dragoljub Pokrajac, Longin Jan Latecki
Spherical co v erage v erification Mark o D. P etk ovi ´ c 1 , Dragoljub P okra jac 2 , Longin Jan Latecki 3 1 F aculty of Scienc e and Mathematics, University of Ni ˇ s, Vi ˇ se gr adska 33, 18000 Ni ˇ s, Serbia E-mail: dexterofnis@gmail.com 2 Computer and Information Scienc es Dep artment, Applie d Mathematics R ese ar ch Center, CREOSA Center, Delawar e State University, Dover, DE, 19901, USA E-mail: dpokrajac@desu.edu 3 Computer and Information Scienc es Dep artment, T emple University, Philadelphia, P A, USA E-mail: latecki@temple.edu Abstract W e consider the problem of co vering h yp ersphere b y a set of spherical h yp ercaps. This sort of prob- lem has numerous practical applications suc h as error correcting codes and reverse k -nearest neigh b or problem. Using the reduction of non degenerated conca v e quadratic programming (QP) problem, w e demonstrate that spherical cov erage verification is NP hard. W e propose a recursiv e algorithm based on reducing the problem to several low er dimension subproblems. W e test the p erformance of the prop osed algorithm on a num b er of generated constellations. W e demonstrate that the prop osed al- gorithm, in spite of its exp onen tial w orst-case complexity , is applicable in practice. In con trast, our results indicate that spherical cov erage verification using QP solvers that utilize heuristics, due to n umerical instabilit y , may pro duce false positives. Keywor ds: geometrical algorithms, quadratic programming, h yp ersphere, c o v erage, h yp ercaps. 1 In tro duction W e consider the problem of determining whether a hypersphere is completely cov ered b y a set of hyper- spherical caps. An equiv alent problem is whether the given set of hyperspherical cones cen tered in the origin, cov ers the whole space R d . A three dimensional v ersion of the problem can be solved using the approac h from [4], but a solution to a generalized problem in arbitrary dimensional space has not, to the b est of our kno wledge, b een prop osed y et. The generalized problem of hyperspherical cov erage by a set of h yp ercaps arises in areas suc h as co ding theory [12] and multidimensional queries [11]. Co vering problems are important in computational geometry and hav e b een extensively studied re- cen tly . Elbassioni and Tiwary [8] considered the following problem: Giv en a set of p olyhedral hypercones C 1 , C 2 , . . . , C k and a conv ex set D ⊂ R d , chec k whether cones cov er the set D or not. They prov ed NP completeness in several cases and connected it to the problem of determining whether a union of the con vex sets is conv ex. Also in [3] the authors considered cov ering a giv en set of p oints with a giv en p olygon, whether in [5] authors considered co vering the set of p oints with t w o disjoin t disks and t wo disjoin t squares. Papers [7] and [18] consider cov ering a sphere (ball) with other spheres (balls). Recen tly , dev elopment of algorithms for incremen tal density-based outlier detection [13, 16] ha ve motiv ated the need to efficiently apply tec hniques for reverse k -nearest neigh b or searc h. The minimal 1 n umber of hypercaps that can completely cov er a hypersphere is related to the theoretical upp er b ound for the n umber of reverse k -nearest neigh b ors of a given p oin t and hence the complexity of incremen tal outlier detection algorithms [15]. Also, the sets of h yp ercaps that can completely co ver the h yp ersphere are basis for practical algorithms for reverse k -nearest neigh b ors. In this pap er, w e consider h yp erspherical co verage v erification: for a giv en set of h yp ercaps on a h yp ersphere w e determine whether the set completely co vers the h yp ersphere. W e are, how ev er, not concerned how to determine the set of h yp ercaps that completely co vers the hypersphere. Using the reduction of conca ve quadratic programming (QP) problem, we demonstrate that spherical co verage verification is NP hard. As a consequence, an algorithm with non-p olynomial worst-case com- plexit y may still b e viable and practical. W e pro vide a generalized recursive algorithm that can p erform co verage verifi cation task for arbitrary dimension d . The prop osed algorithm is based on reducing the problem to sev e ral low er dimension subproblems. In addition, we provide a method that can iden tify a p oin t on a hypersphere not cov ered b y any hypercap, if suc h a p oint exists. W e test the p erformance of the prop osed algorithm on a n umber of generated constellations with differen t dimensionality . W e demonstrate that the prop osed algorithm, in spite of its exponential worst- case complexity , is applicable in practice, with acceptable av erage-case p erformance. In con trast, our results indicate that spherical cov erage v erification using heuristics-based QP solv ers, may pro duce false p ositiv es and suffer from numerical instability . 2 Spherical co v erage v erification In this section, w e formally define a spherical cov erage verification problem and demonstrate that the considered problem SphCovV er can be described as a system of non-linear equations and inequalities. Subsequen tly , we demonstrate that the problem at hand can b e represen ted as quadratic programming problems with linear constraints. 2.1 Problem form ulation Supp ose that w e ha v e n cones C 1 , . . . , C n in d -dimensional space R d , d ≤ 2 centered at p oint O = (0 , . . . , 0) and defined by C i = C ( t i ; θ i ) = x ∈ R d | ( x, t i ) k x kk t i k ≥ θ i . Note that eac h cone C i is defined b y point t i ∈ R d and real n um b er − 1 < θ i < 1. There holds cos ∠ xO t i ≥ θ i for eac h x ∈ C i and x 6 = O . F or an y tw o given points x, y ∈ R d , with ( x, y ) we denote usual scalar pro duct as ( x, y ) = P d i =1 x i y i and with k x k we denote the Euclidian norm k x k = p ( x, x ). Problem 1. ( SphCovV er ) Che ck if c ones C i c over the whole sp ac e R d . Equivalently, che ck if hyp er c aps K i = C i ∩ S d (1) c over an unit hyp erspher e S d (1) = { x ∈ R d | k x k = 1 } . Without los s of generality , w e can assume that all p oints t i b elong to unit h yp ersphere S d (1), i.e., whic h holds k t i k = 1. Let t i = ( t i 1 , . . . , t id ). If x ∈ S d (1) then: x ∈ K i ⇐ ⇒ d X j =1 x i t ij ≥ θ i Observ e that p oint x on the unit h yp ers phere S d − 1 (1) is not cov ered by an y of the cones C 1 , . . . , C n if and only if ( x, t i ) < θ i for all i = 1 , . . . , n . Therefore cones C 1 , . . . , C n co ver the unit hypersphere S d − 1 (1) 2 if and only if the follo wing non-linear system of equations and inequities does not hav e a solution x 2 1 + x 2 2 + . . . + x 2 d = 1 t 11 x 1 + t 12 x 2 + . . . + t 1 d x d < θ 1 . . . t n 1 x 1 + t n 2 x 2 + . . . + t nd x d < θ n . (1) 2.2 Spherical co verage verification as quadratic optimization problem Denote b y S the solution space of the linear system of inequalities obtained b y dropping the first equation in (1). Ob viously S is the conv ex set. Denote by ¯ S the closure of S i.e., the solution space of the inequalities from (1) when eac h < is replaced b y ≤ . The set ¯ S is also con v ex. Let f ( x ) = k x k 2 , m = f ( ˜ x ) = min { f ( x ) | x ∈ ¯ S } and M = f ( x ∗ ) = max { f ( x ) | x ∈ ¯ S } (if ¯ S is unbounded, then M = + ∞ ). In other w ords, ˜ x and x ∗ (i.e., m and M ) are solutions of the following QP problem: (min / max) f ( x ) = k x k 2 = x 2 1 + x 2 2 + . . . + x 2 d s . t . t 11 x 1 + t 12 x 2 + . . . + t 1 d x d ≤ θ 1 . . . t n 1 x 1 + t n 2 x 2 + . . . + t nd x d ≤ θ n . (2) If ¯ S is un b ounded ( M = + ∞ ) then let x ∗ ∈ ¯ S b e any feasible p oint so that k x k > 1. W e say that set ¯ S sp ecified b y constraints from (QP problem (2)) is de gener ate d if it is con tained in some hyperplane H . Note that hyperplane H has to b e of the form ( x, t i ) = θ i , i.e. there ha ve to exist tw o constraints i and j from eq. (2) where t ik = − t j k , k = 1 , . . . , d and θ i = − θ j . In suc h case, the system of equations (1) ob viously has no solutions. The follo wing lemma shows the connection b et ween problem SphCo vV er (i.e., the system (1)) and the QP problem (2): Lemma 1. The system of e quations and ine quities (1) has solutions if and only if M > 1 , m < 1 and ¯ S is non-de gener ate d. Pro of. ( ⇐ :) Let M > 1, m < 1 and ¯ S be non-degenerated. Assume that x ∗ is a b oundary p oint of ¯ S and consider a ball B d ( x ∗ , ρ ) where ρ < k x ∗ k − 1. Since ¯ S is a non-degenerated p olytop e, there exists an in ternal p oint x ∗ 1 ∈ S ∩ B d ( x ∗ , ρ ). If x ∗ is an internal p oint, we just set x ∗ 1 = x ∗ . Note that ρ < k x ∗ k − 1 implies k x ∗ 1 k > 1. The same w ay , w e consider a b oundary p oin t ˜ x of ¯ S and construct a new point ˜ x 1 ∈ S so that‘ k ˜ x 1 k < 1 (here we take ρ < 1 − k ˜ x k ). Since [ ˜ x 1 , x ∗ 1 ] ⊂ S ( S is con v ex), function f ( x ) = k x k 2 is con tinuous on [ ˜ x 1 , x ∗ 1 ] and f ( ˜ x 1 ) < 1 < f ( x ∗ 1 ), there exists p oin t u ∈ ( ˜ x 1 , x ∗ 1 ) ⊂ S so that f ( u ) = 1. In other w ords, u is the solution of system (1). ( ⇒ :) Let u be one solution of system (1). Hence u ∈ S and f ( u ) = 1. Since S is an op en set, there exists a ball B d ( u ; ρ ) ⊂ S . Denote by u 1 and u 2 the in tersection p oin ts of B d ( u ; ρ ) and line O u , so that u ∈ ( O u 1 ). Ob viously holds f ( u 1 ) = (1 + ρ ) 2 and f ( u 2 ) = (1 − ρ ) 2 . No w u 1 , u 2 ∈ B d ( u ; ρ ) ⊂ ¯ S and f ( u 1 ) > 1 > f ( u 2 ) directly implies M ≥ f ( u 1 ) > 1 > f ( u 2 ) ≥ m . Non-degeneracy of ¯ S follows immediately from the fact that (1) has solutions. 3 3 Spherical co v erage v erification is NP-hard In this section, w e pro ve that the spherical co verage v erification problem ( SphCo vV er ) is an NP hard problem. First, w e demonstrate that the concav e non degenerated quadratic programming decision prob- lem ( ConNDQPd ) defined b elow can b e p olynomially reduced to SphCovV er . Then we demonstrate that the problem ConNDQPd is NP complete. Note that the v ariant of the problem ConNDQPd without the non-degeneracy assumption is con- sidered b y F reund and Orlin in [10] (HB problem) where its NP completeness is prov en. The degeneracy assumption mak es the problem ConNDQPd considered here more restrictiv e than the one considered in [10], implying that we need a different pro of of NP completeness. It will b e given in the next subsection. 3.1 P olynomial reduction of ConNDQPd to SphCovV er Consider the following concav e quadratic programming (QP) problem. Problem 2. ( ConNDQPd ) Che ck whether exist x 1 , x 2 , . . . , x d ∈ R so that x 2 1 + x 2 2 + . . . + x 2 d > c a 11 x 1 + a 12 x 2 + . . . + a 1 d x d ≤ b 1 . . . a n 1 x 1 + a n 2 x 2 + . . . + a nd x d ≤ b n . (3) wher e c, a ij , b i ∈ R ( c > 0 ) for i = 1 , . . . , n, j = 1 , . . . , d and the p olytop e sp e cifie d by ≤ c onstr aints fr om (3) is non-de gener ate d. Also consider the following algorithm: Algorithm 1. QP-SphCovV er Input: A n instanc e of the pr oblem ConNDQPd , i.e., matrix A = [ a ij ] ∈ R n × d and ve ctor b = [ b i ] ∈ R n × 1 . 1. Normalize e ach c onstr aint, i.e c ompute t ij = a ij P d j 0 =1 a 2 ij 0 , θ i = b i c P d j 0 =1 a 2 ij 0 , i = 1 , . . . , n ; j = 1 , . . . , d. (4) 2. Solve the minimization pr oblem (2) (a c onvex optimization pr oblem) in p olynomial time and denote its minimum by m . If m ≥ 1 , then output True . If the pr oblem is infe asible, output False . Otherwise c ontinue. 3. Dr op e ach c onstr aint which satisfies θ i > 1 . 4. F orm the instanc e of the pr oblem SphCo vV er fr om the r emaining c onstr aints and solve it. Output the c omplementary r esult. The follo wing theorem pro v es the correctness of Algorithm QP-SphCo vV er . Theorem 2. Algorithm QP-SphCo vV er p olynomial ly r e duc es the pr oblem ConNDQPd to the pr oblem SphCo vV er . Pro of. Note that by (4) we form the equiv alen t problem of form (2) so that k t i k = 1 and c = 1. First assume that m > 1. Then all the feasible p oints satisfy k x k 2 > 1 and hence the output of problem ConNDQPd is True . Assume that m = 1. By assumption, the feasible set of the problem 4 (3) is a non-degenerated p olytop e. Hence it cannot belong to the unit h yp ersphere (otherwise, it would reduce to the single p oint) and there is a feasible p oint x 0 so that k x 0 k > 1. It implies that the answ er of ConNDQPd is True . No w assume that m < 1. Let ˜ x b e the solution of minimization problem (2). Then the following holds | t i 1 ˜ x 1 + t i 2 ˜ x 2 + . . . + t id ˜ x d | = | ( t i , ˜ x ) | ≤ k ˜ x kk t i k = k ˜ x k = m < 1 . (5) Equation (5) implies θ i > − 1, since ˜ x is a feasible point. Without loss of generality , assume that θ 1 , θ 2 , . . . , θ p ≤ 1 and θ p +1 , θ p +2 , . . . , θ n > 1. W e can consider cones C 1 = C ( t 1 ; θ 1 ) , C 2 = C ( t 2 ; θ 2 ) , . . . , C p = C ( t p ; θ p ) , since − 1 < θ i ≤ 1 and k t i k = 1 for all i = 1 , 2 , . . . , p . Also consider the corresp onding system of first p + 1 equations from (1): x 2 1 + x 2 2 + . . . + x 2 d = 1 t 11 x 1 + t 12 x 2 + . . . + t 1 d x d < θ 1 . . . t p 1 x 1 + t p 2 x 2 + . . . + t pd x d < θ p . (6) If cones C 1 , C 2 , . . . , C p co ver R d , it implies that the system (6) (and also (1)) do es not hav e a solution. According to Lemma 1, there m ust hold M ≤ 1, which implies that the answ er of ConNDQPd is False . Assume that cones do not co ver R d and denote by ¯ x the solution of the system (6). Now, since the follo wing relation holds for i = p + 1 , p + 2 , . . . , n : t i 1 ¯ x 1 + t i 2 ¯ x 2 + . . . + t id ¯ x d = ( t i , ¯ x ) ≤ k ˜ x kk t i k = 1 < θ i , w e conclude that ¯ x is also a solution of (1). According to Lemma 1, there holds M > 1 and the answ er of ConNDQPd is True . This completes the pro of. 3.2 NP completeness of ConNDQPd and SphCo vV er W e pro ve that the problem ConNDQPd is NP hard by reducing it to the k -clique decision problem. Recall that, for a given graph G = ( V , E ), set C l ⊆ V is a clique , if for every u, v ∈ C l holds { u, v } ∈ E . In other words, a clique is ev ery set C l of vertices, so that each t w o vertices from C l are adjacent. Clique C l is called k -clique , if it con tains exactly k vertices. The k -clique decision problem (see e.g., [6]) can b e form ulated as follows: Problem 3. ( k -Clique ) Given a gr aph G = ( V , E ) , che ck if ther e exists k -clique. It is known, ([6]) that the problem k -Clique is NP complete. The following lemma demonstrates that it can be p olynomially reduced to the problem ConNDQPd . Lemma 3. Pr oblem k -Clique c an b e p olynomial ly r e duc e d to pr oblem ConNDQPd . Pro of. F or a given graph G = ( V , E ) with the vertex set V = { 1 , 2 , . . . , n } , consider the following instance of problem ConNDQPd : x 2 1 + x 2 2 + . . . + x 2 n > n − − 1 ≤ x i ≤ 1 x i + x j ≤ 0 , ∀{ i, j } / ∈ E x 1 + x 2 + . . . + x n ≥ 2 k − n (7) 5 Here 0 < < 1 and its v alue will b e determined later. If there exists a clique C l of length k in graph G , then b y setting x ∗ i = ( 1 , i ∈ C l − 1 , i / ∈ C l w e obtain one feasible solution ( x ∗ 1 , x ∗ 2 , . . . , x ∗ n ) of the problem (7) satisfying ( x ∗ 1 ) 2 + ( x ∗ 2 ) 2 + . . . + ( x ∗ n ) 2 = n > n − . No w assume that there is no clique of length k in graph G . W e sho w that the decision problem (7) do es not ha ve the solution. Consider the follo wing auxiliary optimization problem max x 1 + x 2 + . . . + x n s . t . − 1 ≤ x i ≤ 1 , i = 1 , 2 , . . . , n x i + x j ≤ 0 , ∀{ i, j } / ∈ E x 2 1 + x 2 2 + . . . + x 2 n > n − . (8) Let ( x 1 , x 2 , . . . , x n ) b e an arbitrary feasible solution. It can b e easily c hec ked that x 2 i > 1 − (due to the quadratic condition in (8)) for i = 1 , 2 , . . . , n implying that x i ∈ [ − 1 , − √ 1 − ) ∪ ( √ 1 − , 1]. Let x i 1 , x i 2 , . . . , x i p ∈ ( √ 1 − , 1] and x q ∈ [ − 1 , − √ 1 − ) for q / ∈ { i 1 , i 2 , . . . , i l } . F or arbitrary 1 ≤ r < s ≤ p there holds { i r , i s } ∈ E , according to x i r + x i s ≥ 2 √ 1 − > 0. In other words, v ertices i 1 , i 2 , . . . , i p form a clique of length p in graph G . According to our assumption that there is no clique of length k in G , it m ust hold that p < k . F urthermore it holds that f = x 1 + x 2 + . . . + x n = p X j =1 x i j + X q / ∈{ i 1 ,i 2 ,...,i p } x q ≤ p − ( n − p ) √ 1 − < p − ( n − p )(1 − ) = 2 p − n + ( n − p ) . (9) No w by choosing = 2 /n and b y using k > p and (9), w e obtain f < 2 p − n + 2 n ( n − p ) < 2 p − n + 2( k − p ) n − p ( n − p ) = 2 k − n. According to the previous expression, each feasible solution of (8) satisfies f < 2 k − n implying that the system (7) has no solutions. As a direct consequence of the Lemma 3 and the NP completeness of the problem k -Clique , and since a verification of a solution for the eq. (3) is p ossible in p olynomial time, the follo wing corollary holds: Corollary 4. Pr oblem ConNDQPd is NP c omplete. No w Corollary 4 and Theorem 2 directly imply: Theorem 5. Spheric al c over age verific ation, i.e., the pr oblem SphCo vV er , is NP har d. 4 Algorithms for spherical co v erage v erification The simplest metho d for spherical co verage v erification is to apply a non-deterministic Mon te-Carlo approac h. The idea is to generate a large num b er N (for example N = 10 10 ) of pseudo-random p oin ts x distributed uniformly on the sphere. F or each generated p oint x we chec k if there is h yp ercap K i con taining the point. If the answer for any p oin t is negativ e, the algorithm outputs False . If the system of h yp ercaps K i , i = 1 , . . . , n co vers the unit sphere, this metho d alw ays outputs the correct answ er True . 6 If there is no cov erage, this m etho d will output the correct answer False with probability that increases with N . Ho wev er, there is no guaran tee that the algorithm will not return false p ositiv es (thus providing answ er True for a non-co vering system of h yp ercaps). Hence, in the subsequen t subsections, we discuss algorithms for spherical cov erage v erification based on application of quadratic programming, and on the reduction of the given problem to lo w er dimensional subproblems. 4.1 QP-based v erification According to Section 2.1, cones C i , i = 1 , . . . , n cov er the space R d if an only if the system (1) has no solutions. According to Lemma 1, to chec k whether the system (1) has solutions, w e need to solv e the QP problems (2) and to determine whether the minimal v alue m and maximal v alue M of the ob jectiv e function satisfy m < 1 < M . The follo wing algorithm for solving problem SphCovV er arises from the aforemen tioned discussion: Algorithm 2. Cov er-QP Input: Caps K i , i = 1 , . . . , n define d by t i ∈ S d (1) and θ i ∈ ( − 1 , 1) . 1. R eturn True if QP pr oblems (2) ar e de gener ate d. Otherwise, solve b oth pr oblems. 2. If any of the pr oblems is not fe asible, or m < 1 < M do es not hold, r eturn True . Otherwise r eturn False . T o apply algorithm Cov er-QP w e need appropriate QP problem solv ers. Note that the minimization problem can b e solved in p olynomial time, since it is con vex. Since conca ve QP is NP complete (see Section 3.2), one of known heuristics can b e applied ([17]). 4.2 Recursiv e algorithm In this section, w e describ e our recursive algorithm for solving the problem SphCovV er . The main idea is to reduce the initial d -dimensional problem to several d − 1-dimensional problems. More precisely , using d-dimensional in version, hypercaps are mapp ed in to regions consisting of h yp erspheres in d − 1 dimensional plane and their exteriors/interiors dep ending on the p osition of the cen ter of inv ersion w.r.t. a hypercap. As prov en in the next subsection, a d -dimensional sphere is co vered b y the hypercaps if and only if the resulting d − 1 plane is completely cov ered by the d − 1 dimensional regions. In turn, this may b e true if a b oundary of eac h region (whic h is itself a d − 1 dimensional h yp ersphere) is completely cov ered b y d − 1 dimensional hypercaps defined b y corresp onding regions. Th us, the co verage of d -dimensional sphere reduces to cov erage of d − 1 dimensional spheres. The follo wing subsections pro vide the rationale for the algorithm and discuss the algorithm formally . 4.2.1 Rationale of the algorithm W e restate the w ell-known definition of the in v ersion in R d . Definition 1. L et c ∈ R d b e a given p oint, and let R b e a p ositive r e al numb er. Inversion ψ c,R ( x ) is a function ψ c,R : R d → R d that maps every p oint x ∈ R d to a p oint y so that: y = c + R 2 k x − c k 2 ( x − c ) . (10) Point c is c al le d the c enter of inversion ψ c,R and R is the r adius of inversion. Let us apply inv ersion ψ = ψ (1 , 0 ,..., 0) , 1 on caps K i and unit hypersphere S d (1). It is well-kno wn that an image of a hypersphere, b y in version whose cen ter b elongs to the h yp ersphere is a hyperplane. Thus, the image of the unit hypersphere is h yp erplane x 1 = 1 / 2. 7 Denote b y ∂ X the b oundary of a giv en set X , Particularly , w e denote by D i = ∂ K i the boundary of cap K i . Also, for a given hypersphere S denote by in t S its interior and b y ext S its exterior. Images of caps K i , i = 1 , . . . , n consist of d − 1-dimensional hyperspheres, b elonging to h yp erplane x 1 = 1 / 2 and their exteriors or interiors, dep ending whether the center of inv ersion is outside or inside the cap. More precisely , the follo wing Lemma holds: Lemma 6. L et c = (1 , 0 , . . . , 0) ∈ R d and assume that c / ∈ D i for every i = 1 , . . . , n . Image of D i , by an inversion ψ c, 1 is a d − 1 -dimensional hyp erspher e S i with c enter β i = ( β i 1 , . . . , β id ) and r adius r i . Mor e over, the image of K i is R i = S i ∪ ext S i , if c ∈ K i and R i = S i ∪ in t S i , if c / ∈ K i . V alues β i and r i ar e given by the fol lowing expr essions: β i 1 = 1 2 , β ij = t ij 2( θ i − t i 1 ) , r i = q 1 − θ 2 i 2( θ i − t i 1 ) . (11) Pro of. T ranslate the co ordinate system to the center of in version, i.e., to p oint c = (1 , 0 , . . . , 0). The equation describing cap K i b ecomes ( x 0 = x − c ): t i 1 x 0 1 + . . . + t id x 0 d ≥ θ i − t i 1 . (12) Let ψ ( x ) = y . By the in volution prop ert y of inv ersion, we can conclude that ψ ( y ) = x or in other w ords: x 0 i = 1 P d j =1 y 0 2 j y 0 i . By replacing the last expression into (12) w e obtain: d X j =1 t ij y 0 j ≥ ( θ i − t i 1 ) d X j =1 y 0 2 j . (13) Replacing y = y 0 + c finally yields: R i : d X j =2 ( y j − β ij ) 2 ≥ r 2 i , θ i − t i 1 > 0 ≤ r 2 i , θ i − t i 1 < 0 , 2 ≤ i ≤ d, y 1 = 1 2 . (14) Here condition θ i − t i 1 > 0 is equiv alent to c ∈ K i . In (14) w e denote: β i 1 = 1 2 , β ij = t ij 2( θ i − t i 1 ) , r 2 i = 1 2( θ i − t i 1 ) 2 d X j =2 t 2 ij + t i 1 ( y 1 − 1) θ i − t i 1 − ( y 1 − 1) 2 . (15) The last expression (15) can b e further simplified using P d j =2 t 2 ij = 1 − t 2 i 1 and y 1 = 1 2 in to (11). In Fig. 1 we illustrate Lemma 6 for three spherical caps D i , i = 1 , 2 , 3 and their images. Observ e that c ∈ D 3 , hence D 3 maps to the exterior of sphere S 3 . Due to Lemma 6, problem SphCo vV er reduces to the following h yp erplane cov er verification ( Hp- Co vV er ) problem: Problem 4. ( HpCo vV er ) F or given d − 1 -dimensional hyp erspher es S i , with c enter β i and r adius r i b elonging to the hyp erplane x 1 = 1 / 2 , and sets R i so that R i = S i ∪ int S i or R i = S i ∪ ext S i ( i = 1 , . . . , n ), che ck if sets R i c over the whole hyp erplane. 8 -2 0 2 x 3 c(1,0,0) D 2 D 3 S 1 S 2 S 3 -1 -0.5 0 0.5 1 -2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2.5 x 1 x 2 Figure 1: Three spherical caps and their images under the inv ersion ψ for d = 3. T o simplify further discussion, if R i = S i ∪ int S i set R i is called internal ; otherwise we call it external . If all the sets R i are in ternal, then w e can immediately conclude that the answer to the problem HpCo vV er is False . This holds obviously from the fact that all sets R i are b ounded and hence is also their union S n i =1 R i . Supp ose, in contrast, that at least one set R i is external. If all sets R do not co ver the whole h yp erplane, there exists one h yp ersphere S i and p oint x ∈ S i so that it is uncov ered by other sets R j , j 6 = i . In other words, the following theorem holds: Theorem 7. Sets R 1 , . . . , R n with differ ent b oundaries S i = ∂ R i c over the whole hyp erplane x 1 = 1 / 2 if and only if every hyp erspher e S i is c over e d by other sets R j , j 6 = i , i.e., S i = ∂ R i ⊂ [ j 6 = i R j . (16) Pro of. ( ⇐ :) If sets R 1 , . . . , R n co ver the whole h yp erplane x 1 = 1 / 2, each hypersphere, S i as a subset of the h yp erplane will b e co vered. ( ⇒ :) Denote an unco vered region of the h yp erplane x 1 = 1 / 2 with Q . The boundary ∂ Q consists of the union of spherical caps. Denote by A one of those caps and by S i the hypersphere whic h A b elongs to. Since S i 6 = S j for i 6 = j , interior int A cannot be cov ered b y remaining h yp erspheres S j , j 6 = i . In other w ords, there exists a p oint x ∈ in t A \ S i 6 = j S j . P oint x is co vered b y some set R k and since x / ∈ S k , there holds x ∈ int R k . Hence, there exists a ball B d − 1 ( x ; δ ) ⊆ int R k and holds B d − 1 ( x ; δ ) ∩ Q = ∅ . This is the con tradiction with the fact that x ∈ A is a b oundary p oin t of Q . If S i = S j for i 6 = j then either R i and R j co ver the h yp erplane (if one of them is in ternal and the other is external) or R i = R j . In the second case, we can eliminate one of them and contin ue. 9 According to Theorem 7 w e need to chec k whether eac h h yp ersphere S i is co vered b y sets R j , j 6 = i . W e distinguish the following cases, dep ending on whether the pairs of h yp erspheres S i and S j = ∂ R j are disjoin t: Case 1. Hyp erspheres S i and S j ha ve a nonempty in tersection. In such a case S i ∩ R j is a h yp ercap K i j defined as: K i j = ( x ∈ S i | ( x − β i , x i j ) k x − β i kk x i j k ≥ θ i j ) , i 6 = j, (17) where w e define: θ i j , r 2 i + d 2 ij − r 2 j 2 r i d ij , R j is in ternal − r 2 i + d 2 ij − r 2 j 2 r i d ij , R j is external , x i j , ( d − 1 ij ( β j − β i ) , R j is in ternal − d − 1 ij ( β j − β i ) , R j is external , d ij , k β i − β j k , i 6 = j. (18) Observ e that for the points x i j from eq. (18), x i j = (( x i j ) 2 , . . . , ( x i j ) d ) ∈ R d − 1 , since ( x i j ) 1 = 0 for ev ery i and j . Case 1 is illustrated in the Fig. 2. Figure 2: In tersection of d − 1-dimensional h yp erspheres S i and S j . Case 2. Hyperspheres S i and S j are disjoin t. The equiv alent condition is θ i j / ∈ ( − 1 , 1), where θ i j is defined b y eq. (18). In suc h case either S i ⊂ R j ( Case 2a ) or S i ∩ R j = ∅ ( Case 2b ) holds. Whic h of these t w o sub cases holds can b e determined e.g., b y choosing the arbitrary p oint (for example x j + ( r j , 0 , . . . , 0)) on ∂ S j and c hecking the inequality (14) for hypersphere S i . F or fixed i , if for any j the condition S i ⊂ R j ( Case 2a ) is satisfied, then eq. (16) holds and S i is co vered. Therefore, the algorithm may contin ue with another v alue of i . Otherwise, it is sufficient to determine whether the sphere S i is cov ered by those hypercaps K j defined by eq. (17) and corresp onding to the pairs of spheres ( S i , S j ) satisfying Case 1 (note that pairs ( S i , S j ) satisfying Case 2b do not need be considered due to disjoin tness of S i and R j ). This is an instance of problem SphCo vV er , for dimension d − 1. Hence, w e reduce original problem SphCo vV er to at most n − 1 equiv alent problems of dimension d − 1. 10 4.2.2 Base case of the algorithm When d = 2, the inv ersion from Lemma 6 maps a 2 D sphere (a circle) in to a straigh t line and 2 D caps (arcs) degenerate into in terv als. Hence, as the base case, we c ho ose case d = 2 of problem HpCo vV er . W e omit the first co ordinate (which is equal to 1 / 2) of each p oint from sets S i and R i . Hence, we assume that there are given sets R i so that R i = [ c i , d i ] ( R i is in ternal) or R i = R \ ( c i , d i ) ( R i is external) for some real num b ers c i < d i . Here c i = β i 2 − r i and d i = β i 2 + r i . The problem is to chec k if sets R i co ver the whole real line R , i.e., S n i =1 R i = R . Without loss of generalit y w e can assume that R 1 , . . . , R s b e external and R s +1 , . . . , R n in ternal. Le t us define: c 0 = max i =1 ,...,s c i , d 0 = min i =1 ,...,s d i . Ob viously , if c 0 ≥ d 0 , R is cov ered b y the sets R i . Otherwise, ( c 0 , d 0 ) = R \ s [ i =1 R i , and w e need to chec k if the interv al ( c 0 , d 0 ) is cov ered b y [ c s +1 , d s +1 ] , . . . , [ c n , d n ]. This can b e p erformed b y sorting segmen ts [ c s +1 , d s +1 ] , . . . , [ c n , d n ] with resp ect to c i and sequentially shortening the target in terv al ( c 0 , d 0 ). Let us assume that the segmen ts are sorted so that c s +1 ≤ c s +2 ≤ . . . ≤ c n . If c 0 < c s +1 then the interv al ( c 0 , c s +1 ) is uncov ered (since c s +1 is minimal); hence R is also uncov ered. Otherwise, we need to c hec k if [ c 00 , d 0 ], where c 00 = max { c 0 , d s +1 } is co v ered by segmen ts [ c s +2 , d s +2 ] , . . . , [ c n , d n ]. This leads to the following algorithm for solving the base case. Algorithm 3. Cov er2-2D ( d = 2 c ase of pr oblem HpCo vV er ) Input: V alues β i 2 , r i ( r i ≥ 0 ) and out i for i = 1 , . . . , n . We assume that out i = True if R i is external and otherwise out i = False . 1. L et c i = β i 2 − r i and d i = β i 2 + r i 2. R e or der the sets R i so that out i = True for i = 1 , . . . , s and out i = False for i = s + 1 , . . . , n . 3. Sort sets R i , i = 1 + s, . . . , n so that c s +1 ≤ . . . ≤ c n . 4. L et c 0 := max i =1 ,...,s c i and d 0 := min i =1 ,...,s d i 5. F or every i = s + 1 , . . . , n do the fol lowing: 5.1. If c 0 ≥ d 0 r eturn True . Otherwise c ontinue. 5.2. If c i > c 0 then r eturn False . Otherwise set c 0 := max { c 0 , d i } and c ontinue. 6. If c 0 ≥ d 0 r eturn True , otherwise r eturn False . Note that the complexit y of the Algorithm Co ver2-2D is O ( n log n ) if an asymptotically optimal sorting algorithm is used for interv als sorting. 4.2.3 Algorithm outline Next, w e formulate complete recursive Algorithm Co v er for solving the general case of problem Sph- Co vV er . Algorithm 4. Cov er Input: Caps K i , i = 1 , . . . , n define d by t i ∈ S d (1) , ( d ≥ 2) , and θ i ∈ ( − 1 , 1) . 1. If n = 1 then r eturn False . 11 2. L et c = (1 , 0 , . . . , 0) ∈ R d . Che ck if c ∈ ∂ K i for some i . In such a c ase, r otate the whole hyp erspher e in plane x 1 x 2 by a smal l angle δ so that c / ∈ ∂ K i , i = 1 , . . . , n (c ondition of L emma 6). 3. Compute ve ctors β i and values r i using e q. (11). If t i 1 > θ i ( c ∈ K i ) set out i = True , otherwise set out i = False . If ther e holds out i = False for every i = 1 , . . . , n then r eturn False . 4. If d = 2 apply A lgorithm Co ver2-2D for β i 2 , r i and out i , i = 1 , . . . , n and r eturn the obtaine d value. Otherwise c ont inue. 5. F or every i = 1 , . . . , n do the fol lowing: 5.1 Determine x i j and θ i j (r elations (18) ), for every j 6 = i , j = 1 , . . . , n . If θ i j / ∈ ( − 1 , 1) , let A = x i j + ( r i , 0 , . . . , 0) . If p oint A b elongs to R j (che ck the r elation (14) ), set i = i + 1 and go to step 5. Otherwise c ontinue. If θ i j ∈ ( − 1 , 1) , form c ap K i j , e q. (17) . If out j = True then set x i j = − x i j and θ i j = − θ i j . 5.2 Apply A lgorithm Co ver on the set of al l forme d c aps K i j . Step 1 of Algorithm Co v er implicitly cov ers the case when the n umber of caps n is smaller than the n umber of dimensions d . In step 2, we introduce the rotation b y a small angle δ . W e may set δ to an arbitrary v alue, for example δ = 0 . 01. If, after the rotation, p oint c is again on arc C i , angle δ needs to b e changed. One p ossibility is to exp onentially decrease it by setting δ = pδ , where 0 < p < 1 (w e used p = 0 . 9) and to rep eat the same pro cedure until the p oint c is not on any arc C i . Step 3 p erforms in version and c hec ks whether the caps map in to external or internal regions. If all regions are in ternal, as discussed earlier, the co verage is False . After this, we chec k whether the base case of recursion is ac hieved. Otherwise, Step 5 chec ks whether the conditions of the Theorem 16 are satisfied and, when needed, performs recursive calls of Co ver . The w orst-case time complexity of algorithm Co v er is exponential in terms of d and polynomial in terms of n as shown by the following theorem: Theorem 8. The time c omplexity of Algorithm Cov er is T ( n, d ) = O n d − 1 log n . Pro of. Algorithm Co ver , Step 5, reduces d -dimensional problem SphCo vV er to at most n problems of size n − 1 and dimension d − 1. The complexity of r e duc e op erations (including inv ersion, eq. (11), and c hecking conditions from Step 5), is O ( n 2 · d ). Hence, the follo wing recursiv e relationship holds T ( n, d ) = O ( n 2 · d ) + nT ( n − 1 , d − 1). Since time complexity of the base case, Algorithm Co ver2-2D , is O ( n log n ), this leads to the statement of the theorem. 4.3 Lo calization of unco v ered p oin t When the solution to SphCovV er problem is False , it ma y b e of in terest to identify a p oint on a h yp ersphere not co vered by an y of the caps. W e demonstrate how this could b e accomplished using results of Algorithm Cov er . The main idea of the prop osed metho d is as follo ws. If unit hypersphere S d (1), is not completely co vered by caps K 1 , . . . , K n , then angles of the corresponding cones C i can be slightly widened so that the resulting system of cones still do es not cov er the space. F urther, as demonstrated in this section, it is p ossible to find a point on the b oundary of the region cov ered b y the enlarged cones, whic h corresp onds to an uncov ered p oint of the original problem. Namely , when the output of Algorithm Co v er is False , we can determine a point u ∈ S d (1) on the b oundary of the cov ered region. Moreov er, suc h a p oint b elongs to an in tersection of several b oundaries D i = ∂ K i , but u / ∈ int K i for ev ery i = 1 , . . . , n . If suc h a p oint b elongs to exactly l b oundaries D i , w e 12 call it ( l, d ) -b oundary p oint. W e prop ose a metho d for computation of the ( l , d )-b oundary p oin t and then extend it to the computation of unco vered point. An example is illustrated in Fig. 3. The main idea is to iden tify a b oundary p oin t in a low er dimensional space (during recursive steps of Co ver ) and connect it with the b oundary p oint of the original problem. Let the output of Algorithm Cov er b e False . Define recursiv ely the sequence of p oints u m ∈ S m (1) and v m ∈ R m b y u m = ψ c m , 1 ( v m ) , v m = β m i + r m i (0 , u m − 1 ) where m = m 0 , m 0 + 1 , . . . , d , m 0 = max { d − l , 2 } and i is the index of uncov ered h yp ersphere in step 5 of m -th recursion call of Algorithm Co ver . Also we denoted c m = (1 , 0 , . . . , 0) | {z } m and l = ( d − 1 , False is returned by Algorithm Cov er2-2D (step 4 of Co ver ) , h − 1 , False is returned at steps 1 or 3 of h -th recursion call of Co ver . while β m i and r m i are corresponding v alues β i and r i obtained from m -th recursion call of Algorithm Co ver . In the first case ( l = 1), initial p oint u 1 is equal to − 1 or 1, depending of whether c i or d i is unco vered by sets R j . Otherwise it is given by u d − h +1 = ( − t d − h +1 1 , False is returned by Step 1 of Co ver , c d − h +1 , False is returned by Step 3 of Co ver . (19) Lemma 9, stated b elow, prov es that u d is ( l, d )-b oundary p oint of the initial problem. Lemma 9. Every p oint u m , m = m 0 , m 0 + 1 , . . . , d , m 0 = min { d − l , 3 } is a ( m − d + l , m ) -b oundary p oint of the c orr esp onding m -dimensional SphCovV er pr oblem in m -th r e cursion c al l of Algorithm Cov er . Pro of. W e first iden tify a b oundary p oint in the base case when the False answer of Algorithm Co ver is detected. Then we prov e by induction that each recursiv e call (step 5.2 of Cov er ) results in an additional boundary to whic h the point b elongs. There are t w o possibilities for the base case of induction: 1. Let the answ er False b e generated b y Algorithm Co ver2-2D . Define u 1 as follo ws u 1 = ( − 1 , p oin t c i is unco vered 1 , p oin t d i is unco vered It is not difficult to observ e that v 2 = (1 / 2 , c i ) in the first and v 2 = (1 / 2 , d i ) in the second case. This point corresp onds to p oin t u 2 = ψ c 2 , 1 ( v 2 ) ∈ S 2 (1) ( c 2 = (1 , 0)) which is a b oundary p oin t of some arc K i and unco v ered b y other arcs K j , j 6 = i . Hence, p oin t u 2 is a required (1 , 2)-b oundary p oin t for the case d = 2. 2. No w let the answer False b e generated b y steps 1 or 3 of Algorithm Co ver at recursion level h . P oint u d − h +1 defined b y eq. (19) do es not b elong to any h yp ercap (interior or boundary) and hence is unco vered. Therefore, u d − h +1 is a (0 , d − h + 1)-b oundary p oin t. W e pro ve the inductiv e step now. Consider the recursiv e call of Algorithm Cov er on m -th level. If it is not otherwise stated, all notation corresponds to the SphCovV er problem b eing solved on m -th recursiv e level.fs Assume that the recursive call in step 5.2 returned False on the i -th subproblem (i.e., co vering of S m − 1 (1) b y caps K i j , j = 1 , . . . , l ). Also assume (b y induction h yp othesis) that returned p oin t u m − 1 is a ( m − 1 − d + l, m − 1)-b oundary p oint. P oint v m = β i + r i (0 , u m − 1 ) b elongs to h yp ersphere S i . Since 13 u m − 1 is a ( m − 1 − d + l , m − 1)-b oundary point, v m do es not b elong to the interior of an y cap K k i and hence do es not b elong to the in terior of R j (i.e., set R j \ S j ) for an y j . This directly implies that u m = ψ c m , 1 ( v m ) is not contained in the in terior of an y hypercap K j . Since p oin t u m is, by assumption, contained in the b oundary of m − 1 − d + l h yp ercaps K i j , there exist spheres S i 1 , . . . , S i m − 1 − d + l so that v m ∈ S i j for j = 1 , . . . , l − 1. Obviously , p oint u m is con tained within the b oundary of l hypercaps K i , K i 1 , . . . , K i m − 1 − d + l . This completes the pro of by induction that u m is a ( m − d + l , m )-b oundary p oint. a) b) c=(1,0,0) u 3 v 3 S 2 c) v 2 Figure 3: a) Spherical caps corresp onding to cones C i with θ = √ 3 / 2 (b oundaries denoted b y solid line) do not co ver S 3 (1) and sperical caps (denoted by colored dashed lines) corresp onding to cones C α i that also do not co ver S 3 (1). b) Inv ersion ψ (1 , 0 , 0) , 1 of spherical cap b oundaries; c) Inv ersion ψ (1 , 0) , 1 of 2D sphere S 2 . Poin t v 2 is b oundary p oint detected in Algorithm Cov er2-2D corresp onding to (1 , 2) b oundary p oin t u 2 and also to p oint v 3 ∈ S 2 and an uncov ered p oint u 3 of S 3 (1). The follo wing Lemma 10 formalizes the fact that given cones C i , i = 1 , . . . , n which do not completely co ver space R d , can b e enlarged (by increasing the cen tral angle for a sufficien tly small v alue) so that the resulting system of cones still do es not co ver the space. Lemma 10. If solution of SphCo vV er pr oblem is False for instanc e of c ones C i = C ( t i ; θ i ) , i = 1 , . . . , n , then ther e exists sufficiently smal l value α so that‘ the solution of SphCovV er pr oblem for c ones C α i = C ( t i ; θ i − α ) , i = 1 , . . . , n , is also False . Pro of. Note that U = R d \ S n i =1 C i is an op en set and there exists at least one in ternal uncov ered p oint x ∈ U . He nce there exists a ball B d ( x ; ρ ) so that B d ( x ; ρ ) ⊂ U . Le t γ i = arccos θ i . All angles γ i can b e enlarged by v alue ∆ γ = 2 arcsin( ρ/ (2 k x k )) and corresponding enlarged cones will not contain point x . Hence it is sufficient to set α = cos(∆ γ ). Note that ( l, d )-b oundary p oint u d,α of the enlarged co verage problem (from Lemma 10) is an in ternal unco vered point of the original problem. This holds since the distance b et w een p oint u d,α and arbitrary cone C i is at least 2 k u d,α k sin arccos( θ i − α ) − arccos( θ i ) 2 . 14 Hence, to find an unco vered p oin t, it is sufficient to determine α (from Lemma 10), resolve the enlarged cov erage problem and compute the boundary p oin t (from Lemma 9). In practice, v alue α can b e computed similarly as rotation angle δ (see step 2 of Algorithm Cov er ). E.g., w e can set α = 0 . 01 and chec k if cones C α 1 , . . . , C α n co ver the space R d . If the result is True , w e can exponentially decrease α (e.g., b y setting α = pα where 0 < p < 1) and repeat the same pro cedure un til the result is False . In our implemen tation, we choose p = 0 . 9. Due to considerations ab ov e, the complete algorithm for computing an internal unco vered point can b e form ulated as follows: Algorithm 5. FindUncov eredP oint Input: Caps K i , i = 1 , . . . , n define d by t i ∈ S d (1) and θ i ∈ [ − 1 , 1] . 1. Apply A lgorithm Cov er with values t i and θ i . If the r esult is True , r eturn True . Otherwise c ontinue. 2. Set α = 0 . 01 . 3. Set θ α i = θ i − α . Apply Algorithm Cover with values t i and θ α i and c ompute ( l , d ) -b oundary p oint u d . 4. If the r esult is True , set α = 0 . 9 α and go to step 2. Otherwise r eturn False and p oint u d . 5 Numerical examples In this section, w e compare p erformance of prop osed recursive algorithm Co ver for spherical co verage v erification, Section 4.2, with an algorithm based on quadratic programming Co ver-QP , Section 4.1. W e demonstrate that the application of Cov er-QP could lead to false p ositives (cov erage incorrectly v erified) while the prop osed recursive algorithm Cov er do es not suffer from suc h a problem. Moreov er, w e demonstrate that the performance of Co ver is satisfactory in practice, in spite of its w orst-case exp onen tial complexit y . Algorithm Cov er is implemented in programming language C . T o implement Co ver-QP , we utilized the programming pack ages Mathematica and Matlab . T o test the influence of differen t non-con vex QP solv ers on the results of the algorithm, we also created an AMPL mo del [9] for the Algortihm Co ver-QP and tested it using MINOS and FortMp solv ers. Implemen tations are tested on sev eral test examples. In the exp eriments, h yp ercaps are determined b y constellations t i , i = 1 . . . , n of p oints and cones ha ve a constan t angle, i.e., θ i = θ . This stipulation comes from applications in metho ds for finding inv erse k -nearest neighbors. Namely , an algorithm for rev erse k -nearest neighbor problem from [1, 15] requires the cov ering constellation with minimal n and θ = cos( π / 6) = √ 3 / 2. T est constellations are generated by the relaxation algorithm from [14]. This algorithm pro duces near-uniform placement of the p oin ts on d -dimensional h yp ersphere. It starts with a randomly generated set of initial p oints, where p oints interact through generalized electromagnetic in teractions and eac h p oin t has equal c harge. The algorithm seeks the solution of the d -dimensional generalization of the Thomp- son’s problem [2], and searches iterativ ely for the equilibrium state (the state with minimal electrostatic energy). 5.1 Accuracy of algorithms Our exp erimental results indicate that Algorithm Cov er-QP , whic h utilizes solv ers for conca v e QP problems, can b e numerically very unstable. As a consequence, the result is a p oten tially large n umber of false p ositiv es (an algorithm falsely indicates that a sphere is cov ered b y caps). Consider constellation four D 85 obtained by relaxation algorithm for d = 4 and n = 85 (the whole constellation is giv en in the App endix and can b e found at tesla.cis.desu.edu/data/Constellations ). 15 Matlab implementation of the Algorithm Co v er-QP returns True . W e tested Mathematica implemen ta- tions for differen t w orking precisions (double precision, 20 digits, 50 digits and 100 digits). In all testings, result was True , but the corresp onding optimal p oin t and ob jectiv e function v alues were differen t (see the follo wing table). Results of testing AMPL mo del w ere similar. Ob jective Optimal point function v alue x 1 x 2 x 3 x 4 Matlab : double precision 0.9756 0.9477 -0.0063 -0.0004 0.2314 Mathematica : double precision 0.935851 0.172452 0.898297 -0.158036 -0.272394 Mathematica : 20 digits 0.891469 -0.390757 0.136777 0.77625 0.342789 Mathematica : 50 digits 0.946319 0.652345 -0.476475 0.236943 0.487436 Mathematica : 100 digits 0.958959 -0.901428 0.034657 -0.231813 0.302405 MINOS 0.953615 0.628285 0.469643 0.426044 -0.39597 FortMP 0.947382 0.037172 0.869515 0.0528789 0.224396 This example clearly sho ws the instabilit y of Algorithm Cov er-QP . One reason for such instabilit y ma y be the fact that optimization heuristics em b edded in the QP solver traps into the local maxim um and do not achiev e the global maximum. Suc h a conclusion is supp orted by the fact that the optimal solution found b y the solv er changes when simple v ariable transformation x = x 0 + c ( c ∈ R d is a constan t v ector) is applied. F or example, if w e put c = (0 . 1 , 0 . 4 , − 0 . 7 , − 0 . 3) (this choice comes from the unco v ered point sho wn in the next paragraph) all solvers return ob jectiv e function v alue 1 . 000409053 and the optimal p oin t (0 . 134335 , 0 . 4576476 , − 0 . 7915343 , − 0 . 3826162). Since the ob jective function v alue is greater than 1, in such a case algorithm returns False . Note that Algorithm Cov er-QP can hav e only false positives (for a non-cov ering constellation, pro- viding answer True ), not false negatives. This can b e explained by the fact that a p oint returned from the QP solv er is alwa ys a feasible solution and the ob jective function v alue in this p oint cannot b e larger than the global maximum. Unlik e Co v er-QP , Algorithm Cov er is not based on heuristics and randomized initial conditions, but on relatively simple and numerically stable algebraic operations. During our experiments, we could not iden tify an y case where Cov er would return incorrect results. Particularly , on the same test example as ab o v e (the constellation four D 85 ) Co v er correctly returns False , and detects the following unco vered p oin t: x uncov ered = (0 . 134309 , 0 . 457496 , − 0 . 791181 , − 0 . 383002) . Since max 1 ≤ i ≤ 85 ( x uncov ered , t i ) = 0 . 865901 < cos( π / 6) (see the App endix for v ectors t i of the constellation four d 85 ) w e can verify that x uncov ered is indeed uncov ered. Note that the result of the Monte-Carlo based approach metho d (see Section 4) on four D 85 with N = 10 10 pseudo-random p oints w as also incorrect (i.e., True ). This again demonstrates wh y Algorithm Co ver seems to b e the only reliable metho d for spherical co verage verification. 5.2 P erformance of recursive algorithm W e test the w orking time of the Algorithm Co ver . Implementation is compiled b y a GNU C compiler and runs on the mac hine with AMD Phenom I I CPU on 3.0 GHz and CentOS 5.2 (Lin ux) operating system. T esting is p erformed for d = 3 , 4 , 5 and 200 ≤ n ≤ 500. Results are shown in the Fig. 4. Execution times shown are obtained b y av eraging through 20 runs of the program on different constellations of the same dimensionalit y . All testing constellations corresponded to a sphere co v ered b y caps (output of the Algorithm Co ver w as True ). This restriction is added since the w orking time on constellations without co verage is considerably smaller due to the fact that Algorithm Co ver do es not complete recursion. As it can b e seen from the graph, Algorithm Co ver is practically applicable and do es not reach its theoretically obtained complexit y . This can b e explained by the fact that not all d − 1-dimensional 16 2 0 0 2 5 0 3 0 0 3 5 0 4 0 0 4 5 0 5 0 0 0 1 2 3 4 5 6 7 8 R u n n i n g t i m e o f a l g o r i t h m C o v e r [ s e c ] N u m b e r o f h y p e r c a p s i n c o n s t e l l a t i o n ( n ) d = 3 d = 4 d = 5 Figure 4: Av erage running time of Algorithm Cov er for d = 3 , 4 , 5. h yp erspheres S i in tersect (Case 1 of subsection 4.1 holds). Practically , the n umber of hyperspheres in tersecting S i is drastically smaller than n , which implies that the corresp onding subproblems ha ve a lo wer dimensionality . The applicabilit y of Algorithm Co v er made possible to determine the upper bound M u ( d ) of minimal n umber M ( d ) of hypercaps with θ = cos( π / 6), co vering the unit hypersphere. In [15], it is pro v ed that M ( d ) is at least M l ( d ) = Θ(2 d √ d ), i.e., it is exp onential as the function of d . In order to obtain M u ( d ), w e generate sev eral constellations for fixed v alue n . If the cov ering constellation is found, n is dec reased b y 1 and the whole pro cedure is repeated. Otherwise w e put M u ( d ) := n . The follo wing table shows the v alues of M u ( d ) for d = 3 , 4 , 5 , 6: d 3 4 5 6 M u ( d ) 22 81 234 715 Ho wev er, the closed-form expression for upper b ound M u ( d ) is still unknown and its existence is left as an op en question. 6 Conclusion W e ha ve considered the spherical cov erage problem: given a set of hypercaps in d -dimensional space, determine whether a d -dimensional h yp ersphere is completely cov ered b y the h yp ercaps. W e ha ve demon- strated that the considered problem is NP hard by reducing concav e quadratic programming (QP) prob- lem to it. W e ha ve discussed tw o algorithms to resolv e the spherical cov erage problem: the first metho d (Algorithm Co ver-QP ) is based on the utilization of quadratic programming. The second metho d (Algo- rithm Co ver ) is recursiv e and based on the reduction of the main problem on O ( n ) problems of dimension d − 1. The recursiv e algorithm also provides a metho d to determine an unco v ered p oint (if such a p oint exists). While the w orst-case time complexit y of the prop osed recursive metho d is O n d − 1 log n , it is of practical interest, due to NP hardness of the considered problem (note also that the asymptotic complexity could b e impro ved if, as the base case, we c ho ose the case d = 3 of problem HpCo vV er and utilize a metho d prop osed in [4] to resolve it). How ev er, n umerical exp eriments indicate that the recursive algorithm almost never reac hes maximal complexity and hence t ypically do es not require prohibitiv e execution time. In contrast, Algorithm Cov er-QP , with heuristics-based conca ve QP problem solv ers, can b e n umerically unstable resulting in false positive detection (indicating false cov erage of the sphere) and hence having limited practical application. 17 Our results indicate that the recursiv e algorithm ma y b e the b est algorithm for the problems having a relativ ely small dimension d . F or the high v alues of d , where direct application of the recursiv e algorithm ma y b e to o time consuming, quadratic programming metho d still could b e used as a presolve method, since it do es not ha ve false negatives. W e conclude the pap er with the following op en problem: Problem 5. Develop an efficient algorithm for c onstruction of the c overing c onstel lation of minimal size, for a given dimension d and angle θ . Ac kno wledgemen t The study is supp orted b y the US Department of T ransp ortation Office of the Secretary Grant No. DTOS59-08-G-0014. M. Petk o vi´ c is supp orted b y the researc h pro ject 174013 of the Serbian Ministry of Science and DoD/DoA (a ward 45395-MA-ISP). D. Pokra jac has also b een partially supp orted by NIH (gran t #2 P20 RR016472-04), DoD/DoA (aw ards 45395-MA-ISP , P-54412-CI-ISP) and NSF (aw ards #0320991, CREST gran t (#HRD-0630388 CREOSA).). P art of this research w as done while M. P etko vi´ c w as visiting Delaw are State Univ ersity , Departmen t of Computer and Information Sciences. Authors extend their sp ecial thanks to Dr. Stephen A. V av asis, professor and univ ersity research c hair of Universit y of W aterlo o, Ontario, Canada, for help in formulating pro of for NP completeness of the concav e QP problem discussed in the pap e r. Dr. Richard McCallister, an asso ciate professor of Dela ware State Univ ersity English and F oreign Languages Departmen t, help ed us in proofreading the man uscript and impro ving its English. In addition, authors wish to thank to Saman tha McDaniel for help on dra wing the figures. Also, authors are thankful to anon ymous reviewers whose inputs significan tly help ed impro ving the quality of the man uscript. References [1] J. Anderson, B. Tjanden, The inverse ne ar est neighb or pr oblem with astr ophysic al applic ations , Pro- ceedings of the 12th Symp osium of Discrete Algorithms (SOD A), W ashington, DC, January 2001. [2] N. Ashb y , W.E. Brittin, Thomson ’s Pr oblem , Amer. J. Ph ys. 54, pp 776–777, 1986. [3] G. Barequet, M. T. Dic kerson, Y. Sc harf, Covering p oints with a p olygon , Comput. Geom. 39 (2008), 143-162. [4] F. Cazals, S. Loriot, Computing the Arr angements of Cir cles on a Spher e with Applic ations in Structur al Biolo gy , Comput. Geom. 42:6-7 (2009), 551–565. [5] S. Cabello, J. M. Diaz-Banez, C. Seara, J. A. Sellares, J. Urrutia, I. V entura, Covering p oint sets with two disjoint disks or squar es , Comput. Geom. 40 (2008), 195-206. [6] T. H. Cormen, C. E. Leiserson, R. L. Riv est, C. Stein Intr o duction to Algorithms, 2nd e dition , The MIT Press, Boston, 2001. [7] I. Dumer, Covering Spher es with Spher es , Discrete and Computational Geometry 38, 665-679 (2007). [8] K. Elbassioni, H.R. Tiw ary , On a c one c overing pr oblem , Comput. Geom. (2010), doi: 10.1016/j.comgeo.2010.07.004. [9] R. F ourer, D.M. Gay , B.W. Kernighan, AMPL: A Mo deling L anguage for Mathematic al Pr o gr am- ming , Duxbury Press, Bro oks/Cole Publishing Company , 2002. 18 [10] R.M. F reund, J.B. Orlin, On the c omplexity of four p olyhe dr al set c ontainment pr oblems , Mathemat- ical programming 33 (1985), 139–145. [11] J. Han, M. Kam b er, Data Mining: Conc epts and T e chniques , Morgan Kaufmann Publishers, 2001. [12] F.J. MacWilliams, N.J.A. Sloane, The the ory of err or-c orr e cting c o des , North-Holland mathematical library , 2006. [13] D. Pokra jac, A. Lazarevic, L.J. Latecki, Incr emental L o c al Outlier Dete ction for Data Str e ams , Pro ceedings IEEE Symp osium on Computational In telligence and Data Mining, CIDM 2007, pp 504–515. [14] D. Pokra jac, J. Milutino vi´ c, I. Ek anem, Applic ation of ele ctr ostatic metho d for uniform plac ement of p oints on a hyp erspher e , Pro ceedings of the international conference of applied electromagnetics - PES 2007, Ni ˇ s, Serbia. [15] D. Pokra jac, M.D. P etko vi ´ c, L.J. Latec ki, A. Lazarevi ´ c, N. Reljin, J. Milutino vi ´ c, Computational Ge- ometry Issues of R everse Ne ar est Neighb or A lgorithm , Haw aii In ternational Conference on Statistics, Mathematics and Related Fields, Honolulu, HI, Jan uary 2008. [16] D. P okra jac, N. Reljin, N. Pej ˇ ci ´ c, A. Lazarevi´ c, Incr emental Conne ctivity-Base d Outlier F actor A lgorithm , Proc. Visions in Computer Science, 2008. [17] S.A. V a v asis, Nonline ar Optimization: Complexity Issues , Oxford Universit y Press, USA, 1991. [18] J.-L. V erger-Gaugry , Covering a Bal l with Smal ler Equal Bal ls in R n , Discrete and Computational Geometry 33, 143-155 (2005). A Constellation four D 85 t 1 t 2 t 3 t 4 0.911722 0.083517 -0.402106 0.009974 0.581301 -0.242364 0.278539 -0.725096 0.218664 0.316284 0.287611 0.877172 0.117213 0.928431 -0.32201 0.143484 -0.737931 -0.158832 -0.035717 0.654946 0.36765 -0.515669 -0.616088 0.468352 0.76061 0.150628 -0.006887 0.631456 -0.321449 -0.24082 0.870032 0.285869 0.344644 0.293037 0.546223 -0.704976 -0.412761 0.540358 0.013305 -0.73312 0.816541 -0.164327 0.467615 0.295963 -0.245525 -0.585043 -0.116853 0.76406 0.338649 -0.178681 -0.219111 0.89743 -0.69018 -0.518921 0.382654 0.328555 0.740476 0.306078 0.049102 -0.596322 -0.706023 -0.187694 0.659171 -0.178316 -0.055919 -0.589895 -0.805087 -0.027052 0.328685 -0.068095 0.878784 0.339218 0.943282 -0.273977 0.075663 -0.171553 -0.159834 0.605736 0.643967 -0.43914 0.059374 0.684132 0.716063 0.125265 0.562506 0.316235 0.756618 -0.105413 0.526948 -0.692313 -0.153046 -0.468621 0.585505 -0.29157 0.718579 -0.236252 -0.339912 -0.865955 -0.293451 0.220156 19 t 1 t 2 t 3 t 4 0.504902 0.689107 -0.350765 -0.383627 -0.418404 0.437219 0.203124 0.769752 0.498145 -0.17126 -0.337765 -0.780023 -0.19152 -0.885608 -0.106055 -0.409599 -0.216553 0.126561 -0.29997 0.920383 0.060341 0.820685 0.197331 0.532819 -0.736942 -0.41477 -0.460301 0.270194 -0.547833 -0.450474 -0.569484 -0.415501 0.204308 -0.777431 0.433251 -0.407619 -0.983534 -0.127577 0.121635 0.039876 -0.106398 0.780104 -0.465982 -0.403706 -0.824089 0.406231 0.311938 -0.241967 0.925197 0.302318 0.226185 -0.038143 0.582566 0.128159 -0.589217 0.544991 0.622845 -0.704481 0.326492 0.095782 -0.836308 0.452032 -0.011796 0.310027 -0.309656 0.548631 -0.775807 0.035212 -0.225897 -0.143505 0.308569 0.912777 0.656179 -0.438453 -0.607877 -0.087605 -0.12433 -0.638736 0.476572 0.591132 -0.425256 0.298745 0.854262 -0.012046 0.480447 0.819603 0.259428 -0.173544 -0.514513 0.096686 0.503361 -0.687428 0.267199 0.5784 -0.305584 0.707585 -0.018035 -0.020976 -0.880919 -0.472457 -0.736443 -0.653951 -0.063084 -0.161308 0.033759 -0.030823 0.959959 -0.276386 -0.471801 0.343863 -0.601331 -0.545493 -0.349645 0.727398 -0.321563 0.495215 0.004942 -0.96662 0.232848 0.106783 0.111029 -0.625644 0.768139 0.078766 0.060656 -0.490269 -0.013342 -0.869356 0.588121 0.537068 0.416267 0.438626 0.338081 -0.727789 0.016177 0.596458 -0.342943 -0.152519 -0.336818 -0.863529 -0.868041 0.131346 -0.46866 -0.098029 0.084497 -0.533808 -0.58593 -0.603817 0.555933 0.174623 -0.695937 -0.419663 0.409128 -0.265222 0.414691 0.768312 -0.398058 -0.719619 0.526168 -0.216435 -0.230651 -0.326573 -0.707467 0.582787 -0.629986 0.173969 -0.5524 0.517404 -0.389677 -0.043493 -0.919642 0.0228 -0.013034 0.298631 -0.792652 0.531369 0.186482 0.689421 0.091131 -0.693987 0.814852 -0.386301 -0.176673 0.39443 0.275837 -0.088184 -0.946315 0.143616 -0.497203 -0.491744 0.211614 -0.682786 -0.11172 0.959906 0.17654 -0.186902 -0.028234 0.117661 0.115501 -0.98591 0.083739 0.314023 -0.444543 -0.834721 -0.805373 -0.010126 -0.060896 -0.589545 -0.125314 0.270743 0.737091 0.606376 0.690523 0.637669 -0.183595 0.287835 -0.481723 0.736347 0.436683 0.1872 -0.719834 0.10622 0.559014 0.397567 0.38002 0.543697 -0.745668 0.062905 -0.653788 0.707026 -0.192191 -0.189048 20 t 1 t 2 t 3 t 4 0.013797 -0.290263 0.625477 -0.72411 0.303814 -0.890275 -0.328701 0.084045 B Pseudo co de of Algorithm Co v er Here w e give the complete pseudo-co de description of Algorithm Cov er . 21 Algorithm 1 Cov er ( n ; d ; X ; θ ) - Spherical co verage verification Require: Cen ters and angles of cones: X = [ t T 1 · · · t T n ] T ∈ R n × d and θ = ( θ 1 , . . . , θ n ) 1: if d = 2 then 2: return Co v er2 ( n ; X ; θ ) 3: end if 4: S := (1 , 0 , . . . , 0) ∈ R d 5: if ( S, t i ) = θ i for some i = 1 , 2 , . . . , n then 6: X := Rot ( n ; d ; X ; θ ) (rotates b y small angle δ ). 7: end if 8: r i := √ 1 − θ 2 i 2( θ i − t i 1 ) for i = 1 , . . . , n (radii of the spheres in d − 1-dimensional hyperplane) 9: β ij := t ij 2( θ i − t i 1 ) for i = 1 , . . . , n and j = 2 , . . . , d (cen ters of the spheres in d − 1-dimensional hyperplane) 10: out i := (( S, t i ) < θ i ) ( True if S b elongs to the interior of cap K i and False otherwise) 11: for i := 1 to n do 12: k := 0 13: Cov er i := False 14: for j := 1 to n do 15: if i 6 = j then 16: k := k + 1 17: t 1 k,l − 1 := β j l − β il for l = 2 , 3 , . . . , d 18: d ij := k t 1 k k 19: θ 1 k := r 2 i + d 2 ij − r 2 j 2 r i d ij 20: if − 1 ≤ θ 1 k ≤ 1 then 21: if out j then 22: θ 1 k := − θ 1 k 23: t 1 kl := − t 1 kl for l = 1 , 2 , . . . , d − 1 24: end if 25: else 26: s := q ( β i 2 + r i − β j 2 ) 2 + P d l =3 ( β il − β j l ) 2 27: if ( s ≤ r j and not out j ) or ( s ≥ r j and out j ) then 28: Co ver i := True ( j -th h yp ersphere completely cov ers i -th) 29: break 30: else 31: k := k − 1 ( j -th h yp ersphere is disjoin t from i -th) 32: end if 33: end if 34: end if 35: end for 36: if not Cov er i then 37: if k = 0 then 38: return False 39: else 40: X 1 := [( t 1 1 ) T · · · ( t 1 n ) T ] T ∈ R k × d − 1 41: θ 1 := ( θ 1 1 , . . . , θ 1 k ) 42: Co ver i := Co ver ( k ; d − 1; X 1 ; θ 1 ) 43: if not Co ver i then 44: return False 45: end if 46: end if 47: end if 48: end for 49: return True 22 Mark o D. P etk ovi ´ c w as born in Ni ˇ s, Serbia, in 1984. He graduated in mathematics and computer science at the F aculty of Science and Mathematics, Ni ˇ s, Serbia in 2006. He also graduated in telecom- m unications at F acult y of Electronic Engineering, Nis, Serbia in 2007. Also he received his PhD degree in Computer Science from the Univ ersity of Ni ˇ s in 2008. Currently he is an Assistant Professor in Com- puter Science at the F aculty of Science and Mathematics, Universit y of Ni ˇ s, Serbia. His researc h interests include the theory and computation of matrix generalized inv erses, computational geometry , source and c hannel co ding, symbolic computation and optimization metho ds. He is the author of about 50 papers (ab out 30 of them in p eer-reviewed international journals). Dragoljub P okra jac was b orn in ˇ Sib enik, Croatia. He receiv ed his BS in electrical engineering (1993) and MS in telecommunications at Universit y of Ni ˇ s, Serbia (1997). He attended the PhD program at W ashington State Univ ersity 1998-2000 and received his PhD in Computer Science from T emple Univ ersity in 2002. He is curren tly a tenured Asso ciate Professor at Dela ware State Universit y , Dov er, DE. Also, he w orked as a senior scientist at VIP Mobile, Inc. 2008-2010. His research interests include mac hine learning, computer vision, video analytics, computational geometry and theory of algorithms. Dr. P okra jac is mem b er of ACM and F&AM Equit y Lo dge 591, Pennsylv ania. Longin Jan Latecki received the PhD degree in computer science from Ham burg Univ ersity , Germany , in 1992. He is a Professor of Computer Science at T emple Universit y , Philadelphia. His main researc h in terests include shap e represen tation and similarity , ob ject detection and recognition in images, rob ot p erception, data mining, and digital geometry . He has published o ver 175 research pap ers and b o oks. He is an editorial board member of Pattern Recognition and In ternational Journal of Mathematical Imaging. He receiv ed the ann ual P attern Recognition So ciet y Aw ard together with Azriel Rosenfeld for the b est article published in the journal P attern Recognition in 1998. 23
Original Paper
Loading high-quality paper...
Comments & Academic Discussion
Loading comments...
Leave a Comment