ホーム

# 0x01 前言

2022.7.28 去年 5.19 的时候创建的文档，9 月份正式决定主 cry 方向

苦于学校无密码手，多次比赛面对 RSA 的题目都束手无策，决定花小部分时间学习一下 RSA，熟悉一下 crypto 的函数。

# 0x02 简介

RSA 加密演算法是一种非对称加密演算法，用于现代、保密通信。通过公钥与私钥进行加解密，保密程度较高。

RSA 涉及部分参数如下：

参数解释
m明文
c密文
N一个大整数，有 N=p*q，p、q 为两个质数
e可任取，但要求 e 与 φ(N) 互质
φ(N)phi=φ(N)=(p-1)*(q-1)
d为 e 关于 φ(N) 的模逆元，即 ed≡1 (mod φ(N))
代码表示 d=gmpy2.invert (e,phi)

其中（N,e）为公钥，（N,d）为私钥。

m<N 且与 N 互质。

大数分解在线网站 factordb

常引入的库：

1
2
3
from Crypto.Util.number import *
import binascii
import gmpy2

yafu 分解出现 mismatched parens

1
yafu-x64 "factor(@)" -batchfile pcat.txt

# 0x03 过程

# 加密过程：

c≡me(modn)c\equiv m^e \pmod n

# 解密过程：

m≡cd(modn)m\equiv c^d \pmod n

# 正确性证明：

已知c=me(modn)c=m^e \pmod n

则cd≡med≡m1+kϕ(n)≡m×mkϕ(n)(modn)c^{d}\equiv m^{ed}\equiv m^{1+k\phi(n)} \equiv m\times m^{k\phi(n)} \pmod n

由欧拉定理得mϕ(n)≡1(modn)m^{\phi(n)}\equiv 1 \pmod n，则有

mkϕ(n)≡1(modn)m^{k\phi(n)}\equiv 1 \pmod n

故cd≡m(modn)c^{d}\equiv m \pmod n 成立

# 基础解密脚本：
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
import binascii
import gmpy2
from Crypto.Util.number import *
def decrypt(p,q,e,c):
    d = gmpy2.invert(e,(p - 1) * (q - 1))
    m = pow(c,d,p * q)
    return m
#这边我用yafu分解了n
p=780900790334269659443297956843
q=1034526559407993507734818408829
n=p*q
e=0x10001
c=0x534280240c65bb1104ce3000bc8181363806e7173418d15762
m = decrypt(p,q,e,c)
print(long_to_bytes(m))
# openssl 相关：
# 破解公钥：
1
2
3
openssl rsa -pubin -text -modulus -in gy.key
Exponent: 65537 (0x10001)#e
Modulus=A9BD4C7A7763370A042FE6BEC7DDC841602DB942C7A362D1B5D372A4D08912D9#n的16进制形式
# 排列组合迭代器：
迭代器实参结果
product()p, q, ... [repeat=1]笛卡尔积，相当于嵌套的 for 循环
permutations()p[, r]长度 r 元组，所有可能的排列，无重复元素
combinations()p, r长度 r 元组，有序，无重复元素
combinations_with_replacement()p, r长度 r 元组，有序，元素可重复
例子结果
:-----------------------------------------:------------------------------------------------
product('ABCD', repeat=2)AA AB AC AD BA BB BC BD CA CB CC CD DA DB DC DD
permutations('ABCD', 2)AB AC AD BA BC BD CA CB CD DA DB DC
combinations('ABCD', 2)AB AC AD BC BD CD
combinations_with_replacement('ABCD', 2)AA AB AC AD BB BC BD CC CD DD
# ctf 题目过验证脚本
# 自动 pass 验证
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
from pwn import *
from Crypto.Util.number import *
import gmpy2
import string
import itertools
import hashlib
table = string.ascii_letters + string.digits
r = remote('node4.buuoj.cn',29242)
def proof():
    r.recvuntil(b'sha256(XXXX')
    line = r.recvline()[:-1].decode()
    print(line)
    tmp = line[line.find('+') + 1:line.find(')')]
    print(tmp)
    aim = line[line.find('== ') + 3:]
    print(aim)
    for i in itertools.product(table,repeat=4):
        ans = ''.join(i)
        if hashlib.sha256((ans + tmp).encode()).hexdigest() == aim:
            print(ans)
            r.recvuntil(b'Give me XXXX:')
            r.sendline(ans.encode())
            return

proof()
r.interactive()
# 手动 pass 验证
1
2
3
4
5
6
7
8
9
10
11
12
13
14
from pwn import *
import string
import itertools
import hashlib
table = string.ascii_letters + string.digits
def proof():
    tmp = 'w5p2zxChtQ7tcdbL'
    aim = 'e4a9453196c76cf7a820abc66bf35d050de737eee1ca2acc9487094fb904a02b'
    for i in itertools.product(table,repeat=4):
        ans = ''.join(i)
        if hashlib.sha256((ans + tmp).encode()).hexdigest() == aim:
            print(ans)
            return
proof()

# 0x04 常见函数

# gmpy2:
# iroot(x,n):

对整数 x 开 n 次方，返回第一个数为开方结果向下取整数据类型为 mpz；第二个为 bool 变量，若能开整方为 True，否则为 False。

# gcdext(a,b):

扩展欧几里得算法，ax+by=gcd (a,b)(a,b 的公约数)，需注意正负号问题。可用来计算模逆元。返回三个值，分别是：最大公约数、x、y

# mpz(x):
# 初始化一个大整数 x
# sympy:
# prevprime(n):

返回数字 n 的上一个质数

# nextprime(n):

返回数字 n 的下一个质数

# solve():

求解数学方程，返回方程式的根

example:

1
2
3
4
5
6
7
8
9
import sympy
x = sympy.Symbol('x')
y = sympy.Symbol('y')
f1 = x + 3 * y - 7
f2 = x + y - 3
result = sympy.solve([f1,f2],[x,y])
print(result)

#{x: 1, y: 2}
# sage 中多项式求根：roots ()
1
2
3
4
5
6
7
8
9
10
11
12
13
state1 = 664406108637237317109372377546194289077117372665932150107678757303243619963182733408311739663233548725195723884930183442927293177078507
state1 = (state1 << 64) + 4696961136341886199
a = 68823589009570846705623113091430161477799561031575612135516351257937127579444
b = 76549169049080319489163980188997771079750670038002598167961495813084486794567
c = 99215492498952976642031024510129092308042391285395704888838178561670205468882
p = 12775129699668988473759438271274836254349225413222075887429387682336494103348583050672280757042383792640084197832605411237644937815012509935794275313643031
P.<x> = PolynomialRing(Zmod(p))
#x = PolynomialRing(Zmod(p),'x').gen()
#roots = f.small_roots(X=2^281,beta=0.45,epsilon=0.05)
f = a * x * x + b * x + c - state1
f = f.monic()
print(f.roots())#模数为质数
print(f.small_roots(X=2^256))
# discrete_log(n,a,b)

加密算法中若已知 n,m,c，则可求 e

c=me(modn)c=m^{e} \pmod n
dd=discrete_log(n,m,cn,m,c)
e=gmpy2.invert(d,phi)e=gmpy2.invert(d,phi)

