Nが素因数分解可能ということのゼロ知識証明について

長いよ

発端

【問題】 Aさんはある巨大な素数pの素因数分解をすることに成功し、素因数分解に成功したことをBさんに伝えようとしています。 しかし、Bさんが先にpの素因数を世間に発表することは避けたいため、素因数を伝えることなく、pの素因数を知っていることを示したいです。 どうすればいいでしょうか。

ということなので、「巨大な素数p」は「素数の積N=pq」ということにしておきます。

問題の定式化

簡単のため言語を \{N \in \mathbb{Z} \mid N = pq \land (p,q) \in\mathbb{P} \} とNが2つの素数p,qの積からなる言語として、プロトコルに従う検証者だけを考えるゼロ知識証明を考えます。（専門用語でいうところのHVZK PoK, honest-verifier zero-knowledge proof of knowledge)

要件は以下の3つです

  1. 正当性: プロトコルへの共通入力が言語に入っているとき、正直な証明者は正直な検証者を納得させられる
  2. 正直な検証者に対するゼロ知識性 (HVZK): 検証者が正直な場合のプロトコルのやりとりをシミュレーションできる
  3. 知識証明 (PoK): 証明者へブラックボックスアクセスできる期待多項式時間で動き知識抽出機械 E が存在し, そのEがxの証拠wを出力できる確率は、証明者が正直な検証者を納得させる確率から無視できる分だけしか離れていないこと。
    • 具体的に書くと, \Pr({E^{\tilde{P}}(x)} \in W(x)) = \Pr(\langle\tilde{P} \leftrightarrow V\rangle(x) \to 1) - \kappa(|x|) でEの期待動作時間が\mathrm{poly}(|x|) で, \kappa(|x|)は無視できる関数

案1. RSA暗号を使う

「合成数Nが半素数で、その素因数分解N=pqができた」と仮定してみます。 Aはm=φ(N)=(p-1)(q-1)と互いに素なeを探しBにeを伝える。 Bは（Aに秘密に）xを選び、xe≡y (mod N)を計算しyをAに伝える。 Aはed-mv = 1を満たす(d, v)を互除法で計算し、yd≡z(mod N)なるzをBに返送。 Bはx≡zを確認。

Aはp, q、そこから得られるmの情報を持っているが、いずれも外部に漏らしていない。 また、Bは x を鍵(e, N)で暗号化しAに送っているが、Aから x が返ってきたことにより、Aが復号鍵を持つことを確認できる。 復号鍵は N の素因数 p, q がないと得られないため、Aが素因数を持っていると証明できる。

案1の問題点

案1のプロトコルにより、証明者Aはe乗根、すなわちd乗を行うことが証明できます。ここまではOKです。

さて、問題は3の知識証明です。案1のプロトコルで検証者を納得させられる場合に、その証明者からp, qの素因数分解ができるでしょうか?

割と有名な事実として、(N, e, d) が揃えば素因数分解を行うことができることが知られています。*1 詳しくはおまけを見てください。 そのため、案1のプロトコルがdを知っていることの証明になっていれば、それは素因数分解ができることと等価なため、案1のプロトコルは素因数分解の証明になっています。

残念ながら一般的には仮定が証明されていません。例えば、素因数分解を行った別の人がRSA方式のようにeとdを計算し、証明者Aには eの値とd乗を行う謎の機械を証明者に渡すかもしれません。もしこの謎の機械が難読化 *2 されていてdが読みとれないのであれば、証明者にはdが分かりません。そのため、証明者には素因数分解ができない可能性があります。RSA問題（=d乗すること）を解けることと素因数分解の等価性はここ40年ほど未解決問題になっています。

案2. Rabin暗号を使う

RSA暗号の暗号文を復号するという操作は素因数分解と等価ではありません。そこで自然な考えとして、復号する操作が素因数分解と等価な暗号を使えばいいのではないかというアイデアがあります。ということで、そのようなことが知られている平方剰余を使ったRabin暗号を使ってみましょう。

  1. (Bが暗号文を生成する) Bは秘密にxを選び y≡x2(mod N) を計算し、Aにyを伝える
  2. (Aが復号する) Aはyの平方zを計算し、Bにzを伝える
  3. (Bが確認する) Bはz2≡y (mod N) ならば受理、そうでなければ拒否とする

案2の問題点

