Como resolver um quadrado mágico 4x4 utilizando algoritmos?

Question

O problema refere-se ao exercício 10, da página 36, do livro The Art of Computer Programming, Volume 4A: Combinatorial Algorithms – Donald E. Knuth.

Answer 1

Quadrados mágicos são matrizes quadradas, cuja soma de cada uma de suas linhas, colunas ou diagonais resulte sempre em um mesmo número. De origem remota, o quadrado mágico é um problema que aparece em diferentes culturas, como a chinesa, a árabe e a indiana. Ele é também o símbolo adotado pela Editora 34, conhecida por suas traduções de clássicos da literatura europeia.

Curiosamente, estamos interessados na resolução de quadrados mágicos de ordem 4x4 resultando no número 34.

Uma solução bem ineficiente, baseada em busca exaustiva, seria testar todas as permutações possíveis entre os números de 1 a 16 para encontrar quais se enquadram nas regras especificadas. No entanto, teríamos \(P_{16} = 16! \approx 20 \) trilhões de possibilidades para testar, de modo que um computador que realize um milhão de operações por segundo demoraria cerca de 240 dias para encontrar todas as respostas, no pior dos casos! Uma alternativa seria utilizar backtracking.

O que faremos aqui será construir soluções a partir do primeiro dígito, verificando se, até então, a soma de cada linha e coluna não supera o limite máximo estabelecido de 34. A função magic_square cria um vetor que contém todos os números factíveis, bem como uma matriz "possibility", que receberá os valores a serem checados.

import numpy as np


def magic_square(n = 4):
    """ Finds one solution to magic squares. 

    n: square side
    """ 
    n_square = n**2 # Max number
    target = ((n_square)*(n_square + 1)/2) / n # Row's / column's target sum
    possibility = np.zeros((n, n), dtype = int) # Will receive the tested values

    values = [x for x in range(1, n_square + 1)] # List with all the numbers

    # Does the trick
    solve(possibility, values, target, n, level = 0)

    return(None)

A função solve, por sua vez, resolverá o problema de maneira recursiva, interrompendo sua execução quando a variável level atingir o valor 16, ie, quando preenchermos o último espaço disponível. Note que colocaremos um número na matriz somente se sua adição não violar as condições pré-estabelecidas.

def solve(possibility, values, target, n, level):
    """ Solves the problem.

    possibility: each attempt
    values: all possible values
    target: row's / column's / diagonal's target sum.
    level: tree node level (max level == n_square)
    """
    if level == n**2:
        print("We've found a solution!\n", possibility)
        return True
    else:
        for i in values:

            # Gets level's row/column equivalents to fill the "possibility" matrix
            row = level // n
            col = level % n

            # Checks whether certain possibity is valid so far or not
            if not(i in possibility):
                if is_valid(possibility, n, target, row, col, next_attempt = i):                    
                    possibility[row, col] = i 
                    #print(possibility)
                    if solve(possibility, values, target, n, level + 1):
                        return True

            possibility[row, col] = 0 # Backtracking

        return False

Por fim, a função is_valid fará a verificação. O pulo do gato aqui é que ela testará as linhas apenas no momento em que elas forem completamente preenchidas, valendo o mesmo para as colunas e diagonais.

def is_valid(possibility, n, target, row, col, next_attempt):
    """ Checks wheter certain input will be a valid magic square, so far. If so,
    we're allowed to add "next_attempt" to possibility[row, col].

    n: square's side
    row/col: current testing row/col.
    next_attempt: next number we'll input.
    """                    
    possibility[row, col] = next_attempt

    # Checks rows (If you've reached the end of a row and its sum != target, return FALSE)
    if (col == n - 1) and sum(possibility[row,:]) != target:
        return False

    # Checks columns
    if (row == n - 1) and sum(possibility[:,col]) != target:
        return False

    # Checks secondary diagonal (does the checking when you reach 1st entry of last line)
    if (row == n - 1) and (col == 0):
        my_sum = 0 
        for i in range(n):
            my_sum += possibility[i, n - 1 - i]

        if my_sum != target:
            return False


    # Checks main diagonal (only when you've reached the very last entry)
    if (row == n - 1) and (col == n - 1):
        my_sum = 0
        for i in range(n):
            my_sum += possibility[i, i]

        if my_sum != target:
            return(False)


    return(True)


if __name__ == '__main__':    
    n = 4
    magic_square(n)

É possível, ainda, encontrar a solução para qualquer quadrado mágico de tamanho \(n >= 3\) por meio desse código, bastando alterar o valor de \(n\).