Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

模範回答 #42

Open
ryuhei-mori opened this issue Oct 16, 2018 · 6 comments
Open

模範回答 #42

ryuhei-mori opened this issue Oct 16, 2018 · 6 comments

Comments

@ryuhei-mori
Copy link
Contributor

@prg1-2018/students

def fib_itr(n: Int) : BigInt = {
  def aux(n: Int, a: BigInt, b: BigInt) : BigInt = n match {
    case 0 => b
    case _ => aux(n-1, a+b, a)
  }
  aux(n, 1, 0)
}
@ryuhei-mori
Copy link
Contributor Author

def fib_matrix(n: Int) : BigInt = {
  def aux(n: Int, a: BigInt, b: BigInt, c: BigInt, x: BigInt, y: BigInt) : BigInt = n match {
    case 0 => y
    case _ => n % 2 match {
      case 0 => aux(n/2, a*a+b*b, a*b+b*c, b*b+c*c, x, y)
      case 1 => aux(n/2, a*a+b*b, a*b+b*c, b*b+c*c, a*x+b*y, b*x+c*y)
    }
  }
  aux(n, 1, 1, 0, 1, 0)
}

2x2 行列は4つの変数で表現できますが、対称行列であることを使うと3つの変数で表現できます。
A^n ではなく A^n b を計算する(b は長さ2のベクトル)と考えれば、もう一つの行列はベクトルで十分なので2つの変数で十分です。

@ryuhei-mori
Copy link
Contributor Author

def fib_matrix(n: Int) : BigInt = {
  case class SymMatrix2(a: BigInt, b: BigInt, c: BigInt)

  def mul(lhs: SymMatrix2, rhs: SymMatrix2) : SymMatrix2 = {
    SymMatrix2(
      lhs.a * rhs.a + lhs.b * rhs.b,
      lhs.a * rhs.b + lhs.b * rhs.c,
      lhs.b * rhs.b + lhs.c * rhs.c)
  }

  def pow(m: SymMatrix2, n: Int) : SymMatrix2 = {
    def aux(m: SymMatrix2, r: SymMatrix2, n: Int) : SymMatrix2 = n match {
      case 0 => r
      case _ => n % 2 match {
        case 0 => aux(mul(m,m), r, n/2);
        case 1 => aux(mul(m,m), mul(r,m), n/2);
      }
    }
    aux(m, SymMatrix2(1,0,1), n)
  }

  pow(SymMatrix2(1,1,0), n).b
}

もちろん case class を使うのもよいです。

@ryuhei-mori
Copy link
Contributor Author

def fib_polynomial(n: Int): BigInt = {
  def aux(n: Int, p0: BigInt, p1: BigInt, q1: BigInt, q2: BigInt): BigInt = {
    if(n == 0) p0
    else if(n % 2 == 0) aux(n/2, p0, -p1 * q1 + p0 * q2, -q1 * q1 + 2 * q2, q2 * q2)
    else aux(n/2, -p0 * q1 + p1, p1 * q2, -q1 * q1 + 2* q2, q2 * q2)
  }
  aux(n, 0, 1, -1, -1)
}

@ryuhei-mori
Copy link
Contributor Author

行列は
[a b]
[b c]
という形ですが、 a = b+c が常に満たされているので、それを利用すると掛け算の回数を減らせます。
また、
A^{2n} = (A^n)^2
A^{2n+1} = A (A^n)^2
という漸化式を使うと奇数のときに余って現われる A が常に最初の行列
[1 1]
[1 0]
なので乗算回数が減らせます。

  def fib_matrix(n: Int): BigInt = {
    case class SymMatrix2(b: BigInt, c: BigInt)

    def square(m: SymMatrix2): SymMatrix2 = {
      SymMatrix2(
        m.b * (m.b + 2 * m.c),
        m.b * m.b + m.c * m.c)
    }

    def pow(n: Int): SymMatrix2 = {
      if(n == 1) return SymMatrix2(1, 0)
      else if(n % 2 == 0) square(pow(n/2))
      else {
        val SymMatrix2(b, c) = square(pow(n/2))
        SymMatrix2(b+c, b)
      }
    }

    if(n <= 1) return n
    val m = pow(n/2)
    if(n % 2 == 0) m.b * (m.b + 2 * m.c)
    else 2 * m.b * (m.b + m.c) + m.c * m.c
  }

@ryuhei-mori
Copy link
Contributor Author

多項式の方は q2 が最初だけ -1 でその後はずっと 1 であることが利用できます。

  def fib_polynomial(n: Int): BigInt = {
//再帰は1じゃなく0で止めるのが基本だが、少しでも乗算を減らしたい
    def aux(n: Int, p0: BigInt, p1: BigInt, q1: BigInt): BigInt = {
      if(n == 1) p1 + p0 * q1
      else if(n % 2 == 0) aux(n/2, p0, p0 + p1 * q1, q1 * q1 - 2)
      else                aux(n/2, p1 + p0 * q1, p1, q1 * q1 - 2)
    }
    
    if(n <= 1) n
    else if(n == 2) 1
//    else if(n % 2 == 0) aux(n/2, 0,  1, 3)
    else if(n % 2 == 0) aux(n/2-1, 1,  0, 3)
    else                aux(n/2, 1, -1, 3)
  }

@ryuhei-mori
Copy link
Contributor Author

行列の考え方よりも多項式の考え方のアルゴリズムの方が速いことが多いです(もちろん一般のk+1項間線形漸化式の計算は大きな k については多項式の考え方のアルゴリズムの方が速いです)。
また、ちゃんと確認していませんが、この方法も速そうです。
https://gmplib.org/manual/Fibonacci-Numbers-Algorithm.html

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant