ekg    = [1,2]	\\ A064413
inv    = [1,2]	\\ A064664
nb     = 2		\\ number of known terms in ekg
unseen = 3		\\ least unseen term in ekg

\\ A064413(n)
f(n) = {
	while (nb<n,
		for (v=unseen, oo,
			if ((v>#inv || inv[v]==0) && gcd(ekg[nb],v)>1,
				if (nb++>#ekg,
					ekg = concat(ekg, vector(#ekg));
				);
				ekg[nb] = v;
				while (v>#inv,
					inv = concat(inv, vector(#inv));
				);
				inv[v] = n;
				while (unseen<=#inv && inv[unseen],
					unseen++;
				);
				break;
			);
		);
	);
	return (ekg[n]);
}

\\ A064664
g(n) = {
	while (n>#inv || inv[n]==0,
		f(nb+1);
	);
	inv[n];
}

\\ neighbors of number n
\\		  k  --   k+1
\\		f(k) -- f(k+1)
neighbors(n) = {
	my (w=[]);

	if (n>1,
		w = concat(w, n-1);
	);
	w = concat(w, n+1);

	my (m=g(n));
	if (m>1,
		w = concat(w, f(m-1));
	);
	w = concat(w, f(m+1));

	w
}

\\ see / seen lib
s = 0
see(v) = s = bitor(s, 2^v)
seen(v) = bittest(s, v)

{
	p = [];
	w = [1];			\\ nodes at distance d
	for (d=0, 50,
		print (d" "#w);
		apply(see, w);

		w = select(v -> !seen(v), Set(concat(apply(neighbors, w))));
	);
}

quit