X-Git-Url: http://info.iut-bm.univ-fcomte.fr/pub/gitweb/simgrid.git/blobdiff_plain/69659a8281d4aeac3476820299f26a09dd285996..8e70d54a095cd7779a53c7a164f775f7a37b29e5:/src/kernel/lmm/bmf_test.cpp

diff --git a/src/kernel/lmm/bmf_test.cpp b/src/kernel/lmm/bmf_test.cpp
index eb3c824c02..f13e83da62 100644
--- a/src/kernel/lmm/bmf_test.cpp
+++ b/src/kernel/lmm/bmf_test.cpp
@@ -1,11 +1,11 @@
-/* Copyright (c) 2019-2022. The SimGrid Team. All rights reserved.               */
+/* Copyright (c) 2019-2023. The SimGrid Team. All rights reserved.               */
 
 /* This program is free software; you can redistribute it and/or modify it
  * under the terms of the license (GNU LGPL) which comes with this package. */
 
+#include "simgrid/math_utils.h"
 #include "src/include/catch.hpp"
 #include "src/kernel/lmm/bmf.hpp"
-#include "src/surf/surf_interface.hpp"
 #include "xbt/log.h"
 
 namespace lmm = simgrid::kernel::lmm;
@@ -34,7 +34,7 @@ TEST_CASE("kernel::bmf Basic tests", "[kernel-bmf-basic]")
     Sys.expand(sys_cnst, rho_1, 1);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 3, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 3, sg_precision_workamount));
   }
 
   SECTION("Two flows")
@@ -60,8 +60,8 @@ TEST_CASE("kernel::bmf Basic tests", "[kernel-bmf-basic]")
     Sys.expand(sys_cnst, rho_2, 10);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 3.0 / 2.0, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), (3.0 / 2.0) / 10.0, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 3.0 / 2.0, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), (3.0 / 2.0) / 10.0, sg_precision_workamount));
   }
 
   SECTION("Variable penalty/priority")
@@ -87,8 +87,8 @@ TEST_CASE("kernel::bmf Basic tests", "[kernel-bmf-basic]")
     Sys.expand(sys_cnst, rho_2, 1);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 2.0 / 3.0, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), 1.0 / 3.0, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 2.0 / 3.0, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 1.0 / 3.0, sg_precision_workamount));
   }
 
   SECTION("Disable variable doesn't count")
@@ -113,8 +113,8 @@ TEST_CASE("kernel::bmf Basic tests", "[kernel-bmf-basic]")
     Sys.expand(sys_cnst, rho_2, 10);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 1.0, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), 0.0, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 1.0, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 0.0, sg_precision_workamount));
   }
 
   SECTION("No consumption variable")
@@ -138,7 +138,7 @@ TEST_CASE("kernel::bmf Basic tests", "[kernel-bmf-basic]")
     Sys.expand(sys_cnst, rho_2, 10);
     Sys.solve();
 
-    REQUIRE(double_positive(rho_1->get_value(), sg_maxmin_precision));
+    REQUIRE(double_positive(rho_1->get_value(), sg_precision_workamount));
   }
 
   SECTION("Bounded variable")
@@ -163,8 +163,8 @@ TEST_CASE("kernel::bmf Basic tests", "[kernel-bmf-basic]")
     Sys.expand(sys_cnst, rho_1, 2);
     Sys.expand(sys_cnst, rho_2, 1);
     Sys.solve();
-    REQUIRE(double_equals(rho_1->get_value(), .1, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), .8, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), .1, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), .8, sg_precision_workamount));
   }
 
   SECTION("Fatpipe")
@@ -191,8 +191,8 @@ TEST_CASE("kernel::bmf Basic tests", "[kernel-bmf-basic]")
     Sys.expand(sys_cnst, rho_2, 1);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 3.0, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), 3.0, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 3.0, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 3.0, sg_precision_workamount));
   }
 
   SECTION("(un)Bounded variable")
@@ -217,8 +217,8 @@ TEST_CASE("kernel::bmf Basic tests", "[kernel-bmf-basic]")
     Sys.expand(sys_cnst, rho_1, 1);
     Sys.expand(sys_cnst, rho_2, 1);
     Sys.solve();
-    REQUIRE(double_equals(rho_1->get_value(), .5, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), .5, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), .5, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), .5, sg_precision_workamount));
   }
 
   SECTION("Dynamic bounds")
@@ -242,15 +242,15 @@ TEST_CASE("kernel::bmf Basic tests", "[kernel-bmf-basic]")
     lmm::Variable* rho_1 = Sys.variable_new(nullptr, 1);
     Sys.expand(sys_cnst, rho_1, 1);
     Sys.solve();
