| Code: |
function c = solver(a)
iworth=find(a(:,3));
if isempty(iworth), c=[1,1]; return; end
mfreight=mean(a(iworth,3));
ia=find( a(:,4) | a(:,3)>=mfreight/15.5 );
ia=find( a(:,4) | a(:,3)>=mfreight/16);
a=a(ia,:);
xy = (a(:,1)-a(1)) + 1i*(a(:,2)-a(1,2));
n=size(a,1);
if ~any(a(1,4) >= abs(xy(2:n)))
% Nofreight or not enough fuel to go anywhere, stay home (Ave)
c=[1 1];
return
end
da=a(:,4)-a(:,3);
[X,Y]=meshgrid(xy);
dist = abs(X-Y);
dist0=dist;
I=1;
c=ones(1,n+1);
g=0;
for k=2:n
g=g+a(I,4);
dist(:,I)=inf;
[Y,I]=min(dist(I(1),:));
g=g-Y;
if g<0
k=k-1;
break;
end
c(k) = I;
end
c=c(1:k+1);
xy1 = xy(c);
i=sub2ind([n n],c(1:k),c(2:k+1));
g = cumsum(a(c(1:k),4)-dist0(i)');
if g(end)<0
if k<3
c=[1 1];
elseif g(end)<dist0(c(k))
% look for last point from which starting point can be reached
g=g(1:k-2)+a(c(2:k-1),4)-dist0(c(2:k-1),1);
i=find(g>=0);
if isempty(i)
c=[1 1];
else
i=i(end)+2;
c=c(1:i);
c(i)=1;
end
else
c = c(1:k+1);
c(end)=1;
end
end
score1 = sum(da(c(2:end-1)));
c2 = diesel2(a,xy,dist0);
score2 = sum(da(c2(2:end-1)));
if(score2<score1)
c=c2;
score1=score2;
end
c = addGas(a,c,dist0);
fail=0;
while ~fail
[c,fail] = addFreight(a,c,dist0);
if fail
break
end;
[c,fail] = addGas(a,c,dist0);
end
[c,score1]=optim(a,c,score1,xy,da);
if da(1)<score1
% dummy solution is the best(!!)
c=[1 1];
end
c=ia(c);
function [c,fail] = addGas(a,c,dist)
fail=1;
gasv=a(:,4);
gasv(c)=0;
c=c(:);
Nc = length(c);
fIdx = find(gasv>0);
if isempty(fIdx), return; end % NO GAS LEFT TO ADD
N = length(fIdx);
% CALCULATE EXCESS GAS AT EACH POINT IN PATH
i=sub2ind(size(dist),c(1:Nc-1),c(2:Nc));
d = dist(i);
eg = cumsum(a(c(1:Nc-1),4)-d);
dispgas=zeros(Nc-1,1);
dispgas(Nc-1)=eg(Nc-1);
for counter=Nc-2:-1:1
dispgas(counter)=min(dispgas(counter+1),eg(counter));
end
% SORT FREIGHT BY VALUE
% TRY TO FIT BIGGEST FIRST
[ignore,idx] = sort(-a(fIdx,4));
fIdx = fIdx(idx);
for n=1:N,
i=fIdx(n);
% SET OF EXCESS DISTANCES ADDED BY VISITING THIS POINT
ad = dist(c,i);
ed = ad(1:Nc-1) + ad(2:Nc) - d;
h=dispgas-ed;
[ignore,closest] = max(h);
if ed(closest)<a(fIdx(n),4)
if h(closest)>=0
fail=0;
d1=dist(c(closest),i);
d2=dist(i,c(closest+1));
d0=dist(c(closest),c(closest+1));
c = [c(1:closest); i; c(closest+1:Nc)];
d=[d(1:closest-1);d1;d2;d(closest+1:Nc-1)];
eg = cumsum(a(c(1:Nc),4)-d);
dispgas=zeros(Nc,1);
dispgas(Nc)=eg(Nc);
for counter=Nc-1:-1:1
dispgas(counter)=min(dispgas(counter+1),eg(counter));
end
Nc=Nc+1;
end
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% addFreight
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [c,fail] = addFreight(a,c,dist)
fail=1;
c=c(:);
Nc = length(c);
a(c,3)=0;
fIdx = find(a(:,3)>0);
if isempty(fIdx), return; end % NO FREIGHT LEFT TO ADD
N = length(fIdx);
% CALCULATE EXCESS GAS AT EACH POINT IN PATH
i=sub2ind(size(dist),c(1:Nc-1),c(2:Nc));
d = dist(i);
eg = cumsum(a(c(1:Nc-1),4)-d);
dispgas=zeros(Nc-1,1);
dispgas(Nc-1)=eg(Nc-1);
for counter=Nc-2:-1:1
dispgas(counter)=min(dispgas(counter+1),eg(counter));
end
% SORT FREIGHT BY VALUE
% TRY TO FIT BIGGEST FIRST
[ignore,idx] = sort(-a(fIdx,3));
fIdx = fIdx(idx);
for n=1:N
i=fIdx(n);
% SET OF EXCESS DISTANCES ADDED BY VISITING THIS POINT
ad = dist(c,i);
ed = ad(1:Nc-1) + ad(2:Nc) - d;
h=dispgas-ed;
[ignore,closest] = max(h);
if h(closest)>=0
fail=0;
d1=dist(c(closest),i);
d2=dist(i,c(closest+1));
d0=dist(c(closest),c(closest+1));
c = [c(1:closest); i; c(closest+1:Nc)];
d=[d(1:closest-1);d1;d2;d(closest+1:Nc-1)];
eg=[eg(1:closest-1);eg(closest)+d0-d1;eg(closest:Nc-1)-(d1+d2-d0)];
dispgas=zeros(Nc,1);
dispgas(Nc)=eg(Nc);
for counter=Nc-1:-1:1
dispgas(counter)=min(dispgas(counter+1),eg(counter));
end;
Nc=Nc+1;
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function c = diesel2(a,xy,Dmat)
n=size(a,1);
freightv=a(:,3);
petrolv=a(:,4);
Dmat1=Dmat(:,1);
pos=1;
jdest=1;
totalpetrol=sum(petrolv);
unitpetrol = totalpetrol/(n-1);
petroln=petrolv/totalpetrol;
totalfreight=sum(freightv);
freightn=freightv/totalfreight;
petrol = petrolv(1);
petrolv(1)=0; freightv(1)=-1;
notbeenthere=logical(ones(n,1));
notbeenthere(1)=0;
c=ones(n+1,1);
notnearlythere=notbeenthere;
while 1
jdest=jdest+1;
% is it possible to go anywhere?
dpos=Dmat(:,pos);
bpossible=(petrol >= dpos & notbeenthere);
if ~any(bpossible)
break; % fuel low, go home if possible
end
% is it safe to go there?
safe = find(bpossible & (petrol + petrolv >= dpos + Dmat1));
if (isempty(safe))
possible=find(bpossible);
% is the next station after this one safe? (Ave)
safe=zeros(n,1);
for li=1:length(possible)
la=possible(li);
lanotbeenthere=notbeenthere;
lanotbeenthere(la)=0;
ladpos=Dmat(la,pos)+Dmat(:,la);
lapetrol=petrol+petrolv(la);
lapossible=(lapetrol>=ladpos & lanotbeenthere);
if any(lapossible & lapetrol+petrolv>=ladpos+Dmat1)
safe(la)=1;
end
end
safe=find(safe);
end
if isempty(safe)
break;
end
possible=safe;
Dpossible=Dmat(possible,pos);
metric = ( (petroln(possible) - Dpossible/totalpetrol) * 2*(n-jdest+1)/n*(1+unitpetrol/petrol)^5 ...
+ freightn(possible)) ...
./ Dpossible;
[maxm,imax]=max(metric);
next=possible(imax);
notnearlythere(next)=0;
r=Dmat(pos,next);
% find places that are on the way
ontheway = find( (Dmat(:,next) + Dmat(:, pos) < 1.088*r & notnearlythere & bpossible) );
if(~isempty(ontheway))
dp = Dmat(ontheway,pos); dn=Dmat(ontheway,next);
ontheway = ontheway( find ( ...
dp+dn+Dmat(1,next) <= petrol+petrolv(ontheway)+petrolv(next) & ... % can get home from
((dp+dn-r)<petrolv(ontheway) | freightv(ontheway)>0 ) ) ... % net freight or petrol gain
);
[maxf,im]=max(freightv(ontheway));
if (~isempty(ontheway))
next=ontheway(im(1));
end
end
c(jdest)=next;
petrol=petrol + petrolv(next) - Dmat(pos,next);
pos=next;
notbeenthere(pos)=0;
notnearlythere=notbeenthere;
end
c(jdest)=1;
c=c(1:jdest);
function [c,s]=optim(a,c,s,xy,da)
opt=1;
g0 = cumsum(a(c(1:end-1),4)-abs(diff(xy(c))));
while opt
opt=0;
da1=da(c);
i=find(da1(2:end-1)>0)+1;
if isempty(i)
return
end
dg=a(c(i),4)-abs(xy(c(i))-xy(c(i-1)))-abs(xy(c(i+1))-xy(c(i)))+abs(xy(c(i+1))-xy(c(i-1)));
j=find(dg<=g0(end));
if isempty(j)
return
end
i=i(j);
dg=dg(j);
[da11,j]=sort(-da1(i));
i=i(j);
dg=dg(j);
for j=1:length(i)
if all(g0(i(j):end)>=dg(j))
s=s+da11(j);
g0=[g0(1:i(j)-2);g0(i(j):end)-dg(j)];
c(i(j))=[];
opt=1;
break;
end
end
end
|