- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
Rob
Here is my code that I followed from SPARSE modeling for MILP. I have 9 facilities to choose from. I have 1,056 customers. I want the CUSTOMER_FACILITIES to be the place where I can set my radius. Also, I want to be able to limit the number of facilities to "p" that i set initially. I have a constraint (con Have_This_Many_FAC_OPEN)
for that purpose . Is that correct? Right now I have this commented. When I ran my model with distance <=10,000 the model ran. However, i do see that each customer is assigned to all 9 different facilities, which means that the constraint "assign_def" is not set up correctly and/OR the CUSTOMER_FACILITIES is not set up right. If I reduce my radius to 800 for example because technically any customer can be served from a facility within 800 miles within continental US. But it results in in feasibility. Not sure what I am doing wrong. The SiteCapacity for all facility is set to be a very big number so that right now there is no in feasibility because of this. I can change sitecapacity later.
%let p=4; /*Set how many Facilities you want to be OPENED*/proc optmodel;
set DIMS=1..2;
set CUSTOMERS;
set FACILITIES;
num a {CUSTOMERS,DIMS};
num demand {CUSTOMERS};
num f {FACILITIES,DIMS};
num SiteCapacity {FACILITIES};
read data STDOPT.COGcustomers into CUSTOMERS=
[_N_] {d in DIMS} <a[_N_,d]=col('a'||d)> demand;
read data STDOPT.cogExistingFacilities into FACILITIES=
[_N_] {d in DIMS} <f[_N_,d]=col('f'||d)> SiteCapacity;
/* distance from customer i to facility j */
num dist {i in CUSTOMERS, j in FACILITIES}
= GEODIST(a[i,1],a[i,2],f[j,1],f[j,2]);
set CUSTOMERS_FACILITIES = {i in CUSTOMERS, j in FACILITIES: dist[i,j] <= 10000};
var Assign {CUSTOMERS_FACILITIES} binary;
var Build {FACILITIES} binary;
min Cost
= sum {<i,j> in CUSTOMERS_FACILITIES} dist[i,j] * Assign[i,j];
/* each customer assigned to exactly one site */
con assign_def {i in CUSTOMERS}:
sum {<(i),j> in CUSTOMERS_FACILITIES} Assign[i,j] = 1;
/* if customer i assigned to site j, then facility must be built at j */
con link {<i,j> in CUSTOMERS_FACILITIES}:
Assign[i,j] <= Build[j];
/* each site can handle at most Site Capacity demand */
con capacity {j in FACILITIES}:
sum {<i,(j)> in CUSTOMERS_FACILITIES} demand[i] * Assign[i,j] <= SiteCapacity[j] * Build[j];
**con Have_This_Many_FAC_OPEN: Sum{j in FACILITIES} Build[j]=p;
solve with milp/timetype=real;
Accepted Solutions
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
Note that GEODIST defaults to kilometers. If you want miles, include 'M' as a fifth argument.
When I run with distance threshold 10000, I get one facility per customer, as expected from the assign_def constraints.
When I run with distance threshold 800, I get infeasible, as you did. These notes in the log indicate several customers that are not within 800 miles of any facility:
NOTE: The constraint 'assign_def[478]' is empty and infeasible.
NOTE: The constraint 'assign_def[486]' is empty and infeasible.
NOTE: The constraint 'assign_def[488]' is empty and infeasible.
NOTE: The constraint 'assign_def[489]' is empty and infeasible.
NOTE: The constraint 'assign_def[490]' is empty and infeasible.
NOTE: The constraint 'assign_def[491]' is empty and infeasible.
NOTE: The constraint 'assign_def[492]' is empty and infeasible.
NOTE: The constraint 'assign_def[902]' is empty and infeasible.
The minimum feasible threshold turns out to be 926.21776127 miles, the distance from customer 492 to facility 6.
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
Yes, the commented constraint is correct, except that you need &p instead of p unless you also declare p as a num:
num p = &p;
I don't immediately see anything wrong with the rest of the code, so maybe something is wrong with the data. Please provide the input data sets so that I can run this myself to investigate further.
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
ok thanks rob. updated the variable p. here is my data set in csv format.
two files (COGCUSTOMERS) and COGEXISTINGFACILITIES.
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
Note that GEODIST defaults to kilometers. If you want miles, include 'M' as a fifth argument.
When I run with distance threshold 10000, I get one facility per customer, as expected from the assign_def constraints.
When I run with distance threshold 800, I get infeasible, as you did. These notes in the log indicate several customers that are not within 800 miles of any facility:
NOTE: The constraint 'assign_def[478]' is empty and infeasible.
NOTE: The constraint 'assign_def[486]' is empty and infeasible.
NOTE: The constraint 'assign_def[488]' is empty and infeasible.
NOTE: The constraint 'assign_def[489]' is empty and infeasible.
NOTE: The constraint 'assign_def[490]' is empty and infeasible.
NOTE: The constraint 'assign_def[491]' is empty and infeasible.
NOTE: The constraint 'assign_def[492]' is empty and infeasible.
NOTE: The constraint 'assign_def[902]' is empty and infeasible.
The minimum feasible threshold turns out to be 926.21776127 miles, the distance from customer 492 to facility 6.
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
Rob.
Ok. I missed the "m" as fifth argument. ran the model for 1000 miles . But if u look at the attached output ( filter for any one customer, i see 9 possibile (all 9) facilities as options. Should i print for which one is the winner by displaying the winner facility using Assign variable?
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
When i printed the assign variable I do see that one customer is assigned to only one facility. so that is what you mentioned.
Let me see if p number of facilites are open as well
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
I tested the number of facilities and it worked as well. All good. Just two questions on outputs
1) i want to see only those customer-facility pairs where the assign variable is 1. I tried this below by throwing a filter of Build = 1 but it did not work.
create data STDOPT.COG_Output(where=(Build=1)) from [i j]={i in CUSTOMERS,j in FACILITIES} Cust_Lat=a[i,1] Cust_Lon=a[i,2] Fac_Lat=f[j,1] Fac_Lon=f[j,2] Distance=dist[i,j] Assign SiteState=SiteState[i];
2) How do I display the Sitename of the facility (the actual name) instead of the Sitenumber. like "Charleston" as opposed to number 9.
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
oops. sorry. I got it. I got the wrong variable in my create data step . When i changed to Assign, it worked
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
Comparison of Assign = 1 is dangerous because of numerical tolerances. It is safer to do Assign > 0.5 like this:
create data STDOPT.COG_Output(where=(Assign>0.5)) from [i j]={i in CUSTOMERS,j in FACILITIES}
or this:
create data STDOPT.COG_Output from [i j]={i in CUSTOMERS,j in FACILITIES: Assign[i,j].sol > 0.5}
If you want to include SiteName in the output, first declare and read it:
str SiteName {FACILITIES};
read data STDOPT.cogExistingFacilities into FACILITIES=
[_N_] {d in DIMS} <f[_N_,d]=col('f'||d)> SiteCapacity SiteName;
Then add it to the CREATE DATA statement:
create data STDOPT.COG_Output from [i j]={i in CUSTOMERS,j in FACILITIES: Assign[i,j].sol > 0.5}
Cust_Lat=a[i,1] Cust_Lon=a[i,2] Fac_Lat=f[j,1] Fac_Lon=f[j,2] Distance=dist[i,j]
Assign SiteName=SiteName[j];
Note the [j] index because SiteName depends on j only and not on both i and j.
- Mark as New
- Bookmark
- Subscribe
- Mute
- RSS Feed
- Permalink
- Report Inappropriate Content
Rob
Thank you. that makes total sense. All questions are answered perfectly. Thank you again