teshsuite/surf/lmm_usage/lmm_usage.c

   1 /* A few tests for the maxmin library                                       */
   2
   3 /* Copyright (c) 2007-2015. The SimGrid Team.
   4  * All rights reserved.                                                     */
   5
   6 /* This program is free software; you can redistribute it and/or modify it
   7  * under the terms of the license (GNU LGPL) which comes with this package. */
   8
   9 #include "xbt/sysdep.h"
  10 #include "surf/maxmin.h"
  11 #include "xbt/log.h"
  12 #include "xbt/module.h"
  13 #include <math.h>
  14
  15 XBT_LOG_NEW_DEFAULT_CATEGORY(surf_test,
  16                              "Messages specific for surf example");
  17
  18 #define PRINT_VAR(var) XBT_DEBUG(#var " = %g",lmm_variable_getvalue(var));
  19 #define SHOW_EXPR(expr) XBT_DEBUG(#expr " = %g",expr);
  20
  21 /*                               */
  22 /*        ______                 */
  23 /*  ==l1==  L2  ==L3==           */
  24 /*        ------                 */
  25 /*                               */
  26
  27 typedef enum {
  28   MAXMIN,
  29   LAGRANGE_RENO,
  30   LAGRANGE_VEGAS
  31 } method_t;
  32
  33 static double dichotomy(double func(double), double min, double max,
  34                         double min_error)
  35 {
  36   double middle;
  37   double min_func, max_func, middle_func;
  38
  39   double overall_error = 2 * min_error;
  40
  41   min_func = func(min);
  42   max_func = func(max);
  43
  44   if ((min_func > 0 && max_func > 0))
  45     return min - 1.0;
  46   if ((min_func < 0 && max_func < 0))
  47     return max + 1.0;
  48   if ((min_func > 0 && max_func < 0))
  49     abort();
  50
  51   SHOW_EXPR(min_error);
  52
  53   while (overall_error > min_error) {
  54     SHOW_EXPR(overall_error);
  55     if ((min_func > 0 && max_func > 0) ||
  56         (min_func < 0 && max_func < 0) || (min_func > 0 && max_func < 0)) {
  57       abort();
  58     }
  59
  60     SHOW_EXPR(min);
  61     SHOW_EXPR(min_func);
  62     SHOW_EXPR(max);
  63     SHOW_EXPR(max_func);
  64
  65     middle = (max + min) / 2.0;
  66     if ((min == middle) || (max == middle)) {
  67       break;
  68     }
  69     middle_func = func(middle);
  70     SHOW_EXPR(middle);
  71     SHOW_EXPR(middle_func);
  72
  73     if (middle_func < 0) {
  74       min = middle;
  75       min_func = middle_func;
  76       overall_error = max_func - middle_func;
  77     } else if (middle_func > 0) {
  78       max = middle;
  79       max_func = middle_func;
  80       overall_error = middle_func - min_func;
  81     } else {
  82       overall_error = 0;
  83     }
  84   }
  85   return ((min + max) / 2.0);
  86 }
  87
  88 double a_test_1 = 0;
  89 double b_test_1 = 0;
  90 static double diff_lagrange_test_1(double x)
  91 {
  92   return -(3 / (1 + 3 * x * x / 2) -
  93            3 / (2 * (3 * (a_test_1 - x) * (a_test_1 - x) / 2 + 1)) +
  94            3 / (2 *
  95                 (3 * (b_test_1 - a_test_1 + x) *
  96                  (b_test_1 - a_test_1 + x) / 2 + 1)));
  97 }
  98
  99 void test1(method_t method);
 100 void test1(method_t method)
 101 {
 102   lmm_system_t Sys = NULL;
 103   lmm_constraint_t L1 = NULL;
 104   lmm_constraint_t L2 = NULL;
 105   lmm_constraint_t L3 = NULL;
 106
 107   lmm_variable_t R_1_2_3 = NULL;
 108   lmm_variable_t R_1 = NULL;
 109   lmm_variable_t R_2 = NULL;
 110   lmm_variable_t R_3 = NULL;
 111
 112   double a = 1.0, b = 10.0;
 113
 114   if (method == LAGRANGE_VEGAS)
 115     lmm_set_default_protocol_function(func_vegas_f, func_vegas_fp,
 116                                       func_vegas_fpi);
 117   else if (method == LAGRANGE_RENO)
 118     lmm_set_default_protocol_function(func_reno_f, func_reno_fpi,
 119                                       func_reno_fpi);
 120
 121   Sys = lmm_system_new(1);
 122   L1 = lmm_constraint_new(Sys, (void *) "L1", a);
 123   L2 = lmm_constraint_new(Sys, (void *) "L2", b);
 124   L3 = lmm_constraint_new(Sys, (void *) "L3", a);
 125
 126   R_1_2_3 = lmm_variable_new(Sys, (void *) "R 1->2->3", 1.0, -1.0, 3);
 127   R_1 = lmm_variable_new(Sys, (void *) "R 1", 1.0, -1.0, 1);
 128   R_2 = lmm_variable_new(Sys, (void *) "R 2", 1.0, -1.0, 1);
 129   R_3 = lmm_variable_new(Sys, (void *) "R 3", 1.0, -1.0, 1);
 130
 131   lmm_update_variable_weight(Sys, R_1_2_3, 1.0);
 132   lmm_update_variable_weight(Sys, R_1, 1.0);
 133   lmm_update_variable_weight(Sys, R_2, 1.0);
 134   lmm_update_variable_weight(Sys, R_3, 1.0);
 135
 136   lmm_expand(Sys, L1, R_1_2_3, 1.0);
 137   lmm_expand(Sys, L2, R_1_2_3, 1.0);
 138   lmm_expand(Sys, L3, R_1_2_3, 1.0);
 139
 140   lmm_expand(Sys, L1, R_1, 1.0);
 141
 142   lmm_expand(Sys, L2, R_2, 1.0);
 143
 144   lmm_expand(Sys, L3, R_3, 1.0);
 145
 146
 147   if (method == MAXMIN) {
 148     lmm_solve(Sys);
 149   } else if (method == LAGRANGE_VEGAS) {
 150     double x = 3 * a / 4 - 3 * b / 8 +
 151         sqrt(9 * b * b + 4 * a * a - 4 * a * b) / 8;
 152     /* Computed with mupad and D_f=1.0 */
 153     double max_deviation = 0.0;
 154     if (x > a) {
 155       x = a;
 156     }
 157     if (x < 0) {
 158       x = 0;
 159     }
 160
 161     lagrange_solve(Sys);
 162
 163     max_deviation =
 164         MAX(max_deviation, fabs(lmm_variable_getvalue(R_1) - x));
 165     max_deviation =
 166         MAX(max_deviation, fabs(lmm_variable_getvalue(R_3) - x));
 167     max_deviation =
 168         MAX(max_deviation, fabs(lmm_variable_getvalue(R_2) - (b - a + x)));
 169     max_deviation =
 170         MAX(max_deviation, fabs(lmm_variable_getvalue(R_1_2_3) - (a - x)));
 171
 172     if (max_deviation > 0.00001) { // Legacy value used in lagrange.c
 173       XBT_WARN("Max Deviation from optimal solution : %g", max_deviation);
 174       XBT_WARN("Found x = %1.20f", x);
 175       XBT_WARN("Deviation from optimal solution (R_1 = %g): %1.20f", x,
 176             lmm_variable_getvalue(R_1) - x);
 177       XBT_WARN("Deviation from optimal solution (R_2 = %g): %1.20f",
 178             b - a + x, lmm_variable_getvalue(R_2) - (b - a + x));
 179       XBT_WARN("Deviation from optimal solution (R_3 = %g): %1.20f", x,
 180             lmm_variable_getvalue(R_3) - x);
 181       XBT_WARN("Deviation from optimal solution (R_1_2_3 = %g): %1.20f",
 182             a - x, lmm_variable_getvalue(R_1_2_3) - (a - x));
 183     }
 184   } else if (method == LAGRANGE_RENO) {
 185     double x;
 186     double max_deviation = 0.0;
 187
 188     a_test_1 = a;
 189     b_test_1 = b;
 190     x = dichotomy(diff_lagrange_test_1, 0, a, 1e-13);
 191
 192     if (x < 0)
 193       x = 0;
 194     if (x > a)
 195       x = a;
 196     lagrange_solve(Sys);
 197
 198     max_deviation =
 199         MAX(max_deviation, fabs(lmm_variable_getvalue(R_1) - x));
 200     max_deviation =
 201         MAX(max_deviation, fabs(lmm_variable_getvalue(R_3) - x));
 202     max_deviation =
 203         MAX(max_deviation, fabs(lmm_variable_getvalue(R_2) - (b - a + x)));
 204     max_deviation =
 205         MAX(max_deviation, fabs(lmm_variable_getvalue(R_1_2_3) - (a - x)));
 206
 207     if (max_deviation > 0.00001) { // Legacy value used in lagrange.c
 208       XBT_WARN("Max Deviation from optimal solution : %g", max_deviation);
 209       XBT_WARN("Found x = %1.20f", x);
 210       XBT_WARN("Deviation from optimal solution (R_1 = %g): %1.20f", x,
 211             lmm_variable_getvalue(R_1) - x);
 212       XBT_WARN("Deviation from optimal solution (R_2 = %g): %1.20f",
 213             b - a + x, lmm_variable_getvalue(R_2) - (b - a + x));
 214       XBT_WARN("Deviation from optimal solution (R_3 = %g): %1.20f", x,
 215             lmm_variable_getvalue(R_3) - x);
 216       XBT_WARN("Deviation from optimal solution (R_1_2_3 = %g): %1.20f",
 217             a - x, lmm_variable_getvalue(R_1_2_3) - (a - x));
 218     }
 219   } else {
 220     xbt_die( "Invalid method");
 221   }
 222
 223   PRINT_VAR(R_1_2_3);
 224   PRINT_VAR(R_1);
 225   PRINT_VAR(R_2);
 226   PRINT_VAR(R_3);
 227
 228   lmm_variable_free(Sys, R_1_2_3);
 229   lmm_variable_free(Sys, R_1);
 230   lmm_variable_free(Sys, R_2);
 231   lmm_variable_free(Sys, R_3);
 232   lmm_system_free(Sys);
 233 }
 234
 235 void test2(method_t method);
 236 void test2(method_t method)
 237 {
 238   lmm_system_t Sys = NULL;
 239   lmm_constraint_t CPU1 = NULL;
 240   lmm_constraint_t CPU2 = NULL;
 241
 242   lmm_variable_t T1 = NULL;
 243   lmm_variable_t T2 = NULL;
 244
 245
 246   if (method == LAGRANGE_VEGAS)
 247     lmm_set_default_protocol_function(func_vegas_f, func_vegas_fp,
 248                                       func_vegas_fpi);
 249   else if (method == LAGRANGE_RENO)
 250     lmm_set_default_protocol_function(func_reno_f, func_reno_fp,
 251                                       func_reno_fpi);
 252
 253   Sys = lmm_system_new(1);
 254   CPU1 = lmm_constraint_new(Sys, (void *) "CPU1", 200.0);
 255   CPU2 = lmm_constraint_new(Sys, (void *) "CPU2", 100.0);
 256
 257   T1 = lmm_variable_new(Sys, (void *) "T1", 1.0, -1.0, 1);
 258   T2 = lmm_variable_new(Sys, (void *) "T2", 1.0, -1.0, 1);
 259
 260   lmm_update_variable_weight(Sys, T1, 1.0);
 261   lmm_update_variable_weight(Sys, T2, 1.0);
 262
 263
 264   lmm_expand(Sys, CPU1, T1, 1.0);
 265   lmm_expand(Sys, CPU2, T2, 1.0);
 266
 267
 268
 269   if (method == MAXMIN) {
 270     lmm_solve(Sys);
 271   } else if (method == LAGRANGE_VEGAS) {
 272     lagrange_solve(Sys);
 273   } else if (method == LAGRANGE_RENO) {
 274     lagrange_solve(Sys);
 275   } else {
 276     xbt_die("Invalid method");
 277   }
 278
 279   PRINT_VAR(T1);
 280   PRINT_VAR(T2);
 281
 282   lmm_variable_free(Sys, T1);
 283   lmm_variable_free(Sys, T2);
 284   lmm_system_free(Sys);
 285 }
 286
 287
 288
 289 void test3(method_t method);
 290 void test3(method_t method)
 291 {
 292   int flows = 11;
 293   int links = 10;
 294
 295   int i = 0;
 296   int j = 0;
 297
 298   double **A;
 299
 300   lmm_system_t Sys = NULL;
 301   lmm_constraint_t *tmp_cnst = NULL;
 302   lmm_variable_t *tmp_var = NULL;
 303   char **tmp_name;
 304
 305
 306   /*array to add the the constraints of fictiv variables */
 307   double B[15] = { 10, 10, 10, 10, 10, 10, 10, 10, 10, 10,
 308
 309     1, 1, 1, 1, 1
 310   };
 311
 312   A = xbt_new0(double *, links + 5);
 313
 314   for (i = 0; i < links + 5; i++) {
 315     A[i] = xbt_new0(double, flows + 5);
 316     for (j = 0; j < flows + 5; j++) {
 317       A[i][j] = 0.0;
 318
 319       if (i >= links || j >= flows) {
 320         A[i][j] = 0.0;
 321       }
 322     }
 323   }
 324
 325   /*matrix that store the constraints/topollogy */
 326   /*double A[15][16]=
 327      {{0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0,    0, 0, 0, 0, 0},
 328      {0, 1, 0, 0, 0, 0, 0, 1, 1, 0, 0,    0, 0, 0, 0, 0},
 329      {0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0,    0, 0, 0, 0, 0},
 330      {0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0,    0, 0, 0, 0, 0},
 331      {1, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0,    0, 0, 0, 0, 0},
 332      {1, 0, 0, 1, 1, 0, 0, 0, 0, 1, 0,    0, 0, 0, 0, 0},
 333      {1, 0, 0, 0, 1, 0, 0, 0, 0, 1, 1,    0, 0, 0, 0, 0},
 334      {0, 0, 1, 0, 1, 0, 1, 0, 0, 1, 1,    0, 0, 0, 0, 0},
 335      {0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1,    0, 0, 0, 0, 0},
 336      {0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 0,    0, 0, 0, 0, 0},
 337
 338      {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,    1, 0, 0, 0, 0},
 339      {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,    0, 1, 0, 0, 0},
 340      {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,    0, 0, 1, 0, 0},
 341      {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,    0, 0, 0, 1, 0},
 342      {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,    0, 0, 0, 0, 1}
 343      }; */
 344
 345   A[0][1] = 1.0;
 346   A[0][7] = 1.0;
 347
 348   A[1][1] = 1.0;
 349   A[1][7] = 1.0;
 350   A[1][8] = 1.0;
 351
 352   A[2][1] = 1.0;
 353   A[2][8] = 1.0;
 354
 355   A[2][1] = 1.0;
 356   A[2][8] = 1.0;
 357
 358   A[3][8] = 1.0;
 359
 360   A[4][0] = 1.0;
 361   A[4][3] = 1.0;
 362   A[4][9] = 1.0;
 363
 364   A[5][0] = 1.0;
 365   A[5][3] = 1.0;
 366   A[5][4] = 1.0;
 367   A[5][9] = 1.0;
 368
 369   A[6][0] = 1.0;
 370   A[6][4] = 1.0;
 371   A[6][9] = 1.0;
 372   A[6][10] = 1.0;
 373
 374   A[7][2] = 1.0;
 375   A[7][4] = 1.0;
 376   A[7][6] = 1.0;
 377   A[7][9] = 1.0;
 378   A[7][10] = 1.0;
 379
 380   A[8][2] = 1.0;
 381   A[8][10] = 1.0;
 382
 383   A[9][5] = 1.0;
 384   A[9][6] = 1.0;
 385   A[9][9] = 1.0;
 386
 387
 388   A[10][11] = 1.0;
 389   A[11][12] = 1.0;
 390   A[12][13] = 1.0;
 391   A[13][14] = 1.0;
 392   A[14][15] = 1.0;
 393
 394
 395   if (method == LAGRANGE_VEGAS)
 396     lmm_set_default_protocol_function(func_vegas_f, func_vegas_fp,
 397                                       func_vegas_fpi);
 398   else if (method == LAGRANGE_RENO)
 399     lmm_set_default_protocol_function(func_reno_f, func_reno_fp,
 400                                       func_reno_fpi);
 401
 402   Sys = lmm_system_new(1);
 403
 404
 405
 406   tmp_name = xbt_new0(char *, 31);
 407
 408   /*
 409    * Creates the constraints
 410    */
 411   tmp_cnst = xbt_new0(lmm_constraint_t, 15);
 412   for (i = 0; i < 15; i++) {
 413     tmp_name[i] = bprintf("C_%03d", i);
 414     tmp_cnst[i] = lmm_constraint_new(Sys, (void *) tmp_name[i], B[i]);
 415   }
 416
 417
 418   /*
 419    * Creates the variables
 420    */
 421   tmp_var = xbt_new0(lmm_variable_t, 16);
 422   for (j = 0; j < 16; j++) {
 423     tmp_name[i + j] = bprintf("X_%03d", j);
 424     tmp_var[j] =
 425         lmm_variable_new(Sys, (void *) tmp_name[i + j], 1.0, -1.0, 15);
 426     lmm_update_variable_weight(Sys, tmp_var[j], 1.0);
 427   }
 428
 429   /*
 430    * Link constraints and variables
 431    */
 432   for (i = 0; i < 15; i++) {
 433     for (j = 0; j < 16; j++) {
 434       if (A[i][j]) {
 435         lmm_expand(Sys, tmp_cnst[i], tmp_var[j], 1.0);
 436       }
 437     }
 438   }
 439
 440
 441
 442   if (method == MAXMIN) {
 443     lmm_solve(Sys);
 444   } else if (method == LAGRANGE_VEGAS) {
 445     lagrange_solve(Sys);
 446   } else if (method == LAGRANGE_RENO) {
 447     lagrange_solve(Sys);
 448   } else {
 449     xbt_die("Invalid method");
 450   }
 451
 452   for (j = 0; j < 16; j++) {
 453     PRINT_VAR(tmp_var[j]);
 454   }
 455
 456   for (j = 0; j < 16; j++)
 457     lmm_variable_free(Sys, tmp_var[j]);
 458   xbt_free(tmp_var);
 459   xbt_free(tmp_cnst);
 460   for (i = 0; i < 31; i++)
 461     xbt_free(tmp_name[i]);
 462   xbt_free(tmp_name);
 463   lmm_system_free(Sys);
 464   for (i = 0; i < links + 5; i++)
 465     xbt_free(A[i]);
 466   xbt_free(A);
 467 }
 468
 469 int main(int argc, char **argv)
 470 {
 471   xbt_init(&argc, argv);
 472
 473   XBT_INFO("***** Test 1 (Max-Min)");
 474   test1(MAXMIN);
 475   XBT_INFO("***** Test 1 (Lagrange - Vegas)");
 476   test1(LAGRANGE_VEGAS);
 477   XBT_INFO("***** Test 1 (Lagrange - Reno)");
 478   test1(LAGRANGE_RENO);
 479
 480
 481
 482   XBT_INFO("***** Test 2 (Max-Min)");
 483   test2(MAXMIN);
 484   XBT_INFO("***** Test 2 (Lagrange - Vegas)");
 485   test2(LAGRANGE_VEGAS);
 486   XBT_INFO("***** Test 2 (Lagrange - Reno)");
 487   test2(LAGRANGE_RENO);
 488
 489
 490   XBT_INFO("***** Test 3 (Max-Min)");
 491   test3(MAXMIN);
 492   XBT_INFO("***** Test 3 (Lagrange - Vegas)");
 493   test3(LAGRANGE_VEGAS);
 494   XBT_INFO("***** Test 3 (Lagrange - Reno)");
 495   test3(LAGRANGE_RENO);
 496
 497   return 0;
 498 }