Tengo un problema de combinatoria que me gustaría poner en el OEIS ; el problema es que no tengo suficientes términos. Este desafío de código es para ayudarme a calcular más términos, y el ganador será el usuario con el envío que contenga la mayor cantidad de términos.

El problema

Supongamos que le doy una matriz triangular de bombillas con longitud lateral $n$ :

     o
    o o
   o o o
  o o o o
 o o o o o
o o o o o o
1 2  ...  n

Voy a encender tres bombillas que forman un triángulo equilátero "vertical" como en el siguiente ejemplo:

     o
    o x
   o o o
  o o o o
 o x o o x
o o o o o o

Antes de encender las luces, su trabajo es eliminar la mayor cantidad posible de bombillas de la matriz, sin perder la capacidad de deducir el triángulo de bombillas que se ha encendido. Para ser claros, si se quitó una bombilla, no se enciende cuando se enciende su posición.

Por ejemplo, si eliminó las siguientes bombillas (marcadas con .), solo vería las siguientes dos luces encendidas (marcadas con x), lo que es suficiente para deducir de manera única la tercera posición (apagada):

     .              .
    . o            . x
   . . o          . . o
  o o o .   =>   o o o .
 o o o o .      o x o o . <- the third unlit position
o . . . o o    o . . . o o

Sea a(n)la cantidad máxima de bombillas que se pueden quitar sin introducir ambigüedades.

Ejemplo

Con un algoritmo ingenuo, he comprobado valores hasta un triángulo con longitud de lado 7, como se ve a continuación:

                                                                      .
                                                      .              . o
                                        .            . o            o . o
                           .           . .          . . o          . o o .
              .           . .         . o o        o o o .        o o . o .
 .           . .         . o o       o o . o      o o o o .      o . o . o o
. .         . o o       . o o o     o . . o o    o . . . o o    o . o . o o o

a(2) = 3    a(3) = 4    a(4) = 5    a(5) = 7     a(6) = 9       a(7) = 11

Puntuación

La presentación que calcula la secuencia [a(2), a(3), ..., a(n)]para las n más grandes gana. Si dos presentaciones tienen secuencias idénticas, entonces la que se publicó anteriormente gana.

Aunque no es necesario para el envío, sería instructivo para mí si publica una construcción de las matrices triangulares resultantes, como en el ejemplo anterior.

Python 3 ,n=8

import itertools
from ortools.sat.python import cp_model


def solve(n):
    model = cp_model.CpModel()
    solver = cp_model.CpSolver()
    cells = {
        (y, x): model.NewBoolVar(str((y, x)))
        for y in range(n) for x in range(y + 1)}
    triangles = [
            {cells[v] for v in ((y1, x1), (y2, x1), (y2, x1 + y2 - y1))}
            for (y1, x1) in cells.keys() for y2 in range(y1 + 1, n)]
    for t1, t2 in itertools.combinations(triangles, 2):
        model.AddBoolOr(t1.symmetric_difference(t2))
    model.Minimize(sum(cells.values()))
    solver.Solve(model)
    return len(cells) - round(solver.ObjectiveValue())


for n in itertools.count(2):
    print('a(%d) = %d' % (n, solve(n)))

Utiliza el solucionador CP-SAT de Google OR-Tools .

Después de ejecutarse durante ~ 30 segundos, genera lo siguiente:

a(2) = 3
a(3) = 4
a(4) = 5
a(5) = 7
a(6) = 9
a(7) = 11
a(8) = 13

n=9$n^6$a(9)=15n=8

Cómo funciona

$T_1$$T_2$$T_1$$T_2$

Por lo tanto, la pregunta puede reformularse como un problema SAT, con una restricción para cada par de triángulos.

PD: Me gustaría incluir un ejemplo para n=8, pero estoy teniendo problemas con el solucionador SAT que aparentemente quiere mantener las soluciones para sí mismo.

Obteniendo las soluciones del programa @ Delfad0r

Extendí el programa de @ Delfad0r a soluciones de salida. También da resultados intermedios, por lo que obtienes resultados como este:

Solving n = 8:
a(8) >= 9
a(8) >= 10
a(8) >= 11
a(8) >= 12
a(8) >= 13
       o
      . o
     . o o
    . o o .
   o o . o o
  o o o o . .
 o . . o o o .
o . . o . o o o
a(8) = 13

Solving n = 9:
a(9) >= 10
a(9) >= 13
a(9) >= 14
a(9) >= 15
        o
       o o
      o . .
     o . o o
    . o . o o
   . o o o o o
  o o o . o . .
 o o o . . . o o
. o o o o o o . .
a(9) = 15

Este cálculo tomó varias horas.

