src/surf/lagrange.c

   1 /*      $Id$     */
   2
   3 /* Copyright (c) 2007 Arnaud Legrand, Pedro Velho. All rights reserved.     */
   4
   5 /* This program is free software; you can redistribute it and/or modify it
   6  * under the terms of the license (GNU LGPL) which comes with this package. */
   7
   8 /*
   9  * Modelling the proportional fairness using the Lagrange Optimization
  10  * Approach. For a detailed description see:
  11  * "ssh://username@scm.gforge.inria.fr/svn/memo/people/pvelho/lagrange/ppf.ps".
  12  */
  13 #include "xbt/log.h"
  14 #include "xbt/sysdep.h"
  15 #include "xbt/mallocator.h"
  16 #include "maxmin_private.h"
  17
  18 #include <stdlib.h>
  19 #ifndef MATH
  20 #include <math.h>
  21 #endif
  22
  23
  24 XBT_LOG_NEW_DEFAULT_SUBCATEGORY(surf_lagrange, surf,
  25                                 "Logging specific to SURF (lagrange)");
  26
  27 XBT_LOG_NEW_SUBCATEGORY(surf_writelambda, surf,
  28                                 "Generates the lambda.in file. WARNING: the size of this file might be a few GBs.");
  29
  30 void lagrange_solve(lmm_system_t sys);
  31
  32 void lagrange_solve(lmm_system_t sys)
  33 {
  34
  35   /*
  36    * Lagrange Variables.
  37    */
  38   double epsilon_min_error = 1e-6;
  39   double overall_error = 1;
  40   double sigma_step = 1e-3;
  41   double capacity_error=0, bound_error=0;
  42
  43
  44   /*
  45    * Variables to manipulate the data structure proposed to model the maxmin
  46    * fairness. See docummentation for more details.
  47    */
  48   xbt_swag_t elem_list = NULL;
  49   lmm_element_t elem1 = NULL;
  50   lmm_element_t elem = NULL;
  51
  52   xbt_swag_t cnst_list = NULL;
  53   lmm_constraint_t cnst1 = NULL;
  54   lmm_constraint_t cnst2 = NULL;
  55   lmm_constraint_t cnst = NULL;
  56   double sum;
  57   xbt_swag_t var_list = NULL;
  58   lmm_variable_t var1 = NULL;
  59   lmm_variable_t var = NULL;
  60   lmm_variable_t var2 = NULL;
  61
  62
  63   /*
  64    * Auxiliar variables.
  65    */
  66   int iteration=0;
  67   int max_iterations= 1000;
  68   double mu_partial=0;
  69   double lambda_partial=0;
  70   double tmp=0;
  71   int i,j;
  72   FILE *gnuplot_file=NULL;
  73   char print_buf[1024];
  74   char *trace_buf=xbt_malloc0(sizeof(char));
  75
  76
  77
  78
  79
  80
  81   if ( !(sys->modified))
  82     return;
  83
  84   /*
  85    * Initialize the var list variable with only the active variables.
  86    * Associate an index in the swag variables. Initialize mu.
  87    */
  88   var_list = &(sys->variable_set);
  89   i=0;
  90   xbt_swag_foreach(var1, var_list) {
  91     if((var1->bound > 0.0) || (var1->weight <= 0.0)){
  92       DEBUG1("## NOTE var1(%d) is a boundless variable", i);
  93       var1->mu = -1.0;
  94     } else
  95       var1->mu = 1.0;
  96     DEBUG2("## var1(%d)->mu:  %e", i, var1->mu);
  97     DEBUG2("## var1(%d)->weight: %e", i, var1->weight);
  98     i++;
  99   }
 100
 101   /*
 102    * Initialize lambda.
 103    */
 104   cnst_list=&(sys->active_constraint_set);
 105   xbt_swag_foreach(cnst1, cnst_list) {
 106     cnst1->lambda = 1.0;
 107     DEBUG2("## cnst1(%p)->lambda:  %e", cnst1, cnst1->lambda);
 108   }
 109
 110   if(XBT_LOG_ISENABLED(surf_writelambda, xbt_log_priority_debug)) {
 111     gnuplot_file = fopen("lambda.in", "w");
 112     fprintf(gnuplot_file, "# iteration    lambda1  lambda2 lambda3 ... lambdaP");
 113   }
 114
 115
 116   /*
 117    * While doesn't reach a minimun error or a number maximum of iterations.
 118    */
 119   while(overall_error > epsilon_min_error && iteration < max_iterations){
 120     iteration++;
 121     /*                        d Dual
 122      * Compute the value of ----------- (\lambda^k, \mu^k) this portion
 123      *                       d \mu_i^k
 124      * of code depends on function f(x).
 125      */
 126     var_list = &(sys->variable_set);
 127     xbt_swag_foreach(var1, var_list) {
 128       mu_partial = 0;
 129       if((var1->bound > 0) || (var1->weight <=0) ){
 130         //for each link with capacity cnsts[i] that uses flow of variable var1 do
 131         for(i=0; i<var1->cnsts_number; i++)
 132           mu_partial += (var1->cnsts[i].constraint)->lambda;
 133
 134         mu_partial = -1.0 / mu_partial + var1->bound;
 135         var1->new_mu = var1->mu - sigma_step * mu_partial;
 136         /* Assert that var1->new_mu is positive */
 137      }
 138     }
 139
 140     if(XBT_LOG_ISENABLED(surf_writelambda, xbt_log_priority_debug)) {
 141       fprintf(gnuplot_file, "\n%d ", iteration);
 142     }
 143
 144     /*                         d Dual
 145      * Compute the value of ------------- (\lambda^k, \mu^k) this portion
 146      *                      d \lambda_i^k
 147      * of code depends on function f(x).
 148      */
 149
 150     DEBUG1("######Lambda partial at iteration  %d", iteration);
 151     cnst_list=&(sys->active_constraint_set);
 152     j=0;
 153     xbt_swag_foreach(cnst1, cnst_list) {
 154       j++;
 155
 156       lambda_partial = 0;
 157
 158       elem_list = &(cnst1->element_set);
 159       xbt_swag_foreach(elem1, elem_list) {
 160         lambda_partial = 0;
 161
 162         var2 = elem1->variable;
 163
 164         if(var2->weight<=0) continue;
 165
 166         tmp = 0;
 167
 168         //for each link with capacity cnsts[i] that uses flow of variable var1 do
 169         if(var2->bound > 0)
 170           tmp += var2->mu;
 171
 172         for(i=0; i<var2->cnsts_number; i++)
 173           tmp += (var2->cnsts[i].constraint)->lambda;
 174
 175         lambda_partial += -1 / tmp;
 176       }
 177
 178       lambda_partial += cnst1->bound;
 179
 180       DEBUG2("###########Lambda partial %p : %e", cnst1, lambda_partial);
 181
 182       cnst1->new_lambda = cnst1->lambda - sigma_step * lambda_partial;
 183
 184       if(XBT_LOG_ISENABLED(surf_writelambda, xbt_log_priority_debug)) {
 185         fprintf(gnuplot_file, " %f", cnst1->lambda);
 186       }
 187     }
 188
 189     /* Updating lambda's and mu's */
 190     var_list = &(sys->variable_set);
 191     xbt_swag_foreach(var1, var_list)
 192       if(!((var1->bound > 0.0) || (var1->weight <= 0.0)))
 193         var1->mu = var1->new_mu;
 194
 195
 196     cnst_list=&(sys->active_constraint_set);
 197     xbt_swag_foreach(cnst1, cnst_list)
 198       cnst1->lambda = cnst1->new_lambda;
 199
 200     /*
 201      * Now computes the values of each variable (\rho) based on
 202      * the values of \lambda and \mu.
 203      */
 204     var_list = &(sys->variable_set);
 205     xbt_swag_foreach(var1, var_list) {
 206       if(var1->weight <=0)
 207         var1->value = 0.0;
 208       else {
 209         tmp = 0;
 210         if(var1->bound >0)
 211           tmp+=var1->mu;
 212         for(i=0; i<var1->cnsts_number; i++)
 213           tmp += (var1->cnsts[i].constraint)->lambda;
 214
 215         var1->value = 1 / tmp;
 216       }
 217
 218
 219       DEBUG2("var1->value (id=%s) : %e", (char *)var1->id, var1->value);
 220     }
 221
 222   /* Printing Objective */
 223   var_list = &(sys->variable_set);
 224   sprintf(print_buf,"MAX-MIN ( ");
 225   trace_buf = xbt_realloc(trace_buf,strlen(trace_buf)+strlen(print_buf)+1);
 226   strcat(trace_buf, print_buf);
 227   xbt_swag_foreach(var, var_list) {
 228     sprintf(print_buf,"'%p'(%f) ",var,var->weight);
 229     trace_buf = xbt_realloc(trace_buf,strlen(trace_buf)+strlen(print_buf)+1);
 230     strcat(trace_buf, print_buf);
 231   }
 232   sprintf(print_buf,")");
 233   trace_buf = xbt_realloc(trace_buf,strlen(trace_buf)+strlen(print_buf)+1);
 234   strcat(trace_buf, print_buf);
 235   DEBUG1("%s",trace_buf);
 236   trace_buf[0]='\000';
 237
 238   /* Printing Constraints */
 239   cnst_list = &(sys->active_constraint_set);
 240   xbt_swag_foreach(cnst, cnst_list) {
 241     sum=0.0;
 242     elem_list = &(cnst->element_set);
 243     sprintf(print_buf,"\t");
 244     trace_buf = xbt_realloc(trace_buf,strlen(trace_buf)+strlen(print_buf)+1);
 245     strcat(trace_buf, print_buf);
 246     xbt_swag_foreach(elem, elem_list) {
 247       sprintf(print_buf,"%f.'%p'(%f) + ",elem->value,
 248               elem->variable,elem->variable->value);
 249       trace_buf = xbt_realloc(trace_buf,strlen(trace_buf)+strlen(print_buf)+1);
 250       strcat(trace_buf, print_buf);
 251       sum += elem->value * elem->variable->value;
 252     }
 253     sprintf(print_buf,"0 <= %f ('%p')",cnst->bound,cnst);
 254     trace_buf = xbt_realloc(trace_buf,strlen(trace_buf)+strlen(print_buf)+1);
 255     strcat(trace_buf, print_buf);
 256
 257     if(!cnst->shared) {
 258       sprintf(print_buf," [MAX-Constraint]");
 259       trace_buf = xbt_realloc(trace_buf,strlen(trace_buf)+strlen(print_buf)+1);
 260       strcat(trace_buf, print_buf);
 261     }
 262     DEBUG1("%s",trace_buf);
 263     trace_buf[0]='\000';
 264     if(!(sum<=cnst->bound))
 265       DEBUG3("Incorrect value (%f is not smaller than %f): %g",
 266                 sum,cnst->bound,sum-cnst->bound);
 267   }
 268
 269   /* Printing Result */
 270   xbt_swag_foreach(var, var_list) {
 271     if(var->bound>0) {
 272       DEBUG4("'%p'(%f) : %f (<=%f)",var,var->weight,var->value, var->bound);
 273       if(var->value<=var->bound)
 274         DEBUG0("Incorrect value");
 275     }
 276     else
 277       DEBUG3("'%p'(%f) : %f",var,var->weight,var->value);
 278   }
 279
 280
 281     /*
 282      * Verify for each capacity constraint (lambda) the error associated.
 283      */
 284     capacity_error = 0;
 285     cnst_list=&(sys->active_constraint_set);
 286     xbt_swag_foreach(cnst1, cnst_list) {
 287       cnst2 = xbt_swag_getNext(cnst1,(var_list)->offset);
 288       if(cnst2 != NULL){
 289         capacity_error += fabs( cnst1->lambda - cnst2->lambda );
 290       }
 291     }
 292
 293     /*
 294      * Verify for each variable the error of round trip time constraint (mu).
 295      */
 296     bound_error = 0;
 297     var_list = &(sys->variable_set);
 298     xbt_swag_foreach(var1, var_list) {
 299       var2 = xbt_swag_getNext(var1,(var_list)->offset);
 300       if(var2 != NULL){
 301         bound_error += fabs( var2->mu - var1->mu);
 302       }
 303     }
 304
 305     overall_error = capacity_error + bound_error;
 306   }
 307
 308
 309
 310
 311   if(overall_error <= epsilon_min_error){
 312     DEBUG1("The method converge in %d iterations.", iteration);
 313   }else{
 314     WARN1("Method reach %d iterations, which is the maxmimun number of iterations allowed.", iteration);
 315   }
 316
 317
 318   if(XBT_LOG_ISENABLED(surf_writelambda, xbt_log_priority_debug)) {
 319     fclose(gnuplot_file);
 320   }
 321
 322
 323   /*
 324    * Now computes the values of each variable (\rho) based on
 325    * the values of \lambda and \mu.
 326    */
 327   var_list = &(sys->variable_set);
 328   xbt_swag_foreach(var1, var_list) {
 329     tmp = 0;
 330     for(i=0; i<var1->cnsts_number; i++){
 331       elem1 = &(var1->cnsts[i]);
 332       tmp += (elem1->constraint)->lambda + var1->mu;
 333     }
 334     var1->weight = 1 / tmp;
 335
 336     DEBUG2("var1->weight (id=%s) : %e", (char *)var1->id, var1->weight);
 337   }
 338
 339
 340
 341
 342 }