# 欧拉函数通式

ϕ(x)=x∏i=1n(1−1pi)\phi(x)=x \prod _{i=1} ^{n} (1- \frac{1}{p^{i}})

1
2
3
4
5
6
7
8
9
10
def euler_phi(n):
    ans = n
    for i in range(2, int(sqrt(n)) + 1):
        if n % i == 0:
            ans = ans // i * (i - 1)
            while n % i == 0:
                n = n / i
    if n > 1:
        ans = ans // n * (n - 1)
    return ans
# nthroot_mod(a,n,p)

找出 x 满足xn≡a(modp)x^{n} \equiv a \pmod p

1
2
3
4
5
6
from sympy.ntheory.residue_ntheory import nthroot_mod
from Crypto.Util.number import *
p = 6378942048622331926713333230890982896005250422135330937700526241780531631833065236713335883143504963836211037847951164472600586972332152350894996920936093
c = 4320014546555113258165901655495054266050880690276767417007499881878580270123774083341715559106036817554803445013957408506225560898516994410034594453027434
m = nthroot_mod(c,2,p,all_roots=True)
print(long_to_bytes(m))

# 0x05 常见攻击

# 1. 已知 e、p、q、c
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
#!/usr/bin/env python
# -*- coding:utf-8 -*-
from Crypto.Util.number import *
import binascii
import gmpy2
import hashlib
p = 104046835712664064779194734974271185635538927889880611929931939711001301561682270177931622974642789920918902563361293345434055764293612446888383912807143394009019803471816448923969637980671221111117965227402429634935481868701166522350570364727873283332371986860194245739423508566783663380619142431820861051179
q = 140171048074107988605773731671018901813928130582422889797732071529733091703843710859282267763783461738242958098610949120354497987945911021170842457552182880133642711307227072133812253341129830416158450499258216967879857581565380890788395068130033931180395926482431150295880926480086317733457392573931410220501
n = p * q
e = 65537
c = 4772758911204771028049020670778336799568778930072841084057809867608022732611295305096052430641881550781141776498904005589873830973301898523644744951545345404578466176725030290421649344936952480254902939417215148205735730754808467351639943474816280980230447097444682489223054499524197909719857300597157406075069204315022703894466226179507627070835428226086509767746759353822302809385047763292891543697277097068406512924796409393289982738071019047393972959228919115821862868057003145401072581115989680686073663259771587445250687060240991265143919857962047718344017741878925867800431556311785625469001771370852474292194
phi = (p-1)*(q-1)
d = int(gmpy2.invert(e,phi))
m = pow(c,d,n)
print(binascii.unhexlify(hex(m)[2:]))
# 2. 共模攻击 已知公共 n 和两组 e、c
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
#!/usr/bin/env python
# -*- coding:utf-8 -*-
import gmpy2
import libnum
from Crypto.Util.number import *
import binascii
from gmpy2 import *
def commonattack(n,e1,e2,c1,c2):
    ee,r,s = gmpy2.gcdext(e1,e2)
    c = pow(c1,r,n) * pow(c2,s,n) % n
    return c,ee
e1 = 65536
e2 = 270270
c1 = 39449016403735405892343507200740098477581039605979603484774347714381635211925585924812727991400278031892391996192354880233130336052873275920425836986816735715003772614138146640312241166362203750473990403841789871473337067450727600486330723461100602952736232306602481565348834811292749547240619400084712149673
c2 = 43941404835820273964142098782061043522125350280729366116311943171108689108114444447295511969090107129530187119024651382804933594308335681000311125969011096172605146903018110328309963467134604392943061014968838406604211996322468276744714063735786505249416708394394169324315945145477883438003569372460172268277
n = 122031686138696619599914690767764286094562842112088225311503826014006886039069083192974599712685027825111684852235230039182216245029714786480541087105081895339251403738703369399551593882931896392500832061070414483233029067117410952499655482160104027730462740497347212752269589526267504100262707367020244613503
temp = gmpy2.gcdext(e1,e2)
r,s = temp[1],temp[2]
m = (pow(c1,r,n)*pow(c2,s,n))%n
_gcd = gmpy2.gcd(e1,e2)
m = gmpy2.iroot(m,_gcd)[0]
print(long_to_bytes(m))

# 3. 已知 n,e,d：
1
2
3
4
5
6
7
8
9
10
11
12
def getpq(n, e,d):
    # for e in range(0, 100000):
    k = e * d - 1
    g = random.randint(0, n)
    while gmpy2.gcd(g, n) != 1:
        g = random.randint(0, n)
    while k % 2 == 0:
        k >>= 1
        temp = gmpy2.powmod(g, k, n)-1
        if gmpy2.gcd(temp, n) > 1 and temp != 0:
            p = gmpy2.gcd(temp, n)
            return p,n // p
# 4. 已知 p 高位 (2021 CISCN rsa)
1
2
3
4
5
6
7
8
9
p3 = getPrime(512)
q3 = getPrime(512)
N3 = p3*q3
print pow(msg3,e3,N3)
print (e3,N3)
print p3>>200
#59213696442373765895948702611659756779813897653022080905635545636905434038306468935283962686059037461940227618715695875589055593696352594630107082714757036815875497138523738695066811985036315624927897081153190329636864005133757096991035607918106529151451834369442313673849563635248465014289409374291381429646
#(65537, 113432930155033263769270712825121761080813952100666693606866355917116416984149165507231925180593860836255402950358327422447359200689537217528547623691586008952619063846801829802637448874451228957635707553980210685985215887107300416969549087293746310593988908287181025770739538992559714587375763131132963783147L)
#7117286695925472918001071846973900342640107770214858928188419765628151478620236042882657992902

用 sage 算出 p:

1
2
3
4
5
6
7
8
9
10
p4 = 7117286695925472918001071846973900342640107770214858928188419765628151478620236042882657992902
n = 113432930155033263769270712825121761080813952100666693606866355917116416984149165507231925180593860836255402950358327422447359200689537217528547623691586008952619063846801829802637448874451228957635707553980210685985215887107300416969549087293746310593988908287181025770739538992559714587375763131132963783147L
pbits = 512
kbits = 200
p4 = p4 << kbits
PR.<x> = PolynomialRing(Zmod(n))
f = p4 + x
root = f.small_roots(X=2^kbits, beta=0.4)
p = p4+root[0]
print(p)
# 5.dp 泄露 (la 佬的博客)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
import gmpy2 as gp
e = 0x10001
n = 81196282992606113591233615204680597645208562279327854026981376917977843644855180528227037752692498558370026353244981467900057157997462760732019372185955846507977456657760125682125104309241802108853618468491463326268016450119817181368743376919334016359137566652069490881871670703767378496685419790016705210391
dp = 5372007426161196154405640504110736659190183194052966723076041266610893158678092845450232508793279585163304918807656946147575280063208168816457346755227057
c = 61505256223993349534474550877787675500827332878941621261477860880689799960938202020614342208518869582019307850789493701589309453566095881294166336673487909221860641809622524813959284722285069755310890972255545436989082654705098907006694780949725756312169019688455553997031840488852954588581160550377081811151
for x in range(1, e):
	if(e*dp%x==1):
		p=(e*dp-1)//x+1
		if(n%p!=0):
			continue
		q=n//p
		phin=(p-1)*(q-1)
		d=gp.invert(e, phin)
		m=gp.powmod(c, d, n)
		if(len(hex(m)[2:])%2==1):
			continue
		print('--------------')
		print(m)
		print(hex(m)[2:])
		print(bytes.fromhex(hex(m)[2:]))
# 6. 低加密指数广播攻击（给多组 n,c）
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
import gmpy2
from functools import reduce
import binascii
from Crypto.Util.number import *
n1 = int(str(331310324212000030020214312244232222400142410423413104441140203003243002104333214202031202212403400220031202142322434104143104244241214204444443323000244130122022422310201104411044030113302323014101331214303223312402430402404413033243132101010422240133122211400434023222214231402403403200012221023341333340042343122302113410210110221233241303024431330001303404020104442443120130000334110042432010203401440404010003442001223042211442001413004),5)
c1 = int(str(310020004234033304244200421414413320341301002123030311202340222410301423440312412440240244110200112141140201224032402232131204213012303204422003300004011434102141321223311243242010014140422411342304322201241112402132203101131221223004022003120002110230023341143201404311340311134230140231412201333333142402423134333211302102413111111424430032440123340034044314223400401224111323000242234420441240411021023100222003123214343030122032301042243),5)

n2 = int(str(302240000040421410144422133334143140011011044322223144412002220243001141141114123223331331304421113021231204322233120121444434210041232214144413244434424302311222143224402302432102242132244032010020113224011121043232143221203424243134044314022212024343100042342002432331144300214212414033414120004344211330224020301223033334324244031204240122301242232011303211220044222411134403012132420311110302442344021122101224411230002203344140143044114),5)
c2 = int(str(112200203404013430330214124004404423210041321043000303233141423344144222343401042200334033203124030011440014210112103234440312134032123400444344144233020130110134042102220302002413321102022414130443041144240310121020100310104334204234412411424420321211112232031121330310333414423433343322024400121200333330432223421433344122023012440013041401423202210124024431040013414313121123433424113113414422043330422002314144111134142044333404112240344),5)

