a000043 = [ 2, 3, 5, 7, 13, 17, 19, 31, 61, 89, 107, 127, 521, 607, 1279, 2203, 2281, 3217, 4253, 4423, 9689, 9941, 11213, 19937, 21701, 23209, 44497, 86243, 110503, 132049, 216091, 756839, 859433, 1257787, 1398269, 2976221, 3021377, 6972593, 13466917, 20996011, 24036583, 25964951, 30402457, 32582657, 37156667, 42643801, 43112609 ]

explore(v, sig, tau, k, i) = {
	if (sig == tau^k,
		best = max(best, v);
	);
	for (j=i, oo,
		if (j>#a000043,
			\\ at this point, we may miss a number such that sigma <= tau^k, so we must quit
			quit
		);

		my (xsig = sig * 2^a000043[j], xtau = 2*tau);
		if (xsig <= xtau^k,
			explore(v*(2^a000043[j]-1), xsig, xtau, k, j+1),
			\\ adding a bigger Mersenne prime won't give a number such that sigma <= tau^k, so we can stop here
			break
		);
	);
}

a(n) = {
	best = 1;
	explore(1, 1, 1, n, 1);
	return (best);
}

for (n=1, 1000, print (n " " a(n)))

quit