addhelp(upto, "all positive integers k where k/A007954(k) <= f and is integer."); upto(f) = { ints = vector(floor(f)); fracs = vector(floor(f)); maxn = f; for(i = 1, 9, process([[[i], i]]); ); res = vector(vecsum(ints)); t = 0; ind = 0; for(i = 1, #ints, for(j = 1, ints[i], t++; ind++; res[ind] = t; ); t+=fracs[i]; ); res } process(v) = { my(c); for(i = 1, #v, processperms(v[i][1]); c = childrenbounded(v[i][1], v[i][2]); process(c); ) } processperms(v) = { vp = vecprod(v); forperm(v, p, c = fromdigits(Vec(p))/vp; if(c > maxn, return ); if(denominator(c) > 1, fracs[floor(c)]++ , ints[floor(c)]++ ) ) } childrenbounded(v, p) = { my(res = List(), v = concat(0, v)); for(i = 1, v[2], v[1] = i; fr = fromdigits(v); if(fr\(p*i) <= maxn, listput(res, [v, p*i]); ) ); res }