Prime Factorization with SQL (Project Euler #3)

main.sql

SELECT prime_factorization(600851475143);

-- prime_factorization 
---------------------
 {71,839,1471,6857}
(1 row)

prime_factorization.sql

CREATE OR REPLACE FUNCTION prime_factorization (bigint) RETURNS bigint[] AS $$
WITH RECURSIVE prime(ary) AS ( -- 素数列を作成
    SELECT
        array_agg(num) AS ary
    FROM
        sieve_of_eratosthenes(trunc(sqrt($1))::bigint) AS t(num) 
)
,factorization AS( -- 素因数分解の実行
    SELECT
        1::bigint AS prime -- 現在参照している素数
        ,$1 AS quotient -- 素因数分解を行う対象
        ,1 AS idx -- 現在参照しているindex
        ,FALSE AS flag -- 参照している素数が因数かどうかのフラグ
    UNION ALL
    SELECT -- primeで作成した素数集合に対し一つずつ因数かどうかの判定を行う
        ary[idx]
        ,CASE
            WHEN (quotient % ary[idx]) = 0 THEN quotient / ary[idx]
            ELSE quotient
        END
        ,CASE -- 因数の場合は同じ数でもう一度判定を行う
            WHEN (quotient % ary[idx]) = 0 THEN idx
            ELSE idx + 1
        END
        ,CASE
            WHEN (quotient % ary[idx]) = 0 THEN TRUE
            ELSE FALSE
        END
    FROM
        factorization, prime
    WHERE
        idx <= array_length(ary, 1)
)
,agg AS ( -- factorizationから素因数の集合を取り出す（操作しやすくする為に配列に集約）
SELECT
    array_agg(prime) AS result
    ,min(quotient) AS min
FROM
    factorization
WHERE
    flag
)
SELECT -- 端処理など
    CASE
        WHEN result[1] IS NULL THEN ARRAY[$1]
        WHEN min = 1 THEN result
        ELSE array_append(result, min)
    END
FROM
    agg
;
$$ LANGUAGE SQL IMMUTABLE STRICT;

sieve_of_eratosthenes.sql

CREATE OR REPLACE FUNCTION sieve_of_eratosthenes(bigint) RETURNS SETOF bigint AS $$
    WITH RECURSIVE
    --数列を用意
    t1(n) AS (
        SELECT generate_series(2, $1)
    ),
    t2 (n, i) AS (
        --初期化
        SELECT
            n
            ,1::bigint
        FROM
            t1
        --再起:再起集合内の最小値で割り切れる値を除いた集合を新たに再起集合とする
        UNION ALL(
            --再起集合をコピー（再起集合は再起部で一度しか呼び出せないため）
            WITH s1 (n) AS(
                SELECT
                    n
                FROM
                    t2
            )
            --再起集合内の最小値を取得
            ,s2 (k) AS(
                SELECT
                    min(n)
                FROM
                    s1
            )
            --最小値で割り切れる値を除外して集合を作成
            SELECT
                n
                ,k
            FROM
                s1,s2
            WHERE
                (n%k) <> 0
                AND k*k < $1 --√nで再起を終了
        )
    )
    --除数に用いた数は素数(1以外)
    SELECT DISTINCT
        i AS n
    FROM
        t2
    WHERE
        i <> 1
    UNION
    --最後に残った集合も素数
    SELECT
        n
    FROM
        t2
    WHERE
        i = (SELECT max(i) FROM t2)
    ORDER BY
        n
$$ LANGUAGE SQL IMMUTABLE STRICT;