残念ながら案2のプロトコルでは素因数分解の知識が漏れてしまいます。 Bがプロトコルに従っていても、Bは確率的ですが素因数分解できてしまいます。 中国式剰余定理で, Bが秘密に選んだx in Z/NZを, (x mod p, x mod q) in Z/pZ x Z/qZと同一視します。 平方剰余yの平方根は (± x mod p, ± x mod q) の4通りあるんですが、Aから見るとどの組み合わせを自乗したのか識別できません。そのため、Aがどのような方法で4通りを決めたとしてもBが選んだものと一致する確率は1/4です。ということで, 非自明な組み合わせ (x mod p, -x mod q) または (-x mod p, x mod q) を証明者Aが選んでzとして送ってくると, Bは gcd(x-z,N)で素因数分解ができてしまいます。

案2の改良に必要な平方根を知っていることを示すゼロ知識証明

ということでいったん話を変えて、平方根を知っていることを示すゼロ知識証明を考えてみます。

証明者Pは検証者Vと対話して、入力(t, N)について, その平方根 w を知っていることを示したいとします。

プロトコル1:

  1. Pはrをランダムに選んで a = r2 mod Nを計算し, Vにaを送る
  2. Vはcを0/1からランダムに選び, cを送る
  3. Pはz = r wc mod Nを計算し, Vにzを送る

  4. Vは z2 ≡ a tc (mod N) ならば受理、そうでなければ拒否とする

  5. 正当性は自明なので略

  6. HVゼロ知識性: 先にcを選んでzを計算してからaを計算すると辻褄が合います。しかもこれは正しいP, Vのやりとりと見分けがつきません。
  7. 知識証明: \tilde{P}を使って, aを得ます。c=0の場合のz_0を得た後, \tilde{P}を巻き戻します（または同じ乱数テープを使う）。同じaに対してc=1の場合のz_1を得ます。このとき, {z_0^{2}}\equiv a \pmod{N}かつ{z_1^{2}}\equiv a t \pmod{N}なので, {(z_1/z_0)^{2}} \equiv t \pmod{N} となりz_1/z_2がtの平方根になります。（確率の問題がありますが、\tilde{P}の乱数テープを取り換えていろんなaで試すことで確率を増幅でき、κを減らすことができます。詳細略）

これを繰り返し（または並列に）行うとtの平方根wを知っていることを証明できていると考えてください（証明は略）

案2の改良

のリンク先に書いてあるプロトコルです。

プロトコル1を使って案2の問題を解決したものが以下です。

  1. (Bが暗号文を生成し、原像を知っていることを示す):
    1. Bは秘密にxを選び y≡x2(mod N) を計算し、Aにyを伝える
    2. Bが証明者, Aが検証者になって, プロトコル1を使い, Bはyの平方根xを知っていることを証明する
  2. (Aが復号する)
    1. Aはyの平方zを計算する
    2. Aが証明者, Bが検証者になって, プロトコル1を使い, Aはyの平方根zを知っていることを証明する
  3. (Bの出力) 2-2の結果をBは出力する

リンク先の小山（情報処理 vol.32, No.6 http://id.nii.ac.jp/1001/00004658/) に書かれてます。このプロトコル自体はTompa and Woll (FOCS 1987, https://ieeexplore.ieee.org/document/4568301) が初出だと思います。

（補足）検証者が正直であるという条件を使う場合は、1.で単にyを送ってプロトコル1を省いても大丈夫な気がします。

案3

これ以外の方法もあります。

link.springer.com

Poupard and Stern (PKC 2000) では、プロトコル1のような感じでゼロ知識証明を行います。

z_1,\dots,z_Kは1以上N未満でNと互素な整数からランダムに選ばれているとします。またパラメータA, B, K, 繰り返し回数のLなどは気を付けて選ぶ必要があります。

  1. A-1:
    1. rを1以上A未満の整数からランダムに選ぶ.
    2. x_i = {z_i^{r}} \bmod{n}を計算し, Bにx_1,\dots,x_Kを送る
  2. B-1:
    1. eを1以上B未満の整数からランダムに選び, Aに送る
  3. A-2:
    1. y = r + (n-φ(n)) × eを計算し, Bにyを送る
  4. B-2:
    1. 0 <= y < Aをチェック
    2. {z_i^{y-ne}} \equiv x_i \pmod{N}をチェック

