# Python code for OEIS A260550
# Michael S. Branicky, Jan 15 2021

from itertools import product
from numba import njit
@njit
def decomposable(A):
  M = max(A)
  if M == 1: return False
  a, b, c, d = A
  for e in range(1, min(a, b, M)):
    for f in range(1, min(a-e+1, b-e+1, M)):
      for g in range(1, min(c, d, M)):
        for h in range(1, min(c-g+1, d-g+1, M)):
          for i in range(1, min(a//e, c//g)+1):
            for j in range(1, min(b//e, d//g)+1):
              for k in range(1, min((a-e*i)//f, (c-g*i)//h)+1):
                for l in range(1, min((b-e*j)//f, (d-g*j)//h)+1):
                  if a != e*i + f*k: continue
                  if b != e*j + f*l: continue
                  if c != g*i + h*k: continue
                  if d != g*j + h*l: continue
                  return True
  return False

def a(n): 
  c = n**4 - (n - 1)**4 # count matrices with at least one 1
  return c + sum(not decomposable(A) for A in product(range(2, n+1), repeat=4))

def a2(n): # don't sum over matrices with 1's
  return sum(not decomposable(A) for A in product(range(2, n+1), repeat=4))

def a2_incremental(n): # incrementally, only sum matrices with an n entry
  return sum(not decomposable(A) for A in product(range(2, n+1), repeat=4)
             if any(A[i]==n for i in range(4)))


# RUN CODE, PRINT/COMPARE/SAVE RESULTS

from time import time
time0 = time()

data = [1, 15, 75, 237, 559, 1157, 2055, 3471, 5449, 8131, 11633, 16361, 22041, 29349, 38329, 48839, 61325, 76479, 93957, 114717, 138041, 164153, 194505, 229625, 268259, 311031, 359719, 413245, 472145, 537835, 608837, 688121, 774877, 867549, 971403, 1080637, 1198233, 1326059, 1467029, 1617451, 1777881, 1948219, 2132381, 2329081, 2539351]

prevrhomb = 0
an = 0
for n in range(1, 1001):
  # an = a(n)
  rhomb = n**4 - (n-1)**4 # Rhombic dodecahedral numbers, A005917
  # an = a2(n) + rhomb
  an += (rhomb - prevrhomb) + a2_incremental(n)
  prevrhomb = rhomb 
  if n < len(data):
    print(n, an, data[n-1]==an, time()-time0)
  else:
    print(n, an, time()-time0)
  with open('b260550.txt', 'a') as f:
    f.write(f"{n} {an}\n")