Calculating (3 + sqrt(5))^n exactly

Question

Today your goal is to find integers a and b given non-negative integer n such that:

You should write a program or a function that takes parameter n and outputs a and b in a format of your choice.

Standard loopholes apply. Additionally, it's intended that you implement the above problem using basic arithmetic yourself. So you may not use built-in exact algebra functionality, rationals, or functions implementing non-trivial mathematical constructs (for example the Lucas sequence).

Shortest code in bytes wins.

---

Example input/output:

0 → 1, 0
1 → 3, 1
2 → 14, 6
3 → 72, 32
4 → 376, 168
5 → 1968, 880
6 → 10304, 4608
7 → 53952, 24128
8 → 282496, 126336
9 → 1479168, 661504

Answer 1

Dyalog APL, 18 bytes

((3∘×+5 1×⌽)⍣⎕)1 0

This is a program that takes input through ⎕.

 (         )         Monadic train:
  3∘×                3 times argument
     +               Plus
      5 1×⌽          (5 1) times the reverse
(           ⍣⎕)      Apply that function (input) times
               1 0   starting with (1 0)

The features used here were implemented well before April 2015, making this answer valid.

Try it here. Note that tryapl.org is a limited subset of Dyalog and does not support ⎕.

Answer 2

Octave, 26 bytes

[3 5;1 3]**input('')*[1;0]

Because (a+b*sqrt(5)) * (3+sqrt(5)) = (3a+5b) + (a+3b) * sqrt(5),

multiplying input vector

| 1 |    /* a = 1 */
| 0 |    /* b = 0 */

which stands for 1 = (3+sqrt(5))^0 by matrix

| 3 5 |
| 1 3 |

seems natural. Instead of looping n times, we rather raise the matrix to the power of n and then multiply it by input vector.

Answer 3

Python 2, 50

a=1;b=0
exec"a,b=3*a+5*b,3*b+a;"*input()
print a,b

Multiplies by 3+sqrt(5) repeatedly by its action on the pair (a,b) representing a+b*sqrt(5). Equivalent to starting with the column vector [1,0] and left-multiplying n times by the matrix [[3,5],[1,3]].

Answer 4

Julia, 22 20 bytes

n->[3 5;1 3]^n*[1;0]

This creates a lambda function which takes a single integer as input and returns a 2-element vector of integers corresponding to the solution [a, b]. To call it, give it a name, e.g. f=n->....

Start by multiplying

We can then translate the right hand side of this equation into a 2-column matrix, where the first corresponds to the coefficient of a and the second to the coefficient of b:

Multiply this matrix by itself n times, then right multiply by the column vector (1, 0), and POOF! Out pops the solution vector.

Examples:

julia> println(f(0))
[1,0]

julia> println(f(5))
[1968,880]

Answer 5

J, 20 bytes

+/@:*(3 5,.1 3&)&1 0

Multiplicate the vector [1 0] with the matrix [[3 5] [1 3]] n times.

2 bytes saved thanks to @algorithmshark.

Usage and test:

   (+/@:*(3 5,.1 3&)&1 0) 5
1968 880

   (+/@:*(3 5,.1 3&)&1 0) every i.6
   1   0
   3   1
  14   6
  72  32
 376 168
1968 880

Answer 6

Pyth, 20 bytes

u,+*3sGyeG+sGyeGQ,1Z

u which is reduce in general, is used here as an apply repeatedly loop. The updating function is G -> ,+*3sGyeG+sGyeG, where G is a 2 tuple. That function translates to 3*sum(G) + 2*G[1], sum(G) + 2*G[1]. s is sum, y is *2.

Answer 7

APL (22)

{⍵+.×⍨2 2⍴3 5 1}⍣⎕⍨2↑1

Explanation:

• {...}⍣⎕⍨2↑1: read a number, and run the following function that many times, using [1,0] as the initial input.
  • 2 2⍴3 5 1: the matrix [[3,5],[1,3]]
  • ⍵+.×⍨: multiply the first number in ⍵ by 3, the second by 5, and sum them, this is the new first number; then multiply the first number in ⍵ by 1, the second by 3, and sum those, that is the new second number.

Answer 8

Jelly, 13 bytes

5W×U++Ḥ
2Bdz¡

Try it online!

How it works

5W×U++Ḥ    Helper link. Argument: [a, b]

