CodeIQ:「ルーム・アンド・ルーフ」問題 の解答

Kawazoe(@riverplus)さんのCodeIQの問題です。
問題はこちらにお世話になりました。解答例もあります。
->https://blog.goo.ne.jp/r-de-r/e/1dcaa86093bd1be808d7e846faeb269d
こちらにも情報があります。
->https://togetter.com/li/1079332

挑戦した川添さんの問題でこれだけは自力で解けませんでした。
二項定理利用などさんざん失敗した挙句、
心ならずもtogetterにあったコメントを参考に、
wikipediaの「ほとんど整数」を参照
->https://ja.wikipedia.org/wiki/%E3%81%BB%E3%81%A8%E3%82%93%E3%81%A9%E6%95%B4%E6%95%B0

a=2+sqrt(3),b=2-sqrt(3)としますと、
a+b=4、a*b=1なのでa^n+b^nは整数で、b^n <1なのでfloor(a^n)=a^n+b^n-1
またa^(x+y)+b^(x+y)
　　=(a^x+b^x)*(a^y+b^y)-(a^x*b^y+a^y*b^x)
　　=(a^x+b^x)*(a^y+b^y)-(a*b)^x*(a^(y-x)+b^(y-x))　
　　=(a^x+b^x)*(a^y+b^y)-(a^(y-x)+b^(y-x))　　　　　(y>=x)
ですので、指数を小さくすることができます。
（a^(y-x)+b^(y-x)　が　x=yの時2、y=x+1の時4）

#julia ver 1.41,1.53,1.10.3

function solve(n,a)
 m = 10^6
 function fn(n)   # a^n+b^n
  if n == 1
    return 4
  else
    d,r = divrem(n,2)
    if r == 0
       return mod(fn(d)^2 - 2,m)
    else
       return mod(fn(d)*fn(d+1)-4,m)
    end   
  end
 end 
 ans = fn(n) - 1  
 s = a == ans ? "ok " : "no "
 println(s,ans)
end 

solve(1, 3)
solve(2, 13)
solve(3, 51) 
solve(4, 193)
solve(10, 524173)
solve(601, 868003)
solve(749501, 493043)
solve(987654, 378701)
solve(999600, 1)

別解
上の式より
f(1)=4,f(2)=f(1)*4-2,
f(n)=(a^(n-1)+b^(n-1))*(a^1+b^1)-(a^((n-1)-1)+b^((n-1)-1))
　　=f(n-1)*f(1)-f(n-2)
上の解のほうが速いですがこの解でも十分な速さです。
この解をPrologで書きました。

%SWI-Prolog ver 7.42,7.64,9.25
%:-initialization(start).   %ideone,Online Prolog Compiler
%start.

fn(2,_,_,R,R).
fn(N,M,X,Y,R):-Z is (4*Y-X) mod M,N1 is N-1,fn(N1,M,Y,Z,R).

solve(1,_,3).
solve(N,M,R):-Y is 4*4-2,(N =:= 2->R1=Y;fn(N,M,4,Y,R1)),R is R1-1.

start:-M=1000000,f(N,A),(solve(N,M,R)->(R==A->Str=" ok  ";Str=" no  ");
    write(fail)),write(" "),write(Str),writeln(N->R),fail.
start.

f(1, 3).
f(2, 13).
f(3, 51).
f(4, 193).
f(10, 524173).
f(601, 868003).
f(749501, 493043).
f(987654, 378701).
f(999600, 1).