n3 = int(str(332200324410041111434222123043121331442103233332422341041340412034230003314420311333101344231212130200312041044324431141033004333110021013020140020011222012300020041342040004002220210223122111314112124333211132230332124022423141214031303144444134403024420111423244424030030003340213032121303213343020401304243330001314023030121034113334404440421242240113103203013341231330004332040302440011324004130324034323430143102401440130242321424020323),5)
c3 = int(str(10013444120141130322433204124002242224332334011124210012440241402342100410331131441303242011002101323040403311120421304422222200324402244243322422444414043342130111111330022213203030324422101133032212042042243101434342203204121042113212104212423330331134311311114143200011240002111312122234340003403312040401043021433112031334324322123304112340014030132021432101130211241134422413442312013042141212003102211300321404043012124332013240431242),5)
e = 3#通常情况是有n组数据，e=n
n = [n1,n2,n3]
c = [c1,c2,c3]
def CRT(mi, ai):
    assert(reduce(gmpy2.gcd,mi)==1)
    assert (isinstance(mi, list) and isinstance(ai, list))
    M = reduce(lambda x, y: x * y, mi)
    ai_ti_Mi = [a * (M // m) * gmpy2.invert(M // m, m) for (m, a) in zip(mi, ai)]
    return reduce(lambda x, y: x + y, ai_ti_Mi) % M
m=gmpy2.iroot(CRT(n, c), e)[0]
print(binascii.unhexlify(hex(m)[2:].strip("L")))
# 7. 光滑数

Pollard’s p−1 （ python -m primefac -vs -m=p-1 xxxxxxx ）（光滑数）

Williams’s p+1（ python -m primefac -vs -m=p+1 xxxxxxx ）（光滑数）

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
def pollard_p_1(n):
	small_primes = [n for n in range(2,10000) if isPrime(n)]#range范围可变
	x = 2
	for elt in small_primes:#超过三次方则无法计算 但速度较快
		x*=elt
	for elt in small_primes:
		x*=elt
	for elt in small_primes:
		x*=elt

	q = GCD((pow(2,x,n)-1)%n,n)
	p = n//q

	return [p,q]

def Pollards_p_1(N):#计算速度较慢
    a = 2
    n = 2
    while True:
        a = pow(a, n, N)
        res = gmpy2.gcd(a-1, N)
        if res != 1 and res != N:
            print('n =', n)
            print('p =', res)
            return res
        n += 1

生成脚本如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
from Crypto.Util.number import *


def gen_prime(nbits, gamma):
    g = 2 * getPrime(int(nbits * gamma))
    while True:
        a = getRandomNBitInteger(int((0.5 - gamma) * nbits - 1))
        p = g * a + 1
        if isPrime(p):
            b = getRandomNBitInteger(int((0.5 - gamma) * nbits - 1))
            q = g * b + 1
            while not isPrime(q) or GCD(a, b) != 1:
                b = getRandomNBitInteger(int((0.5 - gamma) * nbits - 1))
                q = g * b + 1
            return p, q


def gen_key(nbits, gamma):
    p, q = gen_prime(nbits, gamma)
    n = p * q
    lcm = (p * q) // GCD(p, q)
    e = getPrime(16)
    while GCD(e, lcm) != 1:
        e = getPrime(16)
    d = inverse(e, lcm)
    return (n, e), (p, q, d)


def main():
    (n, e), (p, q, d) = gen_key(2048, 0.485)
    print('n =', n)
    print('p =', p)
    print('q =', q)
    print('e =', e)


if __name__ == '__main__':
    main()
# 8. 连分数维纳攻击 WienerAttack

RSA 中维纳攻击使用条件d<13N14d<\frac{1}{3}N^{\frac{1}{4}}

Legendre 's theorem

ξ\xi 是一个正数，若gcd(a,b)==1\text{gcd}(a,b) == 1 且∣ξ−ab∣<12b2|\xi-\frac{a}{b}|<\frac{1}{2b^2}。

则对ξ\xi 的连分数展开可以得到\frac{a}

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
from Crypto.Util.number import *
def continuedFra(x,y):
    cf = []
    while y != 0:
        cf.append(x // y)
        x,y = y,x % y
    return cf

def exp(x,y):
    cf = continuedFra(x,y)
    fz = [cf[0],cf[0] * cf[1] + 1]
    fm = [1,cf[1]]
    for i in range(2,len(cf)):
        z = fz[i - 1] * cf[i] + fz[i - 2]
        m = fm[i - 1] * cf[i] + fm[i - 2]
        fz.append(z)
        fm.append(m)
    return fz,fm

def get_d(ds,e,n):
    for d in ds:
        cc = pow(2,e,n)
        if pow(cc,d,n) == 2:
            return d

n = 70797287183558503294424688012661311347401562320738662322515503372418303614976509133453889639849861382169710102785409641199175222842024949807174225076183289128512968273464757607640459479228291059466767031257567440876318784090686027065614914055345627566459797405374881579074036798019061589506339186308789749727
e = 57119035388441682129004645726731934057385195386916473587539158764997045366658122483124672260727492209278576478223181874445160685700241429752484148215732176167740916689891387164255144702607704963686393882626927382203639610368394490188027376571344473424985111146271923650642752668682664268394866053697522126039
c = 70025167656914452138602016128903814432733132373586129151156774199862557367143321492623284914821677153933851187313412266040556220938639194083352203836954782252965527841367117158757341972245284383407689753010666644090978876520563247697916457595448376368084117328367741247511656671531963649242054417373789797611
_,ds = exp(e,n)
d = get_d(ds,e,n)
m = pow(c,d,n)
print(long_to_bytes(m))

sage 版

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
from Crypto.Util.number import *
def attack_d(n,e):
    lst = continued_fraction(Integer(e)/Integer(n))
    conv = lst.convergents()
    _m = 2
    _c = pow(_m,e,n)
    for i in conv:
        _,d = i.numerator(),i.denominator()
        if pow(_c,d,n) == _m:
            return d
n = 70797287183558503294424688012661311347401562320738662322515503372418303614976509133453889639849861382169710102785409641199175222842024949807174225076183289128512968273464757607640459479228291059466767031257567440876318784090686027065614914055345627566459797405374881579074036798019061589506339186308789749727
e = 57119035388441682129004645726731934057385195386916473587539158764997045366658122483124672260727492209278576478223181874445160685700241429752484148215732176167740916689891387164255144702607704963686393882626927382203639610368394490188027376571344473424985111146271923650642752668682664268394866053697522126039
c = 70025167656914452138602016128903814432733132373586129151156774199862557367143321492623284914821677153933851187313412266040556220938639194083352203836954782252965527841367117158757341972245284383407689753010666644090978876520563247697916457595448376368084117328367741247511656671531963649242054417373789797611
d = attack_d(n,e)
m = pow(c,d,n)
print(long_to_bytes(int(m)))

扩展维纳，参考论文 https://eprint.iacr.org/2008/459.pdf

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
def extendwiener(e, n):
    m = 12345
    c = pow(m, e, n)
    q0 = 1
    list = continued_fraction(e / n)  # 计算连分数
    conv = list.convergents()  # 计算渐近分数
    for i in tqdm(conv):
        q1 = i.denominator()  # 渐近分母

        for r in range(20):
            for s in range(20):
                d = r * q1 + s * q0
                m1 = pow(c, d, n)
                if m1 == m:
                    return d
        q0 = q1
# 9. 多元 coppersmith 板子
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
import itertools
def small_roots(f, bounds, m=1, d=None):
    if not d:
        d = f.degree()

    R = f.base_ring()
    N = R.cardinality()

    f /= f.coefficients().pop(0)
    f = f.change_ring(ZZ)

    G = Sequence([], f.parent())
    for i in range(m + 1):
        base = N ^ (m - i) * f ^ i
        for shifts in itertools.product(range(d), repeat=f.nvariables()):
            g = base * prod(map(power, f.variables(), shifts))
            G.append(g)

    B, monomials = G.coefficient_matrix()
    monomials = vector(monomials)

    factors = [monomial(*bounds) for monomial in monomials]
    for i, factor in enumerate(factors):
        B.rescale_col(i, factor)

    B = B.dense_matrix().LLL()

    B = B.change_ring(QQ)
    for i, factor in enumerate(factors):
        B.rescale_col(i, 1 / factor)

    H = Sequence([], f.parent().change_ring(QQ))
    for h in filter(None, B * monomials):
        H.append(h)
        I = H.ideal()
        if I.dimension() == -1:
            H.pop()
        elif I.dimension() == 0:
            roots = []
            for root in I.variety(ring=ZZ):
                root = tuple(R(root[var]) for var in f.variables())
                roots.append(root)
            return roots

    return []
PR.<x> = PolynomialRing(Zmod(n))
bounds = (1 << 64,1 << 86,1 << 26)
roots = small_roots(f,bounds,m=3)[0]
# 10.Boneh and Durfee attack

直接将 n 和 e 替换板子中的就行，ee 非常大接近于NN，即 dd 较小时。与低解密指数攻击类似，比低解密指数攻击（Wiener Attack）更强，可以解决13N14≤d≤N0.292\cfrac{1}{3}N^{\frac{1}{4}} \leq d \leq N^{0.292} 的问题。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
from __future__ import print_function
import time

############################################
# Config
##########################################

"""
Setting debug to true will display more informations
about the lattice, the bounds, the vectors...
"""
debug = True

"""
Setting strict to true will stop the algorithm (and
return (-1, -1)) if we don't have a correct
upperbound on the determinant. Note that this
doesn't necesseraly mean that no solutions
will be found since the theoretical upperbound is
usualy far away from actual results. That is why
you should probably use `strict = False`
"""
strict = False

"""
This is experimental, but has provided remarkable results
so far. It tries to reduce the lattice as much as it can
while keeping its efficiency. I see no reason not to use
this option, but if things don't work, you should try
disabling it
"""
helpful_only = True
dimension_min = 7 # stop removing if lattice reaches that dimension

############################################
# Functions
##########################################

# display stats on helpful vectors
def helpful_vectors(BB, modulus):
    nothelpful = 0
    for ii in range(BB.dimensions()[0]):
        if BB[ii,ii] >= modulus:
            nothelpful += 1

    print(nothelpful, "/", BB.dimensions()[0], " vectors are not helpful")

# display matrix picture with 0 and X
def matrix_overview(BB, bound):
    for ii in range(BB.dimensions()[0]):
        a = ('%02d ' % ii)
        for jj in range(BB.dimensions()[1]):
            a += '0' if BB[ii,jj] == 0 else 'X'
            if BB.dimensions()[0] < 60:
                a += ' '
        if BB[ii, ii] >= bound:
            a += '~'
        print(a)

# tries to remove unhelpful vectors
# we start at current = n-1 (last vector)
def remove_unhelpful(BB, monomials, bound, current):
    # end of our recursive function
    if current == -1 or BB.dimensions()[0] <= dimension_min:
        return BB

    # we start by checking from the end
    for ii in range(current, -1, -1):
        # if it is unhelpful:
        if BB[ii, ii] >= bound:
            affected_vectors = 0
            affected_vector_index = 0
            # let's check if it affects other vectors
            for jj in range(ii + 1, BB.dimensions()[0]):
                # if another vector is affected:
                # we increase the count
                if BB[jj, ii] != 0:
                    affected_vectors += 1
                    affected_vector_index = jj

            # level:0
            # if no other vectors end up affected
            # we remove it
            if affected_vectors == 0:
                print("* removing unhelpful vector", ii)
                BB = BB.delete_columns([ii])
                BB = BB.delete_rows([ii])
                monomials.pop(ii)
                BB = remove_unhelpful(BB, monomials, bound, ii-1)
                return BB

            # level:1
            # if just one was affected we check
            # if it is affecting someone else
            elif affected_vectors == 1:
                affected_deeper = True
                for kk in range(affected_vector_index + 1, BB.dimensions()[0]):
                    # if it is affecting even one vector
                    # we give up on this one
                    if BB[kk, affected_vector_index] != 0:
                        affected_deeper = False
                # remove both it if no other vector was affected and
                # this helpful vector is not helpful enough
                # compared to our unhelpful one
                if affected_deeper and abs(bound - BB[affected_vector_index, affected_vector_index]) < abs(bound - BB[ii, ii]):
                    print("* removing unhelpful vectors", ii, "and", affected_vector_index)
                    BB = BB.delete_columns([affected_vector_index, ii])
                    BB = BB.delete_rows([affected_vector_index, ii])
                    monomials.pop(affected_vector_index)
                    monomials.pop(ii)
                    BB = remove_unhelpful(BB, monomials, bound, ii-1)
                    return BB
    # nothing happened
    return BB

""" 
Returns:
* 0,0   if it fails
* -1,-1 if `strict=true`, and determinant doesn't bound
* x0,y0 the solutions of `pol`
"""
def boneh_durfee(pol, modulus, mm, tt, XX, YY):
    """
    Boneh and Durfee revisited by Herrmann and May
    
    finds a solution if:
    * d < N^delta
    * |x| < e^delta
    * |y| < e^0.5
    whenever delta < 1 - sqrt(2)/2 ~ 0.292
    """

    # substitution (Herrman and May)
    PR.<u, x, y> = PolynomialRing(ZZ)
    Q = PR.quotient(x*y + 1 - u) # u = xy + 1
    polZ = Q(pol).lift()

    UU = XX*YY + 1

    # x-shifts
    gg = []
    for kk in range(mm + 1):
        for ii in range(mm - kk + 1):
            xshift = x^ii * modulus^(mm - kk) * polZ(u, x, y)^kk
            gg.append(xshift)
    gg.sort()

    # x-shifts list of monomials
    monomials = []
    for polynomial in gg:
        for monomial in polynomial.monomials():
            if monomial not in monomials:
                monomials.append(monomial)
    monomials.sort()
    
    # y-shifts (selected by Herrman and May)
    for jj in range(1, tt + 1):
        for kk in range(floor(mm/tt) * jj, mm + 1):
            yshift = y^jj * polZ(u, x, y)^kk * modulus^(mm - kk)
            yshift = Q(yshift).lift()
            gg.append(yshift) # substitution
    
    # y-shifts list of monomials
    for jj in range(1, tt + 1):
        for kk in range(floor(mm/tt) * jj, mm + 1):
            monomials.append(u^kk * y^jj)

    # construct lattice B
    nn = len(monomials)
    BB = Matrix(ZZ, nn)
    for ii in range(nn):
        BB[ii, 0] = gg[ii](0, 0, 0)
        for jj in range(1, ii + 1):
            if monomials[jj] in gg[ii].monomials():
                BB[ii, jj] = gg[ii].monomial_coefficient(monomials[jj]) * monomials[jj](UU,XX,YY)

    # Prototype to reduce the lattice
    if helpful_only:
        # automatically remove
        BB = remove_unhelpful(BB, monomials, modulus^mm, nn-1)
        # reset dimension
        nn = BB.dimensions()[0]
        if nn == 0:
            print("failure")
            return 0,0

    # check if vectors are helpful
    if debug:
        helpful_vectors(BB, modulus^mm)
    
    # check if determinant is correctly bounded
    det = BB.det()
    bound = modulus^(mm*nn)
    if det >= bound:
        print("We do not have det < bound. Solutions might not be found.")
        print("Try with highers m and t.")
        if debug:
            diff = (log(det) - log(bound)) / log(2)
            print("size det(L) - size e^(m*n) = ", floor(diff))
        if strict:
            return -1, -1
    else:
        print("det(L) < e^(m*n) (good! If a solution exists < N^delta, it will be found)")

    # display the lattice basis
    if debug:
        matrix_overview(BB, modulus^mm)

    # LLL
    if debug:
        print("optimizing basis of the lattice via LLL, this can take a long time")

    BB = BB.LLL()

    if debug:
        print("LLL is done!")

    # transform vector i & j -> polynomials 1 & 2
    if debug:
        print("looking for independent vectors in the lattice")
    found_polynomials = False
    
    for pol1_idx in range(nn - 1):
        for pol2_idx in range(pol1_idx + 1, nn):
            # for i and j, create the two polynomials
            PR.<w,z> = PolynomialRing(ZZ)
            pol1 = pol2 = 0
            for jj in range(nn):
                pol1 += monomials[jj](w*z+1,w,z) * BB[pol1_idx, jj] / monomials[jj](UU,XX,YY)
                pol2 += monomials[jj](w*z+1,w,z) * BB[pol2_idx, jj] / monomials[jj](UU,XX,YY)

            # resultant
            PR.<q> = PolynomialRing(ZZ)
            rr = pol1.resultant(pol2)

            # are these good polynomials?
            if rr.is_zero() or rr.monomials() == [1]:
                continue
            else:
                print("found them, using vectors", pol1_idx, "and", pol2_idx)
                found_polynomials = True
                break
        if found_polynomials:
            break

    if not found_polynomials:
        print("no independant vectors could be found. This should very rarely happen...")
        return 0, 0
    
    rr = rr(q, q)

    # solutions
    soly = rr.roots()

    if len(soly) == 0:
        print("Your prediction (delta) is too small")
        return 0, 0

    soly = soly[0][0]
    ss = pol1(q, soly)
    solx = ss.roots()[0][0]

    #
    return solx, soly

def example():
    ############################################
    # How To Use This Script
    ##########################################

    #
    # The problem to solve (edit the following values)
    #
    N = 16043350038355705689747210462564029770514189048023686603223668062559298551137941167871392968882385042895904416328284532906698458292029131698076286801441138815205476275923051191994880955569436018845008058066246413766713207496731647699425954184192904658373733886603851370232470587424580640498514158985740438571951114737629232289370996612855323683190398113040474182083212739351537697750571221786906269749610022043089724912993422380286926739373579601138533650370187776457252560937403752826750173443816338476494719557412929553824254765927609535526976489831437884368102114072944348983057935047911546598879651335996331537851
    e = 509400401183993386598745991841695789965803577018496383973610712015390731952573188892607068365045750956220514439318420820414674457717326610021697773025034023029106539445589380601502110309168017135751408755233566808484079679030884201065957896422213752055277177097935656421755231895714928093448868109156498366044486135035000248179690996646026809626974650538004341731057827787632279121306065487494611608226671415435208508448888705165570622606586815577489805929822189113929662282665193744133120319820122535526728275841105116073019206852443361320494448818509280320727111572056547660440166869901028959846637388604855898329

    # the modulus
   # N = 0xc2fd2913bae61f845ac94e4ee1bb10d8531dda830d31bb221dac5f179a8f883f15046d7aa179aff848db2734b8f88cc73d09f35c445c74ee35b01a96eb7b0a6ad9cb9ccd6c02c3f8c55ecabb55501bb2c318a38cac2db69d510e152756054aaed064ac2a454e46d9b3b755b67b46906fbff8dd9aeca6755909333f5f81bf74db
    # the public exponent
    #e = 0x19441f679c9609f2484eb9b2658d7138252b847b2ed8ad182be7976ed57a3e441af14897ce041f3e07916445b88181c22f510150584eee4b0f776a5a487a4472a99f2ddc95efdd2b380ab4480533808b8c92e63ace57fb42bac8315fa487d03bec86d854314bc2ec4f99b192bb98710be151599d60f224114f6b33f47e357517

    # the hypothesis on the private exponent (the theoretical maximum is 0.292)
    delta = 0.18 # this means that d < N^delta

    #
    # Lattice (tweak those values)
    #

    # you should tweak this (after a first run), (e.g. increment it until a solution is found)
    m = 4 # size of the lattice (bigger the better/slower)

    # you need to be a lattice master to tweak these
    t = int((1-2*delta) * m)  # optimization from Herrmann and May
    X = 2*floor(N^delta)  # this _might_ be too much
    Y = floor(N^(1/2))    # correct if p, q are ~ same size

    #
    # Don't touch anything below
    #

    # Problem put in equation
    P.<x,y> = PolynomialRing(ZZ)
    A = int((N+1)/2)
    pol = 1 + x * (A + y)

    #
    # Find the solutions!
    #

    # Checking bounds
    if debug:
        print("=== checking values ===")
        print("* delta:", delta)
        print("* delta < 0.292", delta < 0.292)
        print("* size of e:", int(log(e)/log(2)))
        print("* size of N:", int(log(N)/log(2)))
        print("* m:", m, ", t:", t)

    # boneh_durfee
    if debug:
        print("=== running algorithm ===")
        start_time = time.time()

    solx, soly = boneh_durfee(pol, e, m, t, X, Y)

    # found a solution?
    if solx > 0:
        print("=== solution found ===")
        if False:
            print("x:", solx)
            print("y:", soly)

        d = int(pol(solx, soly) / e)
        print("private key found:", d)
    else:
        print("=== no solution was found ===")

    if debug:
        print(("=== %s seconds ===" % (time.time() - start_time)))

if __name__ == "__main__":
    example()
# 11.Fermat’s attack

适用条件p−q<n14p-q<n^{\frac{1}{4}}

1
2
3
4
5
6
7
8
9
10
11
12
13
def fetmat(n):
    a = gmpy2.iroot(n,2)[0] + 1
    while True:
        tmp = pow(a,2) - n
        root = gmpy2.iroot(tmp,2)
        if root[1]:
            b = root[0]
            break
        else:
            a += 1
    p = a + b
    assert n % p == 0
    return p,n // p
# 12. 基础背包
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
import math
message = '001100100011000000110111001100100011000100110011'
publickey = [55183912298421517303955421431369579814, 160385599210839494939371828807545852824, 224502758444343716706249613105679862494, 56428275085389527471043040001075739102, 70825066194175985810709736604999644238, 182189572491352797860720727162090526101, 338808370032950618286884981637428243862, 95482209268790288990116600790422572183, 117198917151126812922246482823806475660, 152299506060707455972784040811410726996, 154044680544950085988334360689074120935, 304432433176792575290325776333694930900, 323321897034769947979527706761809614101, 162919704896506359672489539282825803246, 283265796501157821745353825756916883662, 159057562768209085817668041407274035690, 25404229021501562583392877126234942629, 86336557368883012766272473285361999830, 164629030810267942190431863590923485564, 14694575998113417487143606901052265614, 234394532097626625519703654893008523222, 34922472142514721731573590584419709269, 283354938027977015105304329971606310415, 34121879528119165110069595705057525230, 225528146387697418270401363626446851375, 8124474304656160160595696111736910969, 261961711875980080339266126302090089963, 176245649059964468420774795572610547522, 250731469880971850115028008435779803278, 176810419824277118673516120250781139117, 41282715322172470612304293285435894169, 305280523055695643367680491444375899115, 27346235874550002930164039480991650389, 262961464749769618014465643542268136714, 80402108399625659340698685158610227482, 73773525391843782933279077249038915524, 77272645948529523093402915335991715431, 273256036139070966394874649078342509231, 122974872807900030220556030866672041555, 276859666973499125187599029446746454273, 199843018374344446802949929932827676053, 84208090264720231847218786058514379999, 90689500345153371482995174535171832327, 145922619278379570909170168194194919808, 105295125279809362321044947457314181709, 338976660314009313790261272741754167466, 117454690989761407304143687491085793971, 231795219003920197907651849533524370917]
enc = 3106326820431718919805243645288521239341
n = len(publickey)
d = n / math.log(max(publickey),2)
print(d)
mat = [[0 for _ in range(n + 1)] for _ in range(n + 1)]
for i in range(n):
    mat[i][i] = 2
    mat[i][n] = publickey[i]
    mat[n][i] = 1
mat[n][n] = enc
mat = matrix(ZZ,mat)
tmp = mat.LLL()
ans = ''
for i in tmp[0][:-1]:
    ans += str((-i + 1) // 2)
print(ans == message)
print(tmp)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
#密度小于0.9408
#https://github.com/jvdsn/crypto-attacks/blob/master/attacks/knapsack/low_density.py
import math
import logging
# from shared.lattice import shortest_vectors
import os
import sys
from math import ceil
from math import log2
from math import sqrt

from sage.all import QQ
from sage.all import matrix
def shortest_vectors(B):
    """
    Computes the shortest non-zero vectors in a lattice.
    :param B: the basis of the lattice
    :return: a generator generating the shortest non-zero vectors
    """
    logging.debug(f"Computing shortest vectors in {B.nrows()} x {B.ncols()} matrix...")
    B = B.LLL()

    for row in B.rows():
        if not row.is_zero():
            yield row
            
def attack(a, s):
    """
    Tries to find e_i values such that sum(e_i * a_i) = s.
    This attack only works if the density of the a_i values is < 0.9048.
    More information: Coster M. J. et al., "Improved low-density subset sum algorithms"
    :param a: the a_i values
    :param s: the s value
    :return: the e_i values, or None if the e_i values were not found
    """
    n = len(a)
    d = n / log2(max(a))
    N = ceil(1 / 2 * sqrt(n))
    assert d < 0.9408, f"Density should be less than 0.9408 but was {d}."

    L = matrix(QQ, n + 1, n + 1)
    for i in range(n):
        L[i, i] = 1
        L[i, n] = N * a[i]

    L[n] = [1 / 2] * n + [N * s]
    for v in shortest_vectors(L):
        s_ = 0
        e = []
        for i in range(n):
            ei = 1 - (v[i] + 1 / 2)
            if ei != 0 and ei != 1:
                break

            ei = int(ei)
            s_ += ei * a[i]
            e.append(ei)

        if s_ == s:
            return e
s = #总和
m = #待取量
key = attack(m,s)
print(key)
# 13. 多组低解密指数攻击

https://eprint.iacr.org/2009/037.pdf

给了多组 e，且 d 较小

给了两组 e 的

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
def attack_small_d(N,e,e1,e2,bits):
    alpha2 = bits / int(N).bit_length()
    M1 = int(N ** 0.5)
    M2 = int(N ** (1 + alpha2))
    D = diagonal_matrix(ZZ, [N, M1, M2, 1])
    B = [
        [1, -N,  0,  N ** 2],
        [0, e1,-e1, -e1 * N],
        [0, 0,  e2, -e2 * N],
        [0, 0,   0,  e1 * e2]
    ]
    B = matrix(ZZ,B) * D
    mat = B.LLL()
    v = matrix(ZZ,mat[0])
    x = list(v * B ** (-1))
    phi = int(abs(e1 * x[0][1]) // abs(x[0][0]))
    d = inverse_mod(e,phi)
    return d
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
from Crypto.Util.number import *
import gmpy2
e1 = 1038188773022222237625162518466985515806685046439847462572843423800303080199415368325579099819361640945202226526678764311585592296933622966635366454130900252466567292666094830865094694540899938932650663694321540899409821558619513870070621174837528024755540950294728078653453000484865860920060355130142874882872471337494879275434166435493265146752291857135290428750750609423353579700984426964475724965353873095813181244686536072523222027403912142730610262067287620007571352094447066062529895627497159337248165671672168914514241613626520037706745398642583257608070477729851466078618962204332539106519787878047712382699
e2 = 2837849440271663829778449470456059993823700375016504578318494102782617180188657051885856137280051100635878402423110369686929227684421486015532180997159960436120141492683886715611434986294622600612428406093623932339780091710632795226634412256078611259843109876301975664056868908063144172636320692414857287719870275516722663234436495523740203192523105607062687910252368627072074836944313105637959564954309098651598325997792496430340003856687190484681832529188281328826421428597879086043647647886763379182416419551074016810300511817626177321217978912504879476086100668005286481779806010131350674761039391612993646202901
c = 5973798238952580291825915383143493132916118834759984908567429997405141389115327100612059752092101975323145558282778289524466024564450720131251849100687215493221989801105144437981685382023973692198113306045957788268110316519461230170693204752380105917975206409994893101671098451678847638671373239757408532363808224681853024689663345258120864348816343897379881239786554998688501997609152329902187048422237325117741778968505252184157273467466011959504548459297647302026380076579903441434135973514451254950835559924204821846949520738057940287763572642367638668413987340659205489659594044022422368411980101640782079189025
N = 26901814699902439156457451193693740730489294959491270367027927283506475930489639407729426818974347303153364758700002407059993182986763909124690390655890031474097185414651218374672254140022392199647526025638012909369532528422355530044873378287920255523382224453173638818751280227521077881224963029942704252587893395262633450759457753054490886171089835324182422639138198164026845488515879253564971977801724349440235209377091735281830263780308625603392942624306475075157394231585266792247387837984357822842056801420064918953837917678662504712605611080802179768683537742095990507008809197788025847612652983474906829809607

for i in range(1000):
    alpha2 = i/1000
    M1 = int(gmpy2.mpz(N)**0.5)
    M2 = int( gmpy2.mpz(N)**(1+alpha2) )
    D = diagonal_matrix(ZZ, [N, M1, M2, 1])
    B = Matrix(ZZ, [ [1, -N,   0,  N**2],
                 [0, e1, -e1, -e1*N],
                 [0,  0,  e2, -e2*N],
                 [0,  0,   0, e1*e2] ]) * D
    L = B.LLL()
    v = Matrix(ZZ, L[0])
    x = v * B**(-1)
    phi = (x[0,1]/x[0,0]*e1).floor()
    try:
        d = inverse_mod( 65537, phi)
        m = pow(c,d,N)
        flag = long_to_bytes(m)
        if b'flag' in flag:
            print(long_to_bytes(m))
    except:
        pass

给了三组 e 的

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
import gmpy2
from Crypto.Util.number import *
N = 
e1 = 
e2 = 
e3 = 
c = 

for i in range(1000):
    alpha2 = i/1000
    M1 = int(gmpy2.mpz(N)**(3./2))
    M2 = int( gmpy2.mpz(N) )
    M3 = int(gmpy2.mpz(N)**(3./2 + alpha2))
    M4 = int( gmpy2.mpz(N)**(0.5) )
    M5 = int( gmpy2.mpz(N)**(3./2 + alpha2) )
    M6 = int( gmpy2.mpz(N)**(1.+alpha2) )
    M7 = int( gmpy2.mpz(N)**(1.+alpha2) )
    D = diagonal_matrix(ZZ, [M1, M2, M3, M4, M5, M6, M7, 1])
    B = Matrix(ZZ, [ [1, -N,   0,  N**2,   0,      0,      0,    -N**3],
                 [0, e1, -e1, -e1*N, -e1,      0,   e1*N,  e1*N**2],
                 [0,  0,  e2, -e2*N,   0,   e2*N,      0,  e2*N**2],
                 [0,  0,   0, e1*e2,   0, -e1*e2, -e1*e2, -e1*e2*N],
                 [0,  0,   0,     0,  e3,  -e3*N,  -e3*N,  e3*N**2],
                 [0,  0,   0,     0,   0,  e1*e3,      0, -e1*e3*N],
                 [0,  0,   0,     0,   0,      0,  e2*e3, -e2*e3*N],
                 [0,  0,   0,     0,   0,      0,      0, e1*e2*e3] ]) * D

    L = B.LLL()

    v = Matrix(ZZ, L[0])
    x = v * B**(-1)
    phi_ = (e1*x[0,1]/x[0,0]).floor()
    try:
        d = inverse_mod( 65537, phi_)
        m = pow(c,d,N)
        flag = long_to_bytes(m)
        if b'flag' in flag:
            print(long_to_bytes(m))
    except:
        pass
# 14.Common Private Exponent（共私钥指数攻击，d 相同）

Br=[Me1e2⋯er0−N10⋯000−N2⋯0⋮⋮⋮⋱⋮000⋯−Nr]B_r= \begin{bmatrix} {M}&{e_1}&{e_2}&{\cdots}&{e_{r}}\\ {0}&{-N_1}&{0}&{\cdots}&{0}\\ {0}&{0}&{-N_2}&{\cdots}&{0}\\ {\vdots}&{\vdots}&{\vdots}&{\ddots}&{\vdots}\\ {0}&{0}&{0}&{\cdots}&{-N_r}\\ \end{bmatrix}

规约得到∣b1∣=Md\vert b_1\vert=Md

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
import gmpy2
from Crypto.Util.number import *
es = [e_0,e_1,e_2]
ns = [n_0,n_1,n_2]
l = len(es)
mat = [[0 for _ in range(l + 1)] for _ in range(l + 1)]
print(mat)
M = gmpy2.iroot(ns[-1],2)[0]
mat[0][0] = M
for i in range(1,l + 1):
    
    mat[0][i] = es[i - 1]
    mat[i][i] = -ns[i - 1]


mat = matrix(ZZ,mat)
d = abs(mat.LLL()[0][0]) // M
m = pow(c_0,d,n_0)
print(long_to_bytes(int(m)))
# 15.DLP 解离散对数
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
from Crypto.Util.number import *
# Baby-step Giant-step法
def babystep_giantstep(g, y, p, q=None):
    if q is None:
        q = p - 1
    m = int(q**0.5 + 0.5)
    # Baby step
    table = {}
    gr = 1  # g^r
    for r in range(m):
        table[gr] = r
        gr = (gr * g) % p
    # Giant step
    try:
        gm = pow(g, -m, p)  # gm = g^{-m}
    except:
        return None
    ygqm = y                # ygqm = y * g^{-qm}
    for q in range(m):
        if ygqm in table:
            return q * m + table[ygqm]
        ygqm = (ygqm * gm) % p
    return None

# Pohlig–Hellman法
def pohlig_hellman_DLP(g, y, p,ppp):# pow(g,x,p) == y
    crt_moduli = []
    crt_remain = []
    for q,_ in ppp:
        x = babystep_giantstep(pow(g,(p-1)//q,p), pow(y,(p-1)//q,p), p, q)
        if (x is None) or (x <= 1):
            continue
        crt_moduli.append(q)
        crt_remain.append(x)
    x = crt(crt_remain, crt_moduli)
    return x
g = 696376465415968446607383675953857997
y = 75351884040606337127662457946455960228423443937677603718170904462378938882502061014476822055783421908392386804380503123596242003758891619926133807099465797120624009076182390781918339985157326114840926784410018674639537246981505937318380179042568501449024366208980139650052067021073343322300422190243015076307
p = 135413548968824157679549005083702144352234347621794899960854103942091470496598900341162814164511690126111049767046340801124977369460415208157716471020260549912068072662740722359869775486486528791641600354017790255320219623493658736576842207668208174964413000049133934516641398518703502709055912644416582457721
ppp = [263,587,28142457071,395710839697]
x = pohlig_hellman_DLP(g, y, p,ppp)
assert pow(g, x, p) == y
# 16. 中国剩余定理 CRT
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
import gmpy2
from functools import reduce
def union(x1, x2):
    a1, m1 = x1
    a2, m2 = x2
    d = gmpy2.gcd(m1, m2)
    assert (a2 - a1) % d == 0
    p1,p2 = m1 // d,m2 // d
    _,l1,l2 = gmpy2.gcdext(p1,p2)
    k = -((a1 - a2) // d) * l1
    lcm = gmpy2.lcm(m1,m2)
    ans = (a1 + k * m1) % lcm
    return ans,lcm

def excrt(ai,mi):
    tmp = zip(ai,mi)
    return reduce(union, tmp)
# 17.LCG
1
2
3
4
5
6
7
8
9
10
11
12
13
14
from gmpy2 import *
from Crypto.Util.number import *
from functools import reduce
def crack_unknown_increment(states, modulus, multiplier):
    increment = (states[1] - states[0]*multiplier) % modulus
    return modulus, multiplier, increment
def crack_unknown_multiplier(states, modulus):
    multiplier = (states[2] - states[1]) * invert(states[1] - states[0], modulus) % modulus
    return crack_unknown_increment(states, modulus, multiplier)
def crack_unknown_modulus(states):
    diffs = [s1 - s0 for s0, s1 in zip(states, states[1:])]
    zeroes = [t2*t0 - t1*t1 for t0, t1, t2 in zip(diffs, diffs[1:], diffs[2:])]
    modulus = abs(reduce(gcd, zeroes))
    return crack_unknown_multiplier(states, modulus)
# 18.AMM
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
from Crypto.Util.number import *
import gmpy2
import time
import random
from tqdm import tqdm
def AMM(o, r, q):
    start = time.time()
    print('\n----------------------------------------------------------------------------------')
    print('Start to run Adleman-Manders-Miller Root Extraction Method')
    print('Try to find one {:#x}th root of {} modulo {}'.format(r, o, q))
    g = GF(q)
    o = g(o)
    p = g(random.randint(1, q))
    while p ^ ((q-1) // r) == 1:
        p = g(random.randint(1, q))
    print('[+] Find p:{}'.format(p))
    t = 0
    s = q - 1
    while s % r == 0:
        t += 1
        s = s // r
    print('[+] Find s:{}, t:{}'.format(s, t))
    k = 1
    while (k * s + 1) % r != 0:
        k += 1
    alp = (k * s + 1) // r
    print('[+] Find alp:{}'.format(alp))
    a = p ^ (r**(t-1) * s)
    b = o ^ (r*alp - 1)
    c = p ^ s
    h = 1
    for i in range(1, t):
        d = b ^ (r^(t-1-i))
        if d == 1:
            j = 0
        else:
            print('[+] Calculating DLP...')
            j = - discrete_log(d, a)
            print('[+] Finish DLP...')
        b = b * (c^r)^j
        h = h * c^j
        c = c^r
    result = o^alp * h
    end = time.time()
    print("Finished in {} seconds.".format(end - start))
    print('Find one solution: {}'.format(result))
    return result

def onemod(p,r): 
    t=p-2 
    while pow(t,(p-1) // r,p)==1: 
        t -= 1 
    return pow(t,(p-1) // r,p) 
 
def solution(p,root,e):  
    g = onemod(p,e) 
    may = set() 
    for i in range(e): 
        may.add(root * pow(g,i,p)%p) 
    return may

def find_roots(c,e,p):
    mp = AMM(c,e,p)
    mps = solution(p,mp,e)
    return mps
# 19.known_phi
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
from math import gcd
from math import isqrt
from random import randrange

from sage.all import is_prime


def factorize(N, phi):
    """
    Recovers the prime factors from a modulus if Euler's totient is known.
    This method only works for a modulus consisting of 2 primes!
    :param N: the modulus
    :param phi: Euler's totient, the order of the multiplicative group modulo N
    :return: a tuple containing the prime factors, or None if the factors were not found
    """
    s = N + 1 - phi
    d = s ** 2 - 4 * N
    p = int(s - isqrt(d)) // 2
    q = int(s + isqrt(d)) // 2
    return p, q


def factorize_multi_prime(N, phi):
    """
    Recovers the prime factors from a modulus if Euler's totient is known.
    This method works for a modulus consisting of any number of primes, but is considerably be slower than factorize.
    More information: Hinek M. J., Low M. K., Teske E., "On Some Attacks on Multi-prime RSA" (Section 3)
    :param N: the modulus
    :param phi: Euler's totient, the order of the multiplicative group modulo N
    :return: a tuple containing the prime factors
    """
    prime_factors = set()
    factors = [N]
    while len(factors) > 0:
        # Element to factorize.
        N = factors[0]

        w = randrange(2, N - 1)
        i = 1
        while phi % (2 ** i) == 0:
            sqrt_1 = pow(w, phi // (2 ** i), N)
            if sqrt_1 > 1 and sqrt_1 != N - 1:
                # We can remove the element to factorize now, because we have a factorization.
                factors = factors[1:]

                p = gcd(N, sqrt_1 + 1)
                q = N // p

                if is_prime(p):
                    prime_factors.add(p)
                elif p > 1:
                    factors.append(p)

                if is_prime(q):
                    prime_factors.add(q)
                elif q > 1:
                    factors.append(q)

                # Continue in the outer loop
                break

            i += 1

    return tuple(prime_factors)
編集日 閲覧数 回

*~(￣▽￣)~[お茶]を一杯ください

tsuppari Alipay

Alipay

tsuppari PayPal

PayPal