\\ Kevin Ryde, July 2023 \\ \\ This is some PARI/GP code to calculate terms of A364133, \\ which are those n where A000127(n) has a new record high \\ 2-adic valuation (number of low 0 bits). \\ \\ Terms n < 2^num_bits are found for given "num_bits". \\ The strategy is a depth-first search through the bits of n \\ from low to high and pruning the search at valuations found. \\ Speed is good for a few thousand bits, though no final \\ results are known until the end of the whole search. \\ \\ Requires PARI/GP 2.9.x or higher for "oo" infinity. default(strictargs,1); default(recover,0); \\ 24*A000127(n) A364133_poly = x^4 - 6*x^3 + 23*x^2 - 18*x + 24; \\ Return a vector of all terms of A364133 < 2^num_bits. A364133_vector_to_bits(num_bits) = { my(best_n = vector(4*num_bits+3,v, oo), inc_limit = num_bits - 1, n=0); \\ n = candidate term of A364133 while(1, my(t = subst(A364133_poly,'x,n), v = valuation(t,2)); if(n < best_n[v] && n > 0, best_n[v] = n); v = min(v, inc_limit); while(bittest(n,v), v--); if(v < 0, break); n = bitand(n, bitneg(0,v)) + 1<<v); \\ find records my(ret=List([]), m=oo); forstep(v=#best_n,1,-1, if(best_n[v] < m, m=best_n[v]; listput(ret,m))); Vecrev(ret); } { my(v=A364133_vector_to_bits(20)); print1("A364133 = "); for(i=1,#v, print1(v[i],",")); print("..."); } \\----------------------------------------------------------------------------- \\ Implementation Notes \\ -------------------- \\ \\ P(n) = 24*A000127(n) = subst(A364133_poly,'x,n) ignores the \\ /24 part of A000127 since \\ \\ valuation(24*k) = valuation(k) + 3 \\ \\ which is +3 from factor 8, and no change from factor 3. \\ \\ n values are considered by a depth-first search through \\ bit values 0 or 1 at bit positions from low to high. \\ For example, num_bits = 3 would consider \\ \\ n = 000 n values considered in \\ 100 order of their bit reversal \\ 010 \\ 110 \\ 001 \\ 101 \\ 011 \\ 111 \\ \\ Bit positions 0..i of n (position 0 the least significant) \\ determine the corresponding bit positions 0..i in P(n), \\ \\ v = valuation(P(n)) means P(n) ends ..100..00 in bit \\ positions 0..v. n is a residue class in the sense that \\ any with the same bits 0..v has the same valuation, \\ \\ valuation(P(m)) = v \\ for all m in residue class n mod 2^(v+1) \\ \\ This means bit positions v+1 and above in n can be skipped, \\ as any bits there don't change the valuation. \\ \\ The search order and skipping means n < 2^(v+1), ie. n is \\ always the smallest of the residue class. This is since \\ any reduction "n mod 2^(v+1)" would mean clearing some high \\ 1-bits of n and that would be earlier in depth-first order. \\ \\ A given valuation v might occur at various n. The best, \\ meaning the smallest, is held \\ \\ best_n[v] = n smallest having valuation(P(n))=v \\ \\ \\ At a given v, skip bits above 0..v means increment the bit \\ at position v in n. If already a 1-bit then it can't be \\ increased, so drop to try an increment at v-1, or v-2, etc. \\ Those drops are backtracking, or back up in the search tree. \\ \\ The code finds the increment position at the bittest() loop. \\ At a 0-bit, the increment means set that bit to 1 and clear \\ to 0s all bits above (bitand() of a mask). \\ \\ If n has no 0 bits among 0..v then the loop reaches v = -1 \\ (by bittest(n,-1) = 0), and there are no more n to consider. \\ \\ \\ The final step is to pick record highs from best_n[], \\ \\ n = best_n[v] is index of a record high valuation \\ when best_n[v'] > n for all v' > v \\ \\ The code does this going down from biggest v and each new \\ smallest best_n[v] is a term. \\ \\ best_n[] might have holes where a given valuation never \\ occurred. Those are "oo" infinity and they have no effect \\ on the records. All n < 2^num_bits have been considered \\ and the valuations occurring among those n are all that \\ matter. \\----------------------------------------------------------------------------- print("end");