Построить матрицу, в которой сумма элементов в каждой строке, столбце и диагонали имеют одинаковую сумму

 n*n(n*n+1) div (2*n)

type
  TMatrix = array of array of integer;
 
  function CreateMagigSquare(N: integer; var M: TMatrix): boolean;
 
  var
    MagicConstant:  integer;
    UsedNumber: array of boolean;
 
    function TestMagicSquare: boolean;
    begin
      ..........................
    end;
 
    function TestRow(Row: integer): boolean;
    var
      j: integer;
      Sum: integer;
    begin
      Sum := 0;
      for j := 0 to N - 1 do
        Sum := Sum + M[Row, j];
      TestRow := Sum = MagicConstant;
    end;
 
    function Backtracking(Cell: integer): boolean;
    var
      i, j: integer;
      k: integer;
    begin
      {если квадрат заполнен, проверяем его свойства}
      if Cell = sqr(N) then
      begin
        Backtracking := TestMagicSquare;
        Exit;
      end;
      i := Cell div N;
      j := Cell mod N;
      {для сокращения перебора, будем проверять каждую заполненную строку}
      if (i > 0) and (j = 0) then
        if not TestRow(i - 1) then
        begin
          Backtracking := False;
          Exit;
        end;
      {перебор возможных кандидатов, на место в данной ячейке}
      for k := 1 to sqr(N) do
        if not UsedNumber[k - 1] then
        begin
          UsedNumber[k - 1] := True;
          M[i, j] := k;
          Backtracking := Backtracking(Cell + 1);
          if Backtracking then
            break;
          UsedNumber[k - 1] := False;
        end;
    end;
 
  var
    i: integer;
  begin
    SetLength(UsedNumber, sqr(N));
    for i := 0 to sqr(N) - 1 do
      UsedNumber[i] := False;
    MagicConstant := (n * (sqr(n) + 1)) div 2;
    CreateMagigSquare := BackTracking(0);
  end;