Si se impacienta y presiona Ctrl-Cdespués de encontrar una solución posiblemente no óptima, el programa mostrará esa solución. Por lo tanto, no lleva mucho tiempo obtener esto:

                   .
                  o o
                 . o o
                . o o o
               o o . o o
              o . o o o .
             o . o . o o o
            . o o o o o . o
           o . . o o o o o o
          o o o o o o o o o .
         o o . o o o o . o o o
        o o o o o o . o . o o o
       o . o o . o o o o o o o o
      o o o . o o o o o . o o o o
     o o o . o o o o o o o o . . o
    o o o o o o o o o o o . o . . o
   o o o o . o o o o . o o o o o . o
  o o o o o o o o . o o . . o o o o .
 o o o o . o o . o . o o o o o o . o o
o o . o o . o o o o . o o o . o o o o o
a(20) >= 42

Aquí está el programa extendido:

import itertools
from ortools.sat.python import cp_model

class ReportPrinter(cp_model.CpSolverSolutionCallback):

    def __init__(self, n, total):
        cp_model.CpSolverSolutionCallback.__init__(self)
        self.__n = n
        self.__total = total

    def on_solution_callback(self):
        print('a(%d) >= %d' %
              (self.__n, self.__total-self.ObjectiveValue()) )

def solve(n):
    model = cp_model.CpModel()
    solver = cp_model.CpSolver()
    cells = {
        (y, x): model.NewBoolVar(str((y, x)))
        for y in range(n) for x in range(y + 1)}
    triangles = [
            {cells[v] for v in ((y1, x1), (y2, x1), (y2, x1 + y2 - y1))}
            for (y1, x1) in cells.keys() for y2 in range(y1 + 1, n)]
    for t1, t2 in itertools.combinations(triangles, 2):
        model.AddBoolOr(t1.symmetric_difference(t2))
    model.Minimize(sum(cells.values()))
    print('Solving n = %d:' % n)
    status = solver.SolveWithSolutionCallback(model, ReportPrinter(n,len(cells)))
    if status == cp_model.OPTIMAL:
        rel = '='
    elif status == cp_model.FEASIBLE:
        rel = '>='
    else:
        print('No result for a(%d)\n' % n)
        return
    for y in range(n):
        print(' '*(n-y-1), end='')
        for x in range(y+1):
            print('.o'[solver.Value(cells[(y,x)])],end=' ')
        print()
    print('a(%d) %s %d' % (n, rel, len(cells) - solver.ObjectiveValue()))
    print()

for n in itertools.count(2):
    solve(n)

Python 3

Basado en gran medida en la respuesta de Delfad0r , en su mayoría sigue la misma progresión lógica al verificar los pares de triángulos y validar la configuración si no contiene pares de triángulos que fallen esta validación. Como no utilicé ninguna biblioteca además de itertools y copy, tengo control total sobre guardar los ejemplos encontrados en todo el programa.

examples = dict() # stores examples by key pair of n to a tuple with the triangle and number of lights turned off

for n in range(3, 8):
    tri = [] # model of the triangle, to be filled with booleans representing lights
    tri_points = [] # list of tuples representing points of the triangle
    for i in range(n):
        tri.append([True]*(i + 1))
        for j in range(i+1):
            tri_points.append((i, j))

    t_set = [] # list of all possible triangles from tri, represented by lists of points
    for i in range(n):
        for j in range(len(tri[i])):
            for k in range(1, n - i):
                t_set.append([(i, j), (i + k, j), (i + k, j + k)])

    from itertools import combinations
    import copy

    # validates whether or not a triangle of n lights can have i lights turned off, and saves an example to examples if validated
    def tri_validate(x):
        candidate_list = list(combinations(tri_points, x))
        tri_pairs = list(combinations(t_set, 2))
        for candidate in candidate_list:
            temp_tri = copy.deepcopy(tri)
            valid = False
            for point in candidate:
                (row, col) = point
                temp_tri[row][col] = False
            for pair in tri_pairs:
                valid = False
                (tri1, tri2) = pair
                for point in tri1:
                    if not valid:
                        if point not in tri2:
                            (row, col) = point
                            if temp_tri[row][col]:
                                valid = True
                for point in tri2:
                    if not valid:
                        if point not in tri1:
                            (row, col) = point
                            if temp_tri[row][col]:
                                valid = True
                if not valid:
                    break
            if valid:
                examples[n] = (temp_tri, x)
                return True
        return False

    # iterates up to the point that validation fails, then moves on to the next n
    for i in range(len(tri_points)):
        if tri_validate(i + 1):
            continue
        break

El problema es que no es muy eficiente. Funciona muy rápido n=5, pero comienza a disminuir significativamente más allá de ese punto. En n=6, tarda alrededor de un minuto en correr, y es mucho más lento n=7. Me imagino que se pueden hacer muchas correcciones de eficiencia con este programa, pero es un borrador rápido de una buena solución con mucha más flexibilidad para verificar el funcionamiento interno de este método. Trabajaré progresivamente en esto con el tiempo.