f(n) = {
	my (f=factor(n));
	sum(k=1, #f~,
		my (i=primepi(f[k,1])-1, b=Vecrev(binary(f[k,2])));
		sum (l=1, #b, Mod(b[l],2)*'x^i*'y^(l-1))
	);
}

g(q) = {
	prod(i=0, poldegree(q,'x),
		my (p=prime(1+i),qx=polcoef(q,i,'x));
		prod(j=0, poldegree(qx),
			if (polcoef(qx,j,'y),
				p^(2^j), 1
			);
		);
	);
}

T(n,k)=g(f(n)*f(k))

t=0; for (d=1, 100, for (n=1, d, print (t++ " " T(n,d+1-n))))

quit