geomscale/volesti

Practical volume computation and sampling in high dimensions

volesti's Issues

HPolytope prints -0 when initialised using ine file

Describe the bug
Currently, if you initialise an hpolytope using ine file with 0 entries in it, when printed those entries are printed as -0.

To Reproduce
Steps to reproduce the behavior:

1. Initialize a polytope by reading from the ine file used in the docs.
2. Print the hpolytope read in previous step.
3. See that all 0 entries are displayed as -0.

Expected behavior
All 0 entries should be printed as 0 instead of -0.

vol aborting during computation of 3D prism volume

vol aborts with an error when I try to compute the volume of a prism in R^3 (see .ine file below):

$ ./vol -v -f1 prism.ine 
Verbose mode
Reading input from file...
Input polytope: 3
 8 4 float
6.94419e-310 8.5428e-72 1.63042e-322 <= 0.00850976
6.94419e-310 2.21338e-52 1.63042e-322 <= 0.00850976
6.94419e-310 0 2.42092e-322 <= 0.0134911
6.94419e-310 0 0 <= 0.0134911
6.94419e-310 0 1.63042e-322 <= 0.0153591
6.94419e-310 0 1.63042e-322 <= 0.0153591
6.94419e-310 0 1.63042e-322 <= 1
6.94419e-310 0 4.00193e-322 <= 1
Inner ball time: 0.000223
Inner ball center is: 
0 0 0 
radius is: 1e+30
Experiment 1 
Generate the first random point in P

Compute 1318 random points in P
First random points construction time = 0.001338

Furthest distance from Chebychev point= inf

