addhelp(upto, "searches a(n) for n >= 6. If a(n) = k is found for some n >= 6, it prints [n, k]. Regularily it prints a number m. That number m tells you that a(n) >= (m*inc)^2 for n >= 6 for any n where a(n) isn't found yet. "); upto({n = oo}) = {inc = 10^4; printat = inc; if(n == oo, ulim = oo , ulim = sqrtint(n) ); for(i = 1, ulim, if(i > printat, print1(i \ inc", "); printat += inc; ); c = i*sigma(i); if(c <= n, d = divisors(c); r = 0; if(#d < 12, next(1)); for(j = 1, (#d + 1) >> 1, if(d[j] * sigma(d[j]) == c, r++; ) ); if(r >= 6, print(); print(); print([r, c]); print(); print(); ) ) ) } addhelp(qsols, "The number of solutions x to x*sigma(x) = n."); qsols(n) = { if(n > 10^80, warning("A stopping rule was based on oeis.org/a004394 for input <= 10^80. Please revise over there.") ); my(res = 0); fordiv(n, d, if(sigma(d)*d==n, res++; ); if(d^2 > n, return(res) ) ); res }