Finish 2011-11-09 12:00:00 UTC

1

by amit rajora

Status: Passed
Results: 65797264 (cyc: 33, node: 8094)
CPU Time: 89.249
Score: 658577.0
Submitted at: 2011-11-07 08:02:04 UTC
Scored at: 2011-11-07 08:04:14 UTC

Current Rank: 652nd (Highest: 300th )
Based on: in the spirit(2) (diff)
Basis for: 2 (diff)
Basis for: 1 (diff)

Comments
Please login or create a profile.
Code
function [moves, vine] = solver(board, limit)

worththeeffort=sum(board);
worththeeffort=sum(worththeeffort);


[moves1,vine1,scr1] = nick(board,limit);
if worththeeffort<100000
    moves=moves1;
    vine=vine1;
else
    [moves2,vine2,scr2] = alfi(board,limit);

    if scr1>scr2
        moves = moves1;
        vine = vine1;
    else
        moves = moves2;
        vine = vine2;
    end
end

end

%% nick
function [moves, vine,scr] = nick(board, limit)
% Added filter to crosseyed
%
% N. Howe

[vine1,scr1] = trickle2(board);
moves1 = [];
[moves2,vine2,scr2] = pusher(board,limit);
[nrow,ncol] = size(board);
ntot = nrow*ncol;
[moves3,vine3,scr3] = pusher(board',limit);
[moves4,vine4,scr4] = pusher(rot90(board,2),limit);
[moves5,vine5,scr5] = pusher(fliplr(board'),limit);
[moves6,vine6,scr6] = dofilter(board,limit);
if (scr1==max([scr1 scr2 scr3 scr4 scr5 scr6]))
    moves = moves1;
    vine = vine1;
elseif (scr2==max([scr1 scr2 scr3 scr4 scr5 scr6]))
    moves = moves2;
    vine = vine2;
elseif (scr3==max([scr1 scr2 scr3 scr4 scr5 scr6]))
    id = reshape(1:ntot,nrow,ncol);
    id3 = id';
    moves = id3(moves3);
    vine = id3(vine3);
elseif (scr4==max([scr1 scr2 scr3 scr4 scr5 scr6]))
    id = reshape(1:ntot,nrow,ncol);
    id4 = rot90(id,2);
    moves = id4(moves4);
    vine = id4(vine4);
elseif (scr5==max([scr1 scr2 scr3 scr4 scr5 scr6]))
    id = reshape(1:ntot,nrow,ncol);
    id5 = fliplr(id');
    moves = id5(moves5);
    vine = id5(vine5);
else
    moves = moves6;
    vine = vine6;
end;
scr = max([scr1 scr2 scr3 scr4 scr5 scr6]);
end

function [moves, vine, scr] = pusher(board, limit)

moves = zeros(limit,2);
[nrow,ncol] = size(board);
if (limit < nrow)
    scr = -inf;
    moves = [];
    vine = [];
    return;
end;
ntot = nrow*ncol;
xg = repmat(1:ncol,nrow,1);
yg = repmat((1:nrow)',1,ncol);
pos = reshape(1:ntot,nrow,ncol);
fitloc = reshape(1:ntot,nrow,ncol);
fitloc(:,2:2:end) = flipud(fitloc(:,2:2:end));
[s,r] = sort(board(:),'descend');
rank = zeros(nrow,ncol);
rank(r) = 1:ntot;
remain = true(1,ntot);
nextfit = 1;
steps = 0;
cboard = board;
for i = 1:ntot
    if remain(i)
        [cy,cx] = find(rank==i);
        tx = xg(fitloc(nextfit));
        ty = yg(fitloc(nextfit));
        while (cx~=tx)|(cy~=ty)
            dx = sign(tx-cx);
            dy = sign(ty-cy);
            if (dx~=0)&(dy~=0)
                if cboard(cy,cx+dx)<=cboard(cy+dy,cx)
                    dy = 0;
                else
                    dx = 0;
                end;
            end;
            if (rank(cy+dy,cx+dx)>0)
                vx = cx+dx;
                vy = cy+dy;
                if (dx==0)
                    if (vx>1)&(cboard(vy,vx)>cboard(vy,vx-1))&((vx==ncol)|(cboard(vy,vx+1)>cboard(vy,vx-1)))
                        % move left
                        steps = steps+1;
                        if (steps > limit)
                            break;
                        end;
                        moves(steps,1) = pos(vy,vx);
                        moves(steps,2) = pos(vy,vx-1);
                        cboard(vy,vx-1) = cboard(vy,vx);
                        cboard(vy,vx) = 0;
                        if (rank(vy,vx-1)>0)
                            remain(rank(vy,vx-1)) = 0;
                        end;
                        rank(vy,vx-1) = rank(vy,vx);
                        rank(vy,vx) = 0;
                        %imagesc(cboard);
                        %axis image;
                        %drawnow;
                    elseif (vx < ncol)&(cboard(vy,vx)>cboard(vy,vx+1))
                        % move right
                        steps = steps+1;
                        if (steps > limit)
                            break;
                        end;
                        moves(steps,1) = pos(vy,vx);
                        moves(steps,2) = pos(vy,vx+1);
                        cboard(vy,vx+1) = cboard(vy,vx);
                        cboard(vy,vx) = 0;
                        if (rank(vy,vx+1)>0)
                            remain(rank(vy,vx+1)) = 0;
                        end;
                        rank(vy,vx+1) = rank(vy,vx);
                        rank(vy,vx) = 0;
                        %imagesc(cboard);
                        %axis image;
                        %drawnow;
                    else
                        remain(rank(cy+dy,cx+dx)) = 0;
                    end;
                else
                    if (vy>1)&(cboard(vy,vx)>cboard(vy-1,vx))&((vy==nrow)|(cboard(vy+1,vx)>cboard(vy-1,vx)))
                        % move north
                        steps = steps+1;
                        if (steps > limit)
                            break;
                        end;
                        moves(steps,1) = pos(vy,vx);
                        moves(steps,2) = pos(vy-1,vx);
                        cboard(vy-1,vx) = cboard(vy,vx);
                        cboard(vy,vx) = 0;
                        if (rank(vy-1,vx)>0)
                            remain(rank(vy-1,vx)) = 0;
                        end;
                        rank(vy-1,vx) = rank(vy,vx);
                        rank(vy,vx) = 0;
                        %imagesc(cboard);
                        %axis image;
                        %drawnow;
                    elseif (vy < nrow)&(cboard(vy,vx)>cboard(vy+1,vx))
                        % move right
                        steps = steps+1;
                        if (steps > limit)
                            break;
                        end;
                        moves(steps,1) = pos(vy,vx);
                        moves(steps,2) = pos(vy+1,vx);
                        cboard(vy+1,vx) = cboard(vy,vx);
                        cboard(vy,vx) = 0;
                        if (rank(vy+1,vx)>0)
                            remain(rank(vy+1,vx)) = 0;
                        end;
                        rank(vy+1,vx) = rank(vy,vx);
                        rank(vy,vx) = 0;
                        %imagesc(cboard);
                        %axis image;
                        %drawnow;
                    else
                        remain(rank(cy+dy,cx+dx)) = 0;
                    end;
                end;
            end;
            steps = steps+1;
            if (steps > limit)
                break;
            end;
            moves(steps,1) = pos(cy,cx);
            moves(steps,2) = pos(cy+dy,cx+dx);
            cboard(cy+dy,cx+dx) = s(i); %cboard(cy,cx);
            cboard(cy,cx) = 0;
            %rank(cy+dy,cx+dx) = i; %cboard(cy,cx);
            rank(cy,cx) = 0;
            cx = cx+dx;
            cy = cy+dy;
            %imagesc(cboard);
            %axis image;
            %drawnow;
        end;
        if (steps > limit)
            break;
        end;
        rank(ty,tx) = i;
        nextfit = nextfit+1;
    end;
end;
moves = moves(1:min(steps,limit),:);
[vine,scr] = trickle2(cboard);
end

function [vine,scr] = trickle2(board)
doclusters=any(any(board(:,2:end)==board(:,1:end-1)))||any(any(board(1:end-1,:)==board(2:end,:))); %%%%
if doclusters, 
    [vine,scr] = kurt(board);
else
    [m,n]=size(board);
    board=[zeros(1,n+2);[zeros(m,1),board,zeros(m,1)];zeros(1,n+2)];
    vine=bellman_pixel(board);
    scr=sum(board(vine));
    [iC1,iC2]=ind2sub([m+2,n+2],vine);
    vine=sub2ind([m,n],iC1-1,iC2-1);
end
end
% function [vine,scr] = trickle2(board)
% [vine,scr] = kurt(board);
% end

function [vine,scr] = realtrickle2(board)
% Second effort
% Somewhat better handling of ties.
%
%N. Howe

moves = [];
[nrow,ncol] = size(board);
xg = repmat(1:ncol,nrow,1);
yg = repmat((1:nrow)',1,ncol);
limbo = nrow*ncol+1;
vals = [board(:)' -inf];
%up = [zeros(1,ncol); board(1:end-1,:)];
%down = [board(2:end,:);zeros(1,ncol)];
%left = [zeros(nrow,1); board(:,1:end-1)];
%right = [board(:,2:end);zeros(nrow,1)];
ind = reshape(1:nrow*ncol,nrow,ncol);
up = [repmat(limbo,1,ncol); ind(1:end-1,:)];
down = [ind(2:end,:);repmat(limbo,1,ncol)];
left = [repmat(limbo,nrow,1) ind(:,1:end-1)];
right = [ind(:,2:end) repmat(limbo,nrow,1)];
nbr = [[up(:) down(:) left(:) right(:)]; limbo limbo limbo limbo];
%vnbr = vals(nbr);
%big = repmat(board(:),1,4)<=vnbr;
cum = vals;
next = zeros(1,nrow*ncol+1);
next(limbo) = -1;
levels = unique(vals);
done(nrow*ncol+1) = true;
for lvl = levels(end:-1:2)
    lvli = find(board==lvl);
    lvln = numel(lvli);
    for ii = 1:lvln
        moves = (vals(nbr(lvli,:))>=lvl).*done(nbr(lvli,:));
        [heur,jheur] = max(cum(nbr(lvli,:)).*moves,[],2);
        heur(done(lvli)) = 0;
        [mheur,mheuri] = max(heur);
        if (mheur==0)
            done(lvli) = true;
            break;
        end;
        i = lvli(mheuri);
        cum(i) = cum(i)+mheur;
        next(i) = nbr(i,jheur(mheuri));
        done(i) = true;
    end;
end;
vine = zeros(1,nrow*ncol);
[scr,vine(1)] = max(cum);
vlen = 1;
while (next(vine(vlen))>0)
    vine(vlen+1) = next(vine(vlen));
    vlen = vlen+1;
end;
% add on at end
i = vine(vlen);
done(i) = false;
lvl = board(i);
mv = (vals(nbr(i,:))==lvl).*done(nbr(i,:));
while any(mv)
    vine(vlen+1) = nbr(i,find(mv,1));
    vlen = vlen+1;
    i = vine(vlen);
    done(i) = false;
    lvl = board(i);
    mv = (vals(nbr(i,:))==lvl).*done(nbr(i,:));
end;
vine = vine(1:vlen);

% now do augmentation
j = 1;
used = false(nrow,ncol);
used(vine) = true;
while (j < vlen)
    i1 = vine(j);
    i2 = vine(j+1);
    x1 = xg(i1);
    x2 = xg(i2);
    y1 = yg(i1);
    %y2 = yg(i2);
    if (x1==x2)
        % vertical move
        if (x1>1)&&(~used(i1-nrow))&&(~used(i2-nrow))&&(board(i1)<=board(i1-nrow))&&(board(i1-nrow)<=board(i2-nrow))&&(board(i2-nrow)<=board(i2))
            % augment to the left
            vlen = vlen+2;
            vine = [vine(1:j) i1-nrow i2-nrow vine(j+1:end)];
            used(i1-nrow) = true;
            used(i2-nrow) = true;
        elseif (x1<ncol)&&(~used(i1+nrow))&&(~used(i2+nrow))&&(board(i1)<=board(i1+nrow))&&(board(i1+nrow)<=board(i2+nrow))&&(board(i2+nrow)<=board(i2))
            % augment to the right
            vlen = vlen+2;
            vine = [vine(1:j) i1+nrow i2+nrow vine(j+1:end)];
            used(i1+nrow) = true;
            used(i2+nrow) = true;
        else
            % no augmentation; move on
            j = j+1;
        end;
    else
        % horizontal move
        if (y1>1)&&(~used(i1-1))&&(~used(i2-1))&&(board(i1)<=board(i1-1))&&(board(i1-1)<=board(i2-1))&&(board(i2-1)<=board(i2))
            % augment to the north
            vlen = vlen+2;
            vine = [vine(1:j) i1-1 i2-1 vine(j+1:end)];
            used(i1-1) = true;
            used(i2-1) = true;
        elseif (y1<nrow)&&(~used(i1+1))&&(~used(i2+1))&&(board(i1)<=board(i1+1))&&(board(i1+1)<=board(i2+1))&&(board(i2+1)<=board(i2))
            % augment to the south
            vlen = vlen+2;
            vine = [vine(1:j) i1+1 i2+1 vine(j+1:end)];
            used(i1+1) = true;
            used(i2+1) = true;
        else
            % no augmentation; move on
            j = j+1;
        end;
    end;
end;
scr = sum(board(vine));
end


function [moves,vine,scr] = dofilter(board,limit)
[nrow,ncol] = size(board);
ntot = nrow*ncol;
moves = zeros(limit,2);
nmove = 0;
pos = reshape(1:ntot,nrow,ncol);
m1 = median(board,1);
m2 = median(board,2);
if (sum(sign(diff(m1)))>0.67*ncol)||(sum(sign(diff(m2)))>0.67*nrow)||(sum(sign(diff(m2)))<-0.67*nrow)||(sum(sign(diff(m2)))<-0.67*nrow)
    d = diff(board,1,2);
    mb = median(board,2);
    md = median(d,2);
    rmb = repmat(mb,1,ncol);
    rmd = repmat(md,1,ncol);
    guess = rmb+kron(md,-(ncol-1)/2:(ncol-1)/2);
    bad = (abs(board-guess)>20*abs(rmd));
    for i = nrow:-1:1
        offset = 0;
        for j = ceil(ncol/2):ncol
            if bad(i,j)
                offset = offset+1;
            else
                moves(nmove+1:nmove+offset,:) = [pos(i,j:-1:j-offset+1)' pos(i,j-1:-1:j-offset)'];
                nmove = nmove+offset;
                if (nmove>=limit)
                    break;
                end;
            end;
        end;
        offset = 0;
        for j = ceil(ncol/2):-1:1
            if bad(i,j)
                offset = offset+1;
            else
                moves(nmove+1:nmove+offset,:) = [pos(i,j:1:j+offset-1)' pos(i,j+1:1:j+offset)'];
                nmove = nmove+offset;
                if (nmove>=limit)
                    break;
                end;
            end;
        end;
        if (nmove>=limit)
            break;
        end;
    end;
    moves = moves(1:min(nmove,limit),:);
    for i = 1:size(moves,1)
        board(moves(i,2)) = board(moves(i,1));
        board(moves(i,1)) = 0;
        %imagesc(board);
        %axis image;
        %drawnow;
    end;
    [vine,scr] = trickle2(board);
else
    scr = -inf;
    moves = [];
    vine = [];
end;
end

%% alfi
function [moves, vine, gain] = alfi(board, limit)
[moves, vine, gain,newboard] = solver1(board, limit);
% disp(gain)
% imagesc(newboard);hold on;[i1,i2]=ind2sub(size(board),vine);plot(i2,i1,'k.-'); hold off;
end

function [moves, vine, gain, board] = solver1(board, limit)

moves=[];
ab=accumarray(1+board(:),1);
cab=flipud(cumsum(flipud(ab)))-ab(end);
[m,n]=size(board);
board=[zeros(1,n+2);[zeros(m,1),board,zeros(m,1)];zeros(1,n+2)];
[m,n]=size(board);
boardab=cab(1+board);
if max(ab(2:end))>1,
    % compute same-valued clusters
    [BoardCluster,ClusterValue,IdxList,IdxSegments,Nclusters]=connected(board);
    ClusterSize=diff(IdxSegments);
    ClusterNeighb=neighb(BoardCluster,board,Nclusters);
    % search between-clusters
    ClustersOrder = bellman(ClusterNeighb,ClusterValue'.*sqrt(ClusterSize'));
    % search within-clusters
    iC=bellman_postprocess(ClustersOrder,IdxList,IdxSegments,BoardCluster);
else
    % search between-pixels
    iC = bellman_pixel(board);
end
% post-processing moves
if limit>0
    [board,iC,moves,limit] = postprocess(board,iC,moves,limit,boardab);
    [iC1,iC2]=ind2sub([m,n],moves);
    moves=sub2ind([m-2,n-2],iC1-1,iC2-1);
end

[iC1,iC2]=ind2sub([m,n],iC);
vine=sub2ind([m-2,n-2],iC1-1,iC2-1);
vine=vine(end:-1:1);
board=board(2:end-1,2:end-1);
gain=sum(board(vine));

end

function [board,tiC,moves,limit] = postprocess(board,tiC,moves,limit,boardab)
[m,n]=size(board);
nC=numel(tiC);
if nC>1&&limit>0, %%%
    % grow laterally
    d=abs(diff(tiC))==1;
    E=[m*~d+d; m*d+~d];
    BorderUp={repmat(1+(0:n-1)*m,[m,1]),repmat((1:m)',[1,n])};
    BorderDown={repmat(m+(0:n-1)*m,[m,1]),repmat(m*(n-1)+(1:m)',[1,n])};
    tP=board>0;
    tP(tiC)=false;
    for n1=1:1e3,
        [board,moves,tiC,tP,E,limit,ok]=postprocess_lateral(board,moves,tiC,tP,E,limit,BorderDown,BorderUp,boardab);
        if ~ok, break; end
    end
end
if limit>0
    % grow end-points
    tP=board>0;
    tP(tiC)=false;
    for n1=1:1e3,
        [board,moves,tiC,tP,limit,ok]=postprocess_endpoint(board,moves,tiC,tP,limit,1);
        if ~ok, break; end
    end
    for n1=1:1e3,
        [board,moves,tiC,tP,limit,ok]=postprocess_endpoint(board,moves,tiC,tP,limit,2);
        if ~ok, break; end
    end
end
end

function [board,moves,tiC,tP,E,limit,ok]=postprocess_lateral(board,moves,tiC,tP,E,limit,BorderDown,BorderUp,boardab)
[m,n]=size(board);
ok=0;
% D=Dijkstra(tP,tiC(1));
% maskout=tP&~isinf(D)&board>board(tiC(1));
% nmaskout=sum(D(maskout));
for ndir=1:4
    K=size(tiC,2);
    switch(ndir)
        case 1, idx1=find(tP(tiC(1:2:K-1)+E(2,1:2:K-1))&tP(tiC(2:2:K)+E(2,1:2:K-1)));
        case 2, idx1=find(tP(tiC(1:2:K-1)-E(2,1:2:K-1))&tP(tiC(2:2:K)-E(2,1:2:K-1)));
        case 3, idx1=find(tP(tiC(2:2:K-1)+E(2,2:2:K-1))&tP(tiC(3:2:K)+E(2,2:2:K-1)));
        case 4, idx1=find(tP(tiC(2:2:K-1)-E(2,2:2:K-1))&tP(tiC(3:2:K)-E(2,2:2:K-1)));
    end
    %for n2=1:numel(idx1),
    for n2=numel(idx1):-1:1, %%%
        switch(ndir)
            case 1,
                tempA=tiC(2*idx1(n2)-1)+E(2,2*idx1(n2)-1):E(2,2*idx1(n2)-1):BorderDown{1+(E(2,2*idx1(n2)-1)==m)}(tiC(2*idx1(n2)-1));
                tempB=tiC(2*idx1(n2))+E(2,2*idx1(n2)-1):E(2,2*idx1(n2)-1):BorderDown{1+(E(2,2*idx1(n2)-1)==m)}(tiC(2*idx1(n2)));
            case 2,
                tempA=tiC(2*idx1(n2)-1)-E(2,2*idx1(n2)-1):-E(2,2*idx1(n2)-1):BorderUp{1+(E(2,2*idx1(n2)-1)==m)}(tiC(2*idx1(n2)-1));
                tempB=tiC(2*idx1(n2))-E(2,2*idx1(n2)-1):-E(2,2*idx1(n2)-1):BorderUp{1+(E(2,2*idx1(n2)-1)==m)}(tiC(2*idx1(n2)));
            case 3,
                tempA=tiC(2*idx1(n2))+E(2,2*idx1(n2)):E(2,2*idx1(n2)):BorderDown{1+(E(2,2*idx1(n2))==m)}(tiC(2*idx1(n2)));
                tempB=tiC(2*idx1(n2)+1)+E(2,2*idx1(n2)):E(2,2*idx1(n2)):BorderDown{1+(E(2,2*idx1(n2))==m)}(tiC(2*idx1(n2)+1));
            case 4,
                tempA=tiC(2*idx1(n2))-E(2,2*idx1(n2)):-E(2,2*idx1(n2)):BorderUp{1+(E(2,2*idx1(n2))==m)}(tiC(2*idx1(n2)));
                tempB=tiC(2*idx1(n2)+1)-E(2,2*idx1(n2)):-E(2,2*idx1(n2)):BorderUp{1+(E(2,2*idx1(n2))==m)}(tiC(2*idx1(n2)+1));
        end
        idxZ=find(~tP(tempA),1,'first');
        %idxZ=find(~tP(tempA)|maskout(tempA),1,'first');
        switch(ndir)
            case {1,2},
                %idxA=find(board(tiC(2*idx1(n2)-1))>=board(tempA)&board(tempA)>=board(tiC(2*idx1(n2))),1,'first');
                [nill,idxA]=max((1./(abs(board(tiC(2*idx1(n2)-1))-board(tempA))+1)).*(board(tiC(2*idx1(n2)-1))>=board(tempA)&board(tempA)>=board(tiC(2*idx1(n2)))));
                if ~nill,idxA=[]; end
            case {3,4},
                %idxA=find(board(tiC(2*idx1(n2)))>=board(tempA)&board(tempA)>=board(tiC(2*idx1(n2)+1)),1,'first');
                [nill,idxA]=max((1./(abs(board(tiC(2*idx1(n2)))-board(tempA))+1)).*(board(tiC(2*idx1(n2)))>=board(tempA)&board(tempA)>=board(tiC(2*idx1(n2)+1))));
                if ~nill,idxA=[]; end
        end
        if ~isempty(idxA)&&idxA<idxZ
            idxZ=find(~tP(tempB),1,'first');
            %idxZ=find(~tP(tempB)|maskout(tempA),1,'first');
            switch(ndir)
                case {1,2},
                    idxB=find(board(tempA(idxA))>=board(tempB)&board(tempB)>=board(tiC(2*idx1(n2))),1,'first');
                case {3,4},
                    idxB=find(board(tempA(idxA))>=board(tempB)&board(tempB)>=board(tiC(2*idx1(n2)+1)),1,'first');
            end
            if ~isempty(idxB)&&idxB<idxZ&&idxA+idxB-2<=limit,%-nmaskout,%.0*boardab(tiC(1))*(n+m)/2, %%%
                ok=1;
                switch(ndir)
                    case 1, tiC=[tiC(1:2*idx1(n2)-1),tiC(2*idx1(n2)-1)+E(2,2*idx1(n2)-1), tiC(2*idx1(n2))+E(2,2*idx1(n2)-1), tiC(2*idx1(n2):end)];
                    case 2, tiC=[tiC(1:2*idx1(n2)-1),tiC(2*idx1(n2)-1)-E(2,2*idx1(n2)-1), tiC(2*idx1(n2))-E(2,2*idx1(n2)-1), tiC(2*idx1(n2):end)];
                    case 3, tiC=[tiC(1:2*idx1(n2)),tiC(2*idx1(n2))+E(2,2*idx1(n2)), tiC(2*idx1(n2)+1)+E(2,2*idx1(n2)), tiC(2*idx1(n2)+1:end)];
                    case 4, tiC=[tiC(1:2*idx1(n2)),tiC(2*idx1(n2))-E(2,2*idx1(n2)), tiC(2*idx1(n2)+1)-E(2,2*idx1(n2)), tiC(2*idx1(n2)+1:end)];
                end
                E=[E(:,1:2*idx1(n2)-2),E([2,1],2*idx1(n2)-1), E([1,2],2*idx1(n2)-1), E([2,1],2*idx1(n2)-1), E(:,2*idx1(n2):end)];
                tP(tiC)=false;
                newmoveA=[tempA(2:idxA)',tempA(1:idxA-1)'];
                newmoveB=[tempB(2:idxB)',tempB(1:idxB-1)'];
                board(newmoveA(:,2))=board(tempA(idxA));
                board(newmoveB(:,2))=board(tempB(idxB));
                board(newmoveA(:,1))=0;
                board(newmoveB(:,1))=0;
                moves=[moves;flipud(newmoveA);flipud(newmoveB)];
                limit=limit-size(newmoveA,1)-size(newmoveB,1);
                %break;
            end
        end
    end
end
end

function [board,moves,tiC,tP,limit,ok]=postprocess_endpoint(board,moves,tiC,tP,limit,type)
[m,n]=size(board);
ok=1;
switch(type)
    case 1
        mask=tP&board<2*board(tiC(1));
        D=Dijkstra(mask,tiC(1));
        idx1=find(mask&D-1<=limit&board>=board(tiC(1)));
        if isempty(idx1),
            mask=tP;
            D=Dijkstra(mask,tiC(1));
            idx1=find(mask&D-1<=limit&board>=board(tiC(1)));
            if isempty(idx1),ok=0;end
        end
        if ok
            [nill,idx2]=min(board(idx1)+1e-3*D(idx1)); %%%t3
            idx0=idx1(idx2);
            limit=limit-(D(idx0)-1);
            for n2=D(idx0)-1:-1:1
                if      D(idx0+1)==n2, tidx0=idx0+1;
                elseif  D(idx0-1)==n2, tidx0=idx0-1;
                elseif  D(idx0+m)==n2, tidx0=idx0+m;
                elseif  D(idx0-m)==n2, tidx0=idx0-m;
                end
                moves=[moves; idx0, tidx0];
                board(tidx0)=board(idx0);
                board(idx0)=0;
                idx0=tidx0;
            end
            tiC=[idx0,tiC];
            tP(idx0)=false;
        end
    case 2
        mask=tP;
        D=Dijkstra(mask,tiC(end));
        idx1=find(mask&D-1<=limit&board<=board(tiC(end))&board>0);
        if isempty(idx1),ok=0;end
        if ok
            [nill,idx2]=min(-board(idx1)+1e-3*D(idx1)); %%%t3
            idx0=idx1(idx2);
            limit=limit-(D(idx0)-1);
            for n2=D(idx0)-1:-1:1
                if      D(idx0+1)==n2, tidx0=idx0+1;
                elseif  D(idx0-1)==n2, tidx0=idx0-1;
                elseif  D(idx0+m)==n2, tidx0=idx0+m;
                elseif  D(idx0-m)==n2, tidx0=idx0-m;
                end
                moves=[moves; idx0, tidx0];
                board(tidx0)=board(idx0);
                board(idx0)=0;
                idx0=tidx0;
            end
            tiC=[tiC,idx0];
            tP(idx0)=false;
        end
end
end

function iC = bellman_pixel(board)
[m,n]=size(board);
cD=board;
IDX=zeros(m,n);
touched=board>0;
valid=touched;
while any(touched(:))
    touch=find(touched);
    touched(touch)=false;
    for dx=[1,-1,m,-m],
        cDnew=board(touch+dx)+cD(touch);
        idx=find(board(touch+dx)>board(touch)&cDnew>cD(touch+dx));
        touched(touch(idx)+dx)=valid(touch(idx)+dx);
        cD(touch(idx)+dx)=cDnew(idx);
        IDX(touch(idx)+dx)=touch(idx);
    end
end
[nill,idx]=max(cD(:));
iC=zeros(1,m*n);
for n1=1:m*n,
    if idx>0
        iC(n1)=idx;
        idx=IDX(idx);
    else break;
    end
end
iC=iC(1:n1-1);
end

function iC = bellman(C,D)
N=size(C,1);
c=cell(1,N); for n1=1:N,c{n1}=find(C(:,n1)); end
%C=full(C);
IDX=zeros(N,1);
cD=D;
touched=true(1,N);
while any(touched)
    for touch=find(touched),
        touched(touch)=false;
        idx=c{touch};%find(C(:,touch));
        cDnew=D(idx)+cD(touch);
        idx2=find(cD(idx)<cDnew);
        if ~isempty(idx2)
            cD(idx(idx2))=cDnew(idx2);
            touched(idx(idx2))=true;
            IDX(idx(idx2))=touch;
        end
    end
end
[nill,idx]=max(cD);
iC=zeros(N,1);
for n1=1:N,
    if idx>0
        iC(n1)=idx;
        idx=IDX(idx);
    else break;
    end
end
iC=iC(1:n1-1);
end

function  iC=bellman_postprocess(ClustersOrder,IdxList,IdxSegments,BoardCluster)
[m,n]=size(BoardCluster);
iC=[];
for n1=1:numel(ClustersOrder)
    idx=IdxList(IdxSegments(ClustersOrder(n1)):IdxSegments(ClustersOrder(n1)+1)-1);
    if numel(idx)==1, iC=[iC,idx(1)];
    else
        % longest shortest-path
        ThisCluster=false(m,n);
        ThisCluster(idx)=true;
        D=Dijkstra(ThisCluster,iC);
        if n1==numel(ClustersOrder)
            temp=D(ThisCluster);
        else
            temp=D(ThisCluster).*(BoardCluster([false(1,n);ThisCluster(1:m-1,:)])==ClustersOrder(n1+1)|BoardCluster([ThisCluster(2:m,:);false(1,n)])==ClustersOrder(n1+1)|BoardCluster([false(m,1),ThisCluster(:,1:n-1)])==ClustersOrder(n1+1)|BoardCluster([ThisCluster(:,2:n),false(m,1)])==ClustersOrder(n1+1));
        end
        [nill,idx0]=max(temp(:));
        idx0=idx(idx0);
        E=zeros(2,D(idx0)-1);
        tiC=zeros(1,D(idx0));
        tiC(end)=idx0;
        for n2=D(idx0)-1:-1:1
            if      D(idx0+1)==n2, idx0=idx0+1; E(1,n2)=1;E(2,n2)=m;
            elseif  D(idx0-1)==n2, idx0=idx0-1; E(1,n2)=1;E(2,n2)=m;
            elseif  D(idx0+m)==n2, idx0=idx0+m; E(1,n2)=m;E(2,n2)=1;
            elseif  D(idx0-m)==n2, idx0=idx0-m; E(1,n2)=m;E(2,n2)=1;
            end
            tiC(n2)=idx0;
        end
        % now grow it
        tP=ThisCluster;
        tP(tiC)=false;
        ok=1;
        while ok
            ok=0;
            K=size(tiC,2);
            idx1=find(tP(tiC(1:2:K-1)+E(2,1:2:K-1))&tP(tiC(2:2:K)+E(2,1:2:K-1)));
            for n2=numel(idx1):-1:1
                ok=1;
                tiC=[tiC(1:2*idx1(n2)-1),tiC(2*idx1(n2)-1)+E(2,2*idx1(n2)-1), tiC(2*idx1(n2))+E(2,2*idx1(n2)-1), tiC(2*idx1(n2):end)];
                E=[E(:,1:2*idx1(n2)-2),E([2,1],2*idx1(n2)-1), E([1,2],2*idx1(n2)-1), E([2,1],2*idx1(n2)-1), E(:,2*idx1(n2):end)];
                tP(tiC)=false;
            end
            K=size(tiC,2);
            idx1=find(tP(tiC(1:2:K-1)-E(2,1:2:K-1))&tP(tiC(2:2:K)-E(2,1:2:K-1)));
            for n2=numel(idx1):-1:1
                ok=1;
                tiC=[tiC(1:2*idx1(n2)-1),tiC(2*idx1(n2)-1)-E(2,2*idx1(n2)-1), tiC(2*idx1(n2))-E(2,2*idx1(n2)-1), tiC(2*idx1(n2):end)];
                E=[E(:,1:2*idx1(n2)-2),E([2,1],2*idx1(n2)-1), E([1,2],2*idx1(n2)-1), E([2,1],2*idx1(n2)-1), E(:,2*idx1(n2):end)];
                tP(tiC)=false;
            end
            K=size(tiC,2);
            idx1=find(tP(tiC(2:2:K-1)+E(2,2:2:K-1))&tP(tiC(3:2:K)+E(2,2:2:K-1)));
            for n2=numel(idx1):-1:1
                ok=1;
                tiC=[tiC(1:2*idx1(n2)),tiC(2*idx1(n2))+E(2,2*idx1(n2)), tiC(2*idx1(n2)+1)+E(2,2*idx1(n2)), tiC(2*idx1(n2)+1:end)];
                E=[E(:,1:2*idx1(n2)-1),E([2,1],2*idx1(n2)), E([1,2],2*idx1(n2)), E([2,1],2*idx1(n2)), E(:,2*idx1(n2)+1:end)];
                tP(tiC)=false;
            end
            K=size(tiC,2);
            idx1=find(tP(tiC(2:2:K-1)-E(2,2:2:K-1))&tP(tiC(3:2:K)-E(2,2:2:K-1)));
            for n2=numel(idx1):-1:1
                ok=1;
                tiC=[tiC(1:2*idx1(n2)),tiC(2*idx1(n2))-E(2,2*idx1(n2)), tiC(2*idx1(n2)+1)-E(2,2*idx1(n2)), tiC(2*idx1(n2)+1:end)];
                E=[E(:,1:2*idx1(n2)-1),E([2,1],2*idx1(n2)), E([1,2],2*idx1(n2)), E([2,1],2*idx1(n2)), E(:,2*idx1(n2)+1:end)];
                tP(tiC)=false;
            end
        end
        
        iC=[iC,tiC];
    end
end
end

function B = neighb(A,V,nA)
[m,n] = size(A);
B=sparse(A(:,1:n-1),A(:,2:n),double(V(:,1:n-1)>V(:,2:n)),nA,nA)+sparse(A(:,2:n),A(:,1:n-1),double(V(:,2:n)>V(:,1:n-1)),nA,nA)+sparse(A(1:m-1,:),A(2:m,:),double(V(1:m-1,:)>V(2:m,:)),nA,nA)+sparse(A(2:m,:),A(1:m-1,:),double(V(2:m,:)>V(1:m-1,:)),nA,nA);
B=B>0;
end

function [B,C,p,r,nc] = connected(A)
[m,n] = size(A);
N = m*n ;
K = reshape (1:N, m, n) ;
V = A(:,1:n-1) == A(:,2:n);
H = A(1:m-1,:) == A(2:m,:);
G = sparse(K([V,false(m,1)]),K([false(m,1),V]),1,N,N)+sparse(K([H; false(1,n)]),K([false(1,n); H]),1,N,N);
G = G + G' + speye(N);
[p q r] = dmperm(G);
nc = numel(r) - 1;
C = A(p(r(1:nc)));
B = ones(m, n);
for a = 2:nc
    B(p(r(a):r(a+1)-1)) = a;
end
end

function D = Dijkstra(A,i1)
[m,n]=size(A);
mn=m*n;
D=inf(m,n);
D(~A)=0;
if isempty(i1),
    i1=find(A,1,'first');
    mode=0;
else
    i1=i1(end);
    mode=1;
    % 	if A(i1+1); i1=i1+1;
    %     elseif A(i1-1); i1=i1-1;
    %     elseif A(i1+m); i1=i1+m;
    %     elseif A(i1-m); i1=i1-m;
    %     end
end
D(i1)=1;
for n1=2:mn,
    idx1=D(i1+1)>n1;
    idx2=D(i1-1)>n1;
    idx3=D(i1+m)>n1;
    idx4=D(i1-m)>n1;
    %i1=unique([i1(idx1)+1,i1(idx2)-1,i1(idx3)+m,i1(idx4)-m]);
    X=false(m,n);X(i1(idx1)+1)=true;X(i1(idx2)-1)=true;X(i1(idx3)+m)=true;X(i1(idx4)-m)=true;
    i1=find(X);
    if isempty(i1), break; end
    D(i1)=n1;
end
if mode, D(D>0)=D(D>0)-1; end
end

%% kurt
function [vine, bestsom] = kurt(board)
moves = [];
[bestsom,bestvine] = max(board(:));

A=board;
[m,n]=size(board);

inhoud = unique(A(:))';
heur = sparse(0);
for i=inhoud;
    heur(i+1) = sum(sum(A(A<=i)));
end
H=full(heur(A+1));

if any(size(H) ~= size(A))
    H=H';
end


B=A;


IND =  reshape(1:m*n,[m,n]);
C = num2cell(IND);

Hmod = zeros(m,n);


updated = true(size(A));

iteration = 0;
while 1
    G=B+H;
    G=G.*updated;
    [maxg,startindex] = max(G(:));
    
    [bestsom] = max(B(:));
    
    if maxg<=bestsom
        break
    end
    
    if maxg==0
        break
    end
    
    iteration = iteration+1;
    i = mod(startindex,m);
    j = ceil(startindex/m);
    if i==0
        i=m;
    end
    
    
    if 1<i && A(i-1,j)<=A(i,j)
        verplaatsindex = IND(i-1,j);
        if ~any(verplaatsindex==C{i,j})
            [B(verplaatsindex),which] = max([B(verplaatsindex),B(i,j)+A(verplaatsindex)]);
            if which == 2
                updated(verplaatsindex)=true;
                C{verplaatsindex} = [C{i,j} verplaatsindex];
                %                 Hmod(verplaatsindex) = sum(A(C{verplaatsindex})==A(verplaatsindex));
            end
        end
    end
    
    if i<size(A,1) && A(i+1,j)<=A(i,j)
        verplaatsindex = IND(i+1,j);
        if ~any(verplaatsindex==C{i,j})
            [B(verplaatsindex),which] = max([B(verplaatsindex),B(i,j)+A(verplaatsindex)]);
            if which == 2
                updated(verplaatsindex)=true;
                C{verplaatsindex} = [C{i,j} verplaatsindex];
                %                 Hmod(verplaatsindex) = sum(A(C{verplaatsindex})==A(verplaatsindex));
            end
        end
    end
    if j<size(A,2) && A(i,j+1)<=A(i,j)
        verplaatsindex = IND(i,j+1);
        if ~any(verplaatsindex==C{i,j})
            [B(verplaatsindex),which] = max([B(verplaatsindex),B(i,j)+A(verplaatsindex)]);
            if which == 2
                updated(verplaatsindex)=true;
                C{verplaatsindex} = [C{i,j} verplaatsindex];
                %                 Hmod(verplaatsindex) = sum(A(C{verplaatsindex})==A(verplaatsindex));
            end
        end
    end
    if 1<j && A(i,j-1)<=A(i,j)
        verplaatsindex = IND(i,j-1);
        if ~any(verplaatsindex==C{i,j})
            [B(verplaatsindex),which] = max([B(verplaatsindex),B(i,j)+A(verplaatsindex)]);
            if which == 2
                updated(verplaatsindex)=true;
                C{verplaatsindex} = [C{i,j} verplaatsindex];
                %                 Hmod(verplaatsindex) = sum(A(C{verplaatsindex})==A(verplaatsindex));
            end
        end
    end
    
    updated(i,j)=false;
    
    
    
end

[bestsom,index]=max(B(:));
vine = fliplr(C{index});

end