CRAN version 1.1.2 (#118)

* fix review comments on R package

* update R interface in volume and sample_points

* update Rcpp modules and S4 classes

* replace modules with R classes

* update R documentation

* fix R tests

* fix R bugs

* add hpoly new generator, fix bugs in rounding with biw

* update the Rd files

* fix c++ tests and a bug in volume.coo

* fix c++ test

* update polytope random generators in R, update rd files

* fix bugs in R interface

* seperate sdp functions

R-proj/DESCRIPTION

@@ -11,13 +11,13 @@ Description: Provides an R interface for 'volesti' C++ package. 'volesti' comput
              (a.k.a. zonotopes). There are three algorithms for volume estimation as well as algorithms 
              for sampling, rounding and rotating polytopes. Moreover, 'volesti' provides algorithms for 
              estimating copulas useful in computational finance.
-Version: 1.1.1
-Date: 2020-06-18
+Version: 1.1.2
+Date: 2020-10-18
 Maintainer: Vissarion Fisikopoulos <[email protected]>
 Depends: Rcpp (>= 0.12.17)
 Imports: methods, stats
 LinkingTo: Rcpp, RcppEigen, BH
 Suggests: testthat
 Encoding: UTF-8
-RoxygenNote: 6.1.1
+RoxygenNote: 7.1.1
 BugReports: https://github.com/GeomScale/volume_approximation/issues

R-proj/NAMESPACE

@@ -4,7 +4,6 @@ export(compute_indicators)
 export(copula)
 export(direct_sampling)
 export(exact_vol)
-export(file_to_polytope)
 export(frustum_of_simplex)
 export(gen_cross)
 export(gen_cube)
@@ -16,13 +15,19 @@ export(gen_simplex)
 export(gen_skinny_cube)
 export(inner_ball)
 export(loadSdpaFormatFile)
+export(load_sdpa_format_file)
 export(readSdpaFormatFile)
 export(rotate_polytope)
 export(round_polytope)
 export(sample_points)
 export(volume)
-export(writeSdpaFormatFile)
+export(write_sdpa_format_file)
 export(zonotope_approximation)
+exportClasses(Hpolytope)
+exportClasses(Spectrahedron)
+exportClasses(Vpolytope)
+exportClasses(VpolytopeIntersection)
+exportClasses(Zonotope)
 exportPattern("^[[:alpha:]]+")
 importFrom("methods","new")
 importFrom("stats","cov")

R-proj/R/HpolytopeClass.R

@@ -0,0 +1,42 @@
+#' An R class to represent an H-polytope
+#'
+#' An H-polytope is a convex polytope defined by a set of linear inequalities or equivalently a \eqn{d}-dimensional H-polytope with \eqn{m} facets is defined by a \eqn{m\times d} matrix A and a \eqn{m}-dimensional vector b, s.t.: \eqn{Ax\leq b}.
+#'
+#' \describe{
+#'    \item{A}{An \eqn{m\times d} numerical matrix.}
+#'
+#'    \item{b}{An \eqn{m}-dimensional vector b.}
+#' 
+#'    \item{volume}{The volume of the polytope if it is known, \eqn{NaN} otherwise by default.}
+#'    
+#'    \item{type}{A character with default value 'Hpolytope', to declare the representation of the polytope.}
+#' }
+#' 
+#' @examples
+#' A = matrix(c(-1,0,0,-1,1,1), ncol=2, nrow=3, byrow=TRUE)
+#' b = c(0,0,1)
+#' P = Hpolytope(A = A, b = b)
+#'  
+#' @name Hpolytope-class
+#' @rdname Hpolytope-class
+#' @exportClass Hpolytope
+Hpolytope <- setClass (
+  # Class name
+  "Hpolytope",
+
+  # Defining slot type
+  representation (
+    A = "matrix",
+    b = "numeric",
+    volume = "numeric",
+    type = "character"
+    ),
+
+  # Initializing slots
+  prototype = list(
+    A = as.matrix(0),
+    b = as.numeric(NULL),
+    volume = as.numeric(NaN),
+    type = "Hpolytope"
+  )
+)

R-proj/R/RcppExports.R

@@ -46,9 +46,9 @@ copula <- function(r1, r2 = NULL, sigma = NULL, m = NULL, n = NULL, seed = NULL)
 #' \item{\code{type} }{ A string that declares the type of the body for the exact sampling: a) \code{'unit_simplex'} for the unit simplex, b) \code{'canonical_simplex'} for the canonical simplex, c) \code{'hypersphere'} for the boundary of a hypersphere centered at the origin, d) \code{'ball'} for the interior of a hypersphere centered at the origin.}
 #' \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{seed} }{ A fixed seed for the number generator.}
 #' }
 #' @param n The number of points that the function is going to sample.
-#' @param seed Optional. A fixed seed for the number generator.
 #'
 #' @references \cite{R.Y. Rubinstein and B. Melamed,
 #' \dQuote{Modern simulation and modeling} \emph{ Wiley Series in Probability and Statistics,} 1998.}
@@ -60,8 +60,8 @@ copula <- function(r1, r2 = NULL, sigma = NULL, m = NULL, n = NULL, seed = NULL)
 #' # 100 uniform points from the 2-d unit ball
 #' points = direct_sampling(n = 100, body = list("type" = "ball", "dimension" = 2))
 #' @export
-direct_sampling <- function(body, n, seed = NULL) {
-    .Call(`_volesti_direct_sampling`, body, n, seed)
+direct_sampling <- function(body, n) {
+    .Call(`_volesti_direct_sampling`, body, n)
 }
 
 #' Compute the exact volume of (a) a zonotope (b) an arbitrary simplex in V-representation or (c) if the volume is known and declared by the input object.
@@ -83,7 +83,7 @@ direct_sampling <- function(body, n, seed = NULL) {
 #'
 #' \donttest{# compute the exact volume of a 2-d arbitrary simplex
 #' V = matrix(c(2,3,-1,7,0,0),ncol = 2, nrow = 3, byrow = TRUE)
-#' P = Vpolytope$new(V)
+#' P = Vpolytope(V = V)
 #' vol = exact_vol(P)
 #' }
 #'
@@ -142,6 +142,21 @@ inner_ball <- function(P) {
     .Call(`_volesti_inner_ball`, P)
 }
 
+#' Read a SDPA format file
+#'
+#' @param input_file Name of the input file
+#'
+#' @return A list with two named items: an item "matrices" which is a list of the matrices and an vector "objFunction"
+#'
+#' @examples
+#' path = system.file('extdata', package = 'volesti')
+#' l = load_sdpa_format_file(paste0(path,'/sdpa_n2m3.txt'))
+#' 
+#' @export
+load_sdpa_format_file <- function(input_file = NULL) {
+    .Call(`_volesti_load_sdpa_format_file`, input_file)
+}
+
 #' An internal Rccp function as a polytope generator
 #'
 #' @param kind_gen An integer to declare the type of the polytope.
@@ -174,13 +189,14 @@ rotating <- function(P, T = NULL, seed = NULL) {
 #' Internal rcpp function for the rounding of a convex polytope
 #'
 #' @param P A convex polytope (H- or V-representation or zonotope).
+#' @param settings A list to set the random walk and its walk length
 #' @param seed Optional. A fixed seed for the number generator.
 #'
 #' @keywords internal
 #'
 #' @return A numerical matrix that describes the rounded polytope, a numerical matrix of the inverse linear transofmation that is applied on the input polytope, the numerical vector the the input polytope is shifted and the determinant of the matrix of the linear transformation that is applied on the input polytope.
-rounding <- function(P, seed = NULL) {
-    .Call(`_volesti_rounding`, P, seed)
+rounding <- function(P, settings = NULL, seed = NULL) {
+    .Call(`_volesti_rounding`, P, settings, seed)
 }
 
 #' Sample uniformly or normally distributed points from a convex Polytope (H-polytope, V-polytope, zonotope or intersection of two V-polytopes).
@@ -192,19 +208,19 @@ rounding <- function(P, seed = NULL) {
 #' @param random_walk Optional. A list that declares the random walk and some related parameters as follows:
 #' \itemize{
 #' \item{\code{walk} }{ A string to declare the random walk: i) \code{'CDHR'} for Coordinate Directions Hit-and-Run, ii) \code{'RDHR'} for Random Directions Hit-and-Run, iii) \code{'BaW'} for Ball Walk, iv) \code{'BiW'} for Billiard walk, v) \code{'BCDHR'} boundary sampling by keeping the extreme points of CDHR or vi) \code{'BRDHR'} boundary sampling by keeping the extreme points of RDHR. The default walk is \code{'BiW'} for the uniform distribution or \code{'CDHR'} for the Gaussian distribution.}
-#' \item{\code{walk_length} }{ The number of the steps per generated point for the random walk. The default value is \eqn{5} for \code{'BiW'} and \eqn{\lfloor 10 + d/10\rfloor} otherwise.}
+#' \item{\code{walk_length} }{ The number of the steps per generated point for the random walk. The default value is 1.}
 #' \item{\code{nburns} }{ The number of points to burn before start sampling.}
 #' \item{\code{starting_point} }{ A \eqn{d}-dimensional numerical vector that declares a starting point in the interior of the polytope for the random walk. The default choice is the center of the ball as that one computed by the function \code{inner_ball()}.}
 #' \item{\code{BaW_rad} }{ The radius for the ball walk.}
 #' \item{\code{L} }{ The maximum length of the billiard trajectory.}
+#' \item{\code{seed} }{ A fixed seed for the number generator.}
 #' }
 #' @param distribution Optional. A list that declares the target density and some related parameters as follows:
 #' \itemize{
 #' \item{\code{density} }{ A string: (a) \code{'uniform'} for the uniform distribution or b) \code{'gaussian'} for the multidimensional spherical distribution. The default target distribution is uniform.}
 #' \item{\code{variance} }{ The variance of the multidimensional spherical gaussian. The default value is 1.}
 #'  \item{\code{mode} }{ A \eqn{d}-dimensional numerical vector that declares the mode of the Gaussian distribution. The default choice is the center of the as that one computed by the function \code{inner_ball()}.}
 #' }
-#' @param seed Optional. A fixed seed for the number generator.
 #'
 #' @return A \eqn{d\times n} matrix that contains, column-wise, the sampled points from the convex polytope P.
 #' @examples
@@ -215,52 +231,16 @@ rounding <- function(P, seed = NULL) {
 #' # gaussian distribution from the 2d unit simplex in H-representation with variance = 2
 #' A = matrix(c(-1,0,0,-1,1,1), ncol=2, nrow=3, byrow=TRUE)
 #' b = c(0,0,1)
-#' P = Hpolytope$new(A,b)
+#' P = Hpolytope(A = A, b = b)
 #' points = sample_points(P, n = 100, distribution = list("density" = "gaussian", "variance" = 2))
 #'
 #' # uniform points from the boundary of a 2-dimensional random H-polytope
 #' P = gen_rand_hpoly(2,20)
 #' points = sample_points(P, n = 100, random_walk = list("walk" = "BRDHR"))
 #'
 #' @export
-sample_points <- function(P, n, random_walk = NULL, distribution = NULL, seed = NULL) {
-    .Call(`_volesti_sample_points`, P, n, random_walk, distribution, seed)
-}
-
-#' Write a SDPA format file
-#'
-#' Outputs a spectrahedron (the matrices defining a linear matrix inequality) and a vector (the objective function)
-#' to a SDPA format file.
-#'
-#' @param spectrahedron A spectrahedron in n dimensions; must be an object of class Spectrahedron
-#' @param objectiveFunction A numerical vector of length n
-#' @param outputFile Name of the output file
-#'
-#' @examples
-#' A0 = matrix(c(-1,0,0,0,-2,1,0,1,-2), nrow=3, ncol=3, byrow = TRUE)
-#' A1 = matrix(c(-1,0,0,0,0,1,0,1,0), nrow=3, ncol=3, byrow = TRUE)
-#' A2 = matrix(c(0,0,-1,0,0,0,-1,0,0), nrow=3, ncol=3, byrow = TRUE)
-#' lmi = list(A0, A1, A2)
-#' S = Spectrahedron$new(lmi);
-#' objFunction = c(1,1)
-#' writeSdpaFormatFile(S, objFunction, "output.txt")
-#' @export
-writeSdpaFormatFile <- function(spectrahedron = NULL, objectiveFunction = NULL, outputFile = NULL) {
-    invisible(.Call(`_volesti_writeSdpaFormatFile`, spectrahedron, objectiveFunction, outputFile))
-}
-
-#' Read a SDPA format file
-#'
-#' @param inputFile Name of the input file
-#'
-#' @return A list with two named items: an item "matrices" which is a list of the matrices and an vector "objFunction"
-#'
-#' @examples
-#' path = system.file('extdata', package = 'volesti')
-#' l = loadSdpaFormatFile(paste0(path,'/sdpa_n2m3.txt'))
-#' @export
-loadSdpaFormatFile <- function(inputFile = NULL) {
-    .Call(`_volesti_loadSdpaFormatFile`, inputFile)
+sample_points <- function(P, n, random_walk = NULL, distribution = NULL) {
+    .Call(`_volesti_sample_points`, P, n, random_walk, distribution)
 }
 
 #' The main function for volume approximation of a convex Polytope (H-polytope, V-polytope, zonotope or intersection of two V-polytopes)
@@ -270,15 +250,15 @@ loadSdpaFormatFile <- function(inputFile = NULL) {
 #' @param P A convex polytope. It is an object from class a) Hpolytope or b) Vpolytope or c) Zonotope or d) VpolytopeIntersection.
 #' @param settings Optional. A list that declares which algorithm, random walk and values of parameters to use, as follows:
 #' \itemize{
-#' \item{\code{algorithm} }{ A string to set the algorithm to use: a) \code{'CB'} for CB algorithm, b) \code{'SoB'} for SOB algorithm or b) \code{'CG'} for CG algorithm. The defalut algorithm for H-polytopes is \code{'CB'} when \eqn{d\leq 200} and \code{'CG'} when \eqn{d>200}. For the other representations the default algorithm is \code{'CB'}.}
+#' \item{\code{algorithm} }{ A string to set the algorithm to use: a) \code{'CB'} for CB algorithm, b) \code{'SoB'} for SOB algorithm or b) \code{'CG'} for CG algorithm. The defalut algorithm is \code{'CB'}.}
 #' \item{\code{error} }{ A numeric value to set the upper bound for the approximation error. The default value is \eqn{1} for SOB algorithm and \eqn{0.1} otherwise.}
 #' \item{\code{random_walk} }{ 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, c) \code{'BaW'} for Ball Walk, or \code{'BiW'} for Billiard walk. For CB and SOB algorithms the default walk is \code{'CDHR'} for H-polytopes and \code{'BiW'} for the other representations. For CG algorithm the default walk is \code{'CDHR'} for H-polytopes and \code{'RDHR'} for the other representations.}
 #' \item{\code{walk_length} }{ An integer to set the number of the steps for the random walk. The default value is \eqn{\lfloor 10 + d/10\rfloor} for \code{'SOB'} and \eqn{1} otherwise.}
 #' \item{\code{win_len} }{ The length of the sliding window for CB or CG algorithm. The default value is \eqn{400+3d^2} for CB or \eqn{500+4d^2} for CG.}
 #' \item{\code{hpoly} }{ A boolean parameter to use H-polytopes in MMC of CB algorithm when the input polytope is a zonotope. The default value is \code{TRUE} when the order of the zonotope is \eqn{<5}, otherwise it is \code{FALSE}.}
+#' \item{\code{seed} }{ A fixed seed for the number generator.}
 #' }
-#' @param rounding Optional. A boolean parameter for rounding. The default value is \code{TRUE} for V-polytopes and \code{FALSE} otherwise.
-#' @param seed Optional. A fixed seed for the number generator.
+#' @param rounding A boolean parameter for rounding. The default value is \code{FALSE}.
 #'
 #' @references \cite{I.Z.Emiris and V. Fisikopoulos,
 #' \dQuote{Practical polytope volume approximation,} \emph{ACM Trans. Math. Soft.,} 2018.},
@@ -303,8 +283,47 @@ loadSdpaFormatFile <- function(inputFile = NULL) {
 #' vol = volume(Z, settings = list("random_walk" = "RDHR", "walk_length" = 2))
 #'
 #' @export
-volume <- function(P, settings = NULL, rounding = NULL, seed = NULL) {
-    .Call(`_volesti_volume`, P, settings, rounding, seed)
+volume <- function(P, settings = NULL, rounding = FALSE) {
+    .Call(`_volesti_volume`, P, settings, rounding)
+}
+
+#' Write a SDPA format file
+#'
+#' Outputs a spectrahedron (the matrices defining a linear matrix inequality) and a vector (the objective function)
+#' to a SDPA format file.
+#'
+#' @param spectrahedron A spectrahedron in n dimensions; must be an object of class Spectrahedron
+#' @param objective_function A numerical vector of length n
+#' @param output_file Name of the output file
+#'
+#' @examples
+#' \dontrun{
+#' A0 = matrix(c(-1,0,0,0,-2,1,0,1,-2), nrow=3, ncol=3, byrow = TRUE)
+#' A1 = matrix(c(-1,0,0,0,0,1,0,1,0), nrow=3, ncol=3, byrow = TRUE)
+#' A2 = matrix(c(0,0,-1,0,0,0,-1,0,0), nrow=3, ncol=3, byrow = TRUE)
+#' lmi = list(A0, A1, A2)
+#' S = Spectrahedron(matrices = lmi);
+#' objFunction = c(1,1)
+#' write_sdpa_format_file(S, objFunction, "output.txt")
+#' }
+#' @export
+write_sdpa_format_file <- function(spectrahedron, objective_function, output_file) {
+    invisible(.Call(`_volesti_write_sdpa_format_file`, spectrahedron, objective_function, output_file))
+}
+
+#' Read a SDPA format file
+#'
+#' @param inputFile Name of the input file
+#'
+#' @return A list with two named items: an item "matrices" which is a list of the matrices and an vector "objFunction"
+#'
+#' @examples
+#' path = system.file('extdata', package = 'volesti')
+#' l = loadSdpaFormatFile(paste0(path,'/sdpa_n2m3.txt'))
+#' 
+#' @export
+loadSdpaFormatFile <- function(inputFile = NULL) {
+    .Call(`_volesti_loadSdpaFormatFile`, inputFile)
 }
 
 #' An internal Rccp function for the over-approximation of a zonotope

R-proj/R/SpectrahedronClass.R

@@ -0,0 +1,33 @@
+#' An R class to represent a Spectrahedron
+#'
+#' A spectrahedron is a convex body defined by a linear matrix inequality of the form \eqn{A_0 + x_1 A_1 + ... + x_n A_n \preceq 0}.
+#' The matrices \eqn{A_i} are symmetric \eqn{m \times m} real matrices and \eqn{\preceq 0} denoted negative semidefiniteness.
+#'
+#' \describe{
+#'    \item{matrices}{A List that contains the matrices \eqn{A_0, A_1, ..., A_n}.}
+#' }
+#' 
+#' @examples
+#' A0 = matrix(c(-1,0,0,0,-2,1,0,1,-2), nrow=3, ncol=3, byrow = TRUE)
+#' A1 = matrix(c(-1,0,0,0,0,1,0,1,0), nrow=3, ncol=3, byrow = TRUE)
+#' A2 = matrix(c(0,0,-1,0,0,0,-1,0,0), nrow=3, ncol=3, byrow = TRUE)
+#' lmi = list(A0, A1, A2)
+#' S = Spectrahedron(matrices = lmi);
+#'  
+#' @name Spectrahedron-class
+#' @rdname Spectrahedron-class
+#' @exportClass Spectrahedron
+Spectrahedron <- setClass (
+  # Class name
+  "Spectrahedron",
+
+  # Defining slot type
+  representation (
+    matrices = "list"
+  ),
+
+  # Initializing slots
+  prototype = list(
+    matrices = list()
+  )
+)

R-proj/R/VpolytopeClass.R

@@ -0,0 +1,37 @@
+#' An R class to represent a V-polytope
+#'
+#' A V-polytope is a convex polytope defined by the set of its vertices.
+#'
+#' \describe{
+#'    \item{V}{An \eqn{m\times d} numerical matrix that contains the vertices row-wise.}
+#' 
+#'    \item{volume}{The volume of the polytope if it is known, \eqn{NaN} otherwise by default.}
+#'    
+#'    \item{type}{A character with default value 'Vpolytope', to declare the representation of the polytope.}
+#' }
+#'  
+#' @examples
+#' V = matrix(c(2,3,-1,7,0,0),ncol = 2, nrow = 3, byrow = TRUE)
+#' P = Vpolytope(V = V)
+#' 
+#' @name Vpolytope-class
+#' @rdname Vpolytope-class
+#' @exportClass Vpolytope
+Vpolytope <- setClass (
+  # Class name
+  "Vpolytope",
+
+  # Defining slot type
+  representation (
+    V = "matrix",
+    volume = "numeric",
+    type = "character"
+  ),
+
+  # Initializing slots
+  prototype = list(
+    V = as.matrix(0),
+    volume = as.numeric(NaN),
+    type = "Vpolytope"
+  )
+)