function c = solver(a)
global I J I1 J1
i=sqrt(-1);
al=length(a);
if al==2
c=1; return
else
s=sum(a); hcl=al-1; c=ones(hcl,1);
if s<2, return; end
end
if s==2
j=find(a); brk=floor((j(2)-j(1))/2);
c(j(1)+brk)=-i; c(j(1)+brk+1:end)=-1;
return
end
id = ones(s); [I,J] = find(triu(id,1));
Ia=logical(a); iones=find(Ia); I=iones(I); J=iones(J);
a1=a(iones(1):iones(end));a1l=length(a1);I1=I-iones(1)+1;J1=J-iones(1)+1;
csal=ceil(sqrt(al));
% more spirals
bc=Inf;
for k1=1:floor(al/3)-1,
hc=zeros(hcl,1); k2=1;j=1;v=1;
while j<=al
hc(j:min(j+k1-1,hcl))=v; j=j+k1; hc(j:min(j+k2-1,hcl))=i*v;
j=j+k2; v=-v; k1=k1+1; k2=k2+1;
end
ec=energy(hc);
if ec<bc
bc=ec; c=hc;
end
ec=energy(flipud(hc));
if ec<bc
bc=ec; c=flipud(hc);
end
end
dum=[1 i 1 -i 1 -i -1 -i -1 i -1 -1 i i 1 ...
i 1 1 1 -i 1 -i -i -i -1 -i -1 -1 -1 ...
i -1 -1 i i i i 1 i 1 1 1 1 1 -i ...
1 -i -i -i -i -i -1 -i -1 -1 -1 -1 -1 i ...
-1 -1 i i i i i i 1 i 1 1 1 1 1 ...
1 1 -i 1 -i -i -i -i -i -i -i -1 -i ...
-1 -1 -1 -1 -1 -1 -1 i -1 -1 i i i i i ...
i i i 1 i 1 1 1 1 1 1 1 1 1 -i ...
1 -i -i -i -i -i -i -i -i -i -1 -i -1 ...
-1 -1 -1 -1 -1 -1 -1 -1 i -1 -1 i i i i ...
i i i i i i 1 i 1 1 1 1 1 1 1 ...
1 1 1 1 -i 1 -i -i -i -i -i -i -i ...
-i -i -i -i -1 -i -1 -1 -1 -1 -1 -1 -1 -1 ...
-1 -1 -1 i -1 -1 i i i i i i i i i ...
i i i 1 i 1 1 1 1 1 1 1 1 1 1 ...
1 1 1 -i 1 -i -i -i -i -i -i]';
ec=energy(dum(1:hcl));
if ec<bc
bc=ec; c=dum(1:hcl);
end
ec=energy(dum(hcl:-1:1));
if ec<bc
bc=ec; c=dum(hcl:-1:1);
end
%Hilbert : Added by Eli Horn from Israel
if (hcl<64)
hc=[1 i -1 i i 1 -i 1 i 1 -i -i -1 -i 1 1 1 i -1 i i 1 -i 1 i 1 -i -i -1 -i 1 -i -1 -i 1 -i -i -1 i -1 -i -1 i i 1 i -1 -1 -1 -i 1 -i -i -1 i -1 -i -1 i i 1 i -1]';
hc=hc(1:hcl);
p=[0;cumsum(hc)]; dist=abs(p(I)-p(J));ec=sum(dist);
if ec<bc,
bc=ec; c=hc;
end
end;
% zigzag
iones1=iones(1);
iones2=iones(end)-1;
for n=max(2,min(a1l-1,4)):min(13,floor(a1l/2));
t=[ones(n-1,1);-i; -ones(n-1,1);-i];
t=t(:,ones(floor(a1l/n),1));
if rem(n,2)
hc=[-ones(n-1,1);-i;t(:)];
else
hc=t(:);
end
ec=energy1(hc(1:a1l-1));
if ec<bc
bc=ec;
c=[hc(ones(iones1,1));hc(2:a1l-2);hc(a1l-1+zeros(al-iones2,1))];
end
% zigzag2
hc=hc(end-a1l+2:end);
ec=energy1(hc);
if ec < bc
bc=ec;
c=[hc(ones(iones1,1));hc(2:a1l-2);hc(a1l-1+zeros(al-iones2,1))];
end
hc=zigzag3(a1,n);
ec=energy1(hc);
if ec<bc
bc=ec;
c=[hc(ones(iones1,1));hc(2:a1l-2);hc(a1l-1+zeros(al-iones2,1))];
end
hc=flipud(hc);
ec=energy1(hc);
if ec<bc
bc=ec;
c=[hc(ones(iones1,1));hc(2:a1l-2);hc(a1l-1+zeros(al-iones2,1))];
end
hc=zigzag3(flipud(a1),n);
ec=energy1(hc);
if ec<bc
bc=ec;
c=[hc(ones(iones1,1));hc(2:a1l-2);hc(a1l-1+zeros(al-iones2,1))];
end
hc=flipud(hc); ec=energy1(hc);
if ec<bc
bc=ec;
c=[hc(ones(iones1,1));hc(2:a1l-2);hc(a1l-1+zeros(al-iones2,1))];
end
end
% SHELLS
for shell=1:7
hc=[];
p1 = i; p2 = ones(shell,1); p3 = -i; p4=+1;
p5 = [i; i]; p6 = [-1; -p2; -1]; p7 = [-i; -i]; p8=-1;
while length(hc)<hcl,
hc = [hc; p1; p2; p3; p4; p5; p6; p7; p8];
p1 = [p1; i; i];
p2 = [1; 1; p2; 1; 1];
p3 = [p3; -i; -i];
p5 = [p5; i; i];
p6 = [-1; -1; p6; -1; -1];
p7 = [p7; -i; -i];
end
hc = hc(1:hcl);
ec=energy(hc);
if ec<bc
bc=ec;
c=hc;
end
ec=energy(flipud(hc));
if ec<bc
bc=ec;
c=flipud(hc);
end
end
%zigging zigs
ur=[i;1];dr=[-i;1];ul=[i;-1];dl=[-i;-1];uu=[i;i];
hc=[];
for n=1:ceil(al/4)+1,
for ng=1:floor(csal/2),
hc=[hc;ur;dr];
end
hc=[hc;uu];
for ng=1:floor(csal/2),
hc=[hc;ul;dl];
end
hc=[hc;uu];
end
hc = hc(1:hcl);
ec=energy(hc);
if ec<bc
bc=ec; c=hc;
end
ec=energy(flipud(hc));
if ec<bc
bc=ec; c=flipud(hc);
end
%diamond pattern
hc=[];
for n=1:ceil(sqrt(al/10))+1;
for nn=1:2*n-1,
hc=[hc;dr];
end
for an=1:2*n-1,
hc=[hc;ur];
end
for bn=1:n,
hc=[hc;ul;ul];
end
for cn=1:n,
hc=[hc;dl;dl];
end
end
hc = hc(1:hcl);
ec=energy(hc);
if ec<bc,
bc=ec; c=hc;
end
ec=energy(flipud(hc));
if ec<bc
bc=ec; c=flipud(hc);
end
% SHELLS WITH STEM
for shell=1:8
% hc=[];
% GIVE IT A STEM OF HYDEOPHILICS IF POSSIBLE
phobics = find(a);
hc = 1i*ones(phobics(1),1);
p1 = 1i; p2 = ones(shell,1); p3 = -1i; p4=+1;
p5 = [1i; 1i]; p6 = [-1; -p2; -1]; p7 = [-1i; -1i]; p8=-1;
while length(hc)<hcl,
hc = [hc; p1; p2; p3; p4; p5; p6; p7; p8];
p1 = [p1; 1i; 1i];
p2 = [1; 1; p2; 1; 1];
p3 = [p3; -1i; -1i];
p5 = [p5; 1i; 1i];
p6 = [-1; -1; p6; -1; -1];
p7 = [p7; -1i; -1i];
end
hc = hc(1:hcl);
ec=energy(hc);
if ec<bc
bc=ec;
c=hc;
end
ec=energy(flipud(hc));
if ec<bc
bc=ec;
c=flipud(hc);
end
end
%interleave spiral pattern
tmp = a .* [1:length(a)]';
if (length((tmp(tmp>0)))==0)
hc = ones(al-1,1);
return
end
start = round(mean(tmp(tmp>0)));
d = floor(sqrt(start/2)*2-1);
if(~mod(d,2))
d=d-1;
end
dd=sum((1:2:d)*2);
hc=[];
if (dd < start)
if ( (dd+d+2 <start))
hc=[-1*ones(start-(dd+d+2),1); -i*ones(d+2,1)] ;
else
hc=[-i*ones(start-dd,1)];
end
end
df = 1;
for j = d:-2:1
hc = [hc; df*ones(j,1)];
df = df*i;
hc = [hc; df*ones(j,1)];
df = df*i;
end
df = df*-i;
j=1;
while(1)
if ((length(hc)+j) > al)
hc = [hc; df*ones(al-length(hc),1)];
break
else
hc = [hc; df*ones(j,1)];
df = df*-i;
end
if ((length(hc)+j) > al)
hc = [hc; df*ones(al-length(hc),1)];
break
else
hc = [hc; df*ones(j,1)];
df = df*-i;
end
j=j+2;
end
hc = hc(1:end-1);
ec=energy(hc);
if ec<bc,
bc=ec; c=hc;
end
ec=energy(flipud(hc));
if ec<bc
bc=ec; c=flipud(hc);
end
% snake disc
hc = snakedisc(a);
ec = energy(hc);
if ec<bc; bc=ec; c=hc; end;
hc = flipud(snakedisc(flipud(a)));
ec = energy(hc);
if ec<bc; bc=ec; c=hc; end;
% rle snake
hc = rlesnake(a);
ec = energy(hc);
if ec<bc; bc=ec; c=hc; end;
hc = flipud(rlesnake(flipud(a)));
ec = energy(hc);
if ec<bc; bc=ec; c=hc; end;
%%% Energy
function ec=energy(c)
global I J
p = [0;cumsum(c)];
dist=abs(p(I)-p(J));
ec=sum(dist);
%%% Energy1
function ec=energy1(c)
global I1 J1
p = [0;cumsum(c)];
dist=abs(p(I1)-p(J1));
ec=sum(dist);
%%% Snake Disc. Easy. Dumb.
function c = snakedisc(a)
i = sqrt(-1);
N = length(a);
R = sqrt(N/pi);
c = [];
dd = 1;
ll = 0;
for hh = -R:R
ll = max(round(abs(halfchord(R,hh+0.5))),1);
c = [ c;dd(ones(1,ll-1),1);i;-dd(ones(1,ll),1) ];
dd = -dd;
end
if length(c) < N-1; c(end+1:N-1)=c(end); end;
c = c(1:N-1);
function ch = halfchord(r,h)
ch = sqrt(r*r-h*h);
%%% RLE Snake. Needs improvement.
function c = rlesnake(a)
i = sqrt(-1);
N = length(a);
% run length encode a
edges = [1 1+find(a(2:end) ~= a(1:end-1))'];
val = a(edges)';
run = [ edges(2:end)-edges(1:end-1) N-edges(end)+1 ];
N1 = sum(run(val==1));
S = round(sqrt(N1));
c = [];
dd = 1;
for nn = 1:length(val)
V = val(nn);
R = run(nn);
if V
if R < 1.5*S
c = [c;dd(ones(1,R),1)];
else
if R == 1
c = [c;dd];
elseif R == 2
c = [c;dd;i];
dd = -dd;
else
k = round((R-1)/2);
c = [c;dd(ones(1,k),1);i;-dd(ones(1,R-k-1),1)];
dd = -dd;
end
end
else
if R == 1
c = [c;i];
dd = -dd;
elseif R == 2
c = [c;dd;i];
dd = -dd;
else
k = round((R-1)/2);
c = [c;dd(ones(1,k),1);i;-dd(ones(1,R-k-1),1)];
dd = -dd;
end
end
end
if length(c) > N-1; c=c(1:N-1); end;
function c=zigzag3(a,l)
% zigzag3
i=sqrt(-1);
al=length(a);
k=1;
c=zeros(al-1,1);
d=1;
while 1
c(k:k+l-1)=d;
k=k+l;
if k==al-1
c(k)=i;
break;
elseif k>=al-1
break;
elseif k==al-2
c(k)=i;
c(k+1)=-d;
break;
end
l1=find(a(k+1:end)); % can't be empty
l1=l1(1);
if l1>=3
l1=floor((l1-1)/2);
c(k:k+l1-1)=d;
k=k+l1;
c(k)=i;
c(k+1:k+l1)=-d;
k=k+l1+1;
else
c(k)=i;
k=k+1;
end
d=-d;
end
c=c(1:al-1);
|