以下の内容はhttps://inamori.hateblo.jp/entry/2024/12/05/195450より取得しました。


MojoでProject Euler 76

https://projecteuler.net/problem=76

これは分割数から1引いたものです。
n個をm以下で分割した数を、 D_{n,m}とすると、
 \displaystyle D_{n,m} = \sum_{k=1}^m{D_{n-k,k}}
となります。これで D_{N,N}を求めると、計算量は O(N^3)ですが、これでも十分間に合います。
しかし、母関数を使うともっと速くなります。
 \displaystyle G(x) = \prod_{n=1}^{\infty}{(1+x^n+x^{2n}+\cdots)} = \prod_n{\frac{1}{1-x^n}}
なので、これを x^Nまで計算すればよいです。計算量は O(N^2)です。
さらに五角数定理を使うと速くなります。これは、
 \displaystyle \prod_{n}{(1-x^n)} = \sum_{m \in \mathbb{Z}}{(-1)^mx^{\frac{1}{2}m(3m-1)}}
です。最後の xの指数が五角数です( mが負も取ることに注意が必要です)。なので、元の母関数の分母は最初から計算されていて、1を割ればいいだけです。項数が O(\sqrt{N})しかないので、計算量は O(N^{1.5})となります。
データの持ち方に工夫が必要で、指数と係数のペアのVecを用意すればよいです。このときタプルはKeyElementでないので新たな構造体を定義する必要があるようです。

import sys


#################### Pair ####################

@value
struct Pair(KeyElement, Stringable):
    var x: Int
    var y: Int
    
    fn __eq__(self, other: Pair) -> Bool:
        return self.x == other.x and self.y == other.y
    
    fn __ne__(self, other: Pair) -> Bool:
        return self.x != other.x or self.y != other.y
    
    fn __hash__(self) -> Int:
        return self.x * 10000 + self.y
    
    fn __str__(self) -> String:
        return '(' + str(self.x) + ', ' + str(self.y) + ')'


#################### process ####################

fn create_divisor(N: Int) -> List[Pair]:
    var a = List[Pair](Pair(0, 1))
    for n in range(1, N+1):
        var m1 = n * (3 * n - 1) // 2
        if m1 > N:
            break
        a.append(Pair(m1, (-1)**n))
        var m2 = n * (3 * n + 1) // 2
        if m2 > N:
            break
        a.append(Pair(m2, (-1)**n))
    return a

fn inverse(a: List[Pair], N: Int) -> List[Int]:
    var b = List[Int](1)
    for _ in range(N):
        b.append(0)
    var c = List[Int]()
    for i in range(N+1):
        var be = b[i]
        c.append(be)
        for e in a:
            var j = e[].x + i
            var d = e[].y
            if j > N:
                break
            b[j] -= be * d
    return c

fn f(N: Int) -> Int:
    var a = create_divisor(N)
    var c = inverse(a, N)
    return c[N] - 1

fn main() raises:
    var args = sys.argv()
    var N = atol(args[1])
    print(f(N))



以上の内容はhttps://inamori.hateblo.jp/entry/2024/12/05/195450より取得しました。
このページはhttp://font.textar.tv/のウェブフォントを使用してます

不具合報告/要望等はこちらへお願いします。
モバイルやる夫Viewer Ver0.14