-    REQUIRE(double_equals(rho_1->get_value(), 1, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 1, sg_precision_workamount));
 
     // add another variable, half initial capacity
     lmm::Variable* rho_2 = Sys.variable_new(nullptr, 1);
     Sys.expand(sys_cnst, rho_2, 1);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), .25, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), .25, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), .25, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), .25, sg_precision_workamount));
   }
 
   Sys.variable_free_all();
@@ -291,8 +291,8 @@ TEST_CASE("kernel::bmf Advanced tests", "[kernel-bmf-advanced]")
     Sys.expand(sys_cnst2, rho_2, 1);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 1.0 / 11.0, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), 1.0 / 11.0, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 1.0 / 11.0, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 1.0 / 11.0, sg_precision_workamount));
   }
 
   SECTION("BMF paper example")
@@ -327,9 +327,9 @@ TEST_CASE("kernel::bmf Advanced tests", "[kernel-bmf-advanced]")
     Sys.expand(sys_cnst3, rho_3, 3.0 / 4.0);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 1.0 / 3.0, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), 4.0 / 9.0, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_3->get_value(), 4.0 / 9.0, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 1.0 / 3.0, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 4.0 / 9.0, sg_precision_workamount));
+    REQUIRE(double_equals(rho_3->get_value(), 4.0 / 9.0, sg_precision_workamount));
   }
 
   SECTION("IO - example")
@@ -368,8 +368,8 @@ TEST_CASE("kernel::bmf Advanced tests", "[kernel-bmf-advanced]")
     Sys.expand(sys_cnst3, rho_2, 1);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 1e6 / 2.0, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), 1e6 / 2.0, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 1e6 / 2.0, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 1e6 / 2.0, sg_precision_workamount));
   }
 
   SECTION("Proportional fairness")
@@ -398,9 +398,9 @@ TEST_CASE("kernel::bmf Advanced tests", "[kernel-bmf-advanced]")
     Sys.expand(sys_cnst, rho_3, epsilon);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 1.0 / (2.0 + 2 * epsilon), sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), 1.0 / (2.0 + 2 * epsilon), sg_maxmin_precision));
-    REQUIRE(double_equals(rho_3->get_value(), 1.0 / (1.0 + epsilon), sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 1.0 / (2.0 + 2 * epsilon), sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 1.0 / (2.0 + 2 * epsilon), sg_precision_workamount));
+    REQUIRE(double_equals(rho_3->get_value(), 1.0 / (1.0 + epsilon), sg_precision_workamount));
   }
 
   Sys.variable_free_all();
@@ -434,14 +434,14 @@ TEST_CASE("kernel::bmf Subflows", "[kernel-bmf-subflow]")
     lmm::Variable* rho_1      = Sys.variable_new(nullptr, 1);
     lmm::Variable* rho_2      = Sys.variable_new(nullptr, 1);
 
-    Sys.expand_add(sys_cnst, rho_1, 5);
-    Sys.expand_add(sys_cnst, rho_1, 7);
-    Sys.expand_add(sys_cnst, rho_2, 7);
-    Sys.expand_add(sys_cnst, rho_2, 5);
+    Sys.expand(sys_cnst, rho_1, 5);
+    Sys.expand(sys_cnst, rho_1, 7);
+    Sys.expand(sys_cnst, rho_2, 7);
+    Sys.expand(sys_cnst, rho_2, 5);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), 5.0 / 24.0, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), 5.0 / 24.0, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 5.0 / 24.0, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 5.0 / 24.0, sg_precision_workamount));
   }
 
   SECTION("1 subflows, 1 flow and 1 resource")
@@ -470,14 +470,14 @@ TEST_CASE("kernel::bmf Subflows", "[kernel-bmf-subflow]")
     lmm::Variable* rho_1      = Sys.variable_new(nullptr, 1);
     lmm::Variable* rho_2      = Sys.variable_new(nullptr, 1);
 
-    Sys.expand_add(sys_cnst, rho_1, 10);
-    Sys.expand_add(sys_cnst, rho_1, 5);
+    Sys.expand(sys_cnst, rho_1, 10);
+    Sys.expand(sys_cnst, rho_1, 5);
     Sys.expand(sys_cnst, rho_2, 10);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), (5.0 / 25.0), sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), (5.0 / 25.0), sg_maxmin_precision));
