testsuite/surf/lmm_usage.c

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