とします。

ざっくりいうと、これでカーマイケル数 λ(N) の倍数を知っていることが証明できています。

  1. 正当性: B-2-1のところで正当性が統計的になっている点に注意してください。A, B, Kなどの設定が
  2. HVゼロ知識性: eを先に選んで, yを決めてから, x_iを計算します。
  3. 知識証明 (PoK): 平方根を知っていることを証明するプロトコルの際と同様に、A-1の後、2通りのe, e'でy,y'を得て, λ=|y-y'-n(e-e')|を計算します。{z_i^{y-ne}} \equiv {z_i^{y'-ne'}} \pmod{N}を使って, 色々頑張るとλがカーマイケル数 λ(N) の倍数になっていることが示せます。おまけに書いてある Millerの技法を使ったφ(N)の倍数から素因数分解を行う方法と同様に、λ(N)の倍数から素因数分解ができます。

おまけ (RSA原論文のAppendix Cの方法)

  1. RSAの鍵生成の結果から、 e d \equiv 1 \pmod{\phi(N)}. よって, m = e d - 1 は \phi(N) の倍数
  2. 一方, Millerの論文には, \phi(N) の倍数 mを使って, 確率的に N を素因数分解する方法が書いてある。
    1. ( a \in {1,2,\dots,N-1} をランダムに選ぶ)
    2. a \mid N をチェック
    3. 割り切れないなら, \mathrm{gcd} \big( (a^{m/{2^{k}}} \bmod{N})-1, N\big) を 1 \leq k \leq \lg(N-1)について調べる

ということで1. の拡張リーマン予想 (ERH) を仮定してaを小さいところから選んでいるものを, Miller-Rabinのようにランダムにaを選ぶようにすれば確率的帰着になります。

拡張リーマン予想を仮定せずに決定的に帰着を作りたい場合は、格子を使うA. Mayの手法を使ってください。

*1: RSAの原論文のAppendix C. または, May (CRYPTO 2000) https://link.springer.com/chapter/10.1007/978-3-540-28628-8_13 を参照

*2:難読化自体も暗号学的には色々論点がありますが、仮にできるとしておきます

RTACTF Writeup

speedrun.seccon.jp

🏳 Sexy RSA (目標: 300sec)

844.38 sec

#chall.py
from Crypto.Util.number import getPrime, isPrime
import os

def getSexyPrime(n=512):
    # Sexy prime: https://en.wikipedia.org/wiki/Sexy_prime
    while True:
        p = getPrime(n)
        if isPrime(p+6):
            return p, p+6

if __name__ == '__main__':
    # Plaintext (FLAG)
    m = int.from_bytes(os.getenv("FLAG", "FAKE{sample_flag}").encode(), 'big')

    # Generate key
    p, q = getSexyPrime()
    n = p * q
    e = 65537

    # Encryption
    c = pow(m, e, n)

    # Information disclosure
    print(f"n = 0x{n:x}")
    print(f"e = 0x{e:x}")
    print(f"c = 0x{c:x}")
  • ×tar.gzの展開に詰まる.
  • ×問題を見てから二部探索をする.
  • xintからstrへの変換を忘れる. (Crypto.Util.numberのlong_to_bytesでええやん)
# solve.py
from Crypto.Util.number import getPrime, isPrime
import os
from binascii import unhexlify

def find_approx_sqrt(c):
    ng = -1
    ok = c
    while (abs(ok - ng) > 1):
        mid = (ok + ng)//2
        if mid**2 >= c:
            ok = mid
        else:
            ng = mid
    return ok

if __name__ == '__main__':
    n = 0xe72988e811f04091c3291ac28f1e8332193187f3dc5af01579c36badb06671aa9a9543aa07eba8cdab36d787f1ff98a06db995c43cd5c63581ce050e0b9ba856634dabfaf8c7f271fbd026edd6ea1257b16013a526e0581a688cc6a335e7ee4c1b0633f0532d3d0824824195b6b249c70cf0e458609efc01a6575f084e6de53b
    e = 0x10001
    c = 0x6fadd5d7095bd6f45de69bb4e76080e0ea5f8c5a159de10663133e585b71ae580b99b3e0a8e047a9c51c8091a6b33b01c9ab95668794c3acfb084e939a04cb151757c3b2522da99e03f83e205c7c701066d69b120ca17fcf59061c078d9099e5f4bf6dd6dab206418527035f2c1096861c2896327977ac88c2728faa7504d879
    print(n)

    p = find_approx_sqrt(n)
    while True:
        if p * (p-6) == n:
            q = p-6
            d = pow(e,-1,(p-1)*(q-1))
            m = pow(c,d,n)
            print(m)
            m = unhexlify(hex(int(m)).strip('0xL'))
            print(m)
            break
        else:
            p += 1

