DLP攻击算法

DLP问题。

给定一个循环群G，G的阶是n，生成元是g。给定群中的一个元素h，DLP问题就是找出 x 满足下面的条件：$h=g^x$

Okamoto–Uchiyama 密码系统

由 Tatsuaki Okamoto 和 Shigenori Uchiyama 于 1998 年提出的一种公钥 密码体制 。该系统在 整数模 n 的乘法群中 工作，$(Z/nZ{)}^{\ast }$，其中 n 的形式为 $p^{2}q$并且 p 和 q 是大 素数 。

参考 WIKI： https://en.wikipedia.org/wiki/Okamoto%E2%80%93Uchiyama_cryptosystem

密钥生成

公钥/私钥对生成如下：

1、生成两个大素数 p 和 q

2、计算 $n=p^{2}q$

3、选择一个随机整数 $g\in (2,...,n-1)$ 这样 $g^{p-1}\ne 1mod{p}^2$

4、计算 $h=g^{n}modn$

公钥为 (n,g,h) ，私钥为 (p,q)

加密

一个消息 m < p 可以用公钥 (n,g,h) 加密，如下：

1、选择一个随机整数 $r\in (1,...,n-1)$

2、计算 $c=g^m{h}^{r}modn$

c 为加密后的密文。

解密

1、计算 $a=\frac{(c^{p-1}mod{p}^2)-1}{p}$

2、计算 $b=\frac{(g^{p-1}mod{p}^2)-1}{p}$， a 和 b 将是整数。

3、使用扩展欧几里得算法，计算 b 的逆元模 p

$$b^{\prime }=b^{-1}modp$$

4、计算 $m=a{b}^{\prime }modp$

得到明文 m。

例题：[2022cvctf]Weird dlog

from Crypto.Util.number import *
from secret import flag

NBITS = 1337

def keygen():
    p = getPrime(NBITS)
    q = p + 2
    while not isPrime(q):
        q += 2

    n = p**2 * q
    g = getRandomRange(2, n-1)
    assert pow(g, p-1, p**2) != 1
    return (g, n), p

(g, n), pk = keygen()
print(f"g = {g}")
print(f"n = {n}")

m = pow(g, bytes_to_long(flag), n)
print(f"m = {m}")

'''
g = 4217104162231108571881486363264502722577265670533700464292447122725884060865278499273749984084307890296655750169737247120115556105852793495207776085190125693574834744150776071626953537124287733539294768576947214664831668249537361880607900415980188224790377848428124695064324868240382228326860781143717838988468681317873292597463247890454032395659378742623078198171614526987958808613568346779857292295096270767472963769637419908101399359189794634298434235861743414985078720538901239212554829307311050121532457033249775247040962157091404814597818322567791804549931512078619168741053536240553838132062889112375399237994189289540471931223835109201969517169907064471393467949637030171468739827679615037096191102556235525038398343624846458755117043899170432849676184565935742212246568515988089821245576219250167380483425498097448081250808144546318744227429556017907249040862683999049166654695785359522913386597444234975046933660020802696704230400762773828356573597488744036818115263470798266879540313464252075039418721353291327500707724679404373330401143441889215673146351039551079579127929035980199349135110483207097742313251638620757710743549095431313120163844541459241899262156308740654917968060844284412082299
n = 10941888529432710302272044808172366177817083267353728387778251966266072287607755644721338043351870432706497693313492403987246128369454049222700884094312120262698229189208223218732275272981771624419914106314134685967783268356547329406040929859316975902966386276320517786973713050654580905181729305020144204934138012268416005803949005810888621349553556279807338189673928133725107444318753299258813428869755492871371365509253274275507785302703028229706104025231495516828413453372697736689315270600144033367988185883751955533145131820299821036099425033741978037404521469840343451728996005273312522454265432166678945758946508086710070492821787815923425038513533610548985733666276895713644320042830406173577168099532687994444950160159745610177551080984478286003952657902953582221120083238221655231826293209922243103643973158229947410799996443039134812980230923730059957532891088920291921544489491437892953022797147371056895427933338596342279801238281297736043937512612367418666450917035364442328787819461040287010570790023067874699031894369374666073884261229547448186848837461040401572755725068306333093138881471721393682218994576974011413807781411295540237225784367274045480384276599164320452991405168935726178669
m = 4576653037171661542042510858542458890428902270240861216961609304805830644888835861552943484827465024255946333173438385318547946480432546672274379127945737050221711664915718037389483818788255371177342813451318183796063859501483853616124820170927827054719451326957465962519367227753007179628370407849940782754428678471956932743015842138996456351682565775528904950367274203264025982858945560085867807654673726775851728189907590054081432299552285929407918297497972515287529506153936873120705301250301755547851988767682706931104344998515759782906834209711055747621306735557847407025512085632112612715677993856157667539100717815505576125450494688023610882544096421409195218600313621491935755995519193153556491442386466162185336840903547990417011020008723119132837593921370224096091218428827263318876330392141114596260431968349229888104149062981625468823935145842920387345948673405002062554960084883553497321014737635984993337888457310522533554029754162393609118195463217899952617045335297738291366080224287861881668432617054485772999546446076931624127877484962146375663168327918511556849268478695380248447621861402820926977211453653005530269636900365390125831754366239971553991824331582315041274983141065725144997
'''

