BookmarkSubscribeRSS Feed
🔒 This topic is solved and locked. Need further help from the community? Please sign in and ask a new question.
Santha
Pyrite | Level 9

Rob

I am trying to extend my model to do this p-median. 

What I have so far is model picks "x" number of facilities in the map.

I have a very finite set of pre-determined facilities (no more than 15). From this, I want to pick "x" . and (b) from this 15, I have already picked lets say 2 facilities that I want the model to use. I want the model to give me just one more, there by giving a total number of facilities opened to be 3. 

 

I went through the paper (https://support.sas.com/resources/papers/proceedings/pdfs/sgf2008/203-2008.pdf)

I saw that one of the way to do is to have the facilities as customer. I tried that but did not quite achieve that with my code. I also tried one more approach where I store the list of facilities in an array with their lat and lon and was trying to have a binary variable for the site to be open or close. But realized that the first part is NLP and so i cannot have integer variables. I am a little lost here. I am thinking to solve NLP first but I want to restrict the solve to only the list of candidate facilities and not the entire map. Then do post processing to fine tune the assignment.  But i am not sure. Can you help me with your thoughts?

1 ACCEPTED SOLUTION

Accepted Solutions
RobPratt
SAS Super FREQ

Here are sample data and code for the p-median problem.  If you want to force certain facilities to be open, you can use the FIX statement to fix the corresponding Buiild[j] variables to 1.  Note that this formulation is linear because the nonlinear distances are absorbed into the objective coefficients.

 

%let p            = 5;
%let NumCustomers = 200;
%let xmax         = 100;
%let ymax         = 100;

/* generate random customer locations */
data cdata;
   call streaminit(1);
   do i = 1 to &NumCustomers;
      x = rand('UNIF') * &xmax;
      y = rand('UNIF') * &ymax;
      weight = 1;
      output;
   end;
run;

/* plot customer locations */
proc sgplot data=cdata noautolegend aspect=1; 
   xaxis display=(nolabel); 
   yaxis display=(nolabel); 
   scatter y=y x=x;
run;

/* solve the problem by using MILP */
proc optmodel;
   set CUSTOMERS;
   /* x and y coordinates of customers */
   num x {CUSTOMERS};
   num y {CUSTOMERS};
   num weight {CUSTOMERS};
   read data cdata into CUSTOMERS=[i] x y weight;
   /* distance between customers i and j */
   num distance {i in CUSTOMERS, j in CUSTOMERS}
      = sqrt((x[i] - x[j])^2 + (y[i] - y[j])^2);

   /* declare decision variables */
   var Build {CUSTOMERS} binary;
   var Assign {CUSTOMERS, CUSTOMERS} binary;

   /* declare objective */
   min TotalDistance
      = sum {i in CUSTOMERS, j in CUSTOMERS} weight[i] * distance[i,j] * Assign[i,j];

   /* build exactly p facilities */
   con Cardinality:
      sum {j in CUSTOMERS} Build[j] = &p;
   /* each customer assigned to exactly one site */
   con AssignOnce {i in CUSTOMERS}:
      sum {j in CUSTOMERS} Assign[i,j] = 1;
   /* if customer i assigned to site j, then facility must be built at j */
   con AssignImpliesBuild {i in CUSTOMERS, j in CUSTOMERS}:
      Assign[i,j] <= Build[j];
*   con AssignImpliesBuild {j in CUSTOMERS}:
      sum {i in CUSTOMERS} Assign[i,j] <= card(CUSTOMERS) * Build[j];

   /* solve the MILP */
   solve;

   /* create a data set for use by PROC SGPLOT */
   create data solution from
      [customer site]={i in CUSTOMERS, j in CUSTOMERS: Assign[i,j].sol > 0.5}
      x1=x[i] y1=y[i] x2=x[j] y2=y[j];
quit;

/* plot optimal solution */
proc sgplot data=solution noautolegend aspect=1;
   xaxis display=(nolabel); 
   yaxis display=(nolabel); 
   scatter x=x1 y=y1;
   vector x=x2 y=y2 / xorigin=x1 yorigin=y1 noarrowheads;
run;

 

View solution in original post

1 REPLY 1
RobPratt
SAS Super FREQ

Here are sample data and code for the p-median problem.  If you want to force certain facilities to be open, you can use the FIX statement to fix the corresponding Buiild[j] variables to 1.  Note that this formulation is linear because the nonlinear distances are absorbed into the objective coefficients.

 

%let p            = 5;
%let NumCustomers = 200;
%let xmax         = 100;
%let ymax         = 100;

/* generate random customer locations */
data cdata;
   call streaminit(1);
   do i = 1 to &NumCustomers;
      x = rand('UNIF') * &xmax;
      y = rand('UNIF') * &ymax;
      weight = 1;
      output;
   end;
run;

/* plot customer locations */
proc sgplot data=cdata noautolegend aspect=1; 
   xaxis display=(nolabel); 
   yaxis display=(nolabel); 
   scatter y=y x=x;
run;

/* solve the problem by using MILP */
proc optmodel;
   set CUSTOMERS;
   /* x and y coordinates of customers */
   num x {CUSTOMERS};
   num y {CUSTOMERS};
   num weight {CUSTOMERS};
   read data cdata into CUSTOMERS=[i] x y weight;
   /* distance between customers i and j */
   num distance {i in CUSTOMERS, j in CUSTOMERS}
      = sqrt((x[i] - x[j])^2 + (y[i] - y[j])^2);

   /* declare decision variables */
   var Build {CUSTOMERS} binary;
   var Assign {CUSTOMERS, CUSTOMERS} binary;

   /* declare objective */
   min TotalDistance
      = sum {i in CUSTOMERS, j in CUSTOMERS} weight[i] * distance[i,j] * Assign[i,j];

   /* build exactly p facilities */
   con Cardinality:
      sum {j in CUSTOMERS} Build[j] = &p;
   /* each customer assigned to exactly one site */
   con AssignOnce {i in CUSTOMERS}:
      sum {j in CUSTOMERS} Assign[i,j] = 1;
   /* if customer i assigned to site j, then facility must be built at j */
   con AssignImpliesBuild {i in CUSTOMERS, j in CUSTOMERS}:
      Assign[i,j] <= Build[j];
*   con AssignImpliesBuild {j in CUSTOMERS}:
      sum {i in CUSTOMERS} Assign[i,j] <= card(CUSTOMERS) * Build[j];

   /* solve the MILP */
   solve;

   /* create a data set for use by PROC SGPLOT */
   create data solution from
      [customer site]={i in CUSTOMERS, j in CUSTOMERS: Assign[i,j].sol > 0.5}
      x1=x[i] y1=y[i] x2=x[j] y2=y[j];
quit;

/* plot optimal solution */
proc sgplot data=solution noautolegend aspect=1;
   xaxis display=(nolabel); 
   yaxis display=(nolabel); 
   scatter x=x1 y=y1;
   vector x=x2 y=y2 / xorigin=x1 yorigin=y1 noarrowheads;
run;

 

Ready to join fellow brilliant minds for the SAS Hackathon?

Build your skills. Make connections. Enjoy creative freedom. Maybe change the world. Registration is now open through August 30th. Visit the SAS Hackathon homepage.

Register today!
Multiple Linear Regression in SAS

Learn how to run multiple linear regression models with and without interactions, presented by SAS user Alex Chaplin.

Find more tutorials on the SAS Users YouTube channel.

Discussion stats
  • 1 reply
  • 650 views
  • 0 likes
  • 2 in conversation