Zonotope approximation (#59)

* use cdhr in rounding, improve t-test iterations, change diameter of H-polytopes, minor mpdifications in both interfaces. * improve c++ interface, add c++ test for cooling bodies with billiard walk * improve new c++ test * add random generators * update R random generators * update Rd files * update R volume interface * fix generators in both c++ and R interfaces, improve exact_volume check in c++ interface * fix bug in hpoly zonotope volume approximation * improve c++ documentation (help command) * add zonotope_approximation in R interface. improve R volume and sample_points documentation * fix c++ tests * fix bug in zonotope approximation Co-authored-by: Vissarion Fisikopoulos <[email protected]>
GeomScale · Feb 19, 2020 · 801b2dc · 801b2dc
1 parent e8b5ac7
commit 801b2dc
Show file tree

Hide file tree

Showing 12 changed files with 359 additions and 15 deletions.
diff --git a/R-proj/NAMESPACE b/R-proj/NAMESPACE
@@ -18,6 +18,7 @@ export(rand_rotate)
 export(round_polytope)
 export(sample_points)
 export(volume)
+export(zonotope_approximation)
 exportPattern("^[[:alpha:]]+")
 importFrom("methods","new")
 importFrom("stats", "cov")

diff --git a/R-proj/R/RcppExports.R b/R-proj/R/RcppExports.R
@@ -171,7 +171,7 @@ rounding <- function(P, random_walk = NULL, walk_length = NULL, parameters = NUL
 #' \item{\code{dimension} }{ An integer that declares the dimension when exact sampling is enabled for a simplex or a hypersphere.}
 #' \item{\code{radius} }{ The radius of the \eqn{d}-dimensional hypersphere. The default value is \eqn{1}.}
 #' \item{\code{BW_rad} }{ The radius for the ball walk.}
-#' \item{\code{diameter} }{ The diameter of the polytope. It is used to set the maximum length of the trajectory in billiard walk.}
+#' \item{\code{L} }{The maximum length of the billiard trajectory.}
 #' }
 #' @param InnerPoint A \eqn{d}-dimensional numerical vector that defines a point in the interior of polytope P.
 #'
@@ -233,7 +233,7 @@ sample_points <- function(P = NULL, N = NULL, distribution = NULL, random_walk =
 #'  \item{\code{prob} }{ The probability is used for the empirical confidence interval in ratio estimation of CB algorithm. The default value is \eqn{0.75}.}
 #'  \item{\code{hpoly} }{ A boolean parameter to use H-polytopes in MMC of CB algorithm. The default value is \code{FALSE}.}
 #'  \item{\code{minmaxW} }{ A boolean parameter to use the sliding window with a minmax values stopping criterion.}
-#'  \item{\code{diameter} }{ The diameter of the polytope. It is used to set the maximum length of the trajectory in billiard walk.}
+#'  \item{\code{L} }{The maximum length of the billiard trajectory.}
 #' }
 #'
 #' @references \cite{I.Z.Emiris and V. Fisikopoulos,
@@ -261,3 +261,22 @@ volume <- function(P, walk_length = NULL, error = NULL, inner_ball = NULL, algo
     .Call(`_volesti_volume`, P, walk_length, error, inner_ball, algo, random_walk, rounding, parameters)
 }
 
+#' An internal Rccp function for the over-approximation of a zonotope
+#'
+#' @param Z A zonotope.
+#' @param fit_ratio Optional. A boolean parameter to request the computation of the ratio of fitness.
+#' @param walk_length Optional. The number of the steps for the random walk. The default value is \eqn{\lfloor 10 + d/10\rfloor} for SequenceOfBalls and \eqn{1} for CoolingGaussian.
+#' @param error Optional. Declare the upper bound for the approximation error. The default value is \eqn{1} for SequenceOfBalls and \eqn{0.1} for CoolingGaussian.
+#' @param inner_ball Optional. A \eqn{d+1} vector that contains an inner ball. The first \eqn{d} coordinates corresponds to the center and the last one to the radius of the ball. If it is not given then for H-polytopes the Chebychev ball is computed, for V-polytopes \eqn{d+1} vertices are picked randomly and the Chebychev ball of the defined simplex is computed. For a zonotope that is defined by the Minkowski sum of \eqn{m} segments we compute the maximal \eqn{r} s.t.: \eqn{re_i\in Z} for all \eqn{i=1,\dots ,d}, then the ball centered at the origin with radius \eqn{r/\sqrt{d}} is an inscribed ball.
+#' @param random_walk Optional. A string that declares the random walk method: a) \code{'CDHR'} for Coordinate Directions Hit-and-Run, b) \code{'RDHR'} for Random Directions Hit-and-Run or c) \code{'BW'} for Ball Walk. The default walk is \code{'CDHR'}.
+#' @param rounding Optional. A boolean parameter for rounding. The default value is \code{FALSE}.
+#' @param parameters Optional. A list for the parameters of the volume algorithm
+#'
+#' @section warning:
+#' Do not use this function.
+#'
+#' @return A List that contains a numerical matrix that describes the PCA approximation as a H-polytope and the ratio of fitness.
+zono_approx <- function(Z, fit_ratio = NULL, walk_length = NULL, error = NULL, inner_ball = NULL, random_walk = NULL, rounding = NULL, parameters = NULL) {
+    .Call(`_volesti_zono_approx`, Z, fit_ratio, walk_length, error, inner_ball, random_walk, rounding, parameters)
+}
+
diff --git a/R-proj/R/zonotope_approximation.R b/R-proj/R/zonotope_approximation.R
@@ -0,0 +1,53 @@
+#' A function to over-approximate a zonotope with PCA method and to evaluate the approximation by computing a ratio of fitness.
+#' 
+#' For the evaluation of the PCA method the exact volume of the approximation body is computed and the volume of the input zonotope is computed by CoolingBodies algorithm (BAN). The ratio of fitness is \eqn{R=}.
+#' 
+#' @param Z A zonotope.
+#' @param fit_ratio Optional. A boolean parameter to request the computation of the ratio of fitness.
+#' @param walk_length Optional. The number of the steps for the random walk. The default value is \eqn{\lfloor 10 + d/10\rfloor} for SequenceOfBalls and \eqn{1} for CoolingGaussian.
+#' @param error Optional. Declare the upper bound for the approximation error. The default value is \eqn{1} for SequenceOfBalls and \eqn{0.1} for CoolingGaussian.
+#' @param inner_ball Optional. A \eqn{d+1} vector that contains an inner ball. The first \eqn{d} coordinates corresponds to the center and the last one to the radius of the ball. If it is not given then for H-polytopes the Chebychev ball is computed, for V-polytopes \eqn{d+1} vertices are picked randomly and the Chebychev ball of the defined simplex is computed. For a zonotope that is defined by the Minkowski sum of \eqn{m} segments we compute the maximal \eqn{r} s.t.: \eqn{re_i\in Z} for all \eqn{i=1,\dots ,d}, then the ball centered at the origin with radius \eqn{r/\sqrt{d}} is an inscribed ball.
+#' @param random_walk Optional. A string that declares the random walk method: a) \code{'CDHR'} for Coordinate Directions Hit-and-Run, b) \code{'RDHR'} for Random Directions Hit-and-Run or c) \code{'BW'} for Ball Walk. The default walk is \code{'CDHR'}.
+#' @param rounding Optional. A boolean parameter for rounding. The default value is \code{FALSE}.
+#' @param parameters Optional. A list for the parameters of the algorithms:
+#' \itemize{
+#'  \item{\code{Window} }{ The length of the sliding window for CG algorithm. The default value is \eqn{500+4dimension^2}.}
+#'  \item{\code{BW_rad} }{ The radius for the ball walk. The default value is \eqn{4r/dimension}, where \eqn{r} is the radius of the inscribed ball of the polytope.}
+#'  \item{\code{ub} }{ The lower bound for the ratios in MMC in BAN algorithm. The default value is \eqn{0.1}.}
+#'  \item{\code{lb} }{ The upper bound for the ratios in MMC in BAN algorithm. The default value is \eqn{0.15}.}
+#'  \item{\code{N} }{ An integer that controls the number of points \eqn{\nu N} generated in each convex body in annealing schedule.}
+#'  \item{\code{nu} }{ The degrees of freedom for the t-student distribution in t-tests in BAN algorithm. The default value is \eqn{10}.}
+#'  \item{\code{alpha} }{ The significance level for the t-tests in BAN algorithm. The default values is 0.2.}
+#'  \item{\code{prob} }{ The probability is used for the empirical confidence interval in ratio estimation of BAN algorithm. The default value is \eqn{0.75}.}
+#'  \item{\code{hpoly} }{ A boolean parameter to use H-polytopes in MMC of BAN algorithm. The default value is \code{FALSE}.}
+#'  \item{\code{minmaxW} }{ A boolean parameter to use the sliding window with a minmax values stopping criterion.}
+#'  \item{\code{L} }{The maximum length of the billiard trajectory.}
+#' }
+#' 
+#' @return A list that contains the approximation body in H-representation and the ratio of fitness
+#' 
+#' @examples 
+#' # over-approximate a 2-dimensional zonotope with 10 generators and compute the ratio of fitness
+#' Z = GenZonotope(2,10)
+#' retList = zonotope_approximation(Z = Z, fit_ratio = TRUE)
+#' 
+#' @export
+zonotope_approximation <- function(Z, fit_ratio = NULL, walk_length = NULL, error = NULL,
+                        inner_ball = NULL, random_walk = NULL, rounding = NULL,
+                        parameters = NULL){
+
+  ret_list = zono_approx(Z, fit_ratio, walk_length, error, inner_ball, random_walk, rounding, parameters)
+
+  Mat = ret_list$Mat
+
+  # first column is the vector b
+  b = Mat[,1]
+
+  # remove first column
+  A = Mat[,-c(1)]
+  PP = list("P" = Hpolytope$new(A,b), "fit_ratio" = ret_list$fit_ratio)
+
+  return(PP)
+
+}
+
diff --git a/R-proj/man/rounding.Rd b/R-proj/man/rounding.Rd
diff --git a/R-proj/man/sample_points.Rd b/R-proj/man/sample_points.Rd
diff --git a/R-proj/man/volume.Rd b/R-proj/man/volume.Rd
diff --git a/R-proj/man/zono_approx.Rd b/R-proj/man/zono_approx.Rd
diff --git a/R-proj/man/zonotope_approximation.Rd b/R-proj/man/zonotope_approximation.Rd
diff --git a/R-proj/src/sample_points.cpp b/R-proj/src/sample_points.cpp
@@ -44,7 +44,7 @@
 //' \item{\code{dimension} }{ An integer that declares the dimension when exact sampling is enabled for a simplex or a hypersphere.}
 //' \item{\code{radius} }{ The radius of the \eqn{d}-dimensional hypersphere. The default value is \eqn{1}.}
 //' \item{\code{BW_rad} }{ The radius for the ball walk.}
-//' \item{\code{diameter} }{ The diameter of the polytope. It is used to set the maximum length of the trajectory in billiard walk.}
+//' \item{\code{L} }{The maximum length of the billiard trajectory.}
 //' }
 //' @param InnerPoint A \eqn{d}-dimensional numerical vector that defines a point in the interior of polytope P.
 //'

diff --git a/R-proj/src/volume.cpp b/R-proj/src/volume.cpp
@@ -106,7 +106,7 @@ double generic_volume(Polytope& P, unsigned int walk_step, double e,
         if (!hpoly) {
             vol = vol_cooling_balls(P, var, var_ban, InnerB);
         } else {
-            vars_g <NT, RNGType> varg(n, 1, N, 4 * n * n + 500, 1, e, InnerB.second, rng, C, frac, ratio, delta,
+            vars_g <NT, RNGType> varg(n, 1, N, 6 * n * n + 500, 1, e, InnerB.second, rng, C, frac, ratio, delta,
                                       verbose, rand_only, false, false, birk, false, true, false);
             vol = vol_cooling_hpoly < HPolytope < Point > > (P, var, var_ban, varg, InnerB);
         }
@@ -147,7 +147,7 @@ double generic_volume(Polytope& P, unsigned int walk_step, double e,
 //'  \item{\code{prob} }{ The probability is used for the empirical confidence interval in ratio estimation of CB algorithm. The default value is \eqn{0.75}.}
 //'  \item{\code{hpoly} }{ A boolean parameter to use H-polytopes in MMC of CB algorithm. The default value is \code{FALSE}.}
 //'  \item{\code{minmaxW} }{ A boolean parameter to use the sliding window with a minmax values stopping criterion.}
-//'  \item{\code{diameter} }{ The diameter of the polytope. It is used to set the maximum length of the trajectory in billiard walk.}
+//'  \item{\code{L} }{The maximum length of the billiard trajectory.}
 //' }
 //'
 //' @references \cite{I.Z.Emiris and V. Fisikopoulos,