SUSCTF 2022 Crypto复现

InverseProblem

$Problem$

import numpy as np
from secret import flag

def gravity(n,d=0.25):
    A=np.zeros([n,n])
    for i in range(n):
        for j in range(n):
            A[i,j]=d/n*(d**2+((i-j)/n)**2)**(-1.5)
    return A

n=len(flag)
A=gravity(n)
x=np.array(list(flag))
b=A@x
np.savetxt('b.txt',b)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

$Analysis$

简答来说，题目中gravity()函数生成了一个实数矩阵$AA$，并且大小为len(flag) x len(flag)；由flag以单个字符构成的1 x len(flag)大小的矩阵$XX$；由矩阵乘法$BB=AAXX$，生成结果矩阵$bb$；已知量有实数矩阵$AA$，以及结果矩阵$BB$

我们已知矩阵乘法中的一个乘数矩阵以及结果矩阵，那么可以计算$XX=AABB$，得到flag矩阵$XX$；但是这里我们不能直接这样左乘结果矩阵得$BB$，因为$AA$是实数矩阵，并且几乎是矩阵中的元素几乎是无限小数；那么就涉及到一个问题，浮点数的精度达到了17位（Python默认），到哪一位小数才是实际对矩阵乘法的结果造成了影响的呢？

这里我们不得而知，涉及到85*85个数的相乘问题（len(flag) = 85）过于复杂；但是我们可以先用数学符号表达这种有不确定影响的情况（将不确定部分用向量$s$表示），这实际上就是LWE问题的思路
$$BB+s=AAXX$$
LWE问题定义

随机选取一个矩阵$AA\in {\mathbb{Z}}_q^{m\times n}$，一个随机向量$s\in {\mathbb{Z}}_q^n$，和一个随机的噪音$e\in {ϵ}^m$（也叫做误差向量）

系统输出$g_{AA}(s,e)=AAs+e(\mathrm{mod}q)$

问题在于，给定矩阵$AA$，以及输出$g_{AA}(s,e)$；还原$s$

对比发现我们这里的问题中设置的向量$s$实际上就相当于误差向量$e$

对于LWE问题我们可以将其转换为SVP问题求解，通过构造基向量来将目标向量$s$作为其格的最短向量

可以构造矩阵 ( A 0 − B 1 )

(AA−BB​01​)作为基向量（有多种构造方式以及解决方法）
( X 1 × 85 1 ) ( A 85 × 85 0 85 × 1 − B 1 × 85 1 ) = ( s 1 × 85 1 )

(XX1×851)" role="presentation" style="position: relative;">(XX1×851)$$\left(\begin{array}{c}X{X}_{1\times 85}1\end{array}\right)$$

(AA85×85085×1−BB1×851)" role="presentation" style="position: relative;">(AA85×85−BB1×85085×11)$$\left(\begin{array}{cc}A{A}_{85\times 85}& 0_{85\times 1}\\ -B{B}_{1\times 85}& 1\end{array}\right)$$

=

(s1×851)" role="presentation" style="position: relative;">(s1×851)$$\left(\begin{array}{c}{s}_{1\times 85}1\end{array}\right)$$

(XX1×85​1​)(AA85×85​−BB1×85​​085×1​1​)=(s1×85​1​)
对 ( A 0 − B 1 )

(AA0−BB1)" role="presentation" style="position: relative;">(AA−BB01)$$\left(\begin{array}{cc}AA& 0\\ -BB& 1\end{array}\right)$$

(AA−BB​01​)进行LLL算法得到最短向量$ $；那么对于上面这个矩阵乘法，我们可以
( X 1 × 85 1 ) = ( s 1 × 85 1 ) ( A 85 × 85 0 85 × 1 − B 1 × 85 1 ) − 1

(XX1×851)" role="presentation" style="position: relative;">(XX1×851)$$\left(\begin{array}{c}X{X}_{1\times 85}1\end{array}\right)$$

=

(s1×851)" role="presentation" style="position: relative;">(s1×851)$$\left(\begin{array}{c}{s}_{1\times 85}1\end{array}\right)$$

(AA85×85085×1−BB1×851)" role="presentation" style="position: relative;">(AA85×85−BB1×85085×11)$$\left(\begin{array}{cc}A{A}_{85\times 85}& 0_{85\times 1}\\ -B{B}_{1\times 85}& 1\end{array}\right)$$

^{-1} (XX1×85​1​)=(s1×85​1​)(AA85×85​−BB1×85​​085×1​1​)−1
得到flag矩阵$XX$；而关于为什么这样构造的基向量就可以使目标向量是其格的最短向量，应该是构造类似形式的矩阵（基向量）就可以满足条件（论文出处是这样构造的类似形式的矩阵）

