Add more CPLEX paths

CMakeLists.txt

@@ -7,8 +7,6 @@ include_directories(${CPLEX_INCLUDE_DIR})
 
 find_package(GMP REQUIRED)
 find_package(OpenMP REQUIRED)
-find_package(BLAS REQUIRED)
-find_package(LAPACKE REQUIRED)
 
 set(CMAKE_CXX_STANDARD 11)
 set(CMAKE_CXX_FLAGS  "${CMAKE_CXX_FLAGS} -Werror ${OpenMP_CXX_FLAGS} -O3")

cmake/FindLAPACKE.cmake

@@ -1,11 +0,0 @@
-find_library(LAPACKE_LIBRARIES
-        NAMES lapacke)
-
-find_path(LAPACKE_INCLUDE_DIR NAMES lapacke.h)
-
-include(FindPackageHandleStandardArgs)
-find_package_handle_standard_args(LAPACKE
-        DEFAULT_MSG
-        LAPACKE_LIBRARIES
-        LAPACKE_INCLUDE_DIR)
-mark_as_advanced(LAPACKE_LIBRARIES LAPACKE_INCLUDE_DIR)

infinity/library/src/infinity-nd.c

@@ -22,7 +22,6 @@
 #include <multirow/util.h>
 
 #include <infinity/infinity-nd.h>
-#include <multirow/linalg.h>
 
 static long lp_count = 0;
 static double lp_time = 0;
@@ -42,41 +41,6 @@ static double sfree_mip_time = 0;
 static long scale_ahull_lp_count = 0;
 static double scale_ahull_lp_time = 0;
 