Constructing the sequence of balls
vol: /home/arrigo/volume_approximation-master/test/../include/volume/volume.h:174: NT volume(Polytope&, Parameters&, std::pair<Point, NT>) [with Polytope = HPolytope<point<Cartesian<double> > >; Parameters = vars<double, boost::random::mersenne_twister_engine<unsigned int, 32ul, 624ul, 397ul, 31ul, 2567483615u, 11ul, 4294967295u, 7ul, 2636928640u, 15ul, 4022730752u, 18ul, 1812433253u> >; Point = point<Cartesian<double> >; NT = double]: Assertion `!balls.empty()' failed.
Aborted (core dumped)

prism.ine:

H-representation
begin
8 4 real
0.008509755085097572	0.10187124277134102	0.08149699421707308	-0.0
0.008509755085097572	-0.10187124277134102	-0.08149699421707308	-0.0
0.013491075134910584	-0.05780634428435466	0.15375492457251516	-0.0
0.013491075134910584	0.05780634428435466	-0.15375492457251516	-0.0
0.015359070153590594	-0.15967758705569568	0.07225793035544206	-0.0
0.015359070153590596	0.15967758705569568	-0.07225793035544209	-0.0
1.0	-0.0	-0.0	-1.0
1.0	-0.0	-0.0	1.0
end

Automatic differentiation support

It would be a very nice feature to include automatic differentiation support within the functors, so that the gradients can be computed from the automatic differentiation software.

Automatic differentiation references and libraries:

1. https://mc-stan.org/users/interfaces/math.html
2. https://github.com/autodiff/autodiff
3. https://www.boost.org/doc/libs/master/libs/math/doc/html/math_toolkit/autodiff.html

Error in 3D prism volume computation?

I am trying to compute the volume of a simple 3D prism using vol (the .ine file is below). I know that the volume of the prism is approximately equal to 0.0118 (area of the polygonal base X height=1). I have run vol from the latest development branch with different values of the error tolerance:

$ ./vol --file1 prism.ine -e 0.1
Reading input from file...
Experiment 1 
3 8 1 -1 0.1 [-0.9,-1.1] 131833 10 0.02023386 [0.02023386,0.02023386] 0 1.020234 0 0.419022 0 
$ ./vol --file1 prism.ine -e 0.2
Reading input from file...
Experiment 1 
3 8 1 -1 0.2 [-0.8,-1.2] 32958 10 0.02025999 [0.02025999,0.02025999] 0 1.02026 0 0.129329 0 
$ ./vol --file1 prism.ine -e 0.01
Reading input from file...
Experiment 1 
3 8 1 -1 0.01 [-0.99,-1.01] 13183347 10 0.02038337 [0.02038337,0.02038337] 0 1.020383 0 44.884 0 

and from these runs it seems that the volume computed by vol might be off by a factor 2. Or is vol computing something different, like a signed volume?

H-representation
begin
8 4 real
0.008509755085097572 0.10187124277134102 0.08149699421707308 -0.0
0.008509755085097572 -0.10187124277134102 -0.08149699421707308 -0.0
0.013491075134910584 -0.05780634428435466 0.15375492457251516 -0.0
0.013491075134910584 0.05780634428435466 -0.15375492457251516 -0.0
0.015359070153590594 -0.15967758705569568 0.07225793035544206 -0.0
0.015359070153590596 0.15967758705569568 -0.07225793035544209 -0.0
0.5 -0.0 -0.0 -1.0
0.5 -0.0 -0.0 1.0
end

Inefficient handling of burn in phase in sampling functions

In header include/sampling.hpp the sampling functions do not handle the burn in phase efficiently.
They store the points of burn in phase and then they delete them.

RandomPointGenerator::apply(P, p, c, a, eta, nburns, walk_len, randPoints, push_back_policy, rng); randPoints.clear();

We need a new function to perform the burn in phase efficiently.

Remove boost external libraries

Is your feature request related to a problem? Please describe.
External libraries consume space in the repository and have to be updated manually.

Describe the solution you'd like
Remove boost from external folder and download them in cmake before building the package or whenever needed.

See how this was done for eigen #131

R code for Sampling: volesti vs hitandrun -- S4 classes

Working through the example "Sampling: volesti vs hitandrun"

> d<-10
> P = gen_cube(d, 'H')
> tim = system.time({ sample_points(P, n=N, random_walk = list("walk" = "RDHR", "walk_length"=5)) })
> times_volesti = c(times_volesti, as.numeric(tim)[3])
> constr <- list(constr = P$A, dir=rep('<=', 2*d), rhs=P$b)
Error in P$A : $ operator not defined for this S4 class

P is an S4 object

> isS4(P) 
[1] TRUE

so the code should have "@" rather than "$", like this

N <-  1000
times_volesti  <-  c()
times_hnr <-  c()
for (d in seq(from=10,to=100,by=10)) {
   P <-  gen_cube(d, 'H')
  tim <- system.time({ sample_points(P, n=N, random_walk = list("walk" = "RDHR", "walk_length"=5)) })
  times_volesti <- c(times_volesti, as.numeric(tim)[3])
  
  constr <- list(constr = P@A, dir=rep('<=', 2*d), rhs=P@b)
  tim <- system.time({ hitandrun(constr, n.samples=N, thin=5) })
  times_hnr <-  c(times_hnr, as.numeric(tim)[3])
}

Bye

> sessionInfo()
R version 4.0.2 (2020-06-22)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Linux Mint 19

Matrix products: default
BLAS/LAPACK: /usr/lib/x86_64-linux-gnu/libopenblasp-r0.2.20.so

locale:
 [1] LC_CTYPE=en_AU.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_AU.UTF-8        LC_COLLATE=en_AU.UTF-8    
 [5] LC_MONETARY=en_AU.UTF-8    LC_MESSAGES=en_AU.UTF-8   
 [7] LC_PAPER=en_AU.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_AU.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] ggplot2_3.3.2   hitandrun_0.5-5 volesti_1.1.2   Rcpp_1.0.6     

loaded via a namespace (and not attached):
 [1] withr_2.3.0      dplyr_1.0.2      crayon_1.3.4     grid_4.0.2      
 [5] R6_2.4.1         lifecycle_0.2.0  gtable_0.3.0     magrittr_1.5    
 [9] scales_1.1.1     pillar_1.4.6     rlang_0.4.8      generics_0.0.2  
[13] vctrs_0.3.4      ellipsis_0.3.1   tools_4.0.2      glue_1.4.2      
[17] purrr_0.3.4      munsell_0.5.0    compiler_4.0.2   pkgconfig_2.0.3 
[21] colorspace_1.4-1 tidyselect_1.1.0 tibble_3.0.3    
> 

Where are the samples stored?

I am trying to use vol under Linux to generate uniform random samples from a polytope but I can't figure where the samples are stored. If I run for example

./vol -v -f1 cube10.ine -rand -nsample 1000 

I get

Verbose mode
Reading input from file...
Input polytope: 10
 20 11 float
1 0 0 0 0 0 0 0 0 0 <= 1
0 1 0 0 0 0 0 0 0 0 <= 1
0 0 1 0 0 0 0 0 0 0 <= 1
0 0 0 1 0 0 0 0 0 0 <= 1
0 0 0 0 1 0 0 0 0 0 <= 1
0 0 0 0 0 1 0 0 0 0 <= 1
0 0 0 0 0 0 1 0 0 0 <= 1
0 0 0 0 0 0 0 1 0 0 <= 1
0 0 0 0 0 0 0 0 1 0 <= 1
0 0 0 0 0 0 0 0 0 1 <= 1
-1 0 0 0 0 0 0 0 0 0 <= 1
0 -1 0 0 0 0 0 0 0 0 <= 1
0 0 -1 0 0 0 0 0 0 0 <= 1
0 0 0 -1 0 0 0 0 0 0 <= 1
0 0 0 0 -1 0 0 0 0 0 <= 1
0 0 0 0 0 -1 0 0 0 0 <= 1
0 0 0 0 0 0 -1 0 0 0 <= 1
0 0 0 0 0 0 0 -1 0 0 <= 1
0 0 0 0 0 0 0 0 -1 0 <= 1
0 0 0 0 0 0 0 0 0 -1 <= 1
Generate random points only
Inner ball time: 0.000803
Inner ball center is: 
0 0 0 0 0 0 0 0 0 0 
radius is: 1
Sampling time: 0.001194

Code coverage for volesti

Is your feature request related to a problem? Please describe.
It would be nice to have code coverage for volesti.

Describe the solution you'd like
We could use codecov https://about.codecov.io/

Khachiyan algorithm for Minimum Volume Enclosing Ellipsoids

Is your feature request related to a problem? Please describe.
Reimplement or find a more efficient implementation for Khachiyan algorithm for Minimum Volume Enclosing Ellipsoids.
See https://people.orie.cornell.edu/miketodd/TYKhach.pdf as a reference.
Currently we use this code https://github.com/GeomScale/volume_approximation/tree/develop/external/minimum_ellipsoid
which depends on ublas. Even changing this code to use Eigen or lapack could result in more efficient implementation.
Currently, the implementation is very efficient for d<200.

A CGAL function is also available: https://doc.cgal.org/latest/Bounding_volumes/classCGAL_1_1Approximate__min__ellipsoid__d.html and used by volesti in the past.

@TolisChal feel free to add more information and/or benchmarks if available.

Replace Eigen pointers with iterators

Replace eigen pointers in hpolytope.hpp and in the implementations of largest inscribed ball and elipsoid methods.

Integral Collocation method is numerically unstable

Describe the bug

The integral collocation method implementation of this paper is numerically unstable.
We use the Chebyshev transform class provided by Boost to compute the fixed point operator (Lagrange interpolation at Chebyshev nodes).

The implementation of the ODE solver is located at include/ode_solvers/integral_collocation.hpp.

To Reproduce
Steps to reproduce the behavior:

Simulate the ode x''(t) = -x(t) with x(0) = 0 and x'(0) = 1.
The numerical solver returns a sine function with increasing amplitude.
I believe this has to be with the parameterization of the algorithm.

Expected behavior

The solution is x(t) = sin(t) which is bounded above by 1.

Integrate Codecov's Github app with Geomscale

Is your feature request related to a problem? Please describe.
Codecov has been enabled on the repo by #188. The aim is that whenever a new PR is made, Codecov should report any changes to the coverage as a comment to that PR.

Describe the solution you'd like
This can be done by integrating Codecov's Github app with GeomScale org. I think this can be done by simply installing the Codecov app from here: https://github.com/apps/codecov. More detailed instructions can be followed from https://docs.codecov.com/docs/how-to-create-a-github-app-for-codecov-enterprise

Change representation of multiple domains in ODE solvers interface

Instead of defining std::vector<Polytope*> to represent a list which has elements either a Polytope
or NULL to represent R^d.

We want to change this structure to std::vector<Polytope> and the representation of R^d should be
handled inside the corresponding polytope classes.

Python Interface Installation failed

Describe the bug
When I try to install python interface of the repo, I got this error, shown in the following picture. There is not the boosl/random.hpp file in the repo. I followed the installation guide in the MD file.

Desktop (please complete the following information):
ubuntu 18.04 python 3.7 & 3.8

Hope you can check to see if there are any missing files or code errors.

Fix examples and CI integration

Is your feature request related to a problem? Please describe.
Even after the useful PR #190 some examples are still failing to compile.

Describe the solution you'd like
Some have an easy fix e.g. vpolytope-volume while some other e.g. EnvelopeProblemSOS need some work to fix the cmake (at least in my Ubuntu 20.04 with gcc version 9.3.0 I am getting an error Target "EnvelopeProblem" links to target "Python2::Python" but the target was not found after running cmake .).

Additional context
The ultimate solution would integrate all examples to be build in our CI (i.e. github actions) then we will be sure that future changes will not break them.

Adding License and Contribution Section in README

• Adding License and Contribution Section in README.
• Instead of simple text adding Hyperlink.
  

macOS-latest (devel) check failure

This is general problem as described here: r-lib/actions#140

Remove lpsolve external libraries

Is your feature request related to a problem? Please describe.
External libraries consume space in the repository and have to be updated manually.

Describe the solution you'd like
Remove lpsolve from external folder and download them in cmake before building the package or whenever needed.

See how this was done for eigen #131

Add examples

Is your feature request related to a problem? Please describe.
It would be very useful to add examples of all algorithms such as volume computation of H-polytopes, V-polytopes, Z-polytopes, sampling from those and rounding. Now there are only for spectahedra see https://github.com/GeomScale/volume_approximation/tree/develop/examples

Describe the solution you'd like
An example should be ideally simple, clean and descriptive.

Volume of intersection of two V-polytopes crashes

The following crashes with the R-package version 1.0.2

V1 = matrix(c(2,3,-1,7,0,0),ncol = 2, nrow = 3, byrow = TRUE)
 V2 = matrix(c(1,5,-6,8,2,5),ncol = 2, nrow = 3, byrow = TRUE)
 VP = IntVP$new(V1,V2)
 volume(VP)

No file: convex_bodies/ellipsoids.h

When I tried to install the volesti package using the "build source package" in "Build" of Rstudio, it failed due to the following reason:
In the file "copula.cpp" there is a line: #include "convex_bodies/ellipsoids.h"
However, there is no file named "convex_bodies/ellipsoids.h" .

Hoping that the file "convex_bodies/ellipsoids.h" can be added.

R interface

Fix / improve the followings:

1. storing Polytope objects in one session and loading them in another leads to
   awkward and unexpected errors:
   Error in sprintf("C++ object <%s> of class '%s' <%s>", externalptr_address(pointer),'a0 : 'a0 external pointer is not valid

2. Drop exportPattern("^[[:alpha:]]+") from NAMESPACE file

3. Add flag drop = FALSE to code like A = Mat[,-c(1)]

4. Change input arguments in R functions. Instead of setting the arguments equal to NULL, check internally
   for NULL and replacing with the default.

5. Remove file_to_polytope R function

Volesti 1.1.4 installation error

Hello everyone, I am trying to upgrade Volesti from 1.1.3 to 1.1.4 in R. The operation system is win10, with R 4.1.0 and Rtools properly installed.
The installation instruction in the following link can successfully install volesti 1.1.3, howerver, it cannot install volesti 1.1.4. https://github.com/GeomScale/volesti/blob/2e6572d9ba757de311f2a55b09abf8f9c02e932e/doc/r_interface.md

I tried to install volesti 1.1.4 with instruction in the link, the error information is :
"In file included from ../../include/convex_bodies/hpolytope.h:21,
from ../../include/volume/volume_sequence_of_balls.hpp:21,
from direct_sampling.cpp:19:
../../include/lp_oracles/solve_lp.h:32:10: fatal error: lp_lib.h: No such file or directory
#include "lp_lib.h"
^~~~~~~~~~
compilation terminated."

I think the error is related to the fact that lp_solve is not included in volesti directory root/include/lp_oracles. I understand that lp_solve is designed to be downloaded in the complie process. However, it did not work in this case.

Could you please update the installation instruction applicable to v 1.1.4?
Or is it possible to update the CRAN version to v1.1.4? The link is https://cran.r-project.org/web/packages/volesti/index.html.
I think another possible solution is to provide an additional version with lp_solve included locally, like v1.1.3.

Thanks a lot!

Reach 70% coverage

Since #188 volesti has automatic coverage reports with the help of CircleCI.

However, after merging #188 the coverage of the code is rather low (<60%).

In this issue we could aim for 70% coverage by adding simple and fast unit tests in the test suite.

Documentation

Is your feature request related to a problem? Please describe.
Improve documentation by using some specialized software for this case, e.g. https://www.mkdocs.org/. Other suggestions?

Implement NUTS

Implement NUTS sampling from this paper.
Ideally, we want the NUTS criterion to be modular, i.e. be combined with any HMC-related sampler (e.g. via changing
the ODE solver).

Setting error parameter e to 0

I unintentionally set the value e in volume_example.cpp in the develop branch to 0 (missing out it is an int and setting it to a value smaller 1). This led to consumption of my complete memory and ultimately crashed the code.

Implement boundary oracle for integral collocation method

The Integral Collocation method of this paper is implemented in volesti under
include/ode_solvers/integral_collocation.hpp.

The current state of the solver does not include a boundary oracle (and thus it does not currently support truncation).

The goal of this feature request is to implement a boundary oracle for the method.
More specifically the boundary oracle method has to solve the intersection of the trajectory
(computed via Lagrange interpolation at the Chebyshev nodes) which is defined by its
Chebyshev transform (via using the Chebyshev transform class provided in boost) with
the boundary of the H-polytope.

To solve it effectively, I propose converting the Chebyshev transform to a complex polynomial
of twice the degree such that the real part of its roots represent solutions of the equation.

So, all in all, we need to

1. Compute the polynomial of twice the degree (already implemented)
2. Compute the integral of this polynomial (a[k] x^k -> a[k] / (k + 1) x^{k + 1})
3. Project on each of the facet normals and solve the corresponding polynomial equation numerically (e.g. with MPSolve)
4. Keep the smallest positive of the real parts of the complex roots

Missing flag in a command

Describe the bug
If we execute the command ./generate -cube -h -d 10 from the README, we get the error The number of facets has to be >= d+1.

To Reproduce
Steps to reproduce the behavior:

1. Go to the test folder.
2. Compile with cmake and make.
3. Execute the command ./generate -cube -h -d 10.
4. The error The number of facets has to be >= d+1 appears.

Expected behavior
We expect to generate a 10-dimensional hypercube in H-representation.

Solution
A solution would be either to suggest in the README to run the command with the flag -m, either to set the variable m for the facets to be at least d+1.
For example command ./generate -cube -h -d 10 -m 20 is working fine !

Some headers don't compile independently

Tried executing a basic try code Try1.cpp in vol_test/examples/Sampling_New:

   #include<stdio.h>
   #include "random_walks.hpp"
   using namespace std;
   int main()
   {
        cout<<"Ran"<<endl;
        return 0;
   } 

Returned the following error (sqrt,pow not a member of std):

Resolved that by using #include<cmath> . However cmath ideally should not have been needed to be included.

Bug in vol_Ali

Bug in function vol_Ali (called by function frustum_of_simplex) in the arithmetic operations that adds an error in the computations which vanishes in high dimensions.

Remove macro #VOLESTIPY

The macro VOLESTIPY is used to exclude the header files that call the lpsolve library. See for example ./include/convex_bodies/hpolytope.h, where the lpsolve is used to compute the maximum inscribed ball of an H-polytope.

We should restructure the code such that VOLESTIPY macro is not needed any more.

Remove volestipy directory

If I am not mistaken the volestipy folder
can be now removed.

Put all .cmake files in one directory

Is your feature request related to a problem?
No, it's not related to any problem, it's about the convention that makes our code more structured.

Request Description
there are only 2 .cmake files Eigen.cmake and Boost.cmake, but in the future, there can be too many .cmake files. like.

• requires lp_solver.cmake for #141
• requires CodeCoverage.cmake #139
• requires ClangFormat.cmake for issue #116

and in the future, there can be many .cmake files in volesti.

Solution Description
we can put all these similar kinds of .cmake files in a directory a cmake (we can call any suitable name), and add its path to our root CMakeLists.txt or where ever it's needed. in this way, our code will be more structured.

CMake build error under Windows

When I try to build the solution using CMake under Windows 10 I get this error:

Configuring incomplete, errors occurred!

I have defined a new DLP_SOLVE CMake variable as described in the documentation.

Error in the .ext file representation

Describe the bug
In the Polytope input section of the README, the representation of a V-polytope in a .ext file
contains the line m d numbertype, while it should be m d+1 numbertype.

To Reproduce
Steps to reproduce the behavior:

1. Go to the test folder.
2. Compile with cmake and make.
3. Execute the command ./generate -simplex -v -d 2.
4. Open the file simplex_2.ext. The corresponding line has the value three (d+1) instead of two (d).

max_inscribed_ellipsoid takes too much time

Describe the bug
I am trying to (approximately) compute the number of linear extensions of a poset by computing the volume of it's corresponding order polytope.
My poset instance has 100 dimensions and has the following adjacency matrix: link to the file

The instance is taken from this repo, when run with the ARMC algorithm given in the same repo, it takes around 1.6 secs to compute the approximate answer (3.5e+129).
I tried using the volume algorithms implemented in volesti to compare the results.

• First attempt
  • Rounding: None
  • Algorithm: Volume Cooling Balls
  • Result: 3.2e+129
  • Time taken: 15 secs
• Next attempt: with rounding
  • Rounding: Computed the max inscribed ellipsoid to round the order polytope
  • Algorithm: Volume Cooling Balls
  • Result: 3.15e+129
  • Time taken: 17 mins 30 secs

As, can be seen from the results above, both attempts gave pretty accurate results but the time taken is drastically different ?
What am I doing wrong ?

To Reproduce

• Code snipper for first attempt (without rounding)

#include "Eigen/Eigen"
#include <vector>
#include "cartesian_geom/cartesian_kernel.h"
#include "convex_bodies/hpolytope.h"
#include "convex_bodies/ellipsoid.h"
#include "preprocess/max_inscribed_ellipsoid.hpp"
#include "random_walks/uniform_billiard_walk.hpp"
#include "random_walks/uniform_accelerated_billiard_walk.hpp"
#include "generators/boost_random_number_generator.hpp"
#include "volume/volume_cooling_balls.hpp"

#include <iostream>
#include <fstream>
#include "misc.h"
#include "poset.h"

typedef double NT;
typedef Cartesian <NT> Kernel;
typedef typename Kernel::Point Point;
typedef BoostRandomNumberGenerator<boost::mt19937, NT> RNGType;
typedef HPolytope<Point> HPOLYTOPE;
typedef Ellipsoid<Point> EllipsoidType;
typedef typename HPOLYTOPE::MT MT;
typedef typename HPOLYTOPE::VT VT;

NT calculateLinearExtension(HPOLYTOPE& HP) {
    // Setup parameters for calculating volume and rounding
    unsigned int d = HP.dimension();
    unsigned int walk_len = 1;
    NT e = 0.1;
    NT round_val = 1.0;

    NT volume = volume_cooling_balls<BilliardWalk, RNGType>(HP, e, walk_len).second;
    volume = volume * round_val;

    // multiplying by d factorial, d = number of elements
    for(NT i=(NT)d; i>1; i-=1) {
        volume = volume * i;
    }

    return volume;
}


/**

 Usage: ./cle INSTANCE

 example:
    ./cle instances/bipartite_0.5_008_0.txt

*/
int main(int argc, char* argv[]) {
    std::string filename (argv[1]);
    std::ifstream in;
    in.open(filename);
    std::pair<bool, Poset> read_poset = read_poset_from_file_adj_matrix(in);

    if( !read_poset.first ) {
        std::cerr << "error in reading data from file" << std::endl;
        return -1;
    }

    // ---------------- Create HPOLYTOPE ---------------------
    Poset poset = read_poset.second;
    int num_relations = poset.num_relations();
    int d = poset.num_elem();

    MT A = Eigen::MatrixXd::Zero(2*d + num_relations, d);
    VT b = Eigen::MatrixXd::Zero(2*d + num_relations, 1);

    // first add (ai >= 0) or (-ai <= 0) rows
    A.topLeftCorner(d, d) = -Eigen::MatrixXd::Identity(d, d);

    // next add (ai <= 1) rows
    A.block(d, 0, d, d) = Eigen::MatrixXd::Identity(d, d);
    b.block(d, 0, d, 1) = Eigen::MatrixXd::Constant(d, 1, 1.0);

    // next add the relations
    for(int idx=0; idx<num_relations; ++idx) {
        std::pair<int, int> edge  = poset.get_relation(idx);
        A(2*d + idx, edge.first)  = 1;
        A(2*d + idx, edge.second) = -1;
    }
    HPOLYTOPE HP(d, A, b);
    //--------------------------------------------------------

    std::cout << calculateLinearExtension(HP) << std::endl;
    return 0;
}

• Code snippet for second attempt (with rounding)

#include "Eigen/Eigen"
#include <vector>
#include "cartesian_geom/cartesian_kernel.h"
#include "convex_bodies/hpolytope.h"
#include "convex_bodies/ellipsoid.h"
#include "preprocess/max_inscribed_ellipsoid.hpp"
#include "random_walks/uniform_billiard_walk.hpp"
#include "random_walks/uniform_accelerated_billiard_walk.hpp"
#include "generators/boost_random_number_generator.hpp"
#include "volume/volume_cooling_balls.hpp"

#include <iostream>
#include <fstream>
#include "misc.h"
#include "poset.h"

typedef double NT;
typedef Cartesian <NT> Kernel;
typedef typename Kernel::Point Point;
typedef BoostRandomNumberGenerator<boost::mt19937, NT> RNGType;
typedef HPolytope<Point> HPOLYTOPE;
typedef Ellipsoid<Point> EllipsoidType;
typedef typename HPOLYTOPE::MT MT;
typedef typename HPOLYTOPE::VT VT;

NT calculateLinearExtension(HPOLYTOPE& HP) {
    // Setup parameters for calculating volume and rounding
    unsigned int d = HP.dimension();
    unsigned int walk_len = 1;
    NT e = 0.1;
    NT round_val = 1.0;

    // ------------ START: max inscribed ellipsoid rounding ----------------------------------------------
    unsigned int max_iter = 150;
    NT tol = std::pow(10, -6.0), reg = std::pow(10, -4.0);

    std::pair<Point, NT> innerBall = HP.ComputeInnerBall();
    if (innerBall.second < 0.0)
        throw std::runtime_error("no inscribed ball");

    VT x0 = innerBall.first.getCoefficients();
    Eigen::SparseMatrix<NT> A = HP.get_mat().sparseView();
    std::pair<std::pair<MT, VT>, bool> inscribed_ellipsoid_res = max_inscribed_ellipsoid<MT>(A,
                                                                                             HP.get_vec(),
                                                                                             x0,
                                                                                             max_iter,
                                                                                             tol,
                                                                                             reg);

    if (!inscribed_ellipsoid_res.second) // not converged
        throw std::runtime_error("no inscribed ellipsoid");

    MT A_ell = inscribed_ellipsoid_res.first.first.inverse();
    EllipsoidType inscribed_ellipsoid( A_ell );
    MT L = inscribed_ellipsoid.Lcov();

    HP.shift(inscribed_ellipsoid_res.first.second);
    round_val = L.transpose().determinant();
    HP.linear_transformIt(L);
    // ------------ END: max inscribed ellipsoid rounding ------------------------------------------------

    NT volume = volume_cooling_balls<BilliardWalk, RNGType>(HP, e, walk_len).second;
    volume = volume * round_val;

    // multiplying by d factorial, d = number of elements
    for(NT i=(NT)d; i>1; i-=1) {
        volume = volume * i;
    }

    return volume;
}


/**

 Usage: ./cle INSTANCE

 example:
    ./cle instances/bipartite_0.5_008_0.txt

*/
int main(int argc, char* argv[]) {
    std::string filename (argv[1]);
    std::ifstream in;
    in.open(filename);
    std::pair<bool, Poset> read_poset = read_poset_from_file_adj_matrix(in);

    if( !read_poset.first ) {
        std::cerr << "error in reading data from file" << std::endl;
        return -1;
    }

    // ---------------- Create HPOLYTOPE ---------------------
    Poset poset = read_poset.second;
    int num_relations = poset.num_relations();
    int d = poset.num_elem();

    MT A = Eigen::MatrixXd::Zero(2*d + num_relations, d);
    VT b = Eigen::MatrixXd::Zero(2*d + num_relations, 1);

    // first add (ai >= 0) or (-ai <= 0) rows
    A.topLeftCorner(d, d) = -Eigen::MatrixXd::Identity(d, d);

    // next add (ai <= 1) rows
    A.block(d, 0, d, d) = Eigen::MatrixXd::Identity(d, d);
    b.block(d, 0, d, 1) = Eigen::MatrixXd::Constant(d, 1, 1.0);

    // next add the relations
    for(int idx=0; idx<num_relations; ++idx) {
        std::pair<int, int> edge  = poset.get_relation(idx);
        A(2*d + idx, edge.first)  = 1;
        A(2*d + idx, edge.second) = -1;
    }
    HPOLYTOPE HP(d, A, b);
    //--------------------------------------------------------

    std::cout << calculateLinearExtension(HP) << std::endl;
    return 0;
}

Expected behavior
I expected that both the attempts will have comparable time and possibly rounding will improve the efficiency.

Desktop (please complete the following information):

• OS: Ubuntu, macOS (tested on both)

Remove eigen external libraries

Is your feature request related to a problem? Please describe.
External libraries consume space in the repository and have to be updated manually.

Describe the solution you'd like
Remove them from external folder and download them in cmake before building the package or whenever needed. This should be straightforward for the cases of boost and eigen. For lpsolve removing from external could be more complicated.

EDIT: to make the procedure more clear we stick with eigen in this issue and we will open specified issues for the rest of the libraries

Code formating

Is your feature request related to a problem? Please describe.
The code contains many styles since we have a few contributors now. So a code formatting tool is needed.

Describe the solution you'd like
Clang-format could be used with some well known style, e.g. google format.

Improve documentation

I open this to post a few issues related to the improvement of the R interface documentation.

Generators: explain exactly what you generate, e.g. GenSkinnyCube doesn't say anything about the generated cube.

Rouding:
-- Fix examples in round_polytope (they say "rotate")
-- typos in description
-- references to min ellipsoid method (implementation + paper)

CMakeLists not creating requested executable. Default executable ( Makefile) fails to run

Created a folder called Sampling within vol-test/Examples to test small sampling codes. Modified an appropriate CMakelist file

It does not produce the executable as requested , but simply generates a "makefile". Neither does it generate the sub directory as requested. It does however generate the subdirectories that has been commented out.

The Makefile generated fails to compile, raising an issue with the mkl header, which has been downloaded and installed

Moving to a subdirectory "hpolytope-volume" , it contains an object "hpolytopeVolume" that executes successfully with make

Fatal error: 'ellipsoids.h' file not found

Describe the bug
Problem in installing the R-interface.

To Reproduce
Steps to reproduce the behavior:
I followed the given steps:
-> Install package-dependencies: Rcpp, RcppEigen, BH, lpSolveAPI
-> change working directory to "volume_approximation/R-Proj" using setwd
-> Ran these commands

Rcpp::compileAttributes()  
library(devtools)  
devtools::build()  
devtools::install()  
library(volesti)  

Desktop:

• macOS Catalina

Compiling Example(hpolytope-volume) c++ interface failed

I'm testing the c++ interface of hpolytope-volume.

Under the folder examples/hpolytope-volume, I have tried two types of cmake,
When I run：

cmake -DLP_SOLVE=/usr/lib/lp_solve/liblpsolve55.so

succeed, but after that, when I run

make

Error appeared:

hpolytopeVolume.cpp:10:10: fatal error: Eigen/Eigen: No such file or directory
   10 | #include "Eigen/Eigen"
      |          ^~~~~~~~~~~~~

So I have tried another cmake:

cmake .. -DEXAMPLES=hpolytope-volume

another error repoted:

-- The C compiler identification is GNU 9.3.0
-- The CXX compiler identification is GNU 9.3.0
-- Check for working C compiler: /usr/bin/cc
-- Check for working C compiler: /usr/bin/cc -- works
-- Detecting C compiler ABI info
-- Detecting C compiler ABI info - done
-- Detecting C compile features
-- Detecting C compile features - done
-- Check for working CXX compiler: /usr/bin/c++
-- Check for working CXX compiler: /usr/bin/c++ -- works
-- Detecting CXX compiler ABI info
-- Detecting CXX compiler ABI info - done
-- Detecting CXX compile features
-- Detecting CXX compile features - done
-- Eigen library not found locally, downloading it.
-- Using downloaded Eigen library at:volesti/examples/hpolytope-volume/_deps/eigen-src
-- Boost Library found: /usr/include
-- lp_solve library not found locally, downloading it.
-- Using downloaded lp_solve at: volesti/examples/hpolytope-volume/_deps/lpsolve-src/src
-- Library lp_solve found: /usr/lib/lp_solve/liblpsolve55.so
-- Library lp_solve found: /usr/lib/lp_solve/liblpsolve55.so BLAS-NOTFOUND
-- Using downloaded Eigen library at:volesti/examples/hpolytope-volume/_deps/eigen-src
-- Boost Library found: /usr/include
-- Library lp_solve found: /usr/lib/lp_solve/liblpsolve55.so
-- Found OpenMP_C: -fopenmp (found version "4.5")
-- Found OpenMP_CXX: -fopenmp (found version "4.5")
-- Found OpenMP: TRUE (found version "4.5")
CMake Error: The following variables are used in this project, but they are set to NOTFOUND.
Please set them or make sure they are set and tested correctly in the CMake files:
BLAS
    linked by target "high_dimensional_hmc" in directory volesti/examples/logconcave
    linked by target "high_dimensional_hmc_rand_poly" in directory volesti/examples/logconcave
    linked by target "simple_sde" in directory volesti/examples/logconcave
    linked by target "simple_uld" in directory volesti/examples/logconcave
    linked by target "simple_ode" in directory volesti/examples/logconcave
    linked by target "simple_hmc" in directory volesti/examples/logconcave

-- Configuring incomplete, errors occurred!
See also "volesti/examples/hpolytope-volume/CMakeFiles/CMakeOutput.log".

Could anyone tell me how to solve this problem or compile and test the c++ interface under other folders in volesti?
Thank you.

Instructions for installing volestipy in windows platform

I am trying to install the python package volestipy on a computer with win-10 rather than linux. I hope that the instructions for the command line installation of volestipy can be kindly provided. I think there maybe many others like me with a win-10 laptop and eager to explore volesti using python. Thank you!

Redundant C++ tests

In folder test there are a few C++ tests for the previous code structure and interface, e.g. vol.ccp, that can be safely removed. This will improve the performance of the CI suites.

Sample Points cannot use InnerPoint parameter

Describe the bug
When i use the sample_points function and pass the paramater InnerPoint i get the following
error:
Error in sample_points(P, N = 1, distribution = "gaussian", InnerPoint = point) :
std::bad_alloc

To Reproduce
Code:

library(volesti)
A = rbind(c(1,0), c(0,1),
c(-1,0), c(0,-1))
b = c(1,1,0,0)
P = Hpolytope$new(A, b)
point = c(0.2, 0.2)
x = sample_points(P, distribution="gaussian", InnerPoint=point)

Environment
Ubuntu 18.04 LTS
Rstudio, compiler version: 3.4.4
volesti version: 1.0.3

Add support for Hit and Run for general densities

Describe the solution you'd like

Adapt the existing functionality of volesti for hit-and-run walk (both randomized and coordinate directions versions) to support arbitrary target densities, e.g. log-concave.

The current support include uniform and Gaussian densities.

Volesti v1.1.4 installation failure (R interface)

Hello everyone, I am trying to upgrade Volesti from 1.1.3 to 1.1.4 in R. The operation system is win10, with R 4.1.0 and Rtools properly installed.
The installation instruction in the following link can successfully install volesti 1.1.3, howerver, it cannot install volesti 1.1.4 successfully. https://github.com/GeomScale/volesti/blob/2e6572d9ba757de311f2a55b09abf8f9c02e932e/doc/r_interface.md

I tried to install volesti 1.1.4 with instruction in the link, the error information is :
"In file included from ../../include/convex_bodies/hpolytope.h:21,
from ../../include/volume/volume_sequence_of_balls.hpp:21,
from direct_sampling.cpp:19:
../../include/lp_oracles/solve_lp.h:32:10: fatal error: lp_lib.h: No such file or directory
#include "lp_lib.h"
^~~~~~~~~~
compilation terminated."

I think the error is related to the fact that lp_solve is not included in volesti directory root/include/lp_oracles. I understand that lp_solve is designed to be downloaded in the complie process. However, it did not work in this case.

Could you please update the installation instruction applicable to v 1.1.4?
Or is it possible to update the CRAN version to v1.1.4? The link is https://cran.r-project.org/web/packages/volesti/index.html.
I think another possible solution is to provide an additional version with lp_solve included locally, like v1.1.3.

Thanks a lot!

How to compute the intersection of V-polytopes

Is there a function that does the intersection of V-polytopes, like it is stated for the class IntVP?