因为 $n=p^{2}q$，p,q很接近，所以直接对 n 开三次方根，然后在取值附近搜索符合条件的素数，找到真正的 p,q。

然后用 Okamoto–Uchiyama 密码系统进行解密即可

EXP：

from Crypto.Util.number import long_to_bytes
from gmpy2 import *


n = 10941888529432710302272044808172366177817083267353728387778251966266072287607755644721338043351870432706497693313492403987246128369454049222700884094312120262698229189208223218732275272981771624419914106314134685967783268356547329406040929859316975902966386276320517786973713050654580905181729305020144204934138012268416005803949005810888621349553556279807338189673928133725107444318753299258813428869755492871371365509253274275507785302703028229706104025231495516828413453372697736689315270600144033367988185883751955533145131820299821036099425033741978037404521469840343451728996005273312522454265432166678945758946508086710070492821787815923425038513533610548985733666276895713644320042830406173577168099532687994444950160159745610177551080984478286003952657902953582221120083238221655231826293209922243103643973158229947410799996443039134812980230923730059957532891088920291921544489491437892953022797147371056895427933338596342279801238281297736043937512612367418666450917035364442328787819461040287010570790023067874699031894369374666073884261229547448186848837461040401572755725068306333093138881471721393682218994576974011413807781411295540237225784367274045480384276599164320452991405168935726178669
guess = iroot(n,3)[0]
for i in range(-1000,1000):
    p = i + guess
    if is_prime(p) and n % (p**2) == 0:
        q = n//p**2
        # print('p = ',p)
        # print('q = ',q)

p =  2220056847922888806692961201204181494180468404265144656979415727846610609650612145504228170996668036830178564363339729040733053444504290404813978909079947107025031382736413106168755521395169078323433549912352350586785691382035545583930303726697079796565769117602611817953388750640643945286203371871152311393277382213336946213995315438785489971755082168989679975255419041436378409691325787848675432435671
q =  2220056847922888806692961201204181494180468404265144656979415727846610609650612145504228170996668036830178564363339729040733053444504290404813978909079947107025031382736413106168755521395169078323433549912352350586785691382035545583930303726697079796565769117602611817953388750640643945286203371871152311393277382213336946213995315438785489971755082168989679975255419041436378409691325787848675432436509
assert p **2 * q == n

