\\ * * * * * * * * * * * * * * *
\\ coordination sequence library
\\ * * * * * * * * * * * * * * *

width = -1
height = -1
coord = matrix(0,0)
labels = matrix(0,0)

init(screen) = {
	height = #screen;
	width = #screen[1];
	coord = matrix(width, height, x, y, []);
	labels = matrix(width, height, x, y, "");

	for (s=1, #screen,
		my (y=#screen - s);
		my (row = Vec(screen[s]));
		for (x=0, #row-1,
			my (d = row[x+1]);
			if (    (d >= "A" && d <= "Z")
			     || (d >= "a" && d <= "Z")
			     || (d >= "0" && d <= "9"),
				eval(Str("p" d "=" x+I*y));

				labels[1+x, 1+y] = d;
			);
		);
	);
}

up   (z) = z + I*height
down (z) = z - I*height
right(z) = z + width
left (z) = z - width

\\ unidirectional link
link1(z, zz) = {
	my (x=real(z) % width, y=imag(z) % height);
	coord[1+x, 1+y] = concat(coord[1+x, 1+y], zz-z);
}

\\ bidirectional link
link(a, b) = link1(a, b); link1(b, a);

\\ bidirectional links from same origin
links(a, bb) = for (i=1, #bb, link(a, bb[i]));

\\ neighbours
neighbours(z) = {
	my (x=real(z) % width, y=imag(z) % height);
	my (n = apply (c -> c + z, coord[1+x, 1+y]));
	\\ print ("#" z " -> " n);
	return (n);
}

\\ complex see/seen routines

big limit = 100 000
big values = Set([])
small values = vector(4, i, [])

see index(z) = {
	my (s=1, x=real(z), y=imag(z));
	if (x<1,
		s += 1;
		x = 1-x;
	);
	if (y<0,
		s += 2;
		y = -y;
	);
	return ([s,x,y]);
}

see(z) = {
	my (index = see index(z));
	if (max(index[2],index[3]) > big limit,
		big values = set union(big values, Set([z])),
		if (#small values[index[1]] < index[2],
			small values[index[1]] = concat(small values[index[1]], vector(index[2] - #small values[index[1]]))
		);
		small values[index[1]][index[2]] = bit or(small values[index[1]][index[2]], 2^index[3])
	);
}

seen(z) = {
	my (index = see index(z));
	if (max(index[2],index[3]) > big limit,
		return (set search(big values, z)),
		if (#small values[index[1]] < index[2],
			return (0),
			return (bit test(small values[index[1]][index[2]], index[3]))
		);
	);
}

\\ compute coordination sequence

compute(z, mx) = {
	my (current = [z]);

	for (n=0, mx,
		\\ print ("#" current);
		print (n " " #current);
		if (n==mx, break);
		for (c=1, #current,
			see(current[c]);
		);

		my (others = concat(apply(neighbours, current)));
		current = select(p -> !seen(p), Set(others));
	);
}

{
	init([
		"................",
		"........O.......",
		"..M.......N...*.",
		"....K.......L...",
		"................",
		"H...I...J.......",
		"................",
		"....F.......G...",
		"..D.......E...*.",
		"........C.......",
		"................",
		"....A...B...*..."
	]);

	links(pA, [ pB, pD, down(pM) ]);
    links(pB, [ pC, down(right(pH)), down(pO) ]);
    links(pC, [ pE, pF ]);
    links(pD, [ pF, pH, left(pL) ]);
    links(pE, [ pG, pJ ]);
    links(pF, [ pI ]);
    links(pG, [ pL, down(right(pM)) ]);
    links(pH, [ pM ]);
    links(pI, [ pJ, pK ]);
    links(pJ, [ pN ]);
    links(pK, [ pM, pO ]);
    links(pL, [ pN ]);
    links(pM, [ ]);
    links(pN, [ pO ]);
}

allocatemem(2^30)
compute(pE, 1000)

quit