$Solution$

# type:ignore
import numpy as np

def gravity(n,d=0.25):
    A=np.zeros([n,n])
    for i in range(n):
        for j in range(n):
            A[i,j]=d/n*(d**2+((i-j)/n)**2)**(-1.5)
    return A

enc = []
with open("C:\\Users\\Menglin\\Desktop\\b.txt","rb") as f:
    for line in f.readlines():
        enc.append(float(line.strip().decode()))

b = enc
n = 85
multiple = 10 ^ 20
A = gravity(n)

A = [[int(j * multiple) for j in i] for i in A]
b = [int(i * (-1) * multiple) for i in b]
M = [A[i] + [0] for i in range(n)]
M.append(b + [1])
M = Matrix(ZZ, n + 1, n + 1, M)
ans = M.LLL()[0]
# print(ans)
flag = M.solve_left(ans)
# print(bytes(flag))
print(bytes(flag[:-1]).decode())
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

SpecialCurve3

$Problem$

from Crypto.Util.number import *
from secret import flag,getMyPrime
import hashlib
import random

class SpecialCurve:
    def __init__(self,p,a,b):
        self.p=p
        self.a=a
        self.b=b

    def __str__(self):
        return f'SpecialCurve({self.p},{self.a},{self.b})'

    def add(self,P1,P2):
        x1,y1=P1
        x2,y2=P2
        if x1==0:
            return P2
        elif x2==0:
            return P1
        elif x1==x2 and (y1+y2)%self.p==0:
            return (0,0)
        if P1==P2:
            t=(2*self.a*x1-self.b)*inverse(2*y1,self.p)%self.p
        else:
            t=(y2-y1)*inverse(x2-x1,self.p)%self.p
        x3=self.b*inverse(self.a-t**2,self.p)%self.p
        y3=x3*t%self.p
        return (x3,y3)

    def mul(self,P,k):
        assert k>=0
        Q=(0,0)
        while k>0:
            if k%2:
                k-=1
                Q=self.add(P,Q)
            else:
                k//=2
                P=self.add(P,P)
        return Q

def problem(size,k):
    p=getMyPrime(size)
    x=random.randint(1,p-1)
    y=random.randint(1,p-1)
    e=random.randint(1,p-1)
    a=k*random.randint(1,p-1)**2%p
    b=(a*x**2-y**2)*inverse(x,p)%p
    curve=SpecialCurve(p,a,b)
    G=(x,y)
    Q=curve.mul(G,e)
    print(f'curve={curve}')
    print(f'G={G}')
    print(f'Q={Q}')
    return e

e1=problem(128,1)
e2=problem(256,0)
e3=problem(512,-1)
enc=bytes_to_long(hashlib.sha512(b'%d-%d-%d'%(e1,e2,e3)).digest())^bytes_to_long(flag.encode())
print(f'enc={enc}')
'''
curve=SpecialCurve(233083587295210134948821000868826832947,73126617271517175643081276880688551524,88798574825442191055315385745016140538)
G=(183831340067417420551177442269962013567, 99817328357051895244693615825466756115)
Q=(166671516040968894138381957537903638362, 111895361471674668502480740000666908829)
curve=SpecialCurve(191068609532021291665270648892101370598912795286064024735411416824693692132923,0,58972296113624136043935650439499285317465012097982529049067402580914449774185)
G=(91006613905368145804676933482275735904909223655198185414549961004950981863863, 96989919722797171541882834089135074413922451043302800296198062675754293402989)
Q=(13504049588679281286169164714588439287464466303764421302084687307396426249546, 110661224324697604640962229701359894201176516005657224773855350780007949687952)
curve=SpecialCurve(52373730653143623993722188411805072409768054271090317191163373082830382186155222057388907031638565243831629283127812681929449631957644692314271061305360051,28655236915186704327844312279364325861102737672471191366040478446302230316126579253163690638394777612892597409996413924040027276002261574013341150279408716,42416029226399083779760024372262489355327595236815424404537477696856946194575702884812426801334149232783155054432357826688204061261064100317825443760789993)
G=(15928930551986151950313548861530582114536854007449249930339281771205424453985946290830967245733880747219865184207937142979512907006835750179101295088805979, 29726385672383966862722624018664799344530038744596171136235079529609085682764414035677068447708040589338778102975312549905710028842378574272316925268724240)
Q=(38121552296651560305666865284721153617113944344833289618523344614838728589487183141203437711082603199613749216407692351802119887009907921660398772094998382, 26933444836972639216676645467487306576059428042654421228626400416790420281717654664520663525738892984862698457685902674487454159311739553538883303065780163)
enc=4161358072766336252252471282975567407131586510079023869994510082082055094259455767245295677764252219353961906640516887754903722158044643700643524839069337
'''
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