c = 4576653037171661542042510858542458890428902270240861216961609304805830644888835861552943484827465024255946333173438385318547946480432546672274379127945737050221711664915718037389483818788255371177342813451318183796063859501483853616124820170927827054719451326957465962519367227753007179628370407849940782754428678471956932743015842138996456351682565775528904950367274203264025982858945560085867807654673726775851728189907590054081432299552285929407918297497972515287529506153936873120705301250301755547851988767682706931104344998515759782906834209711055747621306735557847407025512085632112612715677993856157667539100717815505576125450494688023610882544096421409195218600313621491935755995519193153556491442386466162185336840903547990417011020008723119132837593921370224096091218428827263318876330392141114596260431968349229888104149062981625468823935145842920387345948673405002062554960084883553497321014737635984993337888457310522533554029754162393609118195463217899952617045335297738291366080224287861881668432617054485772999546446076931624127877484962146375663168327918511556849268478695380248447621861402820926977211453653005530269636900365390125831754366239971553991824331582315041274983141065725144997
g = 4217104162231108571881486363264502722577265670533700464292447122725884060865278499273749984084307890296655750169737247120115556105852793495207776085190125693574834744150776071626953537124287733539294768576947214664831668249537361880607900415980188224790377848428124695064324868240382228326860781143717838988468681317873292597463247890454032395659378742623078198171614526987958808613568346779857292295096270767472963769637419908101399359189794634298434235861743414985078720538901239212554829307311050121532457033249775247040962157091404814597818322567791804549931512078619168741053536240553838132062889112375399237994189289540471931223835109201969517169907064471393467949637030171468739827679615037096191102556235525038398343624846458755117043899170432849676184565935742212246568515988089821245576219250167380483425498097448081250808144546318744227429556017907249040862683999049166654695785359522913386597444234975046933660020802696704230400762773828356573597488744036818115263470798266879540313464252075039418721353291327500707724679404373330401143441889215673146351039551079579127929035980199349135110483207097742313251638620757710743549095431313120163844541459241899262156308740654917968060844284412082299

a = (pow(c,p-1,p**2)-1) //p
b = (pow(g,p-1,p**2)-1) //p
m = a*(invert(b,p)%p)%p
print(long_to_bytes(m))

亏格为0曲线上的DLP

例题：[2022 SEETF] The True ECC

题目文件地址

关键代码：

"""Represents the curve x^2 + ay^2 = 1 mod p"""
a, p = 376014, (1 << 521) - 1
curve = Ellipse(a, p)

gx = 0x1bcfc82fca1e29598bd932fc4b8c573265e055795ba7d68ca3985a78bb57237b9ca042ab545a66b352655a10b4f60785ba308b060d9b7df2a651ca94eeb63b86fdb
gy = 0xca00d73e3d1570e6c63b506520c4fcc0151130a7f655b0d15ae3227423f304e1f7ffa73198f306d67a24c142b23f72adac5f166da5df68b669bbfda9fb4ef15f8e
G = Point(curve, gx, gy)

if __name__ == "__main__":

    from flag import flag

    alice_pub, alice_priv = gen_secret(G)
    blake_pub, blake_priv = gen_secret(G)

    shared = alice_pub * blake_priv
    ct = encrypt(shared, flag)

    assert shared == blake_pub * alice_priv
    assert decrypt(shared, ct) == flag
   
    print("alice_pub =", alice_pub)
    print("blake_pub =", blake_pub)
    print("ct =", ct)
    
# alice_pub = (2138196312148079184382240325330500803425686967483863166176111074666553873369606997586883533252879522314508512610120185663459330505669976143538280185135503158, 1350098408534349199229781714086607605984656432458083815306756044307591092126215092360795039519565477039721903019874871683998662788499599535946383133093677646)
# blake_pub = (4568773897927114993549462717422746966489956811871132275386853917467440322628373530720948282919382453518184702625364310911519327021756484938266990998628406420, 3649097172525610846241260554130988388479998230434661435854337888813320693155483292808012277418005770847521067027991154263171652473498536483631820529350606213)
# ct = b'q\xfa\xf2\xe5\xe3\xba.H\xa5\x07az\xc0;\xc4%\xdf\xfe\xa0MI>o8\x96M\xb0\xfe]\xb2\xfdi\x8e\x9e\xea\x9f\xca\x98\xf9\x95\xe6&\x1fB\xd5\x0b\xf2\xeb\xac\x18\x82\xdcu\xd5\xd5\x8e<\xb3\xe4\x85e\xddX\xca0;\xe2G\xef7\\uM\x8d0A\xde+\x9fu'

题目分析：

曲线方程为：$x^2+a{y}^2=1(modp)$ ，已知量有 $a,p$ ，曲线基点 $G$ ，曲线上的DH算法得到的 $A=n_a\ast G,B=n_b\ast G$ 的值，以及 AES加密后得到的密文 $C$ 。

AES的密钥为 DH 加密的共享密钥，即 $S=B\ast n_a$

因为 p = (1 << 521) - 1 ，所以 p+1 是光滑数。

