% 初速状态，0表示没填的格子
S = [ 7,4,0,0,8,0,0,1,6,
  9,0,0,0,3,5,0,0,4,
  0,0,0,7,0,0,0,0,0,
  0,7,0,0,0,9,5,0,0,
  6,1,0,0,5,0,0,8,7,
  0,0,2,6,0,0,0,4,0,
  0,0,0,0,0,4,0,0,0,
  3,0,0,5,6,0,0,0,2,
  5,6,0,0,1,0,0,3,9];


% 定义0、1数组 A(i, j, k) = 1，如果方格(i, j)里的数为k；否则为0。
% sudoku问题即求一定假设条件下的解。
p = 3;
A = binvar(p^2,p^2,p^2,'full');

% 以下为限制条件
F = [sum(A,1) == 1]; % 限制每行每个数恰好一个
F = [F, sum(A,2) == 1]; % 限制每列每个数恰好一个
F = [F, sum(A,3) == 1]; % 限制每个单元格子里恰好一个数

for m = 1:p
    for n = 1:p
        for k = 1:p^2
            s = sum(sum(A((m-1)*p+(1:p),(n-1)*p+(1:p),k)));  
            F = [F, s == 1];  % 限制每个3×3的方框里每个数恰好出现一次
        end
    end
end

for i = 1:p^2
    for j = 1:p^2
        if S(i,j)
            F = [F, A(i,j,S(i,j)) == 1]; % 初始给定的数要一直
        end
    end
end

% 直接求解
sol = solvesdp(F); 

Z = 0;
for i = 1:p^2
      Z = Z  + i*double(A(:,:,i)); % 简单相加即可
end
Z % 输出结果