\\ Kevin Ryde, August 2022 \\ This is some PARI/GP code for sequences A120244, A120245, A120246, which \\ are certain recurrences using n = m*(m+1)/2 + k for 1 <= k <= m+1. \\ \\ The individual term functions are successive sqrt reductions and suit \\ large and small n. The vector functions use only additions. \\ \\ Requires PARI/GP 2.13 or higher for sqrtint() "remainder" parameter. default(strictargs,1); default(recover,0); \\----------------------------------------------------------------------------- \\ Individual Terms \\ m_and_k(n) returns a vector of two elements [m,k] which satisfy, for n>=2, \\ n = m*(m+1)/2 + k in range 1 <= k <= m+1 m_and_k(n) = { my(r, s=sqrtint(n<<1,&r)); if(r>s, [s, (r-s)>>1], \\ small k [s-1, (r+s)>>1]); \\ big k } \\ A120244 a((m(m+1)/2 + k) = a(m)+ k A120244(n) = { my(ret=1, k); while(n>1, [n,k]=m_and_k(n); ret+=k); ret; } \\ A120245 a((m(m+1)/2+k) = a(m) + a(k), A120245(n) = { if(n==1,1, my(m,k); [m,k]=m_and_k(n); self()(m) + self()(k)); } \\ A120246 a((m(m+1)/2+k) = m + a(k) A120246(n) = { my(ret=1, m); while(n>1, [m,n]=m_and_k(n); ret+=m); ret; } \\----------------------------------------------------------------------------- \\ Vectors \\ mkrec_vector() returns vector v of the first len many terms from one of \\ the m,k recurrences. \\ m_rec=1 means add a(m), or m_rec=0 means add m. \\ k_rec=1 means add a(k), or k_rec=0 means add k. \\ mkrec_vector(len,m_rec,k_rec) = { my(v=vector(len,n,1)); my(m=1,k=1); for(n=2,len, v[n] = if(m_rec,v[m],m) + if(k_rec,v[k],k); if(k<=m, k++, m++;k=1)); v; } \\ A120244_vector(len) returns the first len many terms of A120244, so \\ v[1..len] = A120244(1..len). Similarly respectively the other functions. \\ A120244_vector(len) = mkrec_vector(len,1,0); A120245_vector(len) = mkrec_vector(len,1,1); A120246_vector(len) = mkrec_vector(len,0,1); \\ checks A120244_vector(1000) == vector(1000,n,A120244(n)) || error(); A120245_vector(1000) == vector(1000,n,A120245(n)) || error(); A120246_vector(1000) == vector(1000,n,A120246(n)) || error(); \\----------------------------------------------------------------------------- \\ Samples { my(len=17); print1("A120244 = "); for(n=1,len, print1(A120244(n),", ")); print("..."); print1("A120245 = "); for(n=1,len, print1(A120245(n),", ")); print("..."); print1("A120246 = "); for(n=1,len, print1(A120246(n),", ")); print("..."); print(); print("m_and_k()"); for(n=2,10, my(m,k); [m,k] = m_and_k(n); if(k==1, print1(" n="n": ")); print1(" ["m","k"]"); if(k==m+1, print())); my(len=17); my(m_of_n(n) = m_and_k(n)[1]); print1("m(n) = A003056(n-1) = "); for(n=2,len, my(m,k); [m,k]=m_and_k(n); print1(m","); if(k == m+1, print1(" "))); print("... for n>=2"); my(k_of_n(n) = m_and_k(n)[2]); print1("k(n) = A002260(n) = "); my(L=2,r=1); for(n=2,len, my(m,k); [m,k]=m_and_k(n); print1(k","); if(k == m+1, print1(" "))); print("... for n>=2"); } \\----------------------------------------------------------------------------- \\ Notes \\ M and K \\ ------- \\ \\ m_and_k(n) calculates together \\ \\ m = A003056(n-1) \\ k = A002260(n) \\ \\ The code is a fairly usual quadratic solution, finding m to go with k>=1, \\ \\ m = floor(x) where solve for \\ x*(x+1)/2 + 1 = n \\ x^2 + x - 2*(n-1) = 0 \\ x = ( -1 + sqrt(1+8*(n-1)) )/2 \\ \\ But floor(x) is the same for any n+frac with 0 <= frac < n so simplify by \\ substituting n -> n+7/8 \\ x = floor( sqrt(2*n) - 1/2 ) \\ \\ The code takes advantage of the sqrtint() "remainder" to apply the -1/2 \\ rounding with a comparison. \\ A120244 M Descents \\ ------------------ \\ The A120244 recurrence descends into m only so the code can be iterative \\ and sum the various k seen. \\ \\ The number of iteration steps is modest since each step roughly halves \\ the bit length. The number of steps grows monotonically with n (seen by \\ induction) and new highs are at \\ \\ n = A006894 = 1, 2, 4, 11, 67, 2279, 2598061, ... \\ \\ since the smallest n' which descends to n is when \\ \\ m=n, k=1 so that n' = n*(n+1)/2 + 1 \\ \\ as per the formula in A006894. \\ \\ A006894 grows rapidly, roughly doubling in bit length each time. For \\ example A006894(30) is over 20 megabytes. \\ A120246 k Descents \\ ------------------ \\ The A120246 recurrence descends into k only so the code can be iterative \\ and sum the various m seen. The number of iteration steps for n is \\ \\ A057945(n-1) = 0, 1, 2, 1, 2, 3, 1, 2, 3, 2 \\ n = 1 2 3 4 ... \\ \\ This is since each step subtracts a triangular number m*(m+1)/2 from n, \\ to leave k, but here with offset 1 throughout so that k>=1. \\ \\ Locations of new highs in A057945 are noted there. Here the offset 1 \\ means record high number of k descents are at \\ \\ n = A006893 + 1 \\ = A007501 \\ = 2, 3, 6, 21, 231, 26796, 359026206, ... \\ \\ This can be seen since the smallest n' which descends to n is when \\ \\ m=n-1, k=n so that n' = (n-1)*n/2 + n \\ = n*(n+1)/2 \\ = Triangular(n) \\ \\ and hence the iterated triangular numbers per the formula in A007501. \\ \\ A007501 grows rapidly, roughly doubling in bit length each time. For \\ example A007501(30) is over 50 megabytes. \\ A120245 M and K Recursions \\ -------------------------- \\ \\ The A120245 recurrence descends to both m and k and is implemented by a \\ simple recursion. \\ \\ There's no attempt to notice m or k values which repeat in the descents. \\ For a "random" n, repeats are rare beyond small n, and small n are fast \\ anyway. \\ \\ The greatest recursion depth in the calculation increases monotonically \\ with n. This is seen by induction, with the number of descents \\ \\ descents \\ a(m) greatest k < m \\ a(k) greatest k = m+1 at "record highs" \\ a(m) equal a(k) otherwise \\ \\ The record highs are the same indices n as for the A120246 k-only \\ descents described above. \\ Vector Functions \\ ---------------- \\ \\ These vector functions share a common loop subroutine with parameters to \\ select how to handle m and k. The loop could easily be written out in \\ full in each function, but a subroutine reduces repetition. \\----------------------------------------------------------------------------- print("end");