5W         Yield [5].
  ×U       Multiply it by the reverse of [a, b]. This yields [5b, a].
    +      Hook; add the argument to the result. This yields [a + 5b, a + b].
     +Ḥ    Fork; add the doubled argument ([2a, 2b]) to the result.
           This yields [3a + 5b, a + 3b].

2Bdz¡      Main link. Argument: n

2B         Convert 2 to binary, yielding [1, 0].
    ¡      Repeat:
  Ç            Apply the helper link...
   ³           n times.

Answer 9

Javascript, 63 61 bytes

I am using a recursive evaluation of the binomial: (x+y)^n = (x+y)(x+y)^{n-1}

New(thanks to @edc65)

F=n=>{for(i=y=0,x=1;i++<n;)[x,y]=[3*x+5*y,x+3*y];return[x,y]}

Old

F=n=>{for(i=y=0,x=1;i<n;i++)[x,y]=[3*x+5*y,x+3*y];return [x,y]}

Answer 10

Mathematica, 31

Nest[{{3,5},{1,3}}.#&,{1,0},#]&

Answer 11

CJam, 21 bytes

0X{_2$3*+@5*@3*+}li*p

Try it online.

How it works

0X       " Stack: [ 0 1 ]                                ";
li{      " Do int(input()) times:                        ";
  _2$    " Stack: [ a b ] -> [ a b b a ]                 ";
  3*+    " Stack: [ a b b a ] -> [ a b (b+3a) ]          ";
  @5*@3* " Stack: [ a b (b+3a) ] -> [ (b+3a) 5a 3b ]     ";
  +      " Stack: [ (b+3a) 5a 3b ] -> [ (b+3a) (5a+3b) ] ";
}*       "                                               ";
p        " Print topmost stack item plus linefeed.       ";
         " Print remaining stack item (implicit).        ";

Answer 12

C, 114 bytes

g(n){int i,a[2]={1,0},b[2];for(i=0;i<n;i++)*b=*a*3+5*a[1],b[1]=*a+3*b[1],*a=*b,a[1]=b[1];printf("%d,%d",*a,a[1]);}

This implements matrix multiplication the boring way. For a more fun (quote: "awesomely horrific") 238 byte solution, look no further!

f(n){int p[2][n+3],i,j,k=0,a[2]={0};for(j=0;j<n+3;j++)p[0][j]=0;*p[1]=0;(*p)[1]=1;for(j=0;j<n;j++,k=!k)for(i=1;i<n+3;i++)p[!k][i]=p[k][i-1]+p[k][i];for(i=1;i<n+2;i++)a[!(i%2)]+=p[k][i]*pow(3,n+1-i)*pow(5,(i-1)/2);printf("%d,%d",*a,a[1]);}

Unraveled:

g(n){
    int i,a[2]={1,0},b[2];
    for(i=0;i<n;i++)
        *b=3**a+5*a[1],b[1]=*a+3*b[1],*a=*b,a[1]=b[1];
    printf("%d,%d",*a,a[1]);
}

This could probably be shortened a bit. Try a test program online!

Answer 13

k2 - 22 char

Function taking one argument.

_mul[(3 5;1 3)]/[;1 0]

_mul is matrix multiplication so we curry it with the matrix (3 5;1 3) and then hit it with the functional power adverb: f/[n;x] applies f to x, n times. Again we curry it, this time with the starting vector 1 0.

  _mul[2 2#3 5 1]/[;1 0] 5
1968 880
  f:_mul[2 2#3 5 1]/[;1 0]
  f'!8  /each result from 0 to 7 inclusive
(1 0
 3 1
 14 6
 72 32
 376 168
 1968 880
 10304 4608
 53952 24128)

This will not work in Kona, because for some reason f/[n;x] isn't implemented correctly. Only the n f/x syntax works, so the shortest fix is {x _mul[(3 5;1 3)]/1 0} at 23 char.

Answer 14

ised, 25 bytes (20 characters)

({:{2,4}·x±Σx:}$1)∘1

I hoped for better, but there are just too many braces needed in ised to make it competent, the operator precedence is not optimal for golfing.

It expects the input to be in $1 memory slot, so this works:

ised '@1{9};' '({:{2,4}·x±Σx:}$1)∘1'

For n=0, the zero is skipped (outputs 1, instead of 1 0). If that's an issue, replace the final 1 with ~[2].

Answer 15

Seriously, 32 bytes, non-competing

,╗43/12`╜";)@4*≈(6*-"£n.X4ì±0`n

Hex Dump:

2cbb34332f313260bd223b2940342af728362a2d229c6e2e58348df130606e7f

Try It Onlline

Obviously not a contender for shortest, but at least the method is original. (Noting that such a problem necessarily indicates a Lucas sequence, as mentioned in the description, this program generates successive terms of the sequences using the recurrence relation

a_n = 6*a_{n-1} - 4*a_{n-2}.)

Answer 16

Haskell, 41 bytes

(iterate(\(a,b)->(3*a+5*b,a+3*b))(1,0)!!)

Usage example: (iterate(\(a,b)->(3*a+5*b,a+3*b))(1,0)!!) 8 -> (282496,126336).

Answer 17

C/C++ 89 bytes

void g(int n,long&a,long&b){if(n){long j,k;g(n-1,j,k);a=3*j+5*k;b=j+3*k;}else{a=1;b=0;}}

Formatted:

    void g(int n, long&a, long&b) {
if (n) {
    long j, k;
    g(n - 1, j, k);
    a = 3 * j + 5 * k;
    b = j + 3 * k;
} else {
    a = 1;
    b = 0;
}}

Same concept:

void get(int n, long &a, long& b) {
    if (n == 0) {
        a = 1;
        b = 0;
        return;
    }
    long j, k;
    get(n - 1, j, k);
    a = 3 * j + 5 * k;
    b = j + 3 * k;
}

The test bench:

#include <iostream>
using namespace std;    
int main() {
    long a, b;
    for (int i = 0; i < 55; i++) {
        g(i, a, b);
        cout << i << "-> " << a << ' ' << b << endl;
    }
    return 0;
}

The output:

0-> 1 0
1-> 3 1
2-> 14 6
3-> 72 32
4-> 376 168
5-> 1968 880
6-> 10304 4608
7-> 53952 24128
8-> 282496 126336
9-> 1479168 661504
10-> 7745024 3463680
11-> 40553472 18136064
12-> 212340736 94961664
13-> 1111830528 497225728
14-> 5821620224 2603507712
15-> 30482399232 13632143360
16-> 159607914496 71378829312
17-> 835717890048 373744402432
18-> 4375875682304 1956951097344
19-> 22912382533632 10246728974336
20-> 119970792472576 53652569456640
21-> 628175224700928 280928500842496
22-> 3289168178315264 1470960727228416
23-> 17222308171087872 7702050360000512
24-> 90177176313266176 40328459251089408
25-> 472173825195245568 211162554066534400
26-> 2472334245918408704 1105661487394848768

Answer 18

K, 37 bytes

f:{:[x;*(1;0)*_mul/x#,2 2#3 1 5;1 0]}

or

f:{:[x;*(1;0)*_mul/x#,(3 1;5 3);1 0]}

They're both the same thing.

Answer 19

Python 3, 49 bytes

w=5**0.5;a=(3+w)**int(input())//2+1;print(a,a//w)

although on my machine, this only gives the correct answer for inputs in the range 0 <= n <= 18.

This implements the closed form formula

w = 5 ** 0.5
u = 3 + w
v = 3 - w
a = (u ** n + v ** n) / 2
b = (u ** n - v ** n) / (2 * w)

and takes advantage of the fact that the v ** n part is small, and can be computed by rounding rather than direct calculation.

Answer 20

Scheme, 97 bytes

(define(r n)(let s([n n][a 1][b 0])(if(= 0 n)(cons a b)(s(- n 1)(+(* a 3)(* b 5))(+ a(* b 3))))))

Answer 21

C 71 bytes (60 with pre-initialised variables)

Scope for golfing yet but just to prove that C does not have to be "awesomely horrific".

f(int n,int*a){for(*a=1,a[1]=0;n--;a[1]=*a+3*a[1],*a=(5*a[1]+4**a)/3);}

If the values in a are initialised to {1,0}, we do better.

f(int n,int*a){for(;n--;a[1]=*a+3*a[1],*a=(5*a[1]+4**a)/3);}

This is iteratively using the mappings a->3a+5b, b->a+3b but avoiding a temporary variable by calculating a from the new value of b instead.

Answer 22

(non-competing) Jelly, 10 bytes

31,53Dæ*³Ḣ

Try it online!

Uses matrix. Computes ([[3,1],[5,3]]**input())[0].