Revision: 4264
Initial Code
Initial URL
Initial Description
Initial Title
Initial Tags
Initial Language
at November 14, 2007 15:15 by mertnuhoglu
Initial Code
--- solve_problems.m function report = solve_problems() % data_files={'Alberta';'Galvao100';'Galvao100';'Galvao150';'Galvao150';'Galvao150'}; % p=[10 10 15 5 15 20]; % data_files={'TestData'}; % p=[2]; data_files={'Alberta'}; p=[10]; for i=1:length(data_files) data_file=char(data_files(i)); dist=distances([data_file '_distances.txt']); demand=load([data_file '_demands.txt']); [bestLB,iterations,debug]=solve_p_median(dist,demand,p(i)); report(i).bestLB=bestLB; report(i).iterations=iterations; report(i).debug=debug; dlmwrite([data_file '_p' int2str(p(i)) '_iterations.txt'],iterations); dlmwrite([data_file '_p' int2str(p(i)) '_debug.txt'],debug); end --- distances.m function distances=distances(data_file) distance_data=load(data_file); distances=zeros(distance_data(end,1)+1); for row=1:length(distance_data) from=distance_data(row,1); to=distance_data(row,2); distance=distance_data(row,3); distances(to,from)=distance; distances(from,to)=distance; end --- solve_p_median.m function [bestLB,iterations,debug]=solve_p_median(dist,demand,p) bestLB=0; bestUB=inf; currentLB=0; currentUB=inf; iterations(1,:)=[0 currentLB bestLB currentUB bestUB]; pi=2; n_c=length(demand); % number of customers n_s=length(dist(:,1)); % number of sites u=ones(n_c,1); % LR (Lagrangean Relaxation) multipliers zero=zeros(n_c,1); x=zeros(n_s,n_c); % site/customer assignments i=1; piUpdateTime=1; improvementOccurred=0; debug=[0 0 2 zeros(1,p) u']; % parameters stored for debugging (iteration, step_size, pi, open facilities, u) while(~stoppingCondition(pi,bestLB,bestUB,i)) for s=1:n_s cost=dist(:,s).*demand; newCost=cost-u; z_LR(s)=sum(min(zero,newCost)); x(s,find(newCost<0))=1; % x_{s,c} = 1 if cost negative end [z_sorted,order]=sort(z_LR); currentLB=sum(z_sorted(1:p))+sum(u); % z_{LR}(u) value is current LB facilities=order(1:p); % open p facilities where z is smallest x(setdiff(1:n_s,facilities),:)=0; % x_{s,c} \leq y_s \forall s,c constraint if(currentLB>bestLB) bestLB=currentLB; end currentUB=findUB(facilities,dist,demand); if(currentUB<bestUB) bestUB=currentUB; end iterations(i+1,:)=[i currentLB bestLB currentUB bestUB]; normOfRelaxedCsts=sum((1-sum(x)).^2); if(normOfRelaxedCsts == 0) % hit the lower bound break end step=pi*(bestUB-currentLB)/normOfRelaxedCsts; % s^t = {\pi (UB* - z_{LR}(u^t)) \over \sum_c (1-\sum_s x_{sc})^2} u=u+step*(1-sum(x))'; % u_c^{t+1} = u_c^t + s^t (1-\sum_s x_{sc}^t) if(~improvementsOccur(iterations,piUpdateTime)) pi=pi/2; piUpdateTime=i; end debug(i+1,:)=[i step pi facilities u']; i=i+1 end function result=stoppingCondition(pi,bestLB,bestUB,iterationNo) result = (pi < 0.005) | (bestUB == bestLB); % result = iterationNo > 15000 | (bestUB == bestLB); function result=improvementsOccur(iterations,piUpdateTime) n=30; currentTime=length(iterations(:,1)); timeSinceLastPiUpdate=currentTime-piUpdateTime; if(currentTime <= n | timeSinceLastPiUpdate <= n) result = 1; else lastImpForLB=whenDidLastImprovementOccur(iterations(end-n:end,3)); lastImpForUB=whenDidLastImprovementOccur(iterations(end-n:end,5)); if(lastImpForLB <= n | lastImpForUB <= n) result = 1; else result = 0; end end function lastImp=whenDidLastImprovementOccur(iterations) lastValue=iterations(end); ixLastImp=length(find(iterations~=lastValue)); lastImp=length(iterations)-ixLastImp; function currentUB=findUB(facilities,dist,demand) feasibleAssignments=assignCustomers(facilities,dist); currentUB=sum(feasibleAssignments.*dist)*demand; function customerAssignments=assignCustomers(facilities,dist) n_s=length(dist(:,1)); % number of sites n_c=length(dist(1,:)); % number of customers customerAssignments=zeros(n_s,n_c); distOpenFacilities=dist(facilities,:); [minDist,order]=min(distOpenFacilities); for c=1:n_c customerAssignments(facilities(order(c)),c)=1; end
Initial URL
Initial Description
Initial Title
Solving p-median location allocation problem with Lagrangian Relaxation heuristics
Initial Tags
Initial Language
MatLab