function moves = solver(sequence,target,budget)
moves = ones(0,4);
s=sequence;
t=target;
score=sum(abs(s-t));
moves=algo_sortdiff(s,t,budget,2,1);
score1=algo_findscore(s,t,moves,budget);
if score1/score>0.2
moves2=asolver(s, t, budget);
score2=algo_findscore(s,t,moves2,budget);
if score2<score1
moves=moves2;
end
end
%-------------------------------------------------------------------
function moves=algo_sortdiff(s,t,budget,I,J)
moves = ones(0,4);
for iter=1:budget
value=[];
diff=abs(s-t);
[val,index]=sort(diff,'descend');
% first do simple
for i=1:I
diff2=abs(s(index(i))-t);
[val2,index2]=sort(diff2,'ascend');
for j=1:J
localmoves=algo_removefromto(s,t,index(i),index2(j));
if numel(localmoves)>0
value(end+1,[1:5])=[localmoves(1) localmoves(2) localmoves(3) localmoves(4) algo_findscore(s,t,localmoves,budget)];
end
end
end
if ~isempty(value)
[mv,mi]=min(value(:,5));
localmoves=value(mi,1:4);
moves=[moves;localmoves];
if ~isempty(localmoves)
s=algo_doMoves(s,localmoves,size(localmoves,1));
end
end
end
%-----------------------------------------------------------------
function moves=algo_removefromto(seq2,tar2,index,index2)
L = length(seq2);
tempmoves=algo_removefromtosub1(seq2,index,index2);
if isempty(tempmoves)
moves=ones(0,4);
else
n1 = size(tempmoves,1);
sn1=zeros(n1,L);
tscore1=zeros(n1,L);
for i=1:n1
sn1(i,:)=algo_splice(seq2,tempmoves(i,1),tempmoves(i,2),tempmoves(i,3),tempmoves(i,4));
tscore1(i,:)=sn1(i,:)-tar2;
end
tscore1=sum(abs(tscore1),2);
[minval,minindex]=min(tscore1);
moves=tempmoves(minindex,:);
end
%-----------------------------------------------------------------
function tempmoves=algo_removefromtosub1(seq2,index,index2)
if index<index2
N=abs(index2-index)+1;
tempmoves = zeros(N,4);
for i=1:N
tempmoves(i,:)=[i index index2-i+1 1];
end
N=min([index-1,numel(seq2)-index2])+1;
tempmoves2 = zeros(N-1,4);
for i=2:N
tempmoves2(i-1,:)=[i index-i+1 index2 1];
end
N=index;
tempmoves3 = zeros(N-1,4);
for i=2:N
tempmoves3(i-1,:)=[i index-i+1 index2-i+1 0];
end
N=numel(seq2)-index2+1;
tempmoves4 = zeros(N-1,4);
for i=2:N
tempmoves4(i-1,:)=[i index index2 0];
end
elseif index>=index2
N=abs(index2-index)+1;
tempmoves = zeros(N,4);
for i=1:N
tempmoves(i,:)=[i index-i+1 index2 1];
end
N=min([numel(seq2)-index,index2-1])+1;
tempmoves2 = zeros(N-1,4);
for i=2:N
tempmoves2(i-1,:)=[i index index2-i+1 1];
end
N=numel(seq2)-index+1;
tempmoves3 = zeros(N-1,4);
for i=2:N
tempmoves3(i-1,:)=[i index index2 0];
end
N=index2;
tempmoves4 = zeros(N-1,4);
for i=2:N
tempmoves4(i-1,:)=[i index-i+1 index2-i+1 0];
end
end
tempmoves=[tempmoves;tempmoves2;tempmoves3;tempmoves4];
%----------------------------------------------------------------
function b = algo_splice(a, len, ai, bi, flip)
b = zeros(size(a));
aSpan = ai + (1:len) - 1;
bSpan = aSpan - ai + bi;
if flip
b(bSpan) = fliplr(a(aSpan));
else
b(bSpan) = a(aSpan);
end
a(aSpan) = 0;
b(~b) = a(~~a);
%----------------------------------------------------------------
function seq = algo_doMoves(seq,moves,budget)
if size(moves,1)>budget
moves = moves(1:budget,:);
end
n = numel(seq);
moves(:,1) = max(1,min(n,moves(:,1)));
moves(:,2) = max(1,min(n-moves(:,1)+1,moves(:,2)));
moves(:,3) = max(1,min(n-moves(:,1)+1,moves(:,3)));
moves(:,4) = moves(:,4) ~= 0;
for i =1:size(moves,1)
le = moves(i,1);
ai = moves(i,2);
bi = moves(i,3);
f = moves(i,4);
if f
if ai>bi
seq = seq([1:bi-1 ai+le-1:-1:ai bi:ai-1 ai+le:end]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai+le-1:-1:ai bi+le:end]);
end
else
if ai>bi
seq = seq([1:bi-1 ai:ai+le-1 bi:ai-1 ai+le:end]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai:ai+le-1 bi+le:end]);
end
end
end
%-----------------------------------------------------------
function score=algo_findscore(s,t,moves,budget)
sn=s;
if ~isempty(moves)
sn=algo_doMoves(s,moves,budget);
end
score=sum(abs(sn-t));
%%###################################################################
function moves = asolver(sequence,target,budget)
last = 2;
moves = [];
for b = 1:budget
moves3 = oldconvolutedsolver(sequence,target,1);
score3=findscore(doMoves(sequence,moves,length(sequence)),target,1,moves3);
moves2 = new2solver(doMoves(sequence,moves,length(sequence)),target,1);
score2 = findscore(doMoves(sequence,moves,length(sequence)),target,1,moves2);
if score3<score2,
moves2 = moves3;
score2=score3;
end
moves3 = anothersolver(doMoves(sequence,moves,length(sequence)),target,1);
score3 = findscore(doMoves(sequence,moves,length(sequence)),target,1,moves3);
if score3<score2,
moves2 = moves3;
score3=score2;
end
moves3 = HannesNaudiersolver(doMoves(sequence,moves,length(sequence)),target,1);
score3 = findscore(doMoves(sequence,moves,length(sequence)),target,1,moves3);
if score3<score2,
moves2 = moves3;
score3=score2;
end
moves = [moves; moves2];
sequence = doMoves(sequence,moves2,length(sequence))
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [moves,score] = swapper(sequence,target,budget,rfac)
moves = zeros(budget,4);
n=numel(sequence);
runningScore=zeros(budget,1);
targrows = repmat(target,n,1);
targcols=targrows';
seqrows = repmat(sequence,n,1);
altered=false(1,n);
i=1;
while i<budget
seqrows(:,altered) = repmat(sequence(altered),n,1);
seqcols = seqrows';
% Delta(i,j) = score increase resulting from swap of i with j
Delta = abs(targcols-seqrows) + abs(seqcols-targrows) - abs(targrows-seqrows) - abs(targcols-seqcols) ;
if rfac>0
Delta = Delta.*(1 + rfac*rand(n));
end
[dmin,i1] = min(Delta);
[dmin,j1] = min(dmin);
i1=i1(j1);
ij =[i1 j1]; i1 = min(ij); j1=max(ij);
% now see if it pays to include neighbours so that we swap a band i1 to i2 with a band j1 to j2
dcum1 = Delta(i1,j1);
i2=i1;j2=j1;
% go R (i2>i1) first:
keepgoing = (j2<n) * (i2<j1-1);
while keepgoing
d = Delta(i2+1,j2+1);
if d<0
i2=i2+1;
j2=j2+1;
dcum1 = dcum1+d; % increment cumulative delta
end
keepgoing = ~( (d>=0) + (j2==n) + (i2>=j1-1) );
end
% go L (decreasing i1) :
keepgoing = (i1>1 )*(i2<j1-1);
while keepgoing
d = Delta(i1-1,j1-1);
if d<0
i1=i1-1;
j1=j1-1;
dcum1 = dcum1+d; % increment cumulative delta
end
keepgoing= ~((d>=0) + (i1==1) + (i2>=j1-1) );
end
L = i2-i1+1;
moves1 = [L j1 i1 0; L i1+L j1 0];
% now repeat the exercise for a flipped sequence
dcum2 = Delta(i1,j1);
i2=i1;j2=j1;
% go R (increasing i2) first:
keepgoing = i2<j1-1;
while keepgoing
% d = swapdelta(sequence,target,i2+1,j1-1);
d = Delta(i2+1,j1-1);
if d<0
i2=i2+1;
j1=j1-1;
dcum2 = dcum2+d; % increment cumulative delta
end
keepgoing = ~( (d>=0) + (i2>=j1-1) );
end
% go L (decreasing i1) :
keepgoing = (i1>1)*(i2<j1-1)*( j2<n);
while keepgoing
% d = swapdelta(sequence,target,i1-1,j2+1);
d = Delta(i1-1,j2+1);
if d<0
i1=i1-1;
j2=j2+1;
dcum2 = dcum2+d; % increment cumulative delta
end
keepgoing = ~((d>=0) + (i1==1) + (j2==n));
end
L = i2-i1+1;
moves2 = [L j1 i1 1; L i1+L j1 1];
if dcum1<=dcum2
moves(i:i+1,:) = moves1;
else
moves(i:i+1,:) = moves2;
end
L= moves(i,1); j1 = moves(i,2); i1=moves(i,3); % note the affected genes
sequence1(1,:) = doMoves(sequence,moves(i,:),n);
score1(1) = sum(abs(sequence1(1,:)-target));
sequence1(2,:) = doMoves(sequence1(1,:),moves(i+1,:),n);
score1(2) = sum(abs(sequence1(2,:)-target));
sequence1(3,:) = doMoves(sequence,moves(i+1,:),n);
score1(3) = sum(abs(sequence1(3,:)-target));
[minscore1,ims] = min(score1);
sequence=sequence1(ims,:);
runningScore(i:i+1) = score1(1:2);
if ims==3
moves(i,:) = moves(i+1,:);
runningScore(i) = score1(3);
i=i+1;
else
i=i+ims;
end
% flag the altered bands in the sequence:
altered = false(1,n);
altered(i1:i1+L-1) = true;
altered(j1:j1+L-1) = true;
end
if i==budget
s = abs(sequence-target);
[maxs,imaxs]=max(s); % worst current gene
[mins,imins1] = min(abs(target-sequence(imaxs))) ; % best place to put the worst gene
move1 = [1 imaxs imins1 false];
seq1 = doMoves(sequence,move1,n);
s1= sum(abs(seq1-target));
[mins,imins2] = min(abs(target(imaxs)-sequence)) ; % best gene to put in the worst place
move2 = [1 imins2 imaxs false];
seq2 = doMoves(sequence,move2,n);
s2= sum(abs(seq2-target));
if s1<s2
moves(i,:)=move1;
sequence = seq1;
else
moves(i,:)=move2;
sequence = seq2;
end
runningScore(i)= sum(abs(sequence-target));
end
score= runningScore(budget) ;
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function moves=oldconvolutedsolver(sequence,target,budget)
seq1=sequence; score1 = 67108864;
L = length(target);
if L < 87
moves1 = HannesNaudiersolver(sequence, target, budget);
score1 = findscore(sequence, target, budget, moves1);
if score1<=50
moves=moves1;
return
end
end
moves = new2solver(sequence, target, budget);
score = findscore(sequence,target,budget,moves);
if score<=50
return
end
if score1<score,
moves = moves1;
score=score1;
end
moves2 = anothersolver(sequence, target, budget);
score2 = findscore(sequence,target,budget,moves);
if score2<=25
moves=moves2;
return
end
if score2<score,
moves = moves2;
end
seq3 = doMoves(seq1,moves,L);
score3 = sum(abs(seq3-target));
mybudget = budget-size(moves,1);
if mybudget > 0
tempmoves=ukilaccuratesolver(seq3,target,mybudget);
seq4=doMoves(seq3,tempmoves,L);
score4=sum(abs(seq4-target));
if score4<score3
moves=[moves;tempmoves];
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function moves=ukilaccuratesolver(s,t,budget)
moves = zeros(budget,4);
L = length(s);
for iter=1:budget
diff=abs(s-t);
[val,index]=sort(diff,'descend');
% first do simple
value = inf(4,5);
for i=1:2
diff2 = abs(s(index(i))-t);
[val2, index2] = sort(diff2);
for j=1:3
[localmoves, score] = removefromto(s,t,index(i),index2(j));
if ~isempty(localmoves)
value(i,:) = [localmoves score];
end
end
end
[mv,mi]=min(value(:,5));
if ~isempty(mi)
localmoves=value(mi,1:4);
moves(iter,:) = localmoves;
s = doMove(s,localmoves(1), localmoves(2), localmoves(3), localmoves(4), L);
else
moves = moves(1:iter-1,:);
return
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function moves=anothersolver(seq2,MC2,budget)
moves=zeros(budget,4);
L = length(seq2);
for i=1:budget
diff=abs(seq2-MC2);
[Dim,index]=max(diff);
diff2=abs(seq2(index)-MC2);[N0T,J_m]=min(diff2);
newmoves=removefromto(seq2,MC2,index,J_m);
if ~isempty(newmoves)
moves(i,:) = newmoves;
seq2=doMoves(seq2,newmoves,L);
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [moves,INITIAL]=removefromto(seq2,tar2,index,index2)
L = length(seq2);
INITIAL=sum(abs(seq2-tar2));
tempmoves=removefromtosub1(seq2,index,index2);
if isempty(tempmoves)
moves=ones(0,4);
else
n1 = size(tempmoves,1);
sn1=zeros(n1,L);
tscore1=zeros(n1,L);
for i=1:n1
sn1(i,:)=oldsplice(seq2,tempmoves(i,1),tempmoves(i,2),tempmoves(i,3),tempmoves(i,4));
tscore1(i,:)=sn1(i,:)-tar2;
end
tscore1=sum(abs(tscore1),2);
[minval,minindex1]=min(tscore1);
moves=tempmoves(minindex1,:);
INITIAL=minval;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function tempmoves=removefromtosub1(seq2,index,index2)
if index<index2
N=abs(index2-index)+1;
tempmoves = zeros(N,4);
for i=1:N
tempmoves(i,:)=[i index index2-i+1 1];
end
N=min([index-1,numel(seq2)-index2])+1;
tempmoves2 = zeros(N-1,4);
for i=2:N
tempmoves2(i-1,:)=[i index-i+1 index2 1];
end
N=index;
tempmoves3 = zeros(N-1,4);
for i=2:N
tempmoves3(i-1,:)=[i index-i+1 index2-i+1 0];
end
N=numel(seq2)-index2+1;
tempmoves4 = zeros(N-1,4);
for i=2:N
tempmoves4(i-1,:)=[i index index2 0];
end
elseif index>=index2
N=abs(index2-index)+1;
tempmoves = zeros(N,4);
for i=1:N
tempmoves(i,:)=[i index-i+1 index2 1];
end
N=min([numel(seq2)-index,index2-1])+1;
tempmoves2 = zeros(N-1,4);
for i=2:N
tempmoves2(i-1,:)=[i index index2-i+1 1];
end
N=numel(seq2)-index+1;
tempmoves3 = zeros(N-1,4);
for i=2:N
tempmoves3(i-1,:)=[i index index2 0];
end
N=index2;
tempmoves4 = zeros(N-1,4);
for i=2:N
tempmoves4(i-1,:)=[i index-i+1 index2-i+1 0];
end
end
tempmoves=[tempmoves;tempmoves2;tempmoves3;tempmoves4];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function moves = new2solver(sequence, target, budget)
swap_counter=2;
sequence_ini=sequence;
L = length(sequence);cumSumShiftLeftMatrix = zeros(L, L+1);
cumSumShiftRightMatrix = cumSumShiftLeftMatrix;cumSumShiftLeftRevMatrix = cumSumShiftLeftMatrix;cumSumShiftRightRevMatrix = cumSumShiftLeftMatrix;moves = zeros(budget,4);moveNr = 1;
while 1
for k = 1:L
index = 2:L-k+2;
target1 = target(index-1);
target2 = target(k:L);
cumSumShiftLeftMatrix(k,index) = sequence(k:L) - target1;
cumSumShiftRightMatrix(k,index) = sequence(1:L-k+1) - target2;
cumSumShiftLeftRevMatrix(k,index) = sequence(L-k+1:-1:1) - target1;
cumSumShiftRightRevMatrix(k,index) = sequence(L:-1:k) - target2;
end
cumSumShiftLeftMatrix=cumsum(abs(cumSumShiftLeftMatrix),2);
cumSumShiftRightMatrix=cumsum(abs(cumSumShiftRightMatrix),2);
cumSumShiftLeftRevMatrix=cumsum(abs(cumSumShiftLeftRevMatrix),2);
cumSumShiftRightRevMatrix=cumsum(abs(cumSumShiftRightRevMatrix),2);
cumSumFromleft = cumSumShiftLeftMatrix(1,:);cumSumFromRight = cumSumFromleft(L+1) - cumSumFromleft;
bestscore = cumSumFromleft(L+1);
relMove = moveNr/budget;
tmpscore = bestscore;
curMoveLengthVector = [(33:4:L-1), 32, 29:-1:1];
curMoveLengthVector(curMoveLengthVector >= L) = [];
zeroIndex = (abs(sequence - target) < 1e-3)*floor(relMove);
blen=1;bkFrom=1;bkTo=1;bFlip=false;
for len = curMoveLengthVector
cssLm = cumSumShiftLeftMatrix(len+1,:);
cssRm = cumSumShiftRightMatrix(len+1,:);
kMax = L - len + 1;
for kFrom = 1:kMax
for kTo = 1:kMax*~zeroIndex(kFrom)
if kFrom <= kTo
score1 =cumSumFromleft(kFrom) + cssLm(kTo) +cumSumFromRight(kTo+len) - cssLm(kFrom);
else
score1 =cumSumFromleft(kTo) + cssRm(kFrom) + cumSumFromRight(kFrom+len) - cssRm(kTo);
end
if len ~= 1
seqRevIndex = L-kFrom-len+2;
if kTo >= seqRevIndex
curScore = score1 +cumSumShiftRightRevMatrix(kTo - seqRevIndex + 1, seqRevIndex+len) - ...
cumSumShiftRightRevMatrix(kTo - seqRevIndex + 1, seqRevIndex);
else
curScore = score1 + ...
cumSumShiftLeftRevMatrix(seqRevIndex - kTo + 1, kTo+len) -cumSumShiftLeftRevMatrix(seqRevIndex - kTo + 1, kTo);
end
if curScore < bestscore
bestscore = curScore;
blen=len;bkFrom=kFrom;bkTo=kTo; bFlip=true;
end
end
if len<=L/2
if kFrom <= kTo
curScore = score1 +cumSumShiftRightMatrix(kTo-kFrom+1, kFrom+len) -cumSumShiftRightMatrix(kTo-kFrom+1, kFrom);
else
curScore = score1 + cumSumShiftLeftMatrix(kFrom-kTo+1, kTo+len) -cumSumShiftLeftMatrix(kFrom-kTo+1, kTo);
end
if curScore < bestscore
bestscore = curScore;blen=len;bkFrom=kFrom;bkTo=kTo;
bFlip=false;
end
end
end
end
end
moves(moveNr,:)=[blen bkFrom bkTo bFlip];sequence = doMove(sequence, blen, bkFrom, bkTo, bFlip, L);
swap_counter=swap_counter-1;
if swap_counter<=0
% if moveNr>1
for jj = 1:2;
rfac = (jj>1)*0.8;
seq1 = splice(sequence_ini,moves,moveNr-2);
[movesS,score]=swapper(seq1,target,2,rfac);
if score < bestscore
moves(moveNr-1:moveNr,:) = movesS;
bestscore = score;
sequence = splice(seq1,movesS,2);
swap_counter=2;
end
end
end
moveNr = moveNr + 1;
if (moveNr > budget)+(tmpscore==bestscore)
moves = moves(1:moveNr-1,:);
return
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function scorecheck = findscore(sequence,target,budget,moves)
seq = splice(sequence,moves,budget);
scorecheck = sum(abs(seq-target));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function seq = splice(seq,moves,budget)
if size(moves,1)>budget
moves = moves(1:budget,:);
end
Z65 = numel(seq);moves(:,1) = max(1,min(Z65,moves(:,1)));
moves(:,2) = max(1,min(Z65-moves(:,1)+1,moves(:,2)));
moves(:,3) = max(1,min(Z65-moves(:,1)+1,moves(:,3)));
moves(:,4) = moves(:,4) ~= 0;
for i =1:size(moves,1)
le = moves(i,1);ai = moves(i,2);bi = moves(i,3);
flip = moves(i,4);
if flip
if ai>bi
seq = seq([1:bi-1 ai+le-1:-1:ai bi:ai-1 ai+le:end]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai+le-1:-1:ai bi+le:end]);
end
else
if ai>bi
seq = seq([1:bi-1 ai:ai+le-1 bi:ai-1 ai+le:end]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai:ai+le-1 bi+le:end]);
end
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function moves = markussolver(sequence, target, budget)
swap_counter=2;
sequence_ini=sequence;
moveLengthVector = [14 12 9 8 7];
moveLengthThresh = [0.9 0.6 0.3 0.1 -1];
L = length(sequence);
index = moveLengthVector > L;
moveLengthVector(index) = [];
moveLengthThresh(index) = [];
cumSumShiftLeftMatrix = zeros(L, L+1);
cumSumShiftRightMatrix = cumSumShiftLeftMatrix;
cumSumShiftLeftRevMatrix = cumSumShiftLeftMatrix;
cumSumShiftRightRevMatrix = cumSumShiftLeftMatrix;
moves = zeros(budget,4);
moveNr = 1;
while moveNr <= budget
for k = 1:L
index = 2:L-k+2;
target1 = target(index-1);
target2 = target(k:L);
cumSumShiftLeftMatrix(k,index) = sequence(k:L) - target1;
cumSumShiftRightMatrix(k,index) = sequence(1:L-k+1) - target2;
cumSumShiftLeftRevMatrix(k,index) = sequence(L-k+1:-1:1) - target1;
cumSumShiftRightRevMatrix(k,index) = sequence(L:-1:k) - target2;
end
cumSumShiftLeftMatrix=cumsum(abs(cumSumShiftLeftMatrix),2);
cumSumShiftRightMatrix=cumsum(abs(cumSumShiftRightMatrix),2);
cumSumShiftLeftRevMatrix=cumsum(abs(cumSumShiftLeftRevMatrix),2);
cumSumShiftRightRevMatrix=cumsum(abs(cumSumShiftRightRevMatrix),2);
cumSumFromleft = cumSumShiftLeftMatrix(1,:);
cumSumFromRight = cumSumFromleft(L+1) - cumSumFromleft;
bestscore = cumSumFromleft(L+1);
tmpscore=bestscore;
relMove = 1 - moveNr / budget;
moveLengthIndex = find(relMove >= moveLengthThresh, 1, 'first');
curLen = moveLengthVector(moveLengthIndex);
curMoveLengthVector = [31:-2:21 curLen+4 curLen+3 curLen+2 curLen 6 5 4 3 2 1];
curMoveLengthVector(curMoveLengthVector >= L) = [];
curMoveLengthVector = fliplr(unique(curMoveLengthVector));
zeroIndex = (abs(sequence - target) < 1e-3)*round(relMove-0.3);
blen=1;bkFrom=1;bkTo=1;bFlip=false;
for len = curMoveLengthVector
cssLm = cumSumShiftLeftMatrix(len+1,:);
cssRm = cumSumShiftRightMatrix(len+1,:);
kMax = L - len + 1;
for kFrom = 1:kMax
for kTo = 1:(kFrom-1)*~zeroIndex(kFrom)
score1 =cumSumFromleft(kTo) + cssRm(kFrom) + ...
cumSumFromRight(kFrom+len) - cssRm(kTo);
if len == 1
curScore=score1+cumSumShiftLeftMatrix(kFrom-kTo+1, kTo+len) -cumSumShiftLeftMatrix(kFrom-kTo+1, kTo);
else
seqRevIndex = L-kFrom-len+2;
if kTo >= seqRevIndex
curScore = score1 +cumSumShiftRightRevMatrix(kTo - seqRevIndex + 1, seqRevIndex+len) - ...
cumSumShiftRightRevMatrix(kTo - seqRevIndex + 1, seqRevIndex);
else
curScore = score1 + ...
cumSumShiftLeftRevMatrix(seqRevIndex - kTo + 1, kTo+len) -cumSumShiftLeftRevMatrix(seqRevIndex - kTo + 1, kTo);
end
end
if curScore < bestscore
bestscore = curScore;
blen=len;
bkFrom=kFrom;
bkTo=kTo;
bFlip=true;
end
end
for kTo = kFrom:kMax
score1 =cumSumFromleft(kFrom) + cssLm(kTo) +cumSumFromRight(kTo+len) - cssLm(kFrom);
if len == 1
curScore = score1 +cumSumShiftRightMatrix(kTo-kFrom+1, kFrom+len) -cumSumShiftRightMatrix(kTo-kFrom+1, kFrom);
else
seqRevIndex = L-kFrom-len+2;
if kTo >= seqRevIndex
curScore = score1 +cumSumShiftRightRevMatrix(kTo - seqRevIndex + 1, seqRevIndex+len) - ...
cumSumShiftRightRevMatrix(kTo - seqRevIndex + 1, seqRevIndex);
else
curScore = score1 + ...
cumSumShiftLeftRevMatrix(seqRevIndex - kTo + 1, kTo+len) -cumSumShiftLeftRevMatrix(seqRevIndex - kTo + 1, kTo);
end
end
if curScore < bestscore
bestscore = curScore;blen=len;bkFrom=kFrom;bkTo=kTo;
bFlip=true;
end
end
end
end
moves(moveNr,:)=[blen bkFrom bkTo bFlip];
sequence = doMove(sequence, blen, bkFrom, bkTo, bFlip, L);
swap_counter=swap_counter-1;
% if swap_counter<=0
if moveNr>1
for jj = 1:2
rfac = (jj>1)*0.3;
seq1 = splice(sequence_ini,moves,moveNr-2);
[movesS,score]=swapper(seq1,target,2,rfac);
if score < bestscore
moves(moveNr-1:moveNr,:) = movesS;
bestscore = score;
sequence = splice(seq1,movesS,2);
swap_counter=2;
end
end
end
moveNr = moveNr + 1;
if (moveNr > budget)+(tmpscore==bestscore)
moves = moves(1:moveNr-1,:);
return
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function seq = oldsplice(seq, le, ai, bi, flip)
if flip
if ai>bi
seq = seq([1:bi-1 ai+le-1:-1:ai bi:ai-1 ai+le:end]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai+le-1:-1:ai bi+le:end]);
end
else
if ai>bi
seq = seq([1:bi-1 ai:ai+le-1 bi:ai-1 ai+le:end]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai:ai+le-1 bi+le:end]);
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function seq = doMove(seq, le, ai, bi, flip, L)
if flip
if ai>bi
seq = seq([1:bi-1 ai+le-1:-1:ai bi:ai-1 ai+le:L]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai+le-1:-1:ai bi+le:L]);
end
else
if ai>bi
seq = seq([1:bi-1 ai:ai+le-1 bi:ai-1 ai+le:L]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai:ai+le-1 bi+le:L]);
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function seq = doMoves(seq,moves,L)
for i =1:size(moves,1)
le = moves(i,1);
ai = moves(i,2);
bi = moves(i,3);
if moves(i,4)
if ai>bi
seq = seq([1:bi-1 ai+le-1:-1:ai bi:ai-1 ai+le:L]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai+le-1:-1:ai bi+le:L]);
end
else
if ai>bi
seq = seq([1:bi-1 ai:ai+le-1 bi:ai-1 ai+le:L]);
else
seq = seq([1:ai-1 ai+le:bi+le-1 ai:ai+le-1 bi+le:L]);
end
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function moves=HannesNaudiersolver(s,t,budget)
best = inf;
for trynum=1:5
prelimmoves=jacksolver(s,t,budget);
score = grade(s,t,budget,prelimmoves);
if score<best
moves=prelimmoves;
best=score;
end
end
function moves = jacksolver(s,t,budget)
n = length(s);
on = ones(1,n);
trep = t(on,:);
moves = zeros(budget,4);
for b = 1:budget
d = abs(trep - s(on,:)') + std(t)*(budget-b)/125*randn(n);
moves(b,:) = bestmove(d);
if ~moves(b,1)
moves(b:end,:) = [];
return
end
le = moves(b,1);
ai = moves(b,2);
bi = moves(b,3);
if moves(b,4)
if ai>bi
s = s([1:bi-1 ai+le-1:-1:ai bi:ai-1 ai+le:end]);
else
s = s([1:ai-1 ai+le:bi+le-1 ai+le-1:-1:ai bi+le:end]);
end
else
if ai>bi
s = s([1:bi-1 ai:ai+le-1 bi:ai-1 ai+le:end]);
else
s = s([1:ai-1 ai+le:bi+le-1 ai:ai+le-1 bi+le:end]);
end
end
end
function m = bestmove(d)
n = size(d,1);
dp = diagcumsum(d);
dm = flipud(diagcumsum(flipud(d)));
cost = 0;
m = zeros(1,4);
for i = 1:n-1
for k = i+2:n+1
old = dp(k,k)-dp(i,i); % orig diag [i,k)
for j = i+1:k-1
if dp(k,i+k-j)-dp(j,i)+dp(j,k)-dp(i,i+k-j)-cost<old || dp(k,i+k-j)-dp(j,i)+dm(i,k)-dm(j,i+k-j)-cost<old || dm(j,i+k-j)-dm(k,i)+dp(j,k)-dp(i,i+k-j)-cost<old
if dp(k,i+k-j)-dp(j,i)+dp(j,k)-dp(i,i+k-j)-cost<old
cost = dp(k,i+k-j)-dp(j,i)+dp(j,k)-dp(i,i+k-j)-old;
m = [k-j j i 0 ];
end
if dp(k,i+k-j)-dp(j,i)+dm(i,k)-dm(j,i+k-j)-cost<old
cost = dp(k,i+k-j)-dp(j,i)+dm(i,k)-dm(j,i+k-j)-old;
m = [j-i i i+k-j 1];
end
if dm(j,i+k-j)-dm(k,i)+dp(j,k)-dp(i,i+k-j)-cost<old
cost = dm(j,i+k-j)-dm(k,i)+dp(j,k)-dp(i,i+k-j)-old;
m = [k-j j i 1 ];
end
end
end
end
end
function R = diagcumsum(M)
% cumulative sums along diagonal
n = size(M,1);
R = zeros(n+1);
for c = 2:n+1
for r = 2:n+1
R(r,c) = M(r-1,c-1) + R(r-1,c-1);
end
end
|
