Yes they've been done to death, but here are two minor remarks that came up in a discussion with a (talented) youngster: [...] if a sequence obeys a (homogeneous) linear recurrence of order <= n, then the determinant of an n+1 by n+1 band matrix formed from any 2n+1 consecutive elements = 0. E.g., the fibonaccis are a 2nd order recurrence, so c393) genmatrix(lambda([i,j],fib(integer+i+j-2)),3) [ fib(integer) fib(integer + 1) fib(integer + 2) ] [ ] (d393) [ fib(integer + 1) fib(integer + 2) fib(integer + 3) ] [ ] [ fib(integer + 2) fib(integer + 3) fib(integer + 4) ] (c394) ratsimp(fibtophi(det(%))) (d394) 0 Now suppose you know that the square triangular numbers are a 3rd order recurrence, and you only know the first four: 0, 1, 36, and 1225; and you want more. The matrix needs 2n+1=7 values--what to do? Think about what the recurrence could possibly do for negative values. Are the squares of -1, -6, and -35 equal to the triangulars of index -2, -9, and -50? Yes! So we actually know seven sequence elements, and we only need six: c387) genmatrix(lambda([i,j],[36,1,0,1,36,1225,x][i+j-1]),4) [ 36 1 0 1 ] [ ] [ 1 0 1 36 ] (d387) [ ] [ 0 1 36 1225 ] [ ] [ 1 36 1225 x ] (c388) sqrt(solve(det(%))) (d388) [sqrt(x) = 204] so 41616 is the next one. The formula for the nth, showing which square equals which triangular: 2 n 2 n (sqrt(2) + 1) (1 - sqrt(2)) 2 (---------------- - ----------------) = 4 sqrt(2) 4 sqrt(2) n n n n (sqrt(2) + 1) (sqrt(2) - 1) 2 (sqrt(2) + 1) (sqrt(2) - 1) 2 ((-------------- - --------------) + 1) (-------------- - --------------) 2 2 2 2 ---------------------------------------------------------------------------. 2 (Interesting computer algebra challenge: verify this identity noninductively.) I don't recall seeing the simple observation that this sequence is trivially constructed from the successive approximants to sqrt(2): Given a rational approximation a/b, the fraction 2 b/a errs on the opposite side, and their average would be the next Newton approximation. (Starting with a=b=1.) These form an exponentially sparse subsequence of the continued fraction (1+/2+/2+/2+/...) approximants resulting from taking the mediant, (a + 2 b)/(b + a), instead of the average. E.g., a 1 3 7 17 41 ... b 1 2 5 12 29 ... Then the square triangular numbers are just a^2 b^2 = binom(a^2+(n mod 2),2): (c607) define(inversetriang(n),rhs(solve(binomial(x+1,2) = n,x)[2])) sqrt(8 n + 1) - 1 (d607) inversetriang(n) := ----------------- 2 (c608) block([m : matrix([1,1],[1,0]),fancy_display : false], thru 9 do (print(m,(m[1,1]/m[2,1])^2,inversetriang((m[1,1]*m[2,1])^2)),m:m.matrix([2,1],[1,0])))$ [a a'] -1 [ ] a^2/b^2 triang (a^2 b^2) = a^2-(n+1 mod 2) [b b'] -------------------------------------------------------- [ 1 1 ] 1 [ ] - 1 [ 1 0 ] 1 [ 3 1 ] 9 [ ] - 8 [ 2 1 ] 4 [ 7 3 ] 49 [ ] -- 49 [ 5 2 ] 25 [ 17 7 ] 289 [ ] --- 288 [ 12 5 ] 144 [ 41 17 ] 1681 [ ] ---- 1681 [ 29 12 ] 841 [ 99 41 ] 9801 [ ] ---- 9800 [ 70 29 ] 4900 [ 239 99 ] 57121 [ ] ----- 57121 [ 169 70 ] 28561 [ 577 239 ] 332929 [ ] ------ 332928 [ 408 169 ] 166464 [ 1393 577 ] 1940449 [ ] ------- 1940449 [ 985 408 ] 970225 where [a a'] 1 1 2 1 n [ ] = ( )( ) and a'/b' is the previous a/b. [b b'] 1 0 1 0 --rwg