$Analysis$

构造了三条曲线，而曲线方程是
$$y^2\equiv a{x}^2-bx(\mathrm{mod}p)$$
可以看作圆锥曲线（而不是椭圆曲线）

怎么看出来的呢？看类SpecialCurve的点加法运算关键代码

        if P1==P2:
            t=(2*self.a*x1-self.b)*inverse(2*y1,self.p)%self.p
1
2

而对于整个类里面定义的计算是类似有限域椭圆曲线计算关系的，（关于斜率$t$的单参数曲线，其加法原理是过原点作两点割线（两点相同为切线）的平行线与曲线相交得到新的一点，满足$t_3=\frac{t_1{t}_2+a}{t_1+t_2}\mathrm{mod}p$）

所以对照计算关系可知，在这里P1==P2条件下的两点相加，应该是求切线；而t就是切线方程的斜率，斜率的计算可以看作是对y求导，那么反过来对导函数积分
$$\begin{gathered} 2y\cdot y^{\prime }\equiv 2ax-b(\mathrm{mod}p) \\ \Rightarrow y^2\equiv a{x}^2-bx(\mathrm{mod}p) \end{gathered}$$
题目在这种圆锥曲线加密的情况下，设置了不同的方程参数，分别是$a\equiv r^2(\mathrm{mod}p)$，也就是说$a$是$p$的二次剩余；$a=0$，圆锥曲线方程变为$y^2\equiv -bx(\mathrm{mod}p)$；$a\equiv -r^2(\mathrm{mod}p)$，非二次剩余

那么可以根据计算斜率的公式$t_3=\frac{t_1{t}_2+a}{t_1+t_2}\mathrm{mod}p$；将已知初始点$G$以及终点$Q$的问题，进一步细化为为已知这两点所在切线的斜率，求$G\cdot e=Q$中的$e$

根据论文 SUSCTF2022_official_wp/圆锥曲线公钥密码算法的参数选择_徐旭东.pdf at main · susers/SUSCTF2022_official_wp (github.com) 所写的将圆锥曲线群变为普通乘法群的方法，我们可以将ECDLP问题变为DLP问题，进而使用求解DLP问题的方法快速求解第一种情况；第二种情况中，对照$t_3\equiv \frac{t_1{t}_2+a}{t_1+t_2}\mathrm{mod}p$，发现斜率可以写作$\frac{1}{t_3}\equiv \frac{1}{t_1}+\frac{1}{t_2}(\mathrm{mod}p)$，相当于加法群，也就是说可以将$G\cdot e=Q$表示为$\frac{1}{t_Q}\equiv \frac{1}{t_G}\cdot e(\mathrm{mod}p)$，那么直接求逆元即可；第三种情况，所在曲线参数没有问题，但是模数$p$对于$p+1$是光滑数，那么应用Pohlig_hellman算法求解ECDLP问题的方法即可

$Solution$

# type:ignore
from Crypto.Util.number import *
from sympy.polys.galoistools import gf_crt
from sympy.polys.domains import ZZ
from tqdm import tqdm
import requests
import hashlib

class SpecialCurve:
    def __init__(self,p,a,b):
        self.p=p
        self.a=a
        self.b=b

    def __str__(self):
        return f'SpecialCurve({self.p},{self.a},{self.b})'

    def add(self,P1,P2):
        x1,y1=P1
        x2,y2=P2
        if x1==0:
            return P2
        elif x2==0:
            return P1
        elif x1==x2 and (y1+y2)%self.p==0:
            return (0,0)
        if P1==P2:
            t=(2*self.a*x1-self.b)*inverse(2*y1,self.p)%self.p
        else:
            t=(y2-y1)*inverse(x2-x1,self.p)%self.p
        x3=self.b*inverse(self.a-t**2,self.p)%self.p
        y3=x3*t%self.p
        return (x3,y3)

    def mul(self,P,k):
        assert k>=0
        Q=(0,0)
        while k>0:
            if k%2:
                k-=1
                Q=self.add(P,Q)
            else:
                k//=2
                P=self.add(P,P)
        return Q


curve1 = SpecialCurve(233083587295210134948821000868826832947,73126617271517175643081276880688551524,88798574825442191055315385745016140538)
G = (183831340067417420551177442269962013567, 99817328357051895244693615825466756115)
Q = (166671516040968894138381957537903638362, 111895361471674668502480740000666908829)
p = curve1.p
a = curve1.a
theta = sqrt(GF(p)(a))
t_G = GF(p)(G[1]) // GF(p)(G[0])
t_Q = GF(p)(Q[1]) // GF(p)(Q[0])
func_G = (t_G + theta) // (t_G - theta)
func_Q = (t_Q + theta) // (t_Q - theta)
e1 = func_Q.log(func_G)
print(e1)
e1 = 184572164865068633286768057743716588370