-    REQUIRE(double_equals(15 * rho_1->get_value(), 10 * rho_2->get_value() * 3 / 2, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), (5.0 / 25.0), sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), (5.0 / 25.0), sg_precision_workamount));
+    REQUIRE(double_equals(15 * rho_1->get_value(), 10 * rho_2->get_value() * 3 / 2, sg_precision_workamount));
   }
 
   SECTION("1 subflows using 2 resources: different max for each resource")
@@ -505,16 +505,16 @@ TEST_CASE("kernel::bmf Subflows", "[kernel-bmf-subflow]")
     lmm::Variable* rho_1       = Sys.variable_new(nullptr, 1, -1, 2);
     lmm::Variable* rho_2       = Sys.variable_new(nullptr, 1, -1, 2);
 
-    Sys.expand_add(sys_cnst, rho_1, 1.0);
-    Sys.expand_add(sys_cnst, rho_1, 1.0);
+    Sys.expand(sys_cnst, rho_1, 1.0);
+    Sys.expand(sys_cnst, rho_1, 1.0);
     Sys.expand(sys_cnst, rho_2, 1);
-    Sys.expand_add(sys_cnst2, rho_1, 1.0 / 2.0);
-    Sys.expand_add(sys_cnst2, rho_1, 1.0 / 2.0);
+    Sys.expand(sys_cnst2, rho_1, 1.0 / 2.0);
+    Sys.expand(sys_cnst2, rho_1, 1.0 / 2.0);
     Sys.expand(sys_cnst2, rho_2, 3.0 / 2.0);
     Sys.solve();
 
-    REQUIRE(double_equals(rho_1->get_value(), (1.0 / 3.0), sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), (1.0 / 3.0), sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), (1.0 / 3.0), sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), (1.0 / 3.0), sg_precision_workamount));
   }
 
   Sys.variable_free_all();
@@ -550,13 +550,13 @@ TEST_CASE("kernel::bmf Loop", "[kernel-bmf-loop]")
                   [&vars, &Sys, &A](const auto&) { vars.push_back(Sys.variable_new(nullptr, 1, -1, A.size())); });
     for (size_t j = 0; j < A.size(); j++) {
       for (size_t i = 0; i < A[j].size(); i++) {
-        Sys.expand_add(sys_cnst[j], vars[i], A[j][i]);
+        Sys.expand(sys_cnst[j], vars[i], A[j][i]);
       }
     }
     Sys.solve();
 
     for (const auto* rho : vars) {
-      REQUIRE(double_positive(rho->get_value(), sg_maxmin_precision));
+      REQUIRE(double_positive(rho->get_value(), sg_precision_workamount));
     }
   }
 
@@ -581,13 +581,63 @@ TEST_CASE("kernel::bmf Bugs", "[kernel-bmf-bug]")
     lmm::Variable* rho_2       = Sys.variable_new(nullptr, 1, 2.27328e-10, 2);
     lmm::Variable* rho_3       = Sys.variable_new(nullptr, 1);
 
-    Sys.expand_add(sys_cnst, rho_1, 1.84467e+19);
-    Sys.expand_add(sys_cnst2, rho_1, 1.84467e+19);
-    Sys.expand_add(sys_cnst, rho_2, 1.84467e+19);
-    Sys.expand_add(sys_cnst, rho_3, 1.91268e+11);
+    Sys.expand(sys_cnst, rho_1, 1.84467e+19);
+    Sys.expand(sys_cnst2, rho_1, 1.84467e+19);
+    Sys.expand(sys_cnst, rho_2, 1.84467e+19);
+    Sys.expand(sys_cnst, rho_3, 1.91268e+11);
     Sys.solve();
   }
 
+  SECTION("DadOu's bug: complete matrix")
+  {
+    constexpr int cnsts = 71;
+    constexpr int flows = 3;
+
+    Eigen::MatrixXd A(cnsts, flows);
+    A << 0, 0, 0, 0, 0, 0, 0, 0, 1.84467e+19, 1.91268e+11, 1.84467e+19, 1.84467e+19, 0, 0, 1.84467e+19, 0, 1.84467e+19,
+        0, 0, 1.84467e+19, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0,
+        4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11,
+        0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 1.256e+09, 0, 0, 2.2372e+10, 0, 0,
+        2.2684e+10, 0, 0, 2.2588e+10, 0, 0, 2.3188e+10, 0, 0, 2.228e+10, 0, 0, 2.058e+10, 0, 0, 2.2684e+10, 0, 0,
+        2.2824e+10, 0, 0, 2.2976e+10, 0, 0, 2.2632e+10, 0, 0, 2.2584e+10, 0, 0, 2.2736e+10, 0, 0, 2.2616e+10, 0, 0,
+        2.0828e+10, 0, 0, 2.3184e+10, 0, 0, 2.2524e+10, 0, 0, 2.278e+10, 0, 0, 2.2164e+10, 0, 0, 1.26e+09, 0, 0,
+        2.1872e+10, 0, 0, 1.4e+09, 0, 0, 2.3184e+10, 0, 0, 8.52e+08, 0, 0, 2.2268e+10, 0, 0, 1.756e+09, 0, 0,
+        2.0636e+10, 0, 0, 3.4e+09, 0, 0, 2.2576e+10, 0, 0, 1.352e+09, 0, 0, 2.2832e+10, 0, 0, 1.2e+09, 0, 0, 2.3092e+10,
+        0, 0, 9.48e+08, 0, 0, 2.2436e+10, 0, 0, 1.4e+09, 0, 0, 2.2572e+10, 0, 0, 1.452e+09, 0, 0, 2.2692e+10, 0, 0,
+        1.3e+09, 0, 0, 2.2832e+10, 0, 0, 1.2e+09, 0, 0, 2.1232e+10, 0, 0, 2.804e+09, 0, 0, 2.3184e+10, 0, 0,
+        8.56001e+08, 0, 0, 2.2512e+10, 0, 0, 1.5e+09, 0, 0;
+
+    Eigen::MatrixXd maxA(cnsts, flows);
+    maxA << 0, 0, 0, 0, 0, 0, 0, 0, 1.84467e+19, 3.1e+09, 1.84467e+19, 1.84467e+19, 0, 0, 1.84467e+19, 0, 1.84467e+19,
+        0, 0, 1.84467e+19, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 1, 4.975e+11, 0, 0,
+        4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11,
+        0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 4.975e+11, 0, 0, 1.256e+09, 0, 0, 3.504e+09, 0, 0,
+        3.056e+09, 0, 0, 3.1e+09, 0, 0, 2.952e+09, 0, 0, 2.404e+09, 0, 0, 2.304e+09, 0, 0, 2.556e+09, 0, 0, 2.4e+09, 0,
+        0, 2.804e+09, 0, 0, 2.552e+09, 0, 0, 2.408e+09, 0, 0, 2.9e+09, 0, 0, 2.4e+09, 0, 0, 2.256e+09, 0, 0, 3.504e+09,
+        0, 0, 3.244e+09, 0, 0, 2.556e+09, 0, 0, 2.952e+09, 0, 0, 1.26e+09, 0, 0, 2.552e+09, 0, 0, 1.4e+09, 0, 0,
+        3.244e+09, 0, 0, 8.52e+08, 0, 0, 2.556e+09, 0, 0, 1.756e+09, 0, 0, 2.256e+09, 0, 0, 3.4e+09, 0, 0, 2.6e+09, 0,
+        0, 1.352e+09, 0, 0, 2.952e+09, 0, 0, 1.2e+09, 0, 0, 2.452e+09, 0, 0, 9.48e+08, 0, 0, 2.804e+09, 0, 0, 1.4e+09,
+        0, 0, 3.1e+09, 0, 0, 1.452e+09, 0, 0, 2.404e+09, 0, 0, 1.3e+09, 0, 0, 2.952e+09, 0, 0, 1.2e+09, 0, 0, 3.056e+09,
+        0, 0, 2.804e+09, 0, 0, 2.4e+09, 0, 0, 8.56001e+08, 0, 0, 2.9e+09, 0, 0, 1.5e+09, 0, 0;
+
+    Eigen::VectorXd C(cnsts);
+    C << 3.2e+20, 3.2e+20, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20,
+        3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20, 3.2e+20, 1.83484e+10,
+        2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09,
+        2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10,
+        2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10,
+        2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10,
+        2.5e+09, 1.83484e+10, 2.5e+09, 1.83484e+10;
+
+    Eigen::VectorXd phi(flows);
+    phi << 1.35273, 2.27328e-10, 2.27328e-10;
+
+    std::vector<bool> shared(cnsts, true);
+    lmm::BmfSolver solver(A, maxA, C, shared, phi);
+    auto rho = solver.solve();
+    REQUIRE((rho.array() > 0).all());
+  }
+
   SECTION("is_bmf bug: all limited by bound")
   {
     /*
@@ -600,13 +650,66 @@ TEST_CASE("kernel::bmf Bugs", "[kernel-bmf-bug]")
     lmm::Variable* rho_1       = Sys.variable_new(nullptr, 1, 1.5, 2);
     lmm::Variable* rho_2       = Sys.variable_new(nullptr, 1, 3, 2);
 
-    Sys.expand_add(sys_cnst, rho_1, 5.0);
-    Sys.expand_add(sys_cnst2, rho_1, 1.0);
-    Sys.expand_add(sys_cnst, rho_2, 1.0);
-    Sys.expand_add(sys_cnst2, rho_2, 1.0);
+    Sys.expand(sys_cnst, rho_1, 5.0);
+    Sys.expand(sys_cnst2, rho_1, 1.0);
+    Sys.expand(sys_cnst, rho_2, 1.0);
+    Sys.expand(sys_cnst2, rho_2, 1.0);
     Sys.solve();
-    REQUIRE(double_equals(rho_1->get_value(), 1.4, sg_maxmin_precision));
-    REQUIRE(double_equals(rho_2->get_value(), 3, sg_maxmin_precision));
+    REQUIRE(double_equals(rho_1->get_value(), 1.4, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 3, sg_precision_workamount));
+  }
+
+  SECTION("s4u-cloud-capping bug: all limited by bound extra case")
+  {
+    /*
+     * Inspired by s4u-cloud-capping test.
+     * Step 9: "Put an activity on a PM and an activity on the VM capped by 10%."
+     * Extracted from original test.
+     * The sum of bounds for 1 resource is smaller than C.
+     */
+
+    lmm::Constraint* sys_cnst = Sys.constraint_new(nullptr, 7.6296e+07);
+    lmm::Variable* rho_1      = Sys.variable_new(nullptr, 1, 7.6296e+06, 1);
+    lmm::Variable* rho_2      = Sys.variable_new(nullptr, 1, 3.8148e+07, 1);
+
+    Sys.expand(sys_cnst, rho_1, 1);
+    Sys.expand(sys_cnst, rho_2, 1);
+    Sys.solve();
+    REQUIRE(double_equals(rho_1->get_value(), 7.6296e+06, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 3.8148e+07, sg_precision_workamount));
+  }
+
+  SECTION("Variable penalty with bounds: thread bug")
+  {
+    /*
+     * Detected by exec-thread.
+     * Fair-sharing vector depends on the penalty too.
+     */
+
+    /* don't change the order of creation, important to trigger the bug */
+    lmm::Constraint* sys_cnst = Sys.constraint_new(nullptr, 4e8);
+    lmm::Variable* rho_2      = Sys.variable_new(nullptr, .25, 4e8, 1);
+    lmm::Variable* rho_1      = Sys.variable_new(nullptr, 1, 1e8, 1);
+    Sys.expand(sys_cnst, rho_2, 1);
+    Sys.expand(sys_cnst, rho_1, 1);
+    Sys.solve();
+
+    REQUIRE(double_equals(rho_1->get_value(), 8e7, sg_precision_workamount));
+    REQUIRE(double_equals(rho_2->get_value(), 3.2e8, sg_precision_workamount));
+  }
+
+  SECTION("Variable penalty with bounds greater than C")
+  {
+    /*
+     * Detected by exec-thread. 6 thread running in a 4 core CPU.
+     */
+
+    lmm::Constraint* sys_cnst = Sys.constraint_new(nullptr, 4e8);
+    lmm::Variable* rho_1      = Sys.variable_new(nullptr, 1.0 / 6.0, 6e8, 1);
+    Sys.expand(sys_cnst, rho_1, 1);
+    Sys.solve();
+
+    REQUIRE(double_equals(rho_1->get_value(), 4e8, sg_precision_workamount));
   }
 
   Sys.variable_free_all();
@@ -630,7 +733,7 @@ TEST_CASE("kernel::bmf Stress-tests", "[.kernel-bmf-stress]")
     for (int j = 0; j < N; j++) {
       lmm::Variable* rho = Sys.variable_new(nullptr, 1, -1, C);
       for (int i = 0; i < C; i++) {
-        Sys.expand_add(sys_cnst[i], rho, data[i * j + j]);
+        Sys.expand(sys_cnst[i], rho, data[i * j + j]);
       }
     }
     Sys.solve();
@@ -645,7 +748,7 @@ TEST_CASE("kernel::bmf Stress-tests", "[.kernel-bmf-stress]")
     for (int j = 0; j < N; j++) {
       for (int i = 0; i < C; i++) {
         lmm::Variable* rho = Sys.variable_new(nullptr, 1);
-        Sys.expand_add(sys_cnst[i], rho, data[i * j + j]);
+        Sys.expand(sys_cnst[i], rho, data[i * j + j]);
       }
     }
     Sys.solve();