🏳 Proth RSA (目標: 900sec)

911.18 sec

# chall.py
from Crypto.Util.number import getRandomInteger, getPrime, isPrime
import os

def getProthPrime(n=512):
    # Proth prime: https://en.wikipedia.org/wiki/Proth_prime
    while True:
        k = getRandomInteger(n)
        p = (2*k + 1) * (1<<n) + 1
        if isPrime(p):
            return p, k

if __name__ == '__main__':
    # Plaintext (FLAG)
    m = int.from_bytes(os.getenv("FLAG", "FAKE{sample_flag}").encode(), 'big')

    # Generate key
    p, k1 = getProthPrime()
    q, k2 = getProthPrime()
    n = p * q
    e = 65537
    s = (k1 * k2) % n

    # Encryption
    c = pow(m, e, n)

    # Information disclosure
    print(f"n = 0x{n:x}")
    print(f"e = 0x{e:x}")
    print(f"s = 0x{s:x}")
    print(f"c = 0x{c:x}")

p=(2k_1+1) \cdot 2^{512} + 1, q=(2k_2+1) \cdot 2^{512} + 1という条件で, s=k_1 \cdot k_2 \bmod{n}が分かっている. ビット長を考えると, s=k_1 \cdot k_2そのもの. ビット長の関係でn-1 = 2^{1024} \cdot (2k_1+1)(2k_2+1) + 2^{512} \cdot (2k_1+2k_2+2) と nの512ビット～1024ビット当たりでk_1+k_2が大体わかる. k_1+k_2を絞り込んだら, それで(p-1)(q-1)を作ってdを求めればOK. グレブナー基底はなくてもいけるようになっていて親切.

from Crypto.Util.number import getRandomInteger, getPrime, isPrime
import os
from binascii import unhexlify

if __name__ == '__main__':
    n = 0xa19028b5c0e77e19fc167374358aa346776e6c20c27499505be59c83ea02014e97af631ba0ccbab881313818fd323c15c82dad8793220ba6679ec4b38787e04d0c1fff0880e04423ea288e443660c63a1607532e47dbaad421723d0546c208447f701cd7e9ee1bb43774d132abbb2e91bf50b67be40ed854dbe6c3071ca3ae3307ac03abd76f74e506594106a22795d4b7938611301248a9957e1a637538a9169cf38daf5d60ffc05ae32ea7e638e16d790ffeebfff655a645c99a513616d3ce00000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001
    e = 0x10001
    s = 0x28640a2d7039df867f059cdd0d62a8d19ddb9b08309d265416f96720fa808053a5ebd8c6e8332eae204c4e063f4c8f05720b6b61e4c882e999e7b12ce1e1f812c11cfed72a5c33cfb8f3d34f650e4c19579cf34745f2588aa2fd08a8746257cb789f23ca232346fcf72468a2b160934911902de3f90620aba5874a2d79a33699
    c = 0x4595c3c923bd191ba07456611f80e656a197ff528a031e2952adedda532b1fa2caef719c929132a3cdf06d0e55e6a00f7eb1f189a614b26759916ec42f83579a75ab5948186769a1a936b019466f918f29e32852675c464b7f0797c6fdc55efcd54fbe2083761b1df3dde0b9a9a35b96e3b216c54770b444b1f02525f0268c44483c6e84a781fe9111e6912130d69f462c519873043d44e4a3f1f938491feeb591b5831d0abe7399bc87244576decaf2925f287d3c2bb4061d560c919d820e364744f2322c7efd37d42563842bcf9b1d6b46218694dcd49758d311c6896e38cf2b55c7114d78cfdfaeba74720ecf30d9133034799b9735e26ec913cc9f26bb0a

    print(f'{((n-1) >> 512) % 2**512:#x}')
    k1_plus_k2 = (((n-1) >> 512) % 2**512)//2 - 1  

    for i in range(10):
        phin = 2**1024 * (4*s + 2 * k1_plus_k2 + i * 2**512 + 1)
        d = pow(e,-1,phin)
        m = pow(c,d,n)
        m = unhexlify(hex(int(m)).strip('0xL'))
        print(m)