通过论文：关于循环群、有限域和离散对数问题的说明

可知上述曲线方程的形式符合论文中的亏格为0的一类特殊曲线：$x^2-D{y}^2=1$

检验 a 是否为模 p 的二次剩余即可判断应该使用的算法：

a, p = -376014, (1 << 521) - 1
assert pow(a, (p - 1) // 2, p) == (-1) % p
print("ok")

可知a 为模 p 的二次非剩余：$a^{(p-1)/2}\equiv -1modp$

所以当 $\lambda (D)=-1$ 时，即这里的 $\lambda (-a)=-1$ ，则有：

1.设 $f(W)=W^2-D\in F_p[W]$ ，$f(W)$ 为有限域 $F_p$ 的一个扩展域。

2.映射公式为：$\varphi :(x,y)⟼x+yW$

理论知识参考文章：情况1：当 $ lambda(D) = -1$时

实现EXP：

from ecc import Point, curve, decrypt

gx = 0x1bcfc82fca1e29598bd932fc4b8c573265e055795ba7d68ca3985a78bb57237b9ca042ab545a66b352655a10b4f60785ba308b060d9b7df2a651ca94eeb63b86fdb
gy = 0xca00d73e3d1570e6c63b506520c4fcc0151130a7f655b0d15ae3227423f304e1f7ffa73198f306d67a24c142b23f72adac5f166da5df68b669bbfda9fb4ef15f8e
G = (gx, gy)
A= (2138196312148079184382240325330500803425686967483863166176111074666553873369606997586883533252879522314508512610120185663459330505669976143538280185135503158, 1350098408534349199229781714086607605984656432458083815306756044307591092126215092360795039519565477039721903019874871683998662788499599535946383133093677646)
B = (4568773897927114993549462717422746966489956811871132275386853917467440322628373530720948282919382453518184702625364310911519327021756484938266990998628406420, 3649097172525610846241260554130988388479998230434661435854337888813320693155483292808012277418005770847521067027991154263171652473498536483631820529350606213)
c = b'q\xfa\xf2\xe5\xe3\xba.H\xa5\x07az\xc0;\xc4%\xdf\xfe\xa0MI>o8\x96M\xb0\xfe]\xb2\xfdi\x8e\x9e\xea\x9f\xca\x98\xf9\x95\xe6&\x1fB\xd5\x0b\xf2\xeb\xac\x18\x82\xdcu\xd5\xd5\x8e<\xb3\xe4\x85e\xddX\xca0;\xe2G\xef7\\uM\x8d0A\xde+\x9fu'

a, p = 376014, (1 << 521) - 1
D = -a
assert pow(D, (p - 1) // 2, p) == (-1) % p
print("a不是模p的二次剩余")

R.<w> = GF(p^2, modulus=x^2 - D)

g = G[0] + w * G[1]
A = A[0] + w * A[1]

# A的私钥
n_a = discrete_log(A, g)
print(n_a)

B = Point(curve, int(B[0]), int(B[1]))

# AB的共享密钥
shared = B * int(n_a)
print(decrypt(shared, c))

SEE

类题：

• CryptoHack - Ellipse Curve Cryptography
• m0leCon CTF 2020 Teaser — King Exchange

Baby-step giant-step(小步大步算法)

在群论中，大步小步算法是 Daniel Shanks 发明的一种中间相遇算法，**用于计算离散对数或者有限阿贝尔群的阶。**其中离散对数问题在公钥加密领域有着非常重要的地位。 许多常用的加密系统都基于离散对数极难计算这一假设——计算越困难，这些系统提供的数据传输就越安全。增加离散对数计算难度的一种方法，是把密码系统建立在更大的群上。

这是一种空间换时间的算法，实质上是求解离散对数的朴素算法（枚举并试乘）的一个相当简单的改进。

算法：

给出一个 $n$ 阶循环群 $G$ 、该群的一个生成元 $\alpha$ 和一个元素 $\beta$ 。试找到一个整数 $x$ 满足 ${\alpha }^x=\beta$

大步小步算法把 $x$ 代换成：$x=i\ast m+j$

$m=\sqrt{n}，0\le i\le m，0\le j\le m$

则有：${\alpha }^x=\beta$ —> ${\alpha }^{i\ast m+j}=\beta$ —> ${\alpha }^j=\beta ({\alpha }^{-m}{)}^i$

该算法先对 $j$ 的不同取值计算出 ${\alpha }^j$ 的值，然后固定一个 $m$ ，并对 $i$ 尝试不同的取值，带入上面同余式的右边，看是否与某个（已经预先算出的）同余式左边的值相匹配。

对于离散对数问题：$y\equiv g^{x}modp$ ，已知 $y,g,p$ 求 $x$

令 $x=i\ast m+j$ ，$m=\sqrt{p}，0\le i\le m，0\le j\le m$

$y\equiv g^x\equiv g^{i\ast m+j}$

$y\ast g^{-i\ast m}\equiv g^j$

时间复杂度为 $O(\sqrt{p})$

python实现代码：

def bsgs(g, y, p):
    """ 计算在模 p 的条件下，以 g 为底 y 的对数 """
    m = int(sqrt(p))
    S = {pow(g, j, p): j for j in range(m)}  # 键值对，匹配键返回值
    gs = pow(g, -m, p)
    for i in range(m):
        if y in S:
            return i * m + S[y]
        y = y * gs % p
    return None

Pollard’s Rho algorithm（基于Miller-Rabin算法，递归求解）

Pollard’s kangaroo algorithm(袋鼠算法)

在计算数论和计算代数中，Pollard 的袋鼠算法（也称为Pollard 的 lambda 算法）是一种求解离散对数问题的算法。该算法由数论学家John M. Pollard于 1978 年在同一篇论文中提出，与他更著名的 Pollard's rho 算法 一起解决了相同的问题。尽管 Pollard 描述了他的算法在单位乘法群中以素数p为模的离散对数问题的应用，但它实际上是一种通用的离散对数算法，它适用于任何有限循环群。

**相关论文：**https://www.ams.org/journals/mcom/1978-32-143/S0025-5718-1978-0491431-9/S0025-5718-1978-0491431-9.pdf

算法：

设 $G$ 为模 $p$ 意义下 $g$ 生成的循环群，我们要做的就是从 $G$ 中找到特定的 $x$ 满足 $\beta \equiv g^{x}modp$ 。$x\in [a,b]$ ，记整数集为 $S$ ，选择一个函数 $f:G\to S$ 。人们可以通过设置搜索整个可能的对数范围 $a=0,b=n-1$

定理 $\left\{x_i\right\},\left\{y_i\right\}$ 满足：

$$\{\begin{array}{l}{x}_0\equiv g^{b}modp\\ x_{i+1}\equiv x_i{g}^{f(x_i)}modp\end{array}$$$$\{\begin{array}{l}{y}_0\equiv \beta modp\\ y_{i+1}\equiv y_i{g}^{f(y_i)}modp\end{array}$$

令 $d_{xn}=\sum _{i=0}^{n-1}f(x_i)$，$d_{yn}=\sum _{i=0}^{n-1}f(y_i)$ ，有：

$x_n\equiv x_0\ast g^{d_{xn}}\equiv g^{b+d_{xn}}modp$

同时：

$y_n\equiv y_0\ast g^{d_{yn}}\equiv \beta \ast g^{d_{yn}}modp$

接下来，我们先计算前 $N$ 项 $x_i$ （一般取 $(b-a{)}^{\frac{1}{2}}$ 项），再不断的计算 $y_i$ ，如果发现：

1 、$x_n=y_j$ ， 对于一些 $j$ 。如果序列 $\left\{x_i\right\},\left\{y_i\right\}$ 以这种方式 “碰撞” ，那么我们有：

$x_n=y_j\to g^{b+d_{xn}}=\beta \ast g^{d_{yj}}$ ，此时有：$\beta =g^{b+d_{xn}-dyj}$ ，即：$x=b+d_{xn}-d_{yj}$

2、 $d_{yj}>b-a+d_{xn}$ ，则需要重新考虑 $S$ 或者 $f$

在知道上下界的情况下，预期的时间复杂度为 $O(\sqrt{b-a})$ ，可以使用sage中的 discrete_log_lambda()函数求解。

例题：

from Crypto.Util.number import *
import random
from gmpy2 import *

flag = b'NSSCTF{%d}' % getPrime(40)

q = getPrime(160)
while True:
    t = 2 * getPrime(1024 - 160) * q
    if isPrime(t + 1):
        p = t + 1
        break

h = random.randint(1, p - 2)
g = powmod(h, (p - 1) // q, p)
x = int(flag[7:-1])
y = powmod(g, x, p)

print(p, g, y, sep=',')

'''
170348437069105675857023420369852707107199919792133184452360056253809337668300908125853674590750392947050710727732463955472261258268602007044222277390154053179972622853725777103969911116086298972690976389030851269618265250128626615573340600153282617025741773802362792814749643535242494474688828746111924564347,165897807691210138242402562954872978913428506002086890183488722074489316352258963927481837058576723028094761779902049658188277585007037781070468378089752847174284810042043509739186585150021863898785139309655899362319277284062985948365352003897500446368706636436322274748243327625738381556542611629507069843590,137261382691605813386002060971037192576990902114543987122319043412750149427553984596309208612435032869944636556014997562080112842853873993200670664292253373789312207761559974016876136547508448082441175904892258008489730428622596221443832423875069612364113749289160065933164504747254869418301359752739411555852
'''

题目分析和EXP：

我们需要解决离散对数问题 ：$y\equiv g^{x}modp$ ，这里 $p$ 的大小大概为 1024 bits 左右，如果采用 BSGS 算法的话，时间复杂度则为： $\sqrt{p}=2^{512}$ ，仍然是一个极大的计算量。

这里有题目条件可知 x=40bits ，那么我们采用Pollard’s kangaroo 算法计算指数$x$的时间复杂度即为：

$O(\sqrt{2^{40}})=O(2^{20})$ ，使用 sage中的 discrete_log_lambda()函数求解。

from sage.groups.generic import bsgs

p, g, y = 170348437069105675857023420369852707107199919792133184452360056253809337668300908125853674590750392947050710727732463955472261258268602007044222277390154053179972622853725777103969911116086298972690976389030851269618265250128626615573340600153282617025741773802362792814749643535242494474688828746111924564347,165897807691210138242402562954872978913428506002086890183488722074489316352258963927481837058576723028094761779902049658188277585007037781070468378089752847174284810042043509739186585150021863898785139309655899362319277284062985948365352003897500446368706636436322274748243327625738381556542611629507069843590,137261382691605813386002060971037192576990902114543987122319043412750149427553984596309208612435032869944636556014997562080112842853873993200670664292253373789312207761559974016876136547508448082441175904892258008489730428622596221443832423875069612364113749289160065933164504747254869418301359752739411555852

F = GF(p)
g = F(g)
y = F(y)

x = discrete_log_lambda(y,g,bounds=(2^39,2^40))
print(x)

Pohlig-Hellman 算法(针对阶是光滑且仅有小素因子)

Pohig-Hellman算法是一种求解群有光滑阶时的离散对数的方法，光滑阶即阶的值可以被分解为多个小素因子。比如求解椭圆曲线上的离散对数问题时，当 n=E.order() 为光滑数，且 p-1 也为光滑数时，使用 Pohig-Hellman 算法能够快速解决曲线上的DLP问题。

注：

• 在椭圆曲线密码学中，通常情况下，椭圆曲线上的阶是一个素数，且它应该是一个大素数，以增加密码学安全性。这样的选择可以增加在椭圆曲线密码学中使用的离散对数问题的难度，从而提高密码学方案的安全性；对于 $Q=nP$ 这个问题，可以理解为 P 点的阶很大且为素数时，计算该离散对数问题求解出 n 的复杂度将会很高。反之当 P 点的阶较小，则计算离散对数问题容易。

• 如果椭圆曲线上的阶不是素数，而是一个光滑数（即可以分解为较小质数的乘积），那么会引发一些安全性问题。在椭圆曲线密码学中，使用具有光滑阶的椭圆曲线可能会导致安全性漏洞，因为攻击者可以利用小质数的结构来降低离散对数问题的复杂度。即我们可以利用 Pohlig-Hellman 算法来有效地解决离散对数问题，从而破坏了椭圆曲线密码系统的安全性。

1、针对椭圆曲线上的离散对数问题Pohlig-Hellman算法实现：

参考论文：Weak Curves In Elliptic Curve Cryptography ，实现将离散对数计算简化为P阶的素数子群，并利用中国剩余定理解决整个阶的离散对数同余系统。

对于椭圆曲线上面的点 $Q,P\in E(F_p)$ ，且$P$为生成点，有 $Q=kP$，计算私钥 $k$ 值。

1. 设 n=P.order() ，即n为生成点P的阶，对其进行因式分解可得：$n=p_1^{e_1}\ast p_2^{e_2}\ast ...\ast p_r^{e_r}$

2. 我们想在因式分解中为每个素数找到满足 $k_i\equiv k(mod{p}_i^{e_i})$ 的值，我们用 $p_i$ 来表示 $k_i$ 的多项式形式：$k_i=z_0+z_1{p}_i+z_2{p}_i^2+...+z_{e-1}{p}_i^{e-1}$ ，并且按顺序计算每个 $z_j$ 的值。

3. 设 $P_0=\frac{n}{p_i}P$ ，$Q_0=\frac{n}{p_i}Q$ ；因此有 $p_i\ast P_0=nP$，所以关于 $Q_0$ 的表达式为：

$$Q_0=\frac{n}{p_i}Q=\frac{n}{p_i}(kP)=k(\frac{n}{p_i}P)=k{P}_0$$

4. 这时我们便将在 $P$ 域上的离散对数分解到了 $P_0$ 域上，因为 $P_0$ 的阶是 $n/p_i$，已经较原本的阶 n 运算的复杂度小了很多，有利于快速求解离散对数问题。

5. 由 $k_i\equiv k(mod{p}_i^{e_i})$ 我们可以得到一个同余组：

   $k\equiv k_1(mod{p}_1^{e_1})$

   $k\equiv k_2(mod{p}_2^{e_2})$

   ......

   $k\equiv k_r(mod{p}_r^{e_r})$

6. 将 $k_i$ 的值都求出来后，即可利用中国剩余定理将私钥 $k$ 的值求出来。

sagemath代码实现：

def Pohlig_Hellman(E,P,G):
    """
    :param E: 有限域Fp上的椭圆曲线 E
    :param P: 公钥 P=k*G
    :param G: 基点 G
    :return: 私钥 k
    """
    # 曲线基点的阶
    n = E.order()
    # 将椭圆曲线的阶分解为若干小素数
    factors, exponents = zip(*factor(n))
    print(factors)
    primes = [factors[i] ^ exponents[i] for i in range(len(factors))][:-2]  # 这一步对于小素数，如果大于我们求解的明文的话，可以省略较大的素数值，加速计算结果。

    dlogs = []
    for fac in primes:
        t = int(int(n) // int(fac))
        # dlog为离散对数分解后对应的私钥值
        dlog = discrete_log(t * P, t * G, operation="+")
        dlogs += [dlog]
        print("factor: " + str(fac) + ", Discrete Log: " + str(dlog))  # calculates discrete logarithm for each prime order
    k = crt(dlogs, primes)  # 中国剩余定理
    return k

例题：

[picoCTF 2017] ECC 2

[NCTF 2022] superecc

(非)超奇异椭圆曲线上的 MOV 攻击

关于超奇异椭圆曲线的知识点：超奇异椭圆曲线

MOV攻击（Menezes–Okamoto–Vanstone attack）是一种将椭圆曲线上的离散对数问题（ECDLP）转换为有限域上的离散对数问题（DLP）的方法。这种攻击方法是由Alfred Menezes、Tatsuaki Okamoto和Scott Vanstone在1993年提出的，它利用了椭圆曲线上的双线性配对（bilinear pairings）。

MOV攻击步骤

1. 选择椭圆曲线和点：

   • 选择一条定义在有限域 $F_q$ 上的椭圆曲线 $E$ 。

   • 选择椭圆曲线上的两个点 $P$ 和 $Q$，其中 $P$ 的阶为 $n$ 是一个大素数，即 $P.order()=n$；$Q=KP$，$K$ 未知。

2. 计算嵌入度：

   • 嵌入度 $k$ 是满足 $n$ 能嵌入到有限域 $F_q$ 中的最小整数。 $k$ 可以通过求解 $q^k\equiv 1(modn)$来找到。如果 $k$ 较小（通常小于20），那么MOV攻击可能适用。

3. 计算Weil配对或Tate配对：

• 使用一个合适的双线性配对函数，如Weil配对或Tate配对，将椭圆曲线上的点映射到有限域$F_{q^k}$ 上。例如，Weil配对 $e_n(P,Q)$ 是一个值在 $F_{q^k}$ 上的 $n$ 次单位根。

4. 映射到有限域：

   • 应用配对函数，比如计算$e_n(P,P)$ 和 $e_n(P,Q)$，得到 $F_{q^k}$ 上的元素。这样，原本在椭圆曲线上的离散对数问题就被转换到了 $F_{q^k}$ 上的离散对数问题。

5. 解决有限域上的DLP：

   • 有限域上的DLP算法（如Pohlig-Hellman算法或Pollard's rho算法）求解。

1、利用 MOV 的关键是嵌入度 k 值较小，才能进行有效的攻击。

2、Weil配对和Tate配对是椭圆曲线密码学中重要的双线性配对，它们在密码学中的应用涉及到将椭圆曲线上的点映射到一个更大的有限域上，从而可以利用更高效的算法解决复杂的离散对数问题。

参考 sagemath教程：Points on elliptic curves

在 sagemath中计算 Weil 配对或 Tate 配对的函数分别为：

weil_pairing(Q,n)

计算该点 $P$ 与同一曲线上另一点 $Q$ 的Weil配对。

参数：

• Q — 同一曲线上的另一点
• $n$ — 点 $P$ 的阶

用法：

P.weil_pairing(Q, n)

tate_pairing(Q,n,k)

Tate配对 $t_n(P,Q)$ 是 Weil 配对的一种特殊情况，通常在超奇异椭圆曲线上使用，计算效率会比 Weil 配对更快一些。

参数：

• $Q$ — 和点 $P$ 在同一曲线 $E$ 上的一点
• $n$ — 点 $P$ 的阶
• k — 嵌入度

用法：

p = 103; A = 1; B = 18; E = EllipticCurve(GF(p), [A, B])
P = E(33, 91); n = P.order(); n
# 嵌入度
k = GF(n)(p).multiplicative_order(); k
# tate配对 Tn(P,Q)
Q = E(87, 51)
P.tate_pairing(Q, n, k)

MOV攻击实现代码

sagemath

def MOV_attack(P, G, n, p):
    """
    求解(非)/超奇异椭圆曲线上的离散对数问题，P = k*G ,返回私钥 k 值
    :param P: 公钥
    :param G: 基点
    :param n: 基点G的阶,n = G.order
    :param p: 有限域 F_p
    :return: P = k*G 中的私钥 k
    """
    # 计算嵌入度
    k = 1
    while pow(p, k, n) != 1:
        k += 1
    if k <= 20:
        print(f'[*] 嵌入度 k 值为: {k} ,这个 DLP 问题能够转移到有限域 GF(p^{k}) 上')

        # 将椭圆曲线上的点映射到更大的有限域 F_p^k 上
        Ee = EllipticCurve(GF(p ^ k), [a, b])
        # 扩展有限域上的同一点
        Pe = Ee(P)
        Ge = Ee(G)

        # 取新曲线上的随机一点
        R = Ee.random_point()
        m = R.order()  # 阶
        d = gcd(m, G.order())
        # 计算Q点，其中Q的阶是M除以d的结果，确保Q的阶与G的阶相兼容。
        Q = (m // d) * R
        # 检查G的阶是否是Q的阶的整数倍，这是为了确保Weil对运算的正确性。
        assert G.order() / Q.order() in ZZ

        # 计算Weil配对
        # P_Q = Pe.weil_pairing(Q, G.order())
        # G_Q = Ge.weil_pairing(Q, G.order())

        # 计算Tate配对 (速度更快)
        P_Q = Pe.tate_pairing(Q, G.order(), k)
        G_Q = Ge.tate_pairing(Q, G.order(), k)
        print("[*] 启动 DLP 求解程序.......")

        na = P_Q.log(G_Q)
        return na
    else:
        return f'嵌入度 k 值为: {k} ,这个 k值太大了！'