curve2 = SpecialCurve(191068609532021291665270648892101370598912795286064024735411416824693692132923,0,58972296113624136043935650439499285317465012097982529049067402580914449774185)
G = (91006613905368145804676933482275735904909223655198185414549961004950981863863, 96989919722797171541882834089135074413922451043302800296198062675754293402989)
Q = (13504049588679281286169164714588439287464466303764421302084687307396426249546, 110661224324697604640962229701359894201176516005657224773855350780007949687952)
p = curve2.p
a = curve2.a
theta = sqrt(GF(p)(a))
t_G = GF(p)(G[1]) // GF(p)(G[0])
t_Q = GF(p)(Q[1]) // GF(p)(Q[0])
e2 = t_G // t_Q
print(e2)
e2 = 131789829046710687154053378348742202935151384644040019239219239301007568911745



curve3 = SpecialCurve(52373730653143623993722188411805072409768054271090317191163373082830382186155222057388907031638565243831629283127812681929449631957644692314271061305360051,28655236915186704327844312279364325861102737672471191366040478446302230316126579253163690638394777612892597409996413924040027276002261574013341150279408716,42416029226399083779760024372262489355327595236815424404537477696856946194575702884812426801334149232783155054432357826688204061261064100317825443760789993)
G = (15928930551986151950313548861530582114536854007449249930339281771205424453985946290830967245733880747219865184207937142979512907006835750179101295088805979, 29726385672383966862722624018664799344530038744596171136235079529609085682764414035677068447708040589338778102975312549905710028842378574272316925268724240)
Q = (38121552296651560305666865284721153617113944344833289618523344614838728589487183141203437711082603199613749216407692351802119887009907921660398772094998382, 26933444836972639216676645467487306576059428042654421228626400416790420281717654664520663525738892984862698457685902674487454159311739553538883303065780163)

def get_factor(n_to_fac: int) -> list:
    response = requests.get(r'http://factordb.com/api', params={"query": str(n_to_fac)})
    facs = []
    for one in response.json().get("factors"):
        facs += [int(one[0])] * one[1]
    facs = [[i, facs.count(i)] for i in set(facs)]
    facs = [facs[i][0] ** facs[i][1] for i in range(len(facs))]
    return facs

def Pohlig_Hellman(G, Q, P, curve):
    p_ = P + 1 # 哪个数是光滑数，就赋值什么数；注意最后结果也要有相应变化
    factors = get_factor(p_)

    INF = (0,0)   #无穷远点
    dlogs = []
    for i in tqdm(factors):
        Now = INF
        tmpG = curve.mul(G, (p_) // i)
        tmpQ = curve.mul(Q, (p_) // i)
        for dlog in range(i):
            Now = curve.add(Now, tmpG)
            if Now == tmpQ:
                dlogs.append(dlog)
                break
    return gf_crt(dlogs, factors, ZZ)

e3 = Pohlig_Hellman(G, Q, curve3.p, curve3) + 1
print(e3)

e3 = 23331486889781766099145299968747599730779731613118514070077298627895623872695507249173953050022392729611030101946661150932813447054695843306184318795467215

enc = 4161358072766336252252471282975567407131586510079023869994510082082055094259455767245295677764252219353961906640516887754903722158044643700643524839069337
flag = bytes_to_long(hashlib.sha512(b'%d-%d-%d'%(e1,e2,e3)).digest()) ^^ enc
print(long_to_bytes(flag))
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

Large case

(10条消息) SUSCTF_Crypto_large case_复现_M3ng@L的博客-CSDN博客

Paper_Tiger

流密码的相关操作，使用的BM算法（Berlekamp-Massey）来解决最后的问题；（关于具体实现原理和过程我没看懂）

详情可以见$Reference$中的参考文章

$Reference$

SUSCTF2022_official_wp/crypto/SpecialCurve3 at main · susers/SUSCTF2022_official_wp (github.com)

SUSCTF Writeup by X1cT34m | 小绿草信息安全实验室 (njupt.edu.cn)

SUSCTF-Crypto方向所有赛题Love_DengFeng’s Blog(love-deng-feng.top)

SUSCTF 2022 - hash_hash - 博客园 (cnblogs.com)

lfsr-berlekamp-massey/berlekamp_massey.py at master · thewhiteninja/lfsr-berlekamp-massey (github.com)