pso2.m

这是从matlab上down的源文件。上面显示调试通过

		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
end
function c = diesel2(a,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;
imax=[0;0];
notbeenthere=logical(ones(n,1));
notbeenthere(1)=0;
c=ones(n+1,1);
notnearlythere=notbeenthere;
while 1
	jdest=jdest+1;
	
	dpos=Dmat(:,pos);
	bpossible=(petrol >= dpos & notbeenthere);
	if ~any(bpossible)
		break;	
	end
	
	safe = find(bpossible & (petrol + petrolv >= dpos + Dmat1));
	if (isempty(safe))
		possible=find(bpossible);
		
		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;
				break;
			end
		end
		safe=find(safe);
	end
	if isempty(safe)
		break;
	end
	possible=safe;
	
    lm=length(possible);
    
    if (lm>1)
        Dpossible=Dmat(possible,pos);
        metric = ( (petroln(possible) - Dpossible/totalpetrol) * 2*(n-jdest+1)/n*(1+unitpetrol/petrol)^5 ...
			+ freightn(possible)) ...
			./ Dpossible;
        
        
        sm = [-realmax; -realmax];
        for imetric=1:lm
            mi=metric(imetric);
            if(mi>sm(2))
                sm(1)=sm(2);  imax(1)=imax(2);
                sm(2)=mi;     imax(2)=imetric;
			elseif mi>sm(1)
                sm(1)=mi;  imax(1)=imetric;
			end 
		end
		rn=0;
        if sm(1)/sm(2)>.995
			rn=fix(rand+.5);
		end
        next=possible(imax(end-rn));
	else
        next=possible(1);
	end
	notnearlythere(next)=0;
	r=Dmat(pos,next);
	
	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)  & ...   
			((dp+dn-r)<petrolv(ontheway) | freightv(ontheway)>0 )   ) ...         
			);
		[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 cb=turbodiesel(a,dist)
freightv=a(:,3);
petrolv=a(:,4);
da=freightv-petrolv;
n=size(a,1);
Nstck=1000;
stck=zeros(Nstck,1);
nstck=1;
c=ones(n+1,1);
cb=[1 1];
avail=logical(ones(n,1));
GAS=ones(n,1);
stck(2)=1;
nstck=2;
totalpetrol=sum(petrolv);
unitpetrol = totalpetrol/(n-1);
petroln=petrolv/totalpetrol;
totalfreight=sum(freightv);
freightn=freightv/totalfreight;
gas=petrolv(1);
petrolv(1)=0; freightv(1)=-1;
dist1=dist(:,1);
if n<14
	IMAX=3;
else
	IMAX=2;
end
fgbest=-inf;
loop=1;
while 1
	current=stck(nstck);
	if current==0
		loop=loop-1;
		if loop<=1
			break;
		end
		avail(c(loop))=1;
		nstck=nstck-1;
		gas=GAS(loop-1);
	else
		c(loop)=current;
		avail(current)=0;
		if loop>1
			gas=gas-dist(c(loop-1),current)+petrolv(current);
			
			if gas>=dist(current)
				fg=sum(da(c(2:loop)));
				if fgbest<fg	
					cb=c(1:loop+1);
					cb(loop+1)=1;
					fgbest=fg;
				end
			end
			if all(da(avail)<=0)
				break;
			end
		end
	
		dpos=dist(:,current);
		bpossible=(gas >= dpos & avail);
		nok=1;
		if any(bpossible)
			
			safe = find(bpossible & (gas + petrolv >= dpos + dist1));
			if isempty(safe)
				reachable=find(bpossible);
				
				bsafe=logical(zeros(n,1));
				for li=1:length(reachable)
					la=reachable(li);
					lanotbeenthere=avail;
					lanotbeenthere(la)=0;
					ladpos=dist(la,current)+dist(:,la);
					lapetrol=gas+petrolv(la);
					lapossible=(lapetrol>=ladpos & lanotbeenthere);
					if any(lapossible & lapetrol+petrolv>=ladpos+dist1)
						bsafe(la)=1;
						break;
					end
				end
				safe=find(bsafe);
				if isempty(safe)
					safe=find(bpossible);
					Dpossible=dist(safe,current);
					metric = ( (petroln(safe) - Dpossible/totalpetrol) * 2*(n-loop+1)/n*(1+unitpetrol/gas)^5 ...
						+ freightn(safe)) ...
						./ Dpossible;
					
					[sm,imax]=sort(-metric);
					if nstck>16
						IMAX2=1;
					else
						IMAX2=IMAX;
					end
					safe=safe(imax(1:min(end,IMAX2)));
				end
			end
			if ~isempty(safe)
				reachable=safe;
				
				Dpossible=dist(reachable,current);
				metric = ( (petroln(reachable) - Dpossible/totalpetrol) * 2*(n-loop+1)/n*(1+unitpetrol/gas)^5 ...
					+ freightn(reachable)) ...
					./ Dpossible;
				
				[sm,imax]=sort(-metric);
				
				if nstck>30
					IMAX2=1;
				else
					IMAX2=IMAX;
				end
				reachable=reachable(imax(1:min(end,IMAX2)));
				
				notnearlythere=avail;
				for i=1:length(reachable)
					next=reachable(i);
					r=dist(current,next);
					
					ontheway = find( (dist(:,next) + dist(:, current) < 1.088*r  & notnearlythere & bpossible)  );
					if(~isempty(ontheway))
						dp = dist(ontheway,current); dn=dist(ontheway,next);
						ontheway = ontheway( find ( ...
							dp+dn+dist(1,next) <= gas+petrolv(ontheway)+petrolv(next)  & ...   
							((dp+dn-r)<petrolv(ontheway) | freightv(ontheway)>0 )   ) ...         
							);
						[maxf,im]=max(freightv(ontheway));
						if (~isempty(ontheway))
							reachable(i)=ontheway(im(1));
						end
					end
					notnearlythere(next)=0;
				end	
			
				if ~isempty(reachable)
					stck(nstck)=0;
					stck(nstck+1:nstck+length(reachable))=reachable;
					nstck=nstck+length(reachable);
					GAS(loop)=gas;
					loop=loop+1;
					nok=0;
				end	
			end	
		end	
		if nok
			if dist(current)<=gas	
				fg=sum(da(c(2:loop)));
				if fgbest<fg	
					cb=c(1:loop+1);
					cb(loop+1)=1;
					fgbest=fg;
				end
			end
			nstck=nstck-1;
			avail(current)=1;
		end
	end	
end	
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
function c = petrol3(a,Dmat)
n=size(a,1);
freightv=a(:,3);
petrolv=a(:,4);
petrol=0;
freight=0;
pos=1;
jdest=1;
totalpetrol=sum(petrolv);
unitpetrol = totalpetrol/(n-1);
petroln=petrolv/totalpetrol;
totalfreight=sum(freightv);
freightn=freightv/totalfreight;
petrol  = petrol +  petrolv(1);
freight = freight + freightv(1);
petrolv(1)=0; freightv(1)=-1;
notbeenthere=logical(ones(n,1));
notbeenthere(1)=0;
c=ones(1,n+1);
Dp = max(Dmat(:,1) - petrolv,0);
while 1
    jdest=jdest+1;
    pmD = petrol - Dmat(:,pos);
    possible = find(  ...
        (pmD >= Dp)  & ...
        notbeenthere  ...
        );
    
    if isempty(possible)
        break
    end
    
    Dpossible=Dmat(possible,pos);
    metric = ( (petroln(possible) - Dpossible/totalpetrol) * 2*(n-jdest+1)/n*(1+unitpetrol/petrol)^5 ...
        + freightn(possible)) ...
        ./ Dpossible;
    imax= 1;maxm = -inf;
    for count = 1:numel(metric)
        if metric(count)>maxm
            imax = count;
            maxm = metric(count);
        end
    end
    
    
    next=possible(imax);
    c(jdest)=next;
    petrol=petrol + petrolv(next) - Dmat(pos,next);
    freight=freight+freightv(next);
    pos=next;
    notbeenthere(pos)=0;
end
c=c(1:jdest);