import numpy as np
from math import floor, sqrt
def A361207_list(maxn):
    #Returns a list of maxn terms
    #Variables
    #k = next number to be added to the array
    #rk = smallest shell of the square spiral that contains maxn'th term
    #maxk = largest number k added to array A
    #A = numpy aray of zeros to which numbers k are added
    #B = list of terms a(n)
    #c = center coordinates of array A
    #h = odd dimensions of array A
    #s = current square to which numbers k are added around
    #x,y = coordinates within array A
    #r1,r2 = length of sides of the square spiral
    
    #First creates array A with dimensions based on maxn 
    rk = floor((1/4)*(sqrt((4*maxn)-3)+3))
    maxk,h = 8*(rk**2)+(3*rk)+1,(rk*4)+1
    A,B,c = np.zeros((h,h),dtype=int),[],int((h-1)/2)
    A[c,c] = 1
    #Adds numbers k through kmax to array A
    k,s = 2,1
    while k <= maxk:
        g = np.where(A == s)
        x,y = int(g[0]),int(g[1])
        if A[x][y+1] == 0:
            A[x][y+1] = k
            k += 1
        if A[x+1][y] == 0:
            A[x+1][y] = k
            k += 1
        if A[x][y-1] == 0:
            A[x][y-1] = k
            k += 1
        if A[x-1][y] == 0:
            A[x-1][y] = k
            k += 1
        else: s += 1
    #Reads the array A as a square spiral and creates a list B, of terms up to maxn
    x,y = c,c
    B.append(int(A[c][c]))
    for i in range(1,rk+1):
        r1,r2 = (2*(i-1))+2,(2*i)+1
        for k in range(1,r1):
            x += 1
            B.append(int(A[y][x]))
        for j in range(1,r1):
            y += 1
            B.append(int(A[y][x]))
        for m in range(1,r2):
            x -= 1
            B.append(int(A[y][x]))
        for n in range(1,r2):
            y -= 1
            B.append(int(A[y][x]))
    for i in range(0,len(B)-maxn): B.pop()
    return B

print(A361207_list(67))