src/surf/lagrange.cpp

   1 /* Copyright (c) 2007-2014. The SimGrid Team.
   2  * All rights reserved.                                                     */
   3
   4 /* This program is free software; you can redistribute it and/or modify it
   5  * under the terms of the license (GNU LGPL) which comes with this package. */
   6
   7 /*
   8  * Modeling the proportional fairness using the Lagrangian Optimization Approach. For a detailed description see:
   9  * "ssh://username@scm.gforge.inria.fr/svn/memo/people/pvelho/lagrange/ppf.ps".
  10  */
  11 #include "xbt/log.h"
  12 #include "xbt/sysdep.h"
  13 #include "maxmin_private.hpp"
  14
  15 #include <stdlib.h>
  16 #ifndef MATH
  17 #include <math.h>
  18 #endif
  19
  20 XBT_LOG_NEW_DEFAULT_SUBCATEGORY(surf_lagrange, surf, "Logging specific to SURF (lagrange)");
  21 XBT_LOG_NEW_SUBCATEGORY(surf_lagrange_dichotomy, surf_lagrange, "Logging specific to SURF (lagrange dichotomy)");
  22
  23 #define SHOW_EXPR(expr) XBT_CDEBUG(surf_lagrange,#expr " = %g",expr);
  24
  25 double (*func_f_def) (lmm_variable_t, double);
  26 double (*func_fp_def) (lmm_variable_t, double);
  27 double (*func_fpi_def) (lmm_variable_t, double);
  28
  29 /*
  30  * Local prototypes to implement the Lagrangian optimization with optimal step, also called dichotomy.
  31  */
  32 //solves the proportional fairness using a Lagrangian optimization with dichotomy step
  33 void lagrange_solve(lmm_system_t sys);
  34 //computes the value of the dichotomy using a initial values, init, with a specific variable or constraint
  35 static double dichotomy(double init, double diff(double, void *), void *var_cnst, double min_error);
  36 //computes the value of the differential of constraint param_cnst applied to lambda
  37 static double partial_diff_lambda(double lambda, void *param_cnst);
  38
  39 static int __check_feasible(xbt_swag_t cnst_list, xbt_swag_t var_list, int warn)
  40 {
  41   void *_cnst, *_elem, *_var;
  42   xbt_swag_t elem_list = nullptr;
  43   lmm_element_t elem = nullptr;
  44   lmm_constraint_t cnst = nullptr;
  45   lmm_variable_t var = nullptr;
  46
  47   double tmp;
  48
  49   xbt_swag_foreach(_cnst, cnst_list) {
  50   cnst = static_cast<lmm_constraint_t>(_cnst);
  51     tmp = 0;
  52     elem_list = &(cnst->enabled_element_set);
  53     xbt_swag_foreach(_elem, elem_list) {
  54       elem = static_cast<lmm_element_t>(_elem);
  55       var = elem->variable;
  56       xbt_assert(var->weight > 0);
  57       tmp += var->value;
  58     }
  59
  60     if (double_positive(tmp - cnst->bound, sg_maxmin_precision)) {
  61       if (warn)
  62         XBT_WARN ("The link (%p) is over-used. Expected less than %f and got %f", cnst, cnst->bound, tmp);
  63       return 0;
  64     }
  65     XBT_DEBUG ("Checking feasability for constraint (%p): sat = %f, lambda = %f ", cnst, tmp - cnst->bound,
  66                cnst->lambda);
  67   }
  68
  69   xbt_swag_foreach(_var, var_list) {
  70     var = static_cast<lmm_variable_t>(_var);
  71     if (!var->weight)
  72       break;
  73     if (var->bound < 0)
  74       continue;
  75     XBT_DEBUG("Checking feasability for variable (%p): sat = %f mu = %f", var, var->value - var->bound, var->mu);
  76
  77     if (double_positive(var->value - var->bound, sg_maxmin_precision)) {
  78       if (warn)
  79         XBT_WARN ("The variable (%p) is too large. Expected less than %f and got %f", var, var->bound, var->value);
  80       return 0;
  81     }
  82   }
  83   return 1;
  84 }
  85
  86 static double new_value(lmm_variable_t var)
  87 {
  88   double tmp = 0;
  89
  90   for (int i = 0; i < var->cnsts_number; i++) {
  91     tmp += (var->cnsts[i].constraint)->lambda;
  92   }
  93   if (var->bound > 0)
  94     tmp += var->mu;
  95   XBT_DEBUG("\t Working on var (%p). cost = %e; Weight = %e", var, tmp, var->weight);
  96   //uses the partial differential inverse function
  97   return var->func_fpi(var, tmp);
  98 }
  99
 100 static double new_mu(lmm_variable_t var)
 101 {
 102   double mu_i = 0.0;
 103   double sigma_i = 0.0;
 104
 105   for (int j = 0; j < var->cnsts_number; j++) {
 106     sigma_i += (var->cnsts[j].constraint)->lambda;
 107   }
 108   mu_i = var->func_fp(var, var->bound) - sigma_i;
 109   if (mu_i < 0.0)
 110     return 0.0;
 111   return mu_i;
 112 }
 113
 114 static double dual_objective(xbt_swag_t var_list, xbt_swag_t cnst_list)
 115 {
 116   void *_cnst;
 117   void *_var;
 118   lmm_constraint_t cnst = nullptr;
 119   lmm_variable_t var = nullptr;
 120
 121   double obj = 0.0;
 122
 123   xbt_swag_foreach(_var, var_list) {
 124     var = static_cast<lmm_variable_t>(_var);
 125     double sigma_i = 0.0;
 126
 127     if (!var->weight)
 128       break;
 129
 130     for (int j = 0; j < var->cnsts_number; j++)
 131       sigma_i += (var->cnsts[j].constraint)->lambda;
 132
 133     if (var->bound > 0)
 134       sigma_i += var->mu;
 135
 136     XBT_DEBUG("var %p : sigma_i = %1.20f", var, sigma_i);
 137
 138     obj += var->func_f(var, var->func_fpi(var, sigma_i)) - sigma_i * var->func_fpi(var, sigma_i);
 139
 140     if (var->bound > 0)
 141       obj += var->mu * var->bound;
 142   }
 143
 144   xbt_swag_foreach(_cnst, cnst_list) {
 145     cnst = static_cast<lmm_constraint_t>(_cnst);
 146     obj += cnst->lambda * cnst->bound;
 147   }
 148
 149   return obj;
 150 }
 151
 152 void lagrange_solve(lmm_system_t sys)
 153 {
 154   /* Lagrange Variables. */
 155   int max_iterations = 100;
 156   double epsilon_min_error = 0.00001; /* this is the precision on the objective function so it's none of the configurable values and this value is the legacy one */
 157   double dichotomy_min_error = 1e-14;
 158   double overall_modification = 1;
 159
 160   /* Variables to manipulate the data structure proposed to model the maxmin fairness. See documentation for details. */
 161   xbt_swag_t cnst_list = nullptr;
 162   void *_cnst;
 163   lmm_constraint_t cnst = nullptr;
 164
 165   xbt_swag_t var_list = nullptr;
 166   void *_var;
 167   lmm_variable_t var = nullptr;
 168
 169   /* Auxiliary variables. */
 170   int iteration = 0;
 171   double tmp = 0;
 172   int i;
 173   double obj;
 174   double new_obj;
 175
 176   XBT_DEBUG("Iterative method configuration snapshot =====>");
 177   XBT_DEBUG("#### Maximum number of iterations        : %d", max_iterations);
 178   XBT_DEBUG("#### Minimum error tolerated             : %e", epsilon_min_error);
 179   XBT_DEBUG("#### Minimum error tolerated (dichotomy) : %e", dichotomy_min_error);
 180
 181   if (XBT_LOG_ISENABLED(surf_lagrange, xbt_log_priority_debug)) {
 182     lmm_print(sys);
 183   }
 184
 185   if (!(sys->modified))
 186     return;
 187
 188   /* Initialize lambda. */
 189   cnst_list = &(sys->active_constraint_set);
 190   xbt_swag_foreach(_cnst, cnst_list) {
 191   cnst = (lmm_constraint_t)_cnst;
 192     cnst->lambda = 1.0;
 193     cnst->new_lambda = 2.0;
 194     XBT_DEBUG("#### cnst(%p)->lambda :  %e", cnst, cnst->lambda);
 195   }
 196
 197   /*
 198    * Initialize the var list variable with only the active variables.
 199    * Associate an index in the swag variables. Initialize mu.
 200    */
 201   var_list = &(sys->variable_set);
 202   i = 0;
 203   xbt_swag_foreach(_var, var_list) {
 204   var = static_cast<lmm_variable_t>(_var);
 205     if (!var->weight)
 206       var->value = 0.0;
 207     else {
 208       int nb = 0;
 209       if (var->bound < 0.0) {
 210         XBT_DEBUG("#### NOTE var(%d) is a boundless variable", i);
 211         var->mu = -1.0;
 212         var->value = new_value(var);
 213       } else {
 214         var->mu = 1.0;
 215         var->new_mu = 2.0;
 216         var->value = new_value(var);
 217       }
 218       XBT_DEBUG("#### var(%p) ->weight :  %e", var, var->weight);
 219       XBT_DEBUG("#### var(%p) ->mu :  %e", var, var->mu);
 220       XBT_DEBUG("#### var(%p) ->weight: %e", var, var->weight);
 221       XBT_DEBUG("#### var(%p) ->bound: %e", var, var->bound);
 222       for (i = 0; i < var->cnsts_number; i++) {
 223         if (var->cnsts[i].value ==   0.0)
 224           nb++;
 225       }
 226       if (nb == var->cnsts_number)
 227         var->value = 1.0;
 228     }
 229   }
 230
 231   /*  Compute dual objective. */
 232   obj = dual_objective(var_list, cnst_list);
 233
 234   /* While doesn't reach a minimum error or a number maximum of iterations. */
 235   while (overall_modification > epsilon_min_error && iteration < max_iterations) {
 236     iteration++;
 237     XBT_DEBUG("************** ITERATION %d **************", iteration);
 238     XBT_DEBUG("-------------- Gradient Descent ----------");
 239
 240     /* Improve the value of mu_i */
 241     xbt_swag_foreach(_var, var_list) {
 242       var = static_cast<lmm_variable_t>(_var);
 243       if (!var->weight)
 244         break;
 245       if (var->bound >= 0) {
 246         XBT_DEBUG("Working on var (%p)", var);
 247         var->new_mu = new_mu(var);
 248 /*   dual_updated += (fabs(var->new_mu-var->mu)>dichotomy_min_error); */
 249 /*   XBT_DEBUG("dual_updated (%d) : %1.20f",dual_updated,fabs(var->new_mu-var->mu)); */
 250         XBT_DEBUG("Updating mu : var->mu (%p) : %1.20f -> %1.20f", var, var->mu, var->new_mu);
 251         var->mu = var->new_mu;
 252
 253         new_obj = dual_objective(var_list, cnst_list);
 254         XBT_DEBUG("Improvement for Objective (%g -> %g) : %g", obj, new_obj, obj - new_obj);
 255         xbt_assert(obj - new_obj >= -epsilon_min_error, "Our gradient sucks! (%1.20f)", obj - new_obj);
 256         obj = new_obj;
 257       }
 258     }
 259
 260     /* Improve the value of lambda_i */
 261     xbt_swag_foreach(_cnst, cnst_list) {
 262       cnst = static_cast<lmm_constraint_t>(_cnst);
 263       XBT_DEBUG("Working on cnst (%p)", cnst);
 264       cnst->new_lambda = dichotomy(cnst->lambda, partial_diff_lambda, cnst, dichotomy_min_error);
 265 /*       dual_updated += (fabs(cnst->new_lambda-cnst->lambda)>dichotomy_min_error); */
 266 /*       XBT_DEBUG("dual_updated (%d) : %1.20f",dual_updated,fabs(cnst->new_lambda-cnst->lambda)); */
 267       XBT_DEBUG("Updating lambda : cnst->lambda (%p) : %1.20f -> %1.20f", cnst, cnst->lambda, cnst->new_lambda);
 268       cnst->lambda = cnst->new_lambda;
 269
 270       new_obj = dual_objective(var_list, cnst_list);
 271       XBT_DEBUG("Improvement for Objective (%g -> %g) : %g", obj, new_obj, obj - new_obj);
 272       xbt_assert(obj - new_obj >= -epsilon_min_error, "Our gradient sucks! (%1.20f)", obj - new_obj);
 273       obj = new_obj;
 274     }
 275
 276     /* Now computes the values of each variable (\rho) based on the values of \lambda and \mu. */
 277     XBT_DEBUG("-------------- Check convergence ----------");
 278     overall_modification = 0;
 279     xbt_swag_foreach(_var, var_list) {
 280       var = static_cast<lmm_variable_t>(_var);
 281       if (var->weight <= 0)
 282         var->value = 0.0;
 283       else {
 284         tmp = new_value(var);
 285
 286         overall_modification = MAX(overall_modification, fabs(var->value - tmp));
 287
 288         var->value = tmp;
 289         XBT_DEBUG("New value of var (%p)  = %e, overall_modification = %e", var, var->value, overall_modification);
 290       }
 291     }
 292
 293     XBT_DEBUG("-------------- Check feasability ----------");
 294     if (!__check_feasible(cnst_list, var_list, 0))
 295       overall_modification = 1.0;
 296     XBT_DEBUG("Iteration %d: overall_modification : %f", iteration, overall_modification);
 297 /*     if(!dual_updated) { */
 298 /*       XBT_WARN("Could not improve the convergence at iteration %d. Drop it!",iteration); */
 299 /*       break; */
 300 /*     } */
 301   }
 302
 303   __check_feasible(cnst_list, var_list, 1);
 304
 305   if (overall_modification <= epsilon_min_error) {
 306     XBT_DEBUG("The method converges in %d iterations.", iteration);
 307   }
 308   if (iteration >= max_iterations) {
 309     XBT_DEBUG ("Method reach %d iterations, which is the maximum number of iterations allowed.", iteration);
 310   }
 311
 312   if (XBT_LOG_ISENABLED(surf_lagrange, xbt_log_priority_debug)) {
 313     lmm_print(sys);
 314   }
 315 }
 316
 317 /*
 318  * Returns a double value corresponding to the result of a dichotomy process with respect to a given
 319  * variable/constraint (\mu in the case of a variable or \lambda in case of a constraint) and a initial value init.
 320  *
 321  * @param init initial value for \mu or \lambda
 322  * @param diff a function that computes the differential of with respect a \mu or \lambda
 323  * @param var_cnst a pointer to a variable or constraint
 324  * @param min_erro a minimum error tolerated
 325  *
 326  * @return a double corresponding to the result of the dichotomy process
 327  */
 328 static double dichotomy(double init, double diff(double, void *), void *var_cnst, double min_error)
 329 {
 330   double min =init;
 331   double max= init;
 332   double overall_error;
 333   double middle;
 334   double middle_diff;
 335   double diff_0 = 0.0;
 336
 337   XBT_IN();
 338
 339   if (fabs(init) < 1e-20) {
 340     min = 0.5;
 341     max = 0.5;
 342   }
 343
 344   overall_error = 1;
 345
 346   diff_0 = diff(1e-16, var_cnst);
 347   if (diff_0 >= 0) {
 348     XBT_CDEBUG(surf_lagrange_dichotomy, "returning 0.0 (diff = %e)", diff_0);
 349     XBT_OUT();
 350     return 0.0;
 351   }
 352
 353   double min_diff = diff(min, var_cnst);
 354   double max_diff = diff(max, var_cnst);
 355
 356   while (overall_error > min_error) {
 357     XBT_CDEBUG(surf_lagrange_dichotomy, "[min, max] = [%1.20f, %1.20f] || diffmin, diffmax = %1.20f, %1.20f",
 358                min, max, min_diff, max_diff);
 359
 360     if (min_diff > 0 && max_diff > 0) {
 361       if (min == max) {
 362         XBT_CDEBUG(surf_lagrange_dichotomy, "Decreasing min");
 363         min = min / 2.0;
 364         min_diff = diff(min, var_cnst);
 365       } else {
 366         XBT_CDEBUG(surf_lagrange_dichotomy, "Decreasing max");
 367         max = min;
 368         max_diff = min_diff;
 369       }
 370     } else if (min_diff < 0 && max_diff < 0) {
 371       if (min == max) {
 372         XBT_CDEBUG(surf_lagrange_dichotomy, "Increasing max");
 373         max = max * 2.0;
 374         max_diff = diff(max, var_cnst);
 375       } else {
 376         XBT_CDEBUG(surf_lagrange_dichotomy, "Increasing min");
 377         min = max;
 378         min_diff = max_diff;
 379       }
 380     } else if (min_diff < 0 && max_diff > 0) {
 381       middle = (max + min) / 2.0;
 382       XBT_CDEBUG(surf_lagrange_dichotomy, "Trying (max+min)/2 : %1.20f", middle);
 383
 384       if ((fabs(min - middle) < 1e-20) || (fabs(max - middle) < 1e-20)){
 385         XBT_CWARN(surf_lagrange_dichotomy, "Cannot improve the convergence! min=max=middle=%1.20f, diff = %1.20f."
 386                   " Reaching the 'double' limits. Maybe scaling your function would help ([%1.20f,%1.20f]).",
 387                   min, max - min, min_diff, max_diff);
 388         break;
 389       }
 390       middle_diff = diff(middle, var_cnst);
 391
 392       if (middle_diff < 0) {
 393         XBT_CDEBUG(surf_lagrange_dichotomy, "Increasing min");
 394         min = middle;
 395         overall_error = max_diff - middle_diff;
 396         min_diff = middle_diff;
 397 /*   SHOW_EXPR(overall_error); */
 398       } else if (middle_diff > 0) {
 399         XBT_CDEBUG(surf_lagrange_dichotomy, "Decreasing max");
 400         max = middle;
 401         overall_error = max_diff - middle_diff;
 402         max_diff = middle_diff;
 403 /*   SHOW_EXPR(overall_error); */
 404       } else {
 405         overall_error = 0;
 406 /*   SHOW_EXPR(overall_error); */
 407       }
 408     } else if (fabs(min_diff) < 1e-20) {
 409       max = min;
 410       overall_error = 0;
 411 /*       SHOW_EXPR(overall_error); */
 412     } else if (fabs(max_diff) < 1e-20) {
 413       min = max;
 414       overall_error = 0;
 415 /*       SHOW_EXPR(overall_error); */
 416     } else if (min_diff > 0 && max_diff < 0) {
 417       XBT_CWARN(surf_lagrange_dichotomy, "The impossible happened, partial_diff(min) > 0 && partial_diff(max) < 0");
 418       xbt_abort();
 419     } else {
 420       XBT_CWARN(surf_lagrange_dichotomy,
 421              "diffmin (%1.20f) or diffmax (%1.20f) are something I don't know, taking no action.",
 422              min_diff, max_diff);
 423       xbt_abort();
 424     }
 425   }
 426
 427   XBT_CDEBUG(surf_lagrange_dichotomy, "returning %e", (min + max) / 2.0);
 428   XBT_OUT();
 429   return ((min + max) / 2.0);
 430 }
 431
 432 static double partial_diff_lambda(double lambda, void *param_cnst)
 433 {
 434   int j;
 435   void *_elem;
 436   xbt_swag_t elem_list = nullptr;
 437   lmm_element_t elem = nullptr;
 438   lmm_variable_t var = nullptr;
 439   lmm_constraint_t cnst = static_cast<lmm_constraint_t>(param_cnst);
 440   double diff = 0.0;
 441   double sigma_i = 0.0;
 442
 443   XBT_IN();
 444   elem_list = &(cnst->enabled_element_set);
 445
 446   XBT_CDEBUG(surf_lagrange_dichotomy, "Computing diff of cnst (%p)", cnst);
 447
 448   xbt_swag_foreach(_elem, elem_list) {
 449   elem = static_cast<lmm_element_t>(_elem);
 450     var = elem->variable;
 451     xbt_assert(var->weight > 0);
 452     XBT_CDEBUG(surf_lagrange_dichotomy, "Computing sigma_i for var (%p)", var);
 453     // Initialize the summation variable
 454     sigma_i = 0.0;
 455
 456     // Compute sigma_i
 457     for (j = 0; j < var->cnsts_number; j++) {
 458       sigma_i += (var->cnsts[j].constraint)->lambda;
 459     }
 460
 461     //add mu_i if this flow has a RTT constraint associated
 462     if (var->bound > 0)
 463       sigma_i += var->mu;
 464
 465     //replace value of cnst->lambda by the value of parameter lambda
 466     sigma_i = (sigma_i - cnst->lambda) + lambda;
 467
 468     diff += -var->func_fpi(var, sigma_i);
 469   }
 470
 471   diff += cnst->bound;
 472
 473   XBT_CDEBUG(surf_lagrange_dichotomy, "d D/d lambda for cnst (%p) at %1.20f = %1.20f", cnst, lambda, diff);
 474   XBT_OUT();
 475   return diff;
 476 }
 477
 478 /** \brief Attribute the value bound to var->bound.
 479  *
 480  *  \param func_fpi  inverse of the partial differential of f (f prime inverse, (f')^{-1})
 481  *
 482  *  Set default functions to the ones passed as parameters. This is a polymorphism in C pure, enjoy the roots of
 483  *  programming.
 484  *
 485  */
 486 void lmm_set_default_protocol_function(double (*func_f) (lmm_variable_t var, double x),
 487                                        double (*func_fp) (lmm_variable_t var, double x),
 488                                        double (*func_fpi) (lmm_variable_t var, double x))
 489 {
 490   func_f_def = func_f;
 491   func_fp_def = func_fp;
 492   func_fpi_def = func_fpi;
 493 }
 494
 495 /**************** Vegas and Reno functions *************************/
 496 /* NOTE for Reno: all functions consider the network coefficient (alpha) equal to 1. */
 497
 498 /*
 499  * For Vegas: $f(x) = \alpha D_f\ln(x)$
 500  * Therefore: $fp(x) = \frac{\alpha D_f}{x}$
 501  * Therefore: $fpi(x) = \frac{\alpha D_f}{x}$
 502  */
 503 #define VEGAS_SCALING 1000.0
 504
 505 double func_vegas_f(lmm_variable_t var, double x)
 506 {
 507   xbt_assert(x > 0.0, "Don't call me with stupid values! (%1.20f)", x);
 508   return VEGAS_SCALING * var->weight * log(x);
 509 }
 510
 511 double func_vegas_fp(lmm_variable_t var, double x)
 512 {
 513   xbt_assert(x > 0.0, "Don't call me with stupid values! (%1.20f)", x);
 514   return VEGAS_SCALING * var->weight / x;
 515 }
 516
 517 double func_vegas_fpi(lmm_variable_t var, double x)
 518 {
 519   xbt_assert(x > 0.0, "Don't call me with stupid values! (%1.20f)", x);
 520   return var->weight / (x / VEGAS_SCALING);
 521 }
 522
 523 /*
 524  * For Reno:  $f(x) = \frac{\sqrt{3/2}}{D_f} atan(\sqrt{3/2}D_f x)$
 525  * Therefore: $fp(x)  = \frac{3}{3 D_f^2 x^2+2}$
 526  * Therefore: $fpi(x)  = \sqrt{\frac{1}{{D_f}^2 x} - \frac{2}{3{D_f}^2}}$
 527  */
 528 #define RENO_SCALING 1.0
 529 double func_reno_f(lmm_variable_t var, double x)
 530 {
 531   xbt_assert(var->weight > 0.0, "Don't call me with stupid values!");
 532
 533   return RENO_SCALING * sqrt(3.0 / 2.0) / var->weight * atan(sqrt(3.0 / 2.0) * var->weight * x);
 534 }
 535
 536 double func_reno_fp(lmm_variable_t var, double x)
 537 {
 538   return RENO_SCALING * 3.0 / (3.0 * var->weight * var->weight * x * x + 2.0);
 539 }
 540
 541 double func_reno_fpi(lmm_variable_t var, double x)
 542 {
 543   double res_fpi;
 544
 545   xbt_assert(var->weight > 0.0, "Don't call me with stupid values!");
 546   xbt_assert(x > 0.0, "Don't call me with stupid values!");
 547
 548   res_fpi = 1.0 / (var->weight * var->weight * (x / RENO_SCALING)) - 2.0 / (3.0 * var->weight * var->weight);
 549   if (res_fpi <= 0.0)
 550     return 0.0;
 551 /*   xbt_assert(res_fpi>0.0,"Don't call me with stupid values!"); */
 552   return sqrt(res_fpi);
 553 }
 554
 555 /* Implementing new Reno-2
 556  * For Reno-2:  $f(x)   = U_f(x_f) = \frac{{2}{D_f}}*ln(2+x*D_f)$
 557  * Therefore:   $fp(x)  = 2/(Weight*x + 2)
 558  * Therefore:   $fpi(x) = (2*Weight)/x - 4
 559  */
 560 #define RENO2_SCALING 1.0
 561 double func_reno2_f(lmm_variable_t var, double x)
 562 {
 563   xbt_assert(var->weight > 0.0, "Don't call me with stupid values!");
 564   return RENO2_SCALING * (1.0 / var->weight) * log((x * var->weight) / (2.0 * x * var->weight + 3.0));
 565 }
 566
 567 double func_reno2_fp(lmm_variable_t var, double x)
 568 {
 569   return RENO2_SCALING * 3.0 / (var->weight * x * (2.0 * var->weight * x + 3.0));
 570 }
 571
 572 double func_reno2_fpi(lmm_variable_t var, double x)
 573 {
 574   xbt_assert(x > 0.0, "Don't call me with stupid values!");
 575   double tmp = x * var->weight * var->weight;
 576   double res_fpi = tmp * (9.0 * x + 24.0);
 577
 578   if (res_fpi <= 0.0)
 579     return 0.0;
 580
 581   res_fpi = RENO2_SCALING * (-3.0 * tmp + sqrt(res_fpi)) / (4.0 * tmp);
 582   return res_fpi;
 583 }