Code covered by the BSD License  

Highlights from
recursive solver to peg solitaire contest

from recursive solver to peg solitaire contest by Yi Cao
a solver to push the 3-minute limit

recsol(board)
function [moves,score,v0,k0] = recsol(board)
%
% Note in the interest of collaboration I'm documenting
% the leading code as much as possible.  Instead of going the
% obscufation route with this, I invite everyone to
% document / take credit for their particular changes
% as the codes get modified.
%
% At a minimum, please don't remove the existing comments
% so someones doesn't have to start from scratch on commenting
% updated code again
%
% Alan Chalker
%

[m,n] = size(board); % find the dims of the board
pegCount = nnz(board>0); % check the number of pegs on the board
rv = 5:m+4; % create a new row index starting at 5th row
cv = 5:n+4; % create a new col index starting at 5th col  

%     i = repmat(rv',[n 1]); % create an index of all the new board row coordinates except for the 1st 4 rows
%     j = reshape(repmat(cv,[m 1]),[m*n 1]); % create an index of all the new board col coordinates except for the 1st 4 cols
    % note [i j] would be a list of all original board coordinates in the
    % new board created below

    mm = m+8; % expand the original height by 8 rows
    nn = n+8;  % expand the original width by 8 cols
    ppBoard = zeros(mm, nn)-1; % create a new board with 4 extra rows / cols of offlimits all the way around the board
    ppBoard(rv,cv) = board; % populate the new board with the original board values
    [i,j]=find(ppBoard>=0);

    I = [i;i;i;i];  % vector of 4 repeats of row coords
    J = [j;j;j;j]; % vector of 4 repeats of col coords
    K = [i;i;i-2;i+2]; % vector of row cords, row cords, 2 rows above, 2 rows below
    L = [j-2;j+2;j;j]; % vector of 2 cols before, 2 cols past, col cords, col cords
    % [I J K L] would be a matrix of ALL potential moves for this board

    K1 = [i;i;i-4;i+4]; % vector of row cords, row cords, 4 rows above, 4 rows below
    L1 = [j-4;j+4;j;j]; % vector of 4 cols before, 4 cols past, col cords, col cords
    % if a move from [I J] to [K L] were done, [K L K1 L1] is the potential next colinear moves

    K2 = [i-2;i+2;i-2;i+2]; % vector of 2 rows above, 2 rows below repeated twice
    L2 = [j-2;j+2;j+2;j-2];  % vector of 2 cols before, 2 cols past, 2 cols past, 2 cols before
    % if a move from [I J] to  [K L] were done, [K L K2 L2] is the half of the potential next orthogonal moves

    K3 = [i+2;i-2;i-2;i+2]; % vector of 2 rows below, 2 rows above, 2 rows above, 2 rows below
    L3 = [j-2;j+2;j-2;j+2]; % vector of 2 cols before, 2 cols past repeated twice
    % if a move from [I J] to [K L] were done, [K L K3 L3] is the other half of the potential next orthogonal moves

    F = I+(J-1)*mm; % convert source spot coordinates [I J] into single index value
    T = K+(L-1)*mm; % convert destination spot coordinates [K L] into single index value
    M = (F+T)*0.5; % calculate the index value of the spot that would be jumped if did move [I J K L]

    % find indexes of moves that don't involve off limits (-1) areas
    bogusMoves = (ppBoard(M) < 0) | (ppBoard(T) < 0);

    % remove moves that involve off limits areas
    I(bogusMoves) = [];
    J(bogusMoves) = [];
    K(bogusMoves) = [];
    L(bogusMoves) = [];
    F(bogusMoves) = [];
    T(bogusMoves) = [];
    M(bogusMoves) = [];
    K1(bogusMoves) = [];
    K2(bogusMoves) = [];
    K3(bogusMoves) = [];
    L1(bogusMoves) = [];
    L2(bogusMoves) = [];
    L3(bogusMoves) = [];


    % convert 2nd jump destination spot coordinates into single index value
    % note invalid jumps are kept in index but not converted
    T1 = K1+(L1-1)*mm; % colinear
    T2 = K2+(L2-1)*mm; % ortho
    T3 = K3+(L3-1)*mm; % ortho
    TT = [T1 T2 T3];

    % calculate the index value of the spot that would be jumped if did move
    M1 = (T+T1)*0.5; % [K L K1 L1]
    M2 = (T+T2)*0.5; % [K L K2 L2]
    M3 = (T+T3)*0.5; % [K L K3 L3]
    MM = [M1 M2 M3];

    % create first possible move list
    MV = [I J K L]; 
    MV1 = [K L K1 L1];  
    MV2 = [K L K2 L2]; 
    MV3 = [K L K3 L3];
    MVV = {MV1, MV2, MV3};
    [rMap,rMap3] = createMoves(M,F,T,M1,M2,M3,T1,T2,T3);

    % run the subsolver function
    R    =[0.12  0.83 0.95 0.19 0.22 0.28];
    dep  =[0     23   9    12   5    18];
    level=[1     1    1    1    2    1];
    score=-Inf;
    k0=0;
    for k=1:6,
        [newMoves,newScore,newV0] = subsoltweak( ...
            ppBoard, ...
            R(k),dep(k),level(k), ...
            F,T,M, ...
            pegCount, ...
            TT,MM,MV,MVV, ...
            rMap,rMap3);
        if newScore>score
            k0=k;
            v0=newV0;
            score=newScore;
            moves=newMoves;
        end
        if size(moves,1)<5
            moves=moves-4;
            return
        end
    end
%     moves=moves-4;
%     return
    R=[0.65 0.85 0.71 0.6 0.81 0.7  0.79 0.79 0.8];
    dd=[0.7 1.35 1.31 1.1 1.45 1.49 1.35 1.34 0.9];
    dep=[0   0   15   23  6    10   5    4    36];
    level=[1 1 1 1 2 1 1 1 1];
    for k=1:8
        [newMoves,newScore,newV0] = subsol( ...
            ppBoard, ...
            dd(k),R(k),dep(k),level(k), ...
            F,T,M, ...
            pegCount,  ...
            TT,MM,rMap); %,  ...
        if (newScore > score)
            k0=k+6;
            v0=newV0;
            score = newScore;
            moves = MV(newMoves,:);
        end
    end
%     disp(k0);

    % correct moves due to board padding
    moves = moves - 4;
end   % close function solver

%======================================================================
function [rMap,rMap3] = createMoves(M,F,T,M1,M2,M3,T1,T2,T3)

ni = max([M;F;T;M1;M2;M3;T1;T2;T3]); % find the index of the lower right most possible move spots

% create three copies of a vector long enough to contain all possible move
% position indexes
nmove1 = zeros(ni,1);
nmove2 = nmove1;
nmove3 = nmove1;

% create three copies of a matrix long enough to contain all possible moves
moveid1 = zeros(ni,4);
moveid2 = moveid1;
moveid3 = moveid1;

N=numel(F);
for k = 1:N % repeat for the number of starting positions
    nmove1(F(k)) = nmove1(F(k))+1; %populate the 1st move position vector with the # of times that source spot is in the possible moves
    moveid1(F(k),nmove1(F(k))) = k; % populate the 1st move list vector with the index value of the corresponding source spots in each col
    % Note: if the 10th position on the board is a peg that could make 3
    % possible moves (i.e. up, down, right), that nmove1(10)=3 and
    % nmoveid1(10,:) = [x y z 0] where x,y,x are the index values in
    % vector F where 10 appears

    % Repeat the above for the spots that could be jumped
    nmove2(M(k)) = nmove2(M(k))+1;
    moveid2(M(k),nmove2(M(k))) = k;

    % Repeat the above for the destination spots
    nmove3(T(k)) = nmove3(T(k))+1;
    moveid3(T(k),nmove3(T(k))) = k;
end
rMap=[moveid1(F,:) moveid2(F,:) moveid3(F,:) ...
    moveid1(M,:) moveid2(M,:) moveid3(M,:) ...
    moveid1(T,:) moveid2(T,:) moveid3(T,:)];
rMap3={[moveid1(F,:) moveid2(F,:) moveid3(F,:) ...
    moveid1(M,:) moveid2(M,:) moveid3(M,:) ...
    moveid1(M1,:) moveid2(M1,:) moveid3(M1,:) ...
    moveid1(T1,:) moveid2(T1,:) moveid3(T1,:)],...
   [moveid1(F,:) moveid2(F,:) moveid3(F,:) ...
    moveid1(M,:) moveid2(M,:) moveid3(M,:) ...
    moveid1(M2,:) moveid2(M2,:) moveid3(M2,:) ...
    moveid1(T2,:) moveid2(T2,:) moveid3(T2,:)],...
   [moveid1(F,:) moveid2(F,:) moveid3(F,:) ...
    moveid1(M,:) moveid2(M,:) moveid3(M,:) ...
    moveid1(M3,:) moveid2(M3,:) moveid3(M3,:) ...
    moveid1(T3,:) moveid2(T3,:) moveid3(T3,:)]};
end

%======================================================================
function [moves,v,v0] = subsol( ...
    board, ...
    d,R,depth,level, ...
    F,T,M, ...
    pegCount, ...
    TT,MM,rMap) %,  ...
moves = zeros(pegCount-1,1); % preallocate maximum possible move list based on number of pegs
v0 = zeros(pegCount-1,1); % preallocate maximum size of score list
Bzero = ~board;  % create inverse board where 1 is a hole and every else is a zero, including offlimits
Bmax = max(board, 0); % create board with offlimits as holes

% search for moves where source is peg, destination is hole and jumped spot is peg
validMoves = (Bmax(F) & Bzero(T) & Bmax(M));
T0=find(board<0,1);
t=1;
while any(validMoves)
    % extract indexes of valid moves
    h = find(validMoves);
    if (numel(h) > 1)
        c = (F(h) == T0); % find indexes of jumps with source peg that is same as last one moved             
        if any(c) % if any current 2 jump moves contain the last peg
            C0 = Bmax(M(h(c))); % seed possible 2nd jumps board with jumped peg value for all jumps the match last jump end
            if max(C0)+(1+R)*min(C0)-Bmax(T0)>0
                h = h(c); % extract the jump index
                bFh=Bmax(F(h));
            else
                bFh=Bmax(F(h));
                C0 = Bmax(M(h)) - bFh.*~c;
            end
        else
                bFh=Bmax(F(h));
            C0 = Bmax(M(h)) - bFh; % calculate score for all valid 1st moves
        end
        CV = max(Bmax(MM(h,:)).*Bzero(TT(h,:)),[],2); % check to see if next colinear move is best

        % add jumped peg to 2nd jump and find best 2 move jump
        C1 = (C0+CV*d)./(bFh+R);
        [v,k] = max( C1 );                

        if ~CV(k) && t>(depth<7) && depth
            [C0,pos]=sort(C0,'descend');
            C1=bFh(pos);
            h=h(pos);
            dep=min(depth,numel(h));
            C2=zeros(dep,1);
            for kk=1:dep
                CV(kk)=0;
                [C0D,CL,CVD]=depsol(Bmax,Bzero,F,M,T,MM,TT,rMap,validMoves,h(kk),level);
                if ~isempty(C0D)
                    C2(kk) = max((C0(kk)+C0D+d*CVD)./(R+CL+C1(kk)));
                end
            end
            [c2,k]=max(C2);
        end        
        
        v0(t) = C0(k); % extract score of best 1st jump score
        k = h(k);  % extract index of best 2 move jump
    else
        k = h; % extract the move position
        v0(t) = Bmax(M(k))-Bmax(F(k))*(F(k)~=T0); % calculate the jump score
    end
    moves(t) = k; % add best move (1st jump) to movelist
    T0 = T(k); % extract 1st jump destination spot
    F0 = F(k); % extract source spot
    M0 = M(k); % calculate jumped spot
    t = t+1; % increment move count

    % Update board.
    Bmax(T0) = Bmax(F0); % copy the jumping peg value
    Bmax(F0) = 0; % zero out the jumping spot peg
    Bmax(M0) = 0; % zero out the jumped spot peg
    Bzero(F0)=true;Bzero(M0)=true;Bzero(T0)=false;
    x=rMap(k,:);
    for i=1:36
        j=x(i);
        if j>0
            validMoves(j)=Bmax(F(j)) && Bmax(M(j)) && Bzero(T(j));
        end
    end        
end

v1 = cumsum(v0); % create cumulative sum of scores in movelist
[v,t] = max(v1); % extract location of best cumulative score
moves = moves(1:t); % output moves up to best location
v0 = v0(1:t);
end

%======================================================================
function [moves,v,v0] = subsoltweak( ...
    board, ...
    R,depth,level, ...
    F,T,M, ...
    pegCount, ...
    TT,MM,MV,MVV, ...
    rMap,rMap3)

moves = zeros(pegCount-1,4); % preallocate maximum possible move list based on number of pegs
v0 = zeros(pegCount-1,1); % preallocate maximum size of score list
Bzero = ~board;  % create inverse board where 1 is a hole and every else is a zero, including offlimits
Bmax = max(board, 0); % create board with offlimits as holes
t=1;
hs = (Bmax(F) & Bzero(T) & Bmax(M)); % search for moves where source is peg, destination is hole and jumped spot is peg
T0=find(board<0,1);
while any(hs)
    h = find(hs);
    if numel(h)>1
        c = (F(h) == T0); % find indexes of jumps with source peg that is same as last one moved
        if any(c) % if any current 2 jump moves contain the last peg
            C0 = Bmax(M(h(c))); % seed possible 2nd jumps board with jumped peg value for all jumps the match last jump end
            if max(C0)+(1+R)*min(C0)-Bmax(T0)>0
                h = h(c); % extract the jump index
                C2=C0;
            else
                bFh=Bmax(F(h)).*~c;
                C0 = Bmax(M(h)) - bFh;
                C2 = C0 - R*bFh;
            end
        else
            bFh=Bmax(F(h));
            C0 = Bmax(M(h))-bFh; % calculate score for all valid 1st moves
            C2 = C0 - R*bFh;
        end

        [CV,kc] = max(Bmax(MM(h,:)).*Bzero(TT(h,:)),[],2); % check to see if next colinear move is best
        [v,k] = max(C2+CV); % add jumped peg to 2nd jump and find best 2 move jump
        cv=CV(k);
        if ~cv && depth && t>1
            [C0,pos]=sort(C0,'descend');
            C2=C2(pos);
            h=h(pos);
            dep=min(depth,numel(h));
            for kk=1:dep
                [C0D,CL,CVD]=depsol(Bmax,Bzero,F,M,T,MM,TT,rMap,hs,h(kk),level);
                if ~isempty(C0D)
                    C2(kk) = C2(kk)+max(C0D-R*CL+CVD);
                end
            end
            [c2,k]=max(C2);
        end
    
        v0(t) = C0(k); % extract score of best 1st jump score
        kc=kc(k);
        k = h(k);  % extract index of best 2 move jump
    else
        k=h;
        v0(t) = Bmax(M(k))-Bmax(F(k))*(F(k)~=T0);
        cv=0;
    end

    moves(t,:) = MV(k,:); % add best move (1st jump) to movelist
    F0 = F(k); % extract source spot
    t = t+1; % increment move count
    if ~cv
        T0=T(k); % extract 1st jump destination spot
        x=rMap(k,:);
    else
        T0=TT(k,kc);
        M0=MM(k,kc);
        moves(t,:)=MVV{kc}(k,:);
        x=rMap3{kc}(k,:);
        v0(t)=cv;
        Bmax(M0)=0;Bzero(M0)=true;
        t=t+1;
    end
    M0=M(k);
    Bmax(T0)=Bmax(F0);Bmax(F0)=0;Bmax(M0)=0;
    Bzero(F0)=true;Bzero(M0)=true;Bzero(T0)=false;
    for i=1:length(x)
        j=x(i);
        if j>0
            hs(j)=Bmax(F(j)) && Bmax(M(j)) && Bzero(T(j));
        end
    end
end
v1 = cumsum(v0); % create cumulative sum of scores in movelist
[v,t] = max(v1); % extract location of best cumulative score
v0=v0(1:t);
moves = moves(1:t,:); % output moves up to best location
end

function [C0,CL,CV]=depsol(BmD,BzD,F,M,T,MM,TT,rMap,hsD,jk,level)
    Fk=F(jk);Mk=M(jk);Tk=T(jk);
    BmD(Tk)=BmD(Fk);BmD(Mk)=0;BmD(Fk)=0;
    BzD(Fk)=true;BzD(Mk)=true;BzD(Tk)=false;
    x=rMap(jk,:);
    for i=1:36
        j=x(i);
        if j>0
            hsD(j)=BmD(F(j)) && BmD(M(j)) && BzD(T(j));
        end
    end
    if any(hsD)
        FhD = F(hsD);
        CL = BmD(FhD);
        C0 = BmD(M(hsD))-CL.*(FhD~=Tk);
        CV = max(BmD(MM(hsD,:)).*BzD(TT(hsD,:)),[],2); % check to see if next colinear move is best
        if level>1
            [v,k]=max(C0);
            h=find(hsD);
            [C01,CL1,CV1]=depsol(BmD,BzD,F,M,T,MM,TT,rMap,hsD,h(k),level-1);
            if ~isempty(C01)
                [C01,kk]=max(C01);
                C0(k)=C0(k)+C01;
                CL(k)=CL(k)+CL1(kk);
                CV(k)=CV1(kk);
            end
        end
    else
        C0=[];CL=[];CV=[];
    end
end

Contact us at files@mathworks.com