-static int cone_bound_find_lambda(const struct RayList *rays,
-                                  const double *f,
-                                  const double *x,
-                                  double *lambda)
-{
-    int rval = 0;
-    double *A = 0;
-    double *b = 0;
-    int dim = rays->dim;
-    int nrays = rays->nrays;
-
-    A = (double*) malloc(nrays * dim * sizeof(double));
-    b = (double*) malloc(dim * sizeof(double));
-    abort_if(!A, "could not allocate A");
-    abort_if(!b, "could not allocate b");
-
-    for (int i = 0; i < dim; i++)
-    {
-        b[i] = x[i] - f[i];
-        for (int j = 0; j < nrays; j++)
-            A[i * nrays + j] = rays->values[j * dim + i];
-    }
-
-    rval = LINALG_solve(nrays, dim, A, b, lambda);
-    abort_if(rval, "LINALG_solve failed");
-
-    for (int i = 0; i < dim; i++)
-        abort_iff(lambda[0] < 0, "lambda[i] is negative (%lf)", i, lambda[i]);
-
-CLEANUP:
-    if(A) free(A);
-    if(b) free(b);
-    return rval;
-}
-
 static int create_sfree_mip(int nrows,
                             int nrays,
                             const double *f,
@@ -401,6 +365,10 @@ static int create_tight_rays_lp(int nrows,
     rval = LP_relax(lp);
     abort_if(rval, "LP_relax failed");
 
+    //rval = LP_write(lp, "tight-rays.lp");
+    //abort_if(rval, "LP_write failed");
+
+
 CLEANUP:
     if(rmatind) free(rmatind);
     if(rmatval) free(rmatval);
@@ -462,8 +430,9 @@ static int create_violated_cone_lp(int nrows,
     rval = LP_relax(lp);
     abort_if(rval, "LP_relax failed");
 
-//    rval = LP_write(lp, "violated-cone.lp");
-//    abort_if(rval, "LP_write failed");
+    //rval = LP_write(lp, "violated-cone.lp");
+    //abort_if(rval, "LP_write failed");
+    //UTIL_pause();
 
 CLEANUP:
     if(rmatind) free(rmatind);
@@ -633,7 +602,7 @@ static int find_interior_point_enum(const int nrows,
                 }
             }
 
-    if(best_value < 0.999) *found = 1;
+    if(!DOUBLE_geq(best_value, 1)) *found = 1;
 
 CLEANUP:
     if(beta2) free(beta2);
@@ -694,9 +663,9 @@ static int find_interior_point_cplex(const int nrows,
     rval = LP_get_obj_val(&lp, &objval);
     abort_if(rval, "LP_get_obj_val failed");
 
-    log_debug("        obj: %.8lf\n", objval);
+    log_debug("    obj = %.8lf\n", objval);
 
-    if(objval > 0.999)
+    if(DOUBLE_geq(objval, 1.0))
     {
         log_debug("  set is lattice-free\n");
         *found = 0;
@@ -888,7 +857,7 @@ static int find_violated_cone(int nrows,
                               double *sbar,
                               int *violated_found)
 {
-    log_debug("    finding violated cone:\n");
+    log_verbose("      find_violated_cone\n");
     int rval = 0;
 
     struct LP lp;
@@ -925,30 +894,30 @@ static int find_violated_cone(int nrows,
     rval = LP_get_obj_val(&lp, &obj);
     abort_if(rval, "LP_get_obj_val failed");
 
-    log_debug("      obj: %.8lf\n", obj);
+    log_verbose("  o=%.8lf\n", obj);
 
-    if(obj > 0.999)
+    if(DOUBLE_geq(obj, 0.999))
     {
         *violated_found = 0;
-        log_debug("    no violated cone found\n");
     }
     else
     {
         *violated_found = 1;
 
-        log_debug("    violated cone found:\n");
-        log_debug("      f: %12.8lf %12.8lf\n", f[0], f[1]);
-        log_debug("      x: %12.8lf %12.8lf\n", x[0], x[1]);
+        log_verbose("Violated cone found\n");
+        log_verbose("  f=%.8lf %.8lf\n", f[0], f[1]);
+        log_verbose("  x=%.8lf %.8lf\n", x[0], x[1]);
 
         for(int i = 0; i < nrays; i++)
         {
             rx[i] = (sbar[i] > 1e-9);
 
-            if(rx[i]) if_debug_level
+            if(rx[i]) if_verbose_level
             {
-                const double m = fmin(epsilon, beta[i]);
+                double m = min(epsilon, beta[i]);
                 const double *r = &rays[nrows * i];
-                time_printf("      r[%d]: %12.8lf %12.8lf\n", i, m * sbar[i] * r[0], m * sbar[i] * r[1]);
+                time_printf("  r[%d]=%.8lf %.8lf\n", i, r[0], r[1]);
+                time_printf("  r[%d]=%.8lf %.8lf\n", i, m * r[0], m * r[1]);
             }
         }
     }
@@ -1009,12 +978,15 @@ static int bound(int nrows,
                 fbar[i] -= min(*epsilon, beta[j]) * r[i] * sbar[j];
         }
 
+        log_verbose("%.12lf %.12lf\n", f[0], f[1]);
+        log_verbose("%.12lf %.12lf\n", fbar[0], fbar[1]);
+
         prev_epsilon = *epsilon;
 
         rval = cone_bound(nrows, nrays, fbar, rays, rx, x, beta, epsilon);
         abort_if(rval, "cone_bound failed");
 
-        log_debug("    cone bound: %.20lf\n", epsilon);
+        log_verbose("        e=%.12lf\n", *epsilon);
         abort_if(prev_epsilon < *epsilon, "epsilon should never increase");
     }
 
@@ -1243,6 +1215,8 @@ CLEANUP:
 int INFINITY_ND_generate_lfree(const struct MultiRowModel *model,
                                struct ConvLFreeSet *lfree)
 {
+    log_debug("INFINITY_ND_generate_lfree\n");
+
     int rval = 0;
     int nrows = model->nrows;
     int nrays = model->rays.nrays;
@@ -1275,6 +1249,24 @@ int INFINITY_ND_generate_lfree(const struct MultiRowModel *model,
 
     int it = 0;
 
+    //lp_time = 0;
+    //lp_count = 0;
+
+    //epsilon_lp_count = 0;
+    //epsilon_lp_time = 0;
+    //
+    //sfree_mip_count = 0;
+    //sfree_mip_time = 0;
+
+    //tight_lp_count = 0;
+    //tight_lp_time = 0;
+
+    //violated_lp_count = 0;
+    //violated_lp_time = 0;
+
+    //scale_ahull_lp_time = 0;
+    //scale_ahull_lp_count = 0;
+
     long x_count = 0;
     double epsilon;
 
@@ -1297,7 +1289,7 @@ int INFINITY_ND_generate_lfree(const struct MultiRowModel *model,
 
         while(1)
         {
-            log_debug("    epsilon: %.8lf\n", epsilon);
+            log_debug("  epsilon = %.8lf\n", epsilon);
 
             int found = 0;
 
@@ -1317,14 +1309,8 @@ int INFINITY_ND_generate_lfree(const struct MultiRowModel *model,
                 if(!found) break;
             }
 
-            if_debug_level {
-                time_printf("    found interior point:\n");
-                time_printf("        x: ");
-                for (int i = 0; i < nrows; i++)
-                    printf("%.2lf ", x[i]);
-                printf("\n");
-
-            }
+            log_debug("    found interior point:\n");
+            for(int i = 0; i < nrows; i++) log_debug("        %.2lf\n", x[i]);
 
             x_count++;
             abort_if(x_count > 1000, "infinite loop");
@@ -1334,8 +1320,10 @@ int INFINITY_ND_generate_lfree(const struct MultiRowModel *model,
                     &epsilon_x, tx);
             abort_if(rval, "bound failed");
 
-            abort_if(isinf(epsilon_x), "epsilon_x is infinite");
-            log_debug("    epsilon_x: %.8lf\n", epsilon_x);
+            if(isinf(epsilon_x)) break;
+
+//            epsilon_x = (floor(epsilon_x * 128) / 128);
+            log_debug("    epsilon_x = %.8lf\n", epsilon_x);
 
             if(DOUBLE_eq(epsilon_x, epsilon))
             {
@@ -1355,12 +1343,11 @@ int INFINITY_ND_generate_lfree(const struct MultiRowModel *model,
 
         int skip_ahull = 0;
 
-        log_debug("    updating beta:\n");
         for(int i = 0; i < nrays; i++)
         {
             if(t[i])
             {
-                beta[i] = min(beta[i] * 0.999, epsilon);
+                beta[i] = min(beta[i], epsilon);
             }
             else if(!skip_ahull)
             {
@@ -1377,11 +1364,14 @@ int INFINITY_ND_generate_lfree(const struct MultiRowModel *model,
                     continue;
                 }
 
+//                alpha = (floor(alpha * 128) / 128);
                 beta[i] = min(beta[i], alpha);
             }
 
-            log_debug("        beta[%2d]: %20.12lf\n", i, beta[i]);
+            log_debug("  beta[%2d] = %20.12lf\n", i, beta[i]);
         }
+
+        log_debug("epsilon = %.6lf\n", epsilon);
     }
 
     log_debug("    %6ld lattice points, %ld iterations\n", x_count, it);

infinity/library/src/infinity.c

@@ -108,7 +108,6 @@ static int create_cut_from_lfree(const struct Tableau *tableau,
                                 struct Row *cut)
 {
     int rval = 0;
-    double *ray = 0;
 
     struct LP lp;
     int nvars = map->nvars;
@@ -123,41 +122,40 @@ static int create_cut_from_lfree(const struct Tableau *tableau,
     rval = INFINITY_create_psi_lp(lfree, &lp);
     abort_if(rval, "create_psi_lp failed");
 
-    ray = (double*) malloc(nrows * sizeof(double));
-    abort_if(!ray, "could not allocate ray");
-
     cut->nz = nvars;
     for(int i = 0; i < nvars; i++)
     {
-        double value, norm = 0;
-        double *original_ray = LFREE_get_ray(rays, map->variable_to_ray[i]);
+        double value;
+        const double *q = LFREE_get_ray(rays, map->variable_to_ray[i]);
         char type = tableau->column_types[map->indices[i]];
 
-        for (int j = 0; j < nrows; j++)
-        {
-            ray[j] = original_ray[j];
-            norm += fabs(ray[j]);
-        }
-
-        if (norm < 0.001)
-            for (int j = 0; j < nrows; j++)
-                ray[j] *= 0.001 / norm;
 
         if(ENABLE_LIFTING && type == MILP_INTEGER)
         {
-            rval = INFINITY_pi(nrows, ray, map->ray_scale[i], &lp, &value);
+            rval = INFINITY_pi(nrows, q, map->ray_scale[i], &lp, &value);
             abort_if(rval, "INFINITY_pi failed");
         }
         else
         {
-            rval = INFINITY_psi(nrows, ray, map->ray_scale[i], &lp, &value);
+            rval = INFINITY_psi(nrows, q, map->ray_scale[i], &lp, &value);
             abort_if(rval, "INFINITY_psi failed");
         }
 
-        value *= 1.001;
-        value = fmax(value, 0.001);
         log_verbose("   psi[%4d] = %20.12lf %d\n", map->indices[i], value);
 
+        if_debug_level
+        {
+            time_printf("    q=[");
+            for(int j = 0; j < nrows; j++)
+                printf("%25.20lf ", q[j] * map->ray_scale[i]);
+            printf("] value=%25.20lf\n", value);
+        }
+
+        value = fmax(value, 1 / 1024.0);
+
+//        value *= 1.001;
+//        value = DOUBLE_max(value, 0.001);
+
         cut->indices[i] = map->indices[i];
         cut->pi[i] = -value;
     }
@@ -165,7 +163,6 @@ static int create_cut_from_lfree(const struct Tableau *tableau,
     cut->pi_zero = -1.0;
 
 CLEANUP:
-    if(ray) free(ray);
     LP_free(&lp);
     return rval;
 }
@@ -189,7 +186,12 @@ static int filter_model(const struct MultiRowModel *original_model,
     struct RayList *filtered_rays = &filtered_model->rays;
     const struct RayList *original_rays = &original_model->rays;
 
-    memcpy(f, original_model->f, nrows * sizeof(double));
+    memcpy(f, original_model->f, 2 * sizeof(double));
+    for(int i = 0; i < nrows; i++)
+    {
+        f[i] = (ceil(f[i] * 128) / 128);
+        if(f[i] <= 0.01) f[i] = 0;
+    }
 
     r = (double*) malloc(nrows * sizeof(double));
     abort_if(!r, "could not allocate r");
@@ -266,18 +268,11 @@ static int append_extra_rays(struct MultiRowModel *model)
 
     for(int i = 0; i < nrows; i++)
     {
-        int found, index;
-        double scale;
-
         for(int j = 0; j < nrows; j++) r[j] = (i == j ? e : 0);
-        rval = CG_find_ray(&model->rays, r, &found, &scale, &index);
-        abort_if(rval, "CG_find_ray failed");
-        if(!found) LFREE_push_ray(&model->rays, r);
+        LFREE_push_ray(&model->rays, r);
 
         for(int j = 0; j < nrows; j++) r[j] = (i == j ? -e : 0);
-        rval = CG_find_ray(&model->rays, r, &found, &scale, &index);
-        abort_if(rval, "CG_find_ray failed");
-        if(!found) LFREE_push_ray(&model->rays, r);
+        LFREE_push_ray(&model->rays, r);
     }
 
     if(nrows == 2)
@@ -414,9 +409,6 @@ int INFINITY_generate_cut(const struct Tableau *tableau, struct Row *cut)
     rval = CG_init_map(&original_map, max_nrays, tableau->nrows);
     abort_if(rval, "CG_init_map failed");
 
-    rval = LFREE_init_conv(&lfree, tableau->nrows, max_nrays);
-    abort_if(rval, "LFREE_init_conv failed");
-
     rval = extract_models(tableau, &original_model, &filtered_model,
             &original_map);
     abort_if(rval, "extract_models failed");
@@ -427,18 +419,41 @@ int INFINITY_generate_cut(const struct Tableau *tableau, struct Row *cut)
         goto CLEANUP;
     }
 
+    rval = LFREE_init_conv(&lfree, tableau->nrows, max_nrays);
+    abort_if(rval, "LFREE_init_conv failed");
+
+//    if(tableau->nrows == 2)
+//        rval = INFINITY_2D_generate_lfree(&filtered_model, &lfree);
+//    else
+//        rval = INFINITY_ND_generate_lfree(&filtered_model, &lfree);
+//
+//    if(rval)
+//    {
+//        rval = ERR_NO_CUT;
+//        goto CLEANUP;
+//    }
+//
+//    for(int i = 0; i < filtered_model.rays.nrays; i++)
+//    {
+//        double *r = LFREE_get_ray(&filtered_model.rays, i);
+//        for(int j = 0; j < tableau->nrows; j++)
+//            r[j] *= lfree.beta[j];
+//    }
+
     rval = append_extra_rays(&filtered_model);
     abort_if(rval, "append_extra_rays failed");
 
-    if_debug_level
+//    if(tableau->nrows == 2)
+//        rval = INFINITY_2D_generate_lfree(&filtered_model, &lfree);
+//    else
+        rval = INFINITY_ND_generate_lfree(&filtered_model, &lfree);
+
+    if(rval)
     {
-        rval = CG_print_model(&filtered_model);
-        abort_if(rval, "CG_print_model failed");
+        rval = ERR_NO_CUT;
+        goto CLEANUP;
     }
 
-    rval = INFINITY_ND_generate_lfree(&filtered_model, &lfree);
-    abort_if(rval, "INFINITY_ND_generate_lfree failed");
-
     if(SHOULD_DUMP_CUTS)
     {
         rval = dump_cut(&lfree);

infinity/library/tests/infinity-nd-test.cpp

@@ -27,7 +27,6 @@ extern "C" {
 #include "../src/infinity-nd.c"
 }
 
-const double E = 1e-6;
 
 TEST(InfinityNDTest, find_violated_cone_test)
 {
@@ -165,11 +164,11 @@ TEST(InfinityNDTest, cone_bound_test_1)
 
     rval = cone_bound(2, 6, f, rays, rx1, x, beta, &epsilon);
     abort_if(rval, "cone_bound failed");
-    EXPECT_NEAR(0.5, epsilon, E);
+    EXPECT_NEAR(0.5, epsilon, 1e-6);
 
     rval = cone_bound(2, 6, f, rays, rx2, x, beta, &epsilon);
     abort_if(rval, "cone_bound failed");
-    EXPECT_NEAR(1.0, epsilon, E);
+    EXPECT_NEAR(1.0, epsilon, 1e-6);
 
     CLEANUP:
     if(rval) FAIL();
@@ -198,7 +197,7 @@ TEST(InfinityNDTest, cone_bound_test_2)
 
     rval = cone_bound(2, 2, f, rays, rx, x1, beta1, &epsilon);
     abort_if(rval, "cone_bound failed");
-    EXPECT_NEAR(1.0, epsilon, E);
+    EXPECT_NEAR(1.0, epsilon, 1e-6);
 
     rval = cone_bound(2, 2, f, rays, rx, x1, beta2, &epsilon);
     abort_if(rval, "cone_bound failed");
@@ -206,7 +205,7 @@ TEST(InfinityNDTest, cone_bound_test_2)
 
     rval = cone_bound(2, 2, f, rays, rx, x2, beta2, &epsilon);
     abort_if(rval, "cone_bound failed");
-    EXPECT_NEAR(1.0, epsilon, E);
+    EXPECT_NEAR(1.0, epsilon, 1e-6);
 
     rval = cone_bound(2, 2, f, rays, rx, x2, beta3, &epsilon);
     abort_if(rval, "cone_bound failed");
@@ -240,7 +239,7 @@ TEST(InfinityNDTest, bound_test_1)
 
     rval = bound(2, 6, f, rays, x, beta1, &epsilon, tx);
     abort_if(rval, "bound failed");
-    EXPECT_NEAR(epsilon, 0.5, E);
+    EXPECT_NEAR(epsilon, 0.5, 1e-6);
     EXPECT_TRUE(tx[0]);
     EXPECT_FALSE(tx[1]);
     EXPECT_FALSE(tx[2]);
@@ -250,7 +249,7 @@ TEST(InfinityNDTest, bound_test_1)
 
     rval = bound(2, 6, f, rays, x, beta2, &epsilon, tx);
     abort_if(rval, "bound failed");
-    EXPECT_NEAR(epsilon, 1.0, E);
+    EXPECT_NEAR(epsilon, 1.0, 1e-6);
     EXPECT_TRUE(tx[0]);
     EXPECT_FALSE(tx[1]);
     EXPECT_FALSE(tx[2]);
@@ -309,11 +308,11 @@ TEST(InfinityNDTest, psi_test)
 
     rval = INFINITY_psi(2, q1, 1.0, &lp, &value);
     abort_if(rval, "GREDDY_ND_psi failed");
-    EXPECT_NEAR(value, 2.0, E);
+    EXPECT_NEAR(value, 2.0, 1e-6);
 
     rval = INFINITY_psi(2, q2, 2.0, &lp, &value);
     abort_if(rval, "GREDDY_ND_psi failed");
-    EXPECT_NEAR(value, 8.0, E);
+    EXPECT_NEAR(value, 8.0, 1e-6);
 
 CLEANUP:
     LP_free(&lp);
@@ -361,64 +360,17 @@ TEST(InfinityNDTest, psi_test_2)
 
     rval = INFINITY_psi(3, q1, 1.0, &lp, &value);
     abort_if(rval, "GREDDY_ND_psi failed");
-    EXPECT_NEAR(value, 1.0, E);
+    EXPECT_NEAR(value, 1.0, 1e-6);
 
     rval = INFINITY_psi(3, q2, 1.0, &lp, &value);
     abort_if(rval, "GREDDY_ND_psi failed");
-    EXPECT_NEAR(value, 2.0, E);
+    EXPECT_NEAR(value, 2.0, 1e-6);
 
 CLEANUP:
     if(rval) FAIL();
 }
 
-TEST(InfinityNDTest, psi_test_3)
-{
-    int rval = 0;
-
-    double f[] = { 0.671875, 0.671875 };
-    double rays[] = {
-            -0.007812500000,  0.000000000000,
-            -0.039062500000,  0.046875000000,
-             0.000000000000,  0.046875000000,
-             0.046875000000,  0.000000000000,
-             0.000000000000, -0.039062500000
-    };
-
-    double beta[] = {
-            66.909090909091,
-            29.440000000000,
-            14.000000000000,
-            14.000000000000,
-            29.440000000000,
-    };
-
-    double q[] = { 0 - f[0], 1 - f[1]};
-
-    struct ConvLFreeSet lfree;
-    lfree.f = f;
-    lfree.beta = beta;
-    lfree.rays.nrays = 5;
-    lfree.rays.values = rays;
-    lfree.nrows = lfree.rays.dim = 2;
-
-    double value;
-    struct LP lp;
-
-    rval = LP_open(&lp);
-    abort_if(rval, "LP_open failed");
-
-    rval = INFINITY_create_psi_lp(&lfree, &lp);
-    abort_if(rval, "INFINITY_create_psi_lp failed");
-
-    rval = INFINITY_psi(lfree.nrows, q, 1.0, &lp, &value);
-    abort_if(rval, "GREDDY_ND_psi failed");
-    EXPECT_NEAR(value, 1.0, E);
-
-CLEANUP:
-    if(rval) FAIL();
-}
-
-TEST(DISABLED_InfinityNDTest, generate_cut_test_1)
+TEST(InfinityNDTest, generate_cut_test_1)
 {
     int rval = 0;
 
@@ -446,12 +398,12 @@ TEST(DISABLED_InfinityNDTest, generate_cut_test_1)
     rval = INFINITY_ND_generate_lfree(&model, &lfree);
     abort_if(rval, "INFINITY_ND_generate_lfree failed");
 
-    EXPECT_NEAR(lfree.beta[0], 0.5, E);
-    EXPECT_NEAR(lfree.beta[1], 0.5, E);
-    EXPECT_NEAR(lfree.beta[2], 0.5, E);
-    EXPECT_NEAR(lfree.beta[3], 0.5, E);
-    EXPECT_NEAR(lfree.beta[4], 1.0, E);
-    EXPECT_NEAR(lfree.beta[5], 1.0, E);
+    EXPECT_NEAR(lfree.beta[0], 0.5, 1e-6);
+    EXPECT_NEAR(lfree.beta[1], 0.5, 1e-6);
+    EXPECT_NEAR(lfree.beta[2], 0.5, 1e-6);
+    EXPECT_NEAR(lfree.beta[3], 0.5, 1e-6);
+    EXPECT_NEAR(lfree.beta[4], 1.0, 1e-6);
+    EXPECT_NEAR(lfree.beta[5], 1.0, 1e-6);
 
 CLEANUP:
     if(rval) FAIL();
@@ -486,12 +438,12 @@ TEST(InfinityNDTest, generate_cut_test_2)
     rval = INFINITY_ND_generate_lfree(&model, &lfree);
     abort_if(rval, "INFINITY_ND_generate_lfree failed");
 
-    EXPECT_NEAR(lfree.beta[0], 0.75, E);
-    EXPECT_NEAR(lfree.beta[1], 2.25, E);
-    EXPECT_NEAR(lfree.beta[2], 0.75, E);
-    EXPECT_NEAR(lfree.beta[3], 2.25, E);
-    EXPECT_NEAR(lfree.beta[4], 0.75, E);
-    EXPECT_NEAR(lfree.beta[5], 2.25, E);
+    EXPECT_NEAR(lfree.beta[0], 0.75, 1e-6);
+    EXPECT_NEAR(lfree.beta[1], 2.25, 1e-6);
+    EXPECT_NEAR(lfree.beta[2], 0.75, 1e-6);
+    EXPECT_NEAR(lfree.beta[3], 2.25, 1e-6);
+    EXPECT_NEAR(lfree.beta[4], 0.75, 1e-6);
+    EXPECT_NEAR(lfree.beta[5], 2.25, 1e-6);
 
 CLEANUP:
     CG_free_model(&model);
@@ -535,68 +487,3 @@ TEST(InfinityNDTest, scale_to_ahull_test)
 CLEANUP:
     if(rval) FAIL();
 }
-
-TEST(InfinityNDTest, cone_bound_find_lambda_test)
-{
-    int rval = 0;
-
-    double f[] = { 1/2.0, 1/3.0 };
-    double x[] = { 3.0, 3.0 };
-    double ray_values[] = {
-            1.0, 2.0,
-            3.0, 2.0,
-    };
-
-    struct RayList rays = {.values = ray_values, .nrays = 2, .dim = 2};
-    double lambda[2];
-
-    rval = cone_bound_find_lambda(&rays, f, x, lambda);
-    abort_if(rval, "cone_bound_find_lambda failed");
-
-    EXPECT_NEAR(lambda[0], 3/4.0, E);
-    EXPECT_NEAR(lambda[1], 7/12.0, E);
-
- CLEANUP:
-    if(rval) FAIL();
-}
-
-TEST(InfinityNDTest, cone_bound_find_lambda_test_2)
-{
-    int rval = 0;
-
-    double f[] = { 1/2.0, 1/3.0, 5/6.0 };
-    double x[] = { 3.0, 3.0, 6.0 };
-    double ray_values[] = {
-            1.0, 2.0, 3.0,
-            3.0, 2.0, 5.0,
-    };
-
-    struct RayList rays = {.values = ray_values, .nrays = 2, .dim = 3};
-    double lambda[2];
-
-    rval = cone_bound_find_lambda(&rays, f, x, lambda);
-    abort_if(rval, "cone_bound_find_lambda failed");
-
-    EXPECT_NEAR(lambda[0], 3/4.0, E);
-    EXPECT_NEAR(lambda[1], 7/12.0, E);
-
-CLEANUP:
-    if(rval) FAIL();
-}
-
-TEST(InfinityNDTest, cone_bound_find_lambda_test_3)
-{
-    int rval = 0;
-    double f[] = { 1/2.0, 1/2.0 };
-    double x[] = { 3.0, 3.0 };
-    double ray_values[] = {
-             0.0, -1.0,
-            -1.0,  0.0,
-    };
-
-    struct RayList rays = {.values = ray_values, .nrays = 2, .dim = 2};
-    double lambda[2];
-
-    rval = cone_bound_find_lambda(&rays, f, x, lambda);
-    EXPECT_NE(rval, 0);
-}

lifting/benchmark/run_experiments.sh

@@ -17,6 +17,9 @@ if [ ! -f $RUN ]; then
 fi
 
 INSTANCES="instances/filtered/all.txt"
+SAMPLES_SLOW=10
+SAMPLES_MEDIUM=100
+SAMPLES_FAST=1000
 SEED=1240
 
 # ORIGINAL
@@ -33,29 +36,29 @@ DIR=orig-100
 mkdir -p $DIR; rm -f $DIR/*log $DIR/*yaml
 
 title Heuristic
-$RUN $COMMON_OPTS --samples 100000 --heuristic --log $DIR/heur.log --stats $DIR/heur.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_FAST --heuristic --log $DIR/heur.log --stats $DIR/heur.yaml || exit
 
 title Bound Original
-$RUN $COMMON_OPTS --samples 1000 --bound --log $DIR/bound.log --stats $DIR/bound.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_FAST --bound --log $DIR/bound.log --stats $DIR/bound.yaml || exit
 
 title Bound Pre-processing
-$RUN $COMMON_OPTS --samples 1000 --bound --preprocess --log $DIR/bound-pre.log --stats $DIR/bound-pre.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_FAST --bound --preprocess --log $DIR/bound-pre.log --stats $DIR/bound-pre.yaml || exit
 
 title Naive Bounding-Box
-$RUN $COMMON_OPTS --samples 100 --naive --log $DIR/naive-bbox.log --stats $DIR/naive-bbox.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_MEDIUM --naive --log $DIR/naive-bbox.log --stats $DIR/naive-bbox.yaml || exit
 
 title Naive Fixed-M
 M=50
-$RUN $COMMON_OPTS --samples 100 --naive --fixed-bounds $M --log $DIR/naive-fixed-$M.log --stats $DIR/naive-fixed-$M.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_SLOW --naive --fixed-bounds $M --log $DIR/naive-fixed-$M.log --stats $DIR/naive-fixed-$M.yaml || exit
 
 title Naive Bounding-Box Pre-processing
-$RUN $COMMON_OPTS --samples 100 --naive --preprocess --log $DIR/naive-bbox-pre.log --stats $DIR/naive-bbox-pre.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_SLOW --naive --preprocess --log $DIR/naive-bbox-pre.log --stats $DIR/naive-bbox-pre.yaml || exit
 
 title MIP
-$RUN $COMMON_OPTS --samples 10 --mip --log $DIR/mip.log --stats $DIR/mip.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_SLOW --mip --log $DIR/mip.log --stats $DIR/mip.yaml || exit
 
 title MIP Pre-processing
-$RUN $COMMON_OPTS --samples 10 --mip --preprocess --log $DIR/mip-pre.log --stats $DIR/mip-pre.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_SLOW --mip --preprocess --log $DIR/mip-pre.log --stats $DIR/mip-pre.yaml || exit
 
 # SHEAR
 # ------------------------------------------------------------------------------
@@ -70,28 +73,21 @@ COMMON_OPTS="$COMMON_OPTS --check-answers $ANSWERS"
 DIR=shear-100
 mkdir -p $DIR; rm -f $DIR/*log $DIR/*yaml
 
-title Heuristic + Shear
-$RUN $COMMON_OPTS --samples 100000 --heuristic --log $DIR/heur.log --stats $DIR/heur.yaml || exit
+title Bound Pre-processing + Shear
+$RUN $COMMON_OPTS --samples $SAMPLES_FAST --bound --preprocess --log $DIR/bound-pre.log --stats $DIR/bound-pre.yaml || exit
 
 title Bound Original + Shear
-$RUN $COMMON_OPTS --samples 1000 --bound --log $DIR/bound.log --stats $DIR/bound.yaml || exit
-
-title Bound Pre-processing + Shear
-$RUN $COMMON_OPTS --samples 1000 --bound --preprocess --log $DIR/bound-pre.log --stats $DIR/bound-pre.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_MEDIUM --bound --log $DIR/bound.log --stats $DIR/bound.yaml || exit
 
 title Naive Bounding-Box + Shear
-$RUN $COMMON_OPTS --samples 10 --naive --log $DIR/naive-bbox.log --stats $DIR/naive-bbox.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_SLOW --naive --log $DIR/naive-bbox.log --stats $DIR/naive-bbox.yaml || exit
 
 title Naive Fixed-M + Shear
 M=50
-$RUN $COMMON_OPTS --samples 100 --naive --fixed-bounds $M --log $DIR/naive-fixed-$M.log --stats $DIR/naive-fixed-$M.yaml || exit
-
-title MIP + Shear
-$RUN $COMMON_OPTS --samples 10 --mip --log $DIR/mip.log --stats $DIR/mip.yaml || exit
-
-title MIP Pre-processing + Shear
-$RUN $COMMON_OPTS --samples 10 --mip --preprocess --log $DIR/mip-pre.log --stats $DIR/mip-pre.yaml || exit
+$RUN $COMMON_OPTS --samples $SAMPLES_MEDIUM --naive --fixed-bounds $M --log $DIR/naive-fixed-$M.log --stats $DIR/naive-fixed-$M.yaml || exit
 
+title MIP
+$RUN $COMMON_OPTS --samples $SAMPLES_SLOW --mip --log $DIR/mip.log --stats $DIR/mip.yaml || exit
 
 # TABLES
 # ------------------------------------------------------------------------------

multirow/CMakeLists.txt

@@ -8,7 +8,6 @@ set(COMMON_SOURCES
     src/util.c
     src/rational.c
     src/stats.c
-        src/linalg.c
     include/multirow/cg.h
     include/multirow/double.h
     include/multirow/geometry.h
@@ -18,26 +17,17 @@ set(COMMON_SOURCES
     include/multirow/rational.h
     include/multirow/stats.h
     include/multirow/params.h
-        include/multirow/util.h
-        include/multirow/linalg.h)
+    include/multirow/util.h)
 
 set(TEST_SOURCES
     tests/double-test.cpp
     tests/cg-test.cpp
-        tests/geometry-test.cpp
-        tests/linalg-test.cpp)
+    tests/geometry-test.cpp)
 
 add_library(multirow_static ${COMMON_SOURCES})
 set_target_properties(multirow_static PROPERTIES OUTPUT_NAME lifting)
-target_link_libraries(multirow_static
-        ${CPLEX_LIBRARIES}
-        ${BLAS_LIBRARIES}
-        ${LAPACKE_LIBRARIES})
-target_include_directories(multirow_static PUBLIC ${CMAKE_CURRENT_SOURCE_DIR})
+target_link_libraries(multirow_static ${CPLEX_LIBRARIES})
+target_include_directories (multirow_static PUBLIC ${CMAKE_CURRENT_SOURCE_DIR})
 
 add_executable(multirow-test.run ${COMMON_SOURCES} ${TEST_SOURCES})
-target_link_libraries(multirow-test.run
-        ${CPLEX_LIBRARIES}
-        ${BLAS_LIBRARIES}
-        ${LAPACKE_LIBRARIES}
-        gtest_main)
+target_link_libraries(multirow-test.run ${CPLEX_LIBRARIES} gtest_main)

multirow/include/multirow/cg.h

@@ -116,6 +116,4 @@ int CG_total_nz(const struct Tableau *tableau);
 
 double CG_replace_x(const struct Row *row, const double *x);
 
-int CG_print_model(const struct MultiRowModel *model);
-
 #endif //MULTIROW_CG_H

multirow/include/multirow/lfree2d.h

@@ -40,6 +40,7 @@ struct RayList
     int dim;
 };
 
+
 struct ConvLFreeSet
 {
     double *f;
@@ -81,6 +82,4 @@ int LFREE_init_conv(struct ConvLFreeSet *lfree, int dim, int max_nrays);
 
 void LFREE_free_conv(struct ConvLFreeSet *lfree);
 
-int LFREE_print_set(const struct ConvLFreeSet *lfree);
-
 #endif //LFREE_2D_H

multirow/include/multirow/linalg.h

@@ -1,40 +0,0 @@
-/* Copyright (c) 2015-2017 Alinson Xavier
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program.  If not, see <http://www.gnu.org/licenses/>.
- */
-
-#pragma once
-
-/**
- * Receives an n-dimensional vector x, a scalar alpha and sets x <- alpha * x
- */
-void LINALG_scale(int n, double *x, double alpha);
-
-/**
- * Receives two n-dimensional vectors x and y and returns the dot product of x
- * and y.
- */
-double LINALG_dot(int n, double *x, double *y);
-
-/**
- * Receives an n-dimensional vector x and returns the 1-norm of x.
- */
-double LINALG_norm(int n, double *x);
-
-/**
- * Given a full rank m-by-n matrix A and an m-dimensional vector b, this
- * function finds x such that Ax = b. Returns zero if the operation is
- * successful and non-zero otherwise.
- */
-int LINALG_solve(int n, int m, double *A, double *b, double *x);

multirow/include/multirow/params.h

@@ -41,7 +41,7 @@
 /*
  * Maximum number of sets that should be considered
  */
-#define MAX_N_SETS 1000
+#define MAX_N_SETS 1
 
 /*
  * Number of rays that should be generated per set.

multirow/include/multirow/util.h

@@ -78,7 +78,6 @@ void time_printf(const char *fmt,
 #define UNUSED(x) (void)(x)
 
 #define min(x,y) ((x) < (y) ? (x) : (y))
-#define max(x,y) ((x) > (y) ? (x) : (y))
 #define sign(x) ((x) > 0 ? 1 : ((x) < 0 ? -1 : 0))
 #define swap(x, y, T) do { T SWAP = x; x = y; y = SWAP; } while (0)
 

multirow/src/cg.c

@@ -945,10 +945,7 @@ int CG_add_multirow_cuts(struct CG *cg,
                 LP_free_row(&cut);
                 goto NEXT_COMBINATION;
             }
-            else if(rval) {
-                dump_tableau(&tableau, count);
-                abort_iff(1, "generate failed (cut %d)", count);
-            }
+            else abort_iff(rval, "generate failed (cut %d)", count);
 
             if_verbose_level dump_cut(&cut, count);
 
@@ -1118,32 +1115,4 @@ int CG_total_nz(const struct Tableau *tableau)
     return total_nz;
 }
 
-int CG_print_model(const struct MultiRowModel *model)
-{
-    int rval = 0;
-    int nrows = model->nrows;
-
-    time_printf("       f = [");
-    for (int i = 0; i < nrows; i++)
-        printf("%20.12lf ", model->f[i]);
-    printf("]\n");
-
-    for (int i = 0; i < model->rays.nrays; i++)
-    {
-        double *ray = LFREE_get_ray(&model->rays, i);
-        double norm = 0;
-
-        time_printf("ray[%3d] = [", i);
-        for (int j = 0; j < nrows; j++)
-        {
-            printf("%20.12lf ", ray[j]);
-            norm += fabs(ray[j]);
-        }
-        printf("] norm=%20.12lf \n", norm);
-    }
-
-CLEANUP:
-    return rval;
-}
-
 #endif // TEST_SOURCE

multirow/src/lfree2d.c

@@ -559,41 +559,6 @@ CLEANUP:
     return rval;
 }
 
-int LFREE_print_set(const struct ConvLFreeSet *set)
-{
-    int rval = 0;
-
-    log_debug("f=%12.6lf %12.6lf\n", set->f[0], set->f[1]);
-
-    for (int i = 0; i < set->rays.nrays; i++)
-    {
-        double *ray = LFREE_get_ray(&set->rays, i);
-        log_debug("ray[%-3d] = ", i);
-        for (int j = 0; j < set->nrows; j++)
-            printf("%20.12lf ", ray[j]);
-        printf("\n");
-    }
-
-    for (int i = 0; i < set->rays.nrays; i++)
-    {
-        double beta = set->beta[i];
-        log_debug("beta[%-3d] = %20.12lf\n", i, beta);
-    }
-
-    for (int i = 0; i < set->rays.nrays; i++)
-    {
-        double *ray = LFREE_get_ray(&set->rays, i);
-        double beta = set->beta[i];
-        log_debug("vertex[%-3d] = ", i);
-        for (int j = 0; j < set->nrows; j++)
-            printf("%20.12lf ", ray[j] * beta + set->f[j]);
-        printf("\n");
-    }
-
-CLEANUP:
-    return rval;
-}
-
 void LFREE_free_conv(struct ConvLFreeSet *lfree)
 {
     if(!lfree) return;

multirow/src/linalg.c

@@ -1,52 +0,0 @@
-/*
- * Copyright (C) 2016 Álinson Santos Xavier <isoron@gmail.com>
- *
- * This file is part of Loop Habit Tracker.
- *
- * Loop Habit Tracker is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by the
- * Free Software Foundation, either version 3 of the License, or (at your
- * option) any later version.
- *
- * Loop Habit Tracker is distributed in the hope that it will be useful, but
- * WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
- * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for
- * more details.
- *
- * You should have received a copy of the GNU General Public License along
- * with this program. If not, see <http://www.gnu.org/licenses/>.
- */
-
-#include <math.h>
-#include <cblas.h>
-#include <lapacke.h>
-#include <multirow/util.h>
-
-void LINALG_scale(int n, double *x, double alpha)
-{
-    cblas_dscal(n, alpha, x, 1);
-}
-
-double LINALG_dot(int n, double *x, double *y)
-{
-    return cblas_ddot(n, x, 1, y, 1);
-}
-
-double LINALG_norm(int n, double *x)
-{
-    return cblas_dasum(n, x, 1);
-}
-
-int LINALG_solve(int n, int m, double *A, double *b, double *x)
-{
-    int rval = 0;
-    double A_copy[n * m];
-    memcpy(x, b, m * sizeof(double));
-    memcpy(A_copy, A, n * m * sizeof(double));
-
-    LAPACKE_dgels(LAPACK_ROW_MAJOR, 'N', m, n, 1, A_copy, n, x, 1);
-    abort_if(rval, "LAPACKE_dgesv failed");
-
-CLEANUP:
-    return rval;
-}

multirow/tests/linalg-test.cpp

@@ -1,107 +0,0 @@
-/* Copyright (c) 2015-2017 Alinson Xavier
- *
- * This program is free software: you can redistribute it and/or modify
- * it under the terms of the GNU General Public License as published by
- * the Free Software Foundation, either version 3 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
- * GNU General Public License for more details.
- *
- * You should have received a copy of the GNU General Public License
- * along with this program.  If not, see <http://www.gnu.org/licenses/>.
- */
-
-#include <gtest/gtest.h>
-
-extern "C" {
-#include <multirow/util.h>
-#include <multirow/linalg.h>
-}
-
-const double E = 1e-6;
-
-TEST(LinAlgTest, scale_test)
-{
-    double x[] = {1.0, 2.0, 3.0};
-    LINALG_scale(3, x, 2.0);
-    EXPECT_NEAR(x[0], 2.0, E);
-    EXPECT_NEAR(x[1], 4.0, E);
-    EXPECT_NEAR(x[2], 6.0, E);
-}
-
-TEST(LinAlgTest, dot_test)
-{
-    double x[] = { 1.0, 2.0, 3.0 };
-    double y[] = { 3.0, 4.0, 5.0 };
-    double dot = LINALG_dot(3, x, y);
-    EXPECT_NEAR(dot, 26.0, E);
-}
-
-TEST(LinAlgTest, norm_test)
-{
-    double x[] = { 1.0,  2.0, -3.0 };
-    double y[] = { 3.0, -4.0, -5.0 };
-    double x_norm = LINALG_norm(3, x);
-    double y_norm = LINALG_norm(3, y);
-    EXPECT_NEAR(x_norm, 6.0, E);
-    EXPECT_NEAR(y_norm, 12.0, E);
-}
-
-TEST(LinAlgTest, solve_test)
-{
-    int rval = 0;
-
-    double A[] = {
-            2.0, 1.0,  3.0,
-            2.0, 6.0,  8.0,
-            6.0, 8.0, 18.0,
-    };
-    double b[] = { 1.0, 3.0, 5.0 };
-    double x[] = { 0.0, 0.0, 0.0 };
-
-    rval = LINALG_solve(3, 3, A, b, x);
-    abort_if(rval, "LINALG_solve failed");
-
-    // Should compute x correctly
-    EXPECT_NEAR(x[0], 0.3, E);
-    EXPECT_NEAR(x[1], 0.4, E);
-    EXPECT_NEAR(x[2], 0.0, E);
-
-    // Should not modify A and b
-    EXPECT_EQ(A[0], 2.0);
-    EXPECT_EQ(A[1], 1.0);
-    EXPECT_EQ(A[2], 3.0);
-    EXPECT_EQ(b[0], 1.0);
-    EXPECT_EQ(b[1], 3.0);
-    EXPECT_EQ(b[2], 5.0);
-
-CLEANUP:
-    if(rval) FAIL();
-}
-
-TEST(LinAlgTest, solve_test_2)
-{
-    int rval = 0;
-
-    double A[] = {
-            1.0, 3.0,
-            0.0, 0.0,
-            2.0, 2.0,
-            3.0, 5.0,
-    };
-    double b[] = { 5/2.0, 0.0, 8/3.0, 31/6.0 };
-    double x[2];
-
-    rval = LINALG_solve(2, 4, A, b, x);
-    abort_if(rval, "LINALG_solve failed");
-
-    // Should compute x correctly
-    EXPECT_NEAR(x[0], 3/4.0, E);
-    EXPECT_NEAR(x[1], 7/12.0, E);
-
-CLEANUP:
-    if(rval) FAIL();
-}