🏳 Neighbor RSA (目標: 1800sec)

1836.31 sec

# chall.sage
import os

# Plaintext (FLAG)
plaintext = os.getenv("FLAG", "FAKE{sample_flag}").encode()
plai  = plaintext[:len(plaintext)//2]
ntext = plaintext[len(plaintext)//2:]
m1 = int.from_bytes(plai  + os.urandom(128), 'big')
m2 = int.from_bytes(ntext + os.urandom(128), 'big')

# Generate key
e = 65537
p = random_prime(1<<2048)
q = random_prime(1<<512)
r = random_prime(1<<512)
n1 = p * q
n2 = next_prime(p) * r
assert m1 < n1 and m2 < n2

# Encryption
c1 = pow(m1, e, n1)
c2 = pow(m2, e, n2)

# Information disclosure
print(f"e = {hex(e)}")
print(f"n1 = {hex(n1)}")
print(f"n2 = {hex(n2)}")
print(f"c1 = {hex(c1)}")
print(f"c2 = {hex(c2)}")

n_1=p \cdot q と n_2=(p+\delta) \cdot rでm1とm2が暗号化されている. 条件式を見ると, pが2048ビット素数, p+δはpの次の素数, qとrが512ビット. \frac{n_1}{n_2} = \frac{q}{r} \cdot \frac{p}{p+\delta} \approx \frac{q}{r} となるので ディオファントス近似 - Wikipedia して, 連分数展開してqとrを推測してm1,m2を表示していくとどこかで正解が出る.

(格子でも解けるらしい)

#solve.py
import os
from binascii import unhexlify
import ContinuedFractions

e = 0x10001
n1 = 0xa8ed020c3dd125d503bf124052d643ba1405f2c349244122140e79e7d2244304a1590762c61ac83900c2aced76007b2e3f320464fd51fcfad167ebdc87e69329230869e0a3e153b44ed3b04bfe94174bc8b5ee1a3fa8036b6b9e834666aa07229a431b477e589d94f9a4cfed25b195215b0c694b86e874413b8a00cb064809c8e3677632cde9b43b87a0b812c2024b0c821b5c10764fd4de2d18af55d897d94aeded80b71e36fd73014f75641a8c5b38b36faa020e7cf1327a707bb7d42503bcc28768ef184d66b9ba16efd019b68268885a2da302cd326e78b1d473bcf7cd62442ccd25dc85d23aeb5408922b6b00f13584bea394f1bca4cc431f3c29c5d98ec1683453cc0c526abe4aec08781c7a53f50f2047b4995b9bea6a7a9f6b5425b29be6e867764efaa050799f716af78273041372cfe4f3c88a62329f6f1feff99
n2 = 0x16ea1bde86ea11cdec196a9173258efca235da66f8e3d5437e39e1b2e2574dd3f93d65104ca0225d6119519ae9ea9c035e0f85f02212c0992d0705723fa8b97ed6cff860c4d8fb65f0214a0047feca64e662dcbf025fff47590305e90e5d070d39871880828f5e960ab2ef330129ed5752c3b4debe827376a632b06487740fff4b622a88de23649e3e6993cd332b0284b84eb8765d58527209cc202c89d479421131a2f64ae517ee1e62e6c0f329c306569e427113ec6a8b9d96d73e95580d3a33f6add681f9a9156f0681eb1804183dfa8cebbe921d2fb1d43b256f727d46c5859cc5229f7e555ad25397e5cd14620ebbaefa0a0a520bada3ef8b115481734242af6befbd9b069d4a03281094c0f4aca4e6fdcbe2558b104fc2b383e1c70f0e5a07d1a623f9fc2309ca1d09b69aa1869e280fbc50de2adbada7ea545743b12b
c1 = 0x690e49037fee7649033ffaaa71e4730d2d7143fab97beb22e2afdf6eca449cad3f95b60295f592e7e84833e08b3468d61a34c1d1123f4c683c79d68bbe27dd0af203fc50ef7ebe98b1bc1221918470f058a8fb7645eacb569931835bd7f80494dbb67fbaa592ec19d9b4930c787a2ce1267f8088229b5031e710d6cd5720756923ccb64444939a0f09a51c87488650d4d02551fd4ed7a2fd248825ec34c5df8b6077a6d0d75c5832f9140420c92d3d00cf51e3b0665f5a6d031cb369ddebbb5ce77f2176cd12bb0add5aeda6ae88c4ceade0c1fd0ff3960d3ee36a0c6455ae3027f33e660663d0e2298654e19e8c8a06b4de991fac3b4c1673825b3d9f8f5c675f920a7d137f85ba723bf741321904e0c3c601f5c18d02e1e5b7b118e62e91a7926a9b1eda3cc53e2a6cbc95553e1990ec3f6cceddf283410d6e6849a26f89b
c2 = 0x5c4c7dce82753a68dcdbcdce9af52c9b7af2f561c08b8e23b27c6145d4c3df29d498303bee1bd29829a2e0ae9faaf243b387c39d69daccba07dace7bb420115ffaa69f89a3ea4e1ef0e08eb19043e012a090b79e51d6ae8446ca76e88abe5adbdbe25a731d7ee9aa333a84447edafbc360b505ff293c751571c6bf29dee99fdc443b756f182eb588b4a03de3d35dc4f23736d7239cfbd0ca13fa7b234bc4064a2053ab0045f4833250c8c9de91798502b09d4312ee52f3dc5229dfcb73b42f7c3440932839e6e790bb0db1788fbd7c60365121bbe3858ecedd3d48261d081c380e7ddf6ca570c13cc89c0af2011b4978b22d5456d1122dabd7b2068ab30e301a674809732daede77a27ae13e1bc4779e15d51210f6c10be159907ec1a59bfaf8db6cf290a348f734fd88e3c2b7df6bda84665b810cfe55bc3645d8d118c9172
print(ceil(n1/n2))

# See https://elliptic-shiho.hatenablog.com/entry/2015/12/18/205804
# Use https://github.com/pablocelayes/rsa-wiener-attack
def Wiener(n1,n2):
    frac = ContinuedFractions.rational_to_contfrac(n1,n2)
    convergents = ContinuedFractions.convergents_from_contfrac(frac)
    
    for (q,r) in convergents:        
        print(q,r)
        #check if d is actually the key
        if q != 0 and n1 % q == 0 and n2 % r == 0:

            d1 = pow(e,-1,(n1//q-1)*(q-1))
            m1 = pow(c1,d1,n1)
            m1 = unhexlify(hex(int(m1)).strip('0xL'))
            print(m1)

            d2 = pow(e,-1,(n2//r-1)*(r-1))
            m2 = pow(c2,d2,n2)
            m2 = unhexlify(hex(int(m2)).strip('0xL'))
            print(m2)

            break

Wiener(n1,n2)

Leaky RSA (目標: 1800sec)

Your time is 86428.35 seconds ... 解けず.

# chall.py
from Crypto.Util.number import getPrime, isPrime, inverse
import os

if __name__ == '__main__':
    # Plaintext (FLAG)
    m = int.from_bytes(os.getenv("FLAG", "FAKE{sample_flag}").encode(), 'big')

    # Generate key
    p = getPrime(600)
    q = getPrime(500)
    n = p*p*q*q
    e = 65537
    s = ((p**3 - 20211219*q) * inverse(p*p+q*q,n)) % n # tekitou (ha?)

    # Encryption
    c = pow(m, e, n)

    # Information disclosure
    print(f"n = 0x{n:x}")
    print(f"e = 0x{e:x}")
    print(f"c = 0x{c:x}")
    print(f"s = 0x{s:x}")

s \equiv (p^3 - 20211219 q)/(p^2 + q^2) \pmod{n} なので, p^3 - s \cdot p^2 - s \cdot q^2 - 20211219 q \equiv 0 \pmod{n} . nが大体2200ビットで, p, qが600ビット, 500ビットだからHowgrave-Grahamの補題を使えばいけそう思ってCoppersmith攻撃を試みて失敗. \sqrt{n} = pqが分かるのでp,p+q,qを解とした p^3-20211219 q - s (p+q)^2 + 2spq \equiv 0 \pmod{n}も試して失敗. 他の人のwriteupを見て納得.