Crypto趣题-曲线

该文章主要记录一些曲线相关的趣题

EdRSA

题目来源：暂不明确

题目：

#sagemath
from Crypto.Util.number import *
from secrets import flag

def add(P, Q):
    (x1, y1) = P
    (x2, y2) = Q

    x3 = (x1*y2 + y1*x2) * inverse(1 + d*x1*x2*y1*y2, p) % p
    y3 = (y1*y2 - a*x1*x2) * inverse(1 - d*x1*x2*y1*y2, p) % p
    return (x3, y3)

def mul(x, P):
    Q = (0, 1)
    while x > 0:
        if x % 2 == 1:
            Q = add(Q, P)
        P = add(P, P)
        x = x >> 1
    return Q

p = 64141017538026690847507665744072764126523219720088055136531450296140542176327
a = 362
d = 7
gx=bytes_to_long(flag)
PR.<y>=PolynomialRing(Zmod(p))
f=(d*gx^2-1)*y^2+(1-a*gx^2)
gy=int(f.roots()[0][0])

assert (a*gx^2+gy^2)%p==(1+d*gx^2*gy^2)%p

G=(gx,gy)
e=0x10001
print("eG = {}".format(mul(e, G)))

#eG = (602246821311345089174443402780402388933602410138142480089649941718527311147, 17625197557740535449294773567986004828160284887369041337984750097736030549853)

简单看一下题目流程，题目定义了一个有限域上一条曲线的点加和点乘操作，其中曲线形式为：

assert (a*gx^2+gy^2)%p==(1+d*gx^2*gy^2)%p

写出表达式如下：

$$ax^2+y^2 = 1+dx^2y^2\quad(mod\;p)$$

搜索一下，发现这是标准型的扭曲爱德华曲线：(Twisted Edwards Curves)

曲线 | Lazzaro (lazzzaro.github.io)

而仔细核对一下点加与点乘，发现都是完全对的上的。因此问题就转化为，已知Edcurve上的一个e倍点，求解该e倍点对应的原点G的横坐标，即为flag。

想一想，如果这是一条常见形式的椭圆曲线，求解方式是什么？步骤如下：

• 用sage中的order()函数求解出该椭圆曲线的阶n
• 求出e关于阶n的逆元，记为t
• 求倍点G=t*(eG)，横坐标即为所求

那么再回头，这个求解过程对于Edcurve肯定也是类似的，不过问题就在于，sage中没有办法直接求出Edcurve这种形式的曲线的阶，因此确定思路：

• 将Edcurve通过换元映射，变换为常见的椭圆曲线的形式
• 求解出对应椭圆曲线的阶，记为s
• 求倍点G’ = s*(eG’)
• 将求解出的G’再变换回Edcurve上得到G，其横坐标即为所求

因此难点就在于如何通过换元进行曲线映射，这里陈述一下换元过程：(以下除法均为有限域上除法，即乘逆元)

第一步：转化为蒙哥马利曲线方程(Montgomery)：

参考：Edwards Curves (wordpress.com)

$$x' = \frac{1+y}{1-y}$$ $$y'=\frac{1+y}{x(1-y)}$$ $$B=\frac{4}{a-d}$$ $$A=\frac{2(a+d)}{a-d}$$

代入到Edcurve的曲线方程之后，曲线会转化为蒙哥马利曲线，其方程形式如下：

$$B(y')^2 = (x')^3+A(x')^2+(x')\quad(mod\;p)$$

第二步：转化为椭圆曲线方程(Weierstrass)：

参考：Montgomery curve - Wikipedia

$$x'' = \frac{3x'+A}{3B}$$ $$y''=\frac{y'}{B}$$ $$a=\frac{3-A^2}{3B^2}$$ $$b=\frac{2A^3-9A}{27B^3}$$

此时蒙哥马利曲线就变成了椭圆曲线方程形式：

$$(y'')^2 = (x'')^3+a(x'')+b\quad(mod\;p)$$

然后求该曲线的阶，从而求解出远原点G，并且重新逆变换回Edcurve，得到的横坐标即为flag。

exp.ipynb：

from Crypto.Util.number import *
p = 64141017538026690847507665744072764126523219720088055136531450296140542176327
a = 362
d = 7
c = 1
e = 0x10001
eG = (602246821311345089174443402780402388933602410138142480089649941718527311147, 17625197557740535449294773567986004828160284887369041337984750097736030549853)

#part2 map to ECC
F = GF(p)
dd = F(d*c^4)
A = F(2) * F(a+dd) / F(a-dd)
B = F(4) / F(a-dd)
a = F(3-A^2) / F(3*B^2)
b = F(2*A^3-9*A) / F(27*B^3)

def edwards_to_ECC(x,y):
    x1 = F(x) / F(c)
    y1 = F(y) / F(c)
    #now curve is a*x^2+y^2 = 1+dd*x^2*y^2

    x2 = F(1+y1) / F(1-y1)
    y2 = F(x2) / F(x1)
    #now curve is By^2 = x^3 + Ax^2 + x

    x3 = (F(3*x2) + F(A)) / F(3*B)
    y3 = F(y2) / F(B)
    #now curve is y^2 = x^3 + ax + b

    return (x3,y3)
 
def ECC_to_edwards(x,y):
    x2 = (F(x) * F(3*B) - F(A)) / F(3)
    y2 = F(y) * F(B)
    #now curve is By^2 = x^3 + Ax^2 + x

    x1 = F(x2) / F(y2)
    y1 = F(1) - (F(2) / F(x2+1))
    #now curve is a*x^2+y^2 = 1+dd*x^2*y^2

    x_ = F(x1) * F(c)
    y_ = F(y1) * F(c)
    #now curve is a*x^2+y^2 = c^2(1+d*x^2*y^2)
    
    return (x_,y_)
 
E = EllipticCurve(GF(p), [a, b])
order = E.order()
eG = E(edwards_to_ECC(eG[0],eG[1]))
t = inverse(e,order)
G = t*eG
G = ECC_to_edwards(G[0],G[1])
print(long_to_bytes(int(G[0])))

flag：

DASCTF{y0u_kn0w_edcurv3_w3LL!!}

EC_Party-I

题目来源：“华为杯”第二届中国研究生网络安全创新大赛

题目：

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


assert flag[:5]==b'flag{' and flag[-1:]==b'}'
flag = flag[5:-1]
m = bytes_to_long(flag)


def rabin(m):
    m = m+os.urandom(32)
    p = getPrime(384)
    q = getPrime(384)
    Fp = GF(p)
    Fq = GF(q)
    n = p*q
    e = 2
    a = random.randint(0, p-1)
    b = random.randint(0, p-1)

    Ep = EllipticCurve(Zmod(p), [a, b])
    Eq = EllipticCurve(Zmod(q), [a, b])
    En = EllipticCurve(Zmod(n), [a, b])
    ord_p = Ep.order()
    ord_q = Eq.order()

    xm = bytes_to_long(m)
    while True:
        try:
            Gp = Ep.lift_x(Fp(xm))
            Gq = Eq.lift_x(Fq(xm))
            ym = crt([int(Gp.xy()[1]),int(Gq.xy()[1])],[p,q])
            break
        except :
            xm += 1
            continue

    M = En((xm,ym))
    C = e*M
    pk = [a, b, n, C]
    leak = ord_p*ord_q
    return pk, leak


print(rabin(flag))
"""
[138681122158674534796479818810828100269024674330030901179877002756402543027343312824423418859769980312713625658733, 4989541340743108588577899263469059346332852532421276369038720203527706762720292559751463880310075002363945271507040, 762981334990685089884160169295988791471426441106522959345412318178660817286272606245181160960267776171409174142433857335352402619564485470678152764621235882232914864951345067231483720755544188962798600739631026707678945887174897543, (19591102741441427006422487362547101973286873135330241799412389205281057650306427438686318050682578531286702107543065985988634367524715153650482199099194389191525898366546842016339136884277515665890331906261550080128989942048438965, 728465071542637655949094554469510039681717865811604984652385614821789556549826602178972137405550902004858456181137844771163710123158955524137202319902378503104952106036911634918189377295743976966073577013775200078470659428344462772), 762981334990685089884160169295988791471426441106522959345445792076415993922016249232021560266153453470937452118572318136597282436269660557904217923887981072203978473274822142705255987334355747997513083011853917049784914749699536828]
"""

梳理题目加密流程：

• 题目将flag转化成大整数后，转化为基点M的横坐标
• 生成两个大素数p、q，n=p*q
• 生成随机数a、b，并以此生成三条椭圆曲线Ep、Eq、En
• 求出基点M在Ep、Eq上的倍点，并用中国剩余定理求出组合后的纵坐标
• 给出泄露信息leak=order(Ep)*order(Eq)

那么显然如果能获得n的分解，题目就没有难度了，而泄漏的leak就是n的分解的重要依据。其具体原理可以参考：(因为我其实也并没有理解清楚)

My-CTF-Challenges/HITCON CTF 2022/Chimera/README.md at master · maple3142/My-CTF-Challenges (github.com)

简单来说，就是先获得leak的因式分解：

由于leak=order(Ep)*order(Eq)，因此leak乘任何点都应该是En、Ep、Eq的共同O点(无穷远点)，但是将leak挨个除以其因子，再与C点相乘，就可能会产生倍点在Ep上，而不在Eq上的情况，这会使求解关于n的逆元不存在，sage便会在这时抛出一个报错，而在报错信息中就能看到kp以及n的分解：(这是我的理解，如果有不对的地方欢迎师傅指出)

这个数字就是kp，与n求gcd即可得到p，也可以在最后一行报错信息中直接看到n的分解。

其实这种分解方式就是Lenstra elliptic-curve factorization.的核心原理，但我没有完全理解。

求解出p、q后，就顺势获得了两条曲线，接下来就是如何由曲线上的倍点(2M)求解出原点的问题，一般有两类解法：

1、如果2与曲线阶互素，则可以直接求解2的逆元，将倍点乘上逆元即得原点

2、如果不互素，则可以联立椭圆曲线本身方程及倍点方程，在有限域下求根

而在本题中，2与Eq的阶互素，因此采用第一种解法；与Ep的阶不互素，因此采用第二种解法，第二种解法联立方程过程如下：(记M为(x1,y1),2M为(x2,y2))

$$y^2=x^3+ax+b\quad(mod\;p)$$ $$x_2=k^2−2x_1\quad(mod\;p)$$ $$k = \frac{3x_1^2 + a}{2y_1}\quad(mod\;p)$$

联立上述三式可得方程：

$$\frac{(3x_1^2+a)^2}{(2y_1)^2} - 2x_1 - x_2 = 0\quad(mod\;p)$$

即：

$$(3x_1^2+a)^2 - 2x_1*4(x_1^3+ax_1+b) - x_2*4(x_1^3+ax_1+b) = 0\quad(mod\;p)$$

此时方程中仅有x1一个未知数，在模p下解方程即可，解完后用中国剩余定理将模p与模q下的解组合即得flag。

exp.ipynb：

from Crypto.Util.number import *

e = 2
a,b,n,C,leak = [138681122158674534796479818810828100269024674330030901179877002756402543027343312824423418859769980312713625658733, 4989541340743108588577899263469059346332852532421276369038720203527706762720292559751463880310075002363945271507040, 762981334990685089884160169295988791471426441106522959345412318178660817286272606245181160960267776171409174142433857335352402619564485470678152764621235882232914864951345067231483720755544188962798600739631026707678945887174897543, (19591102741441427006422487362547101973286873135330241799412389205281057650306427438686318050682578531286702107543065985988634367524715153650482199099194389191525898366546842016339136884277515665890331906261550080128989942048438965, 728465071542637655949094554469510039681717865811604984652385614821789556549826602178972137405550902004858456181137844771163710123158955524137202319902378503104952106036911634918189377295743976966073577013775200078470659428344462772), 762981334990685089884160169295988791471426441106522959345445792076415993922016249232021560266153453470937452118572318136597282436269660557904217923887981072203978473274822142705255987334355747997513083011853917049784914749699536828]
E = EllipticCurve(Zmod(n),[a,b])
C = E(C)

#factordb
#2^2,3^4,13,199,307,647,157259,297617,8452217,411927661365999433,1157516701716180046249,1338688620929080207819,31226697952255326809332037614333,581208663471376553417319728009366348095695079579751839149645355600351572890241761173016580183555305805091712621
leak_fac = [2,3,13,199,307,647,157259,297617,8452217,411927661365999433,1157516701716180046249,1338688620929080207819,31226697952255326809332037614333,581208663471376553417319728009366348095695079579751839149645355600351572890241761173016580183555305805091712621]

'''
for i in leak_fac:
    temp = leak // i * C
'''

kp = 422522588482185975632147929645103216180089543839772868484032620301855259044079430236142637458472152787337779544279109415705783248825091269039840404202119229567311048216047356951966653331710686649176005328509793328313264251738045723
p = GCD(kp,n)
q = n//p

Ep = EllipticCurve(Zmod(p), [a, b])
Eq = EllipticCurve(Zmod(q), [a, b])
ord_p = Ep.order()
ord_q = Eq.order()
#print(ord_p)
#print(ord_q)

dq = inverse(e,ord_q)
Q = dq*Eq(C)
mq = int(Q[0])

#print(Ep(C))
PR.<x> = PolynomialRing(GF(p))
f = (3*(x**2)+a)**2 - 2*x*(4*(x**3+a*x+b)) - int(C[0])*4*(x**3+a*x+b)
res = f.roots()
#print(res)

n = [p,q]
if(res):
    for i in res:
        c = [int(i[0]),mq]
        m = crt(c,n)
        print(long_to_bytes(m))

flag：

flag{3crab1n_s0unds_go0d}

strange curve

题目来源：巅峰极客 2022

题目：

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

def add(P,Q):
    (x1,y1)=P
    (x2,y2)=Q


    x3=(x1+x2)*(1+y1*y2)*invert((1+x1*x2)*(1-y1*y2),p)%p
    y3=(y1+y2)*(1+x1*x2)*invert((1-x1*x2)*(1+y1*y2),p)%p

    return (x3,y3)

def mul(e,P):
    Q=(0,0)
    e=e%p
    while e:
        if e&1:
            Q=add(Q,P)
        P=add(P,P)
        e>>=1
    return Q

def Legendre(a,p):
    return (pow((a%p+p)%p,(p-1)//2,p))%p

def get_ts(p):
    p=p-1
    count=0
    while p%2==0:
        count+=1
        p=p//2
    return count,p

def get_nonre(p):
    a=random.randint(1,p)
    while Legendre(a,p)==1:
        a=random.randint(1,p)
    return a

def amm2(a,p):
    t,s=get_ts(p)
    ta=pow(get_nonre(p),s,p)
    tb=pow(a,s,p)
    h=1
    for i in range(1,t):
        d=pow(tb,2**t-1-i,p)
        if d==1:
            k=0
        else:
            k=1
        tb=(tb*pow(ta,2*k,p))%p
        h=(h*pow(ta,k,p))%p
        ta=pow(ta,2,p)
    return h*pow(a,(s+1)//2,p)%p  

def solve(a,b,c,p):
    tmpa=1
    tmpb=b*inverse(a,p)%p
    tmpc=c*inverse(a,p)%p
    assert Legendre(tmpb**2*inverse(4,p)-tmpc,p)==1
    res1=(amm2(tmpb**2*inverse(4,p)-tmpc,p)-tmpb*inverse(2,p))%p
    res2=(-amm2(tmpb**2*inverse(4,p)-tmpc,p)-tmpb*inverse(2,p))%p
    return (res1,res2)

def lift(x,a,b,p):
    tmp=b*(x**2-1)*inverse(a*x,p)%p
    return solve(1,-tmp,-1,p)[0]

p=9410547699903726871336507117271550134683191140146415131394654141737636910570480327296351841515571767317596027931492843621727002889086193529096531342265353
a=54733430689690725746438325219044741824500093621550218736194675295708808435509
b=75237024593957256761258687646797952793573177095902495908321724558796076392871
x=bytes_to_long(flag)

while True:
    try:
        y=lift(x,a,b,p)
        break
    except:
        x+=1
        continue

assert a*x*(y**2-1)%p==b*y*(x**2-1)%p

P=(x,y)
e=65537

eP=mul(e,P)
print(f"P = {P}")
print(f"eP = {eP}")
'''
P = (56006392793427940134514899557008545913996191831278248640996846111183757392968770895731003245209281149, 5533217632352976155681815016236825302418119286774481415122941272968513081846849158651480192550482691343283818244963282636939305751909505213138032238524899)
eP = (mpz(8694229840573103722999959579565187489450818138005222030156495740841851804943200684116883831426548909867463656993852596745698999492932194245562062558787005), mpz(9279986963919197374405152604360936066932975197577643570458423456304679111057526702737279809805694360981565554506626018364382736924914907001214909905449002))
'''

非预期解很容易，flag就是P的横坐标，最多再需要爆破一下就行，这里主要讲一下预期解。

我认为题目预期应该是不给P，而只给eP的，因此我就以只有eP这个条件开始解题，分析过程如下：

首先，前面的很多函数先不看，先关注题目给的曲线方程：

$$ax(y^2-1) \equiv by(x^2-1) \quad(mod\;p)$$

而就在前几天的2023 DASCTF CBCTF中，出了一道huff曲线题，而huff曲线的一般形式为：

$$x(ay^2-1) \equiv y(bx^2-1) \quad(mod\;p)$$

可以发现其实很像，我们只需要做以下映射就可以把题目曲线变成一个标准huff曲线：

$$x' = ax$$ $$y' = by$$

那么题目曲线就变成了：

$$x(\frac{y^2}{b^2}-1) \equiv y(\frac{x^2}{a^2}-1) \quad(mod\;p)$$

可以发现这就是个标准huff曲线：

$$x(a'y^2-1) \equiv y(b'x^2-1) \quad(mod\;p)$$

其中：

$$a' = (b^2)^{-1} \quad(mod\;p)$$ $$b' = (a^2)^{-1} \quad(mod\;p)$$

而huff曲线又可以通过如下方式映射为一条Weiestrass Curve，也就是常见的椭圆曲线：

$$(x,y)→(\frac{b'x-a'y}{y-x},\frac{b'-a'}{y-x})$$

该Weiestrass Curve方程如下：

$$y^2≡x^3+(a'+b')x^2+a'b'x \quad mod \quad p$$

映射为这样的曲线后，sage就可以直接求出阶。而映射由于是双射所以不会改变曲线阶，所以我们其实也就求得了原huff曲线的阶，然后就可以求e关于阶的逆元d，就有：

$$d*eP' = P'$$

然后再将第一次映射逆回去就能得到原P点坐标了。

exp：

from Crypto.Util.number import *

p=9410547699903726871336507117271550134683191140146415131394654141737636910570480327296351841515571767317596027931492843621727002889086193529096531342265353
a=54733430689690725746438325219044741824500093621550218736194675295708808435509
b=75237024593957256761258687646797952793573177095902495908321724558796076392871
e = 65537
eP = (8694229840573103722999959579565187489450818138005222030156495740841851804943200684116883831426548909867463656993852596745698999492932194245562062558787005,9279986963919197374405152604360936066932975197577643570458423456304679111057526702737279809805694360981565554506626018364382736924914907001214909905449002)

def mapping(point):
    x = point[0]
    y = point[1]
    Ex = (b_*x-a_*y) * inverse(y-x,p) % p
    Ey = (b_-a_) * inverse(y-x,p) % p
    return (Ex,Ey)

a_ = inverse(b**2,p)
b_ = inverse(a**2,p)
E = EllipticCurve(GF(p),[0,a_+b_,0,a_*b_,0])
#print(E.order())

eP_ = (eP[0]*a%p,eP[1]*b%p)
order = 9410547699903726871336507117271550134683191140146415131394654141737636910570514004897728229958723858012338384995335419723570802793276851855535834618146832
d = inverse(e,order)

class CB_curve:
    def __init__(self):
        self.p = p
        self.a = a_
        self.b = b_

    def add(self, P, Q):
        if P == -1:
            return Q
        (x1, y1) = P
        (x2, y2) = Q
        x3 =  (x1+x2)*(1+self.a*y1*y2)*inverse((1+self.b*x1*x2)*(1-self.a*y1*y2),self.p)% self.p
        y3 =  (y1+y2)*(1+self.b*x1*x2)*inverse((1-self.b*x1*x2)*(1+self.a*y1*y2),self.p)% self.p
        return (x3, y3)

    def mul(self, x, P):
        Q = -1
        while x > 0:
            if x & 1:
                Q = self.add(Q, P)
            P = self.add(P, P)
            x = x >> 1
        return Q
    
    def negG(self,G):
        return self.mul(order-1,G)

curve = CB_curve()
P_ = curve.mul(d,eP_)
print(long_to_bytes(int(inverse(a,p)*P_[0] % p)))

#flag{b7f209df-1284-4bdf-b030-28197483c47b}

CB_curve其实就是huff曲线的实现，偷懒用了出题师傅的(ᕑᗢᓫ∗)˒

SMM

题目来源：2023福建省数据安全竞赛

题目：

from Crypto.Util.number import long_to_bytes, bytes_to_long
from Crypto.Util.Padding import pad
from Crypto.Cipher import AES
from random import getrandbits
from secret import flag
import hashlib
import base64

msg = pad(flag, 16)
 
ecc_table1 = {
    'n': 'FFFFFFFEFFFFFFFFFFFFFFFFFFFFFFFF7203DF6B21C6052B53BBF40939D54123',
    'p': 'FFFFFFFEFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF00000000FFFFFFFFFFFFFFFF',
    'g': '32c4ae2c1f1981195f9904466a39c9948fe30bbff2660be1715a4589334c74c7'
         'bc3736a2f4f6779c59bdcee36b692153d0a9877cc62a474002df32e52139f0a0',
    'a': 'FFFFFFFEFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF00000000FFFFFFFFFFFFFFFC',
    'b': '28E9FA9E9D9F5E344D5A9E4BCF6509A7F39789F515AB8F92DDBCBD414D940E93',
}

ecc_table2 = {
    'n': '00000000000000000000000000000024c0a1eef669c78f5a60af4b7eece0cec3',
    'p': '00000000000000000000000000000024c0a1eef669c78f5a60af4b7eece0cec3',
    'g': '00000000000000000000000000000001333d7197c5da65ad36ebbb3589634ad6'
         '00000000000000000000000000000007b81a5934e5b39c8d36f449527767209b',
    'a': '00000000000000000000000000000024c0a1eef669c78f5a60af4b7e22c0cec3',
    'b': '00000000000000000000000000000024c0a1eef669c78f5a60af04915ce0cec3',
}
 
class TSM2(object):
    def __init__(self, sk, ecc_table):
        self.ecc_table = ecc_table
        self.n = int(ecc_table['n'], 16)
        self.para_len = len(ecc_table['n'])
        self.ecc_a3 = (int(ecc_table['a'], base=16) +
                       3) % int(ecc_table['p'], base=16)
 
        self.sk = sk
        self.pk = self._kg(self.sk, ecc_table['g'])
 
    def sign(self, data, K):
        e = data
        d = self.sk
        k = K
 
        P1 = self._kg(k, self.ecc_table['g'])
        x = int(P1[0:self.para_len], 16)
        R = ((e + x) % int(self.ecc_table['n'], base=16))
        if R == 0 or R + k == int(self.ecc_table['n'], base=16):
            return None
        d_1 = pow(
            d+1, int(self.ecc_table['n'], base=16) - 2, int(self.ecc_table['n'], base=16))
        S = (d_1*(k + R) - R) % int(self.ecc_table['n'], base=16)
        if S == 0:
            return None
        else:
            return '%064x%064x' % (R, S)
 
    def verify(self, Sign, data):
        r = int(Sign[0:self.para_len], 16)
        s = int(Sign[self.para_len:2 * self.para_len], 16)
        e = int(data.hex(), 16)
        t = (r + s) % self.n
        if t == 0:
            return 0
 
        P1 = self._kg(s, self.ecc_table['g'])
        P2 = self._kg(t, self.pk)
 
        if P1 == P2:
            P1 = '%s%s' % (P1, 1)
            P1 = self._double_point(P1)
        else:
            P1 = '%s%s' % (P1, 1)
            P1 = self._add_point(P1, P2)
            P1 = self._convert_jacb_to_nor(P1)
 
        x = int(P1[0:self.para_len], 16)
        return r == ((e + x) % self.n)
 
    def _kg(self, k, Point):
        if (k % self.n) == 0:
            return '0' * 128
        Point = '%s%s' % (Point, '1')
        mask_str = '8'
        for i in range(self.para_len - 1):
            mask_str += '0'
        mask = int(mask_str, 16)
        Temp = Point
        flag = False
        for n in range(self.para_len * 4):
            if flag:
                Temp = self._double_point(Temp)
            if (k & mask) != 0:
                if flag:
                    Temp = self._add_point(Temp, Point)
                else:
                    flag = True
                    Temp = Point
            k = k << 1
        return self._convert_jacb_to_nor(Temp)
 
    def _double_point(self, Point):
        l = len(Point)
        len_2 = 2 * self.para_len
        if l < self.para_len * 2:
            return None
        else:
            x1 = int(Point[0:self.para_len], 16)
            y1 = int(Point[self.para_len:len_2], 16)
            if l == len_2:
                z1 = 1
            else:
                z1 = int(Point[len_2:], 16)
 
            T6 = (z1 * z1) % int(self.ecc_table['p'], base=16)
            T2 = (y1 * y1) % int(self.ecc_table['p'], base=16)
            T3 = (x1 + T6) % int(self.ecc_table['p'], base=16)
            T4 = (x1 - T6) % int(self.ecc_table['p'], base=16)
            T1 = (T3 * T4) % int(self.ecc_table['p'], base=16)
            T3 = (y1 * z1) % int(self.ecc_table['p'], base=16)
            T4 = (T2 * 8) % int(self.ecc_table['p'], base=16)
            T5 = (x1 * T4) % int(self.ecc_table['p'], base=16)
            T1 = (T1 * 3) % int(self.ecc_table['p'], base=16)
            T6 = (T6 * T6) % int(self.ecc_table['p'], base=16)
            T6 = (self.ecc_a3 * T6) % int(self.ecc_table['p'], base=16)
            T1 = (T1 + T6) % int(self.ecc_table['p'], base=16)
            z3 = (T3 + T3) % int(self.ecc_table['p'], base=16)
            T3 = (T1 * T1) % int(self.ecc_table['p'], base=16)
            T2 = (T2 * T4) % int(self.ecc_table['p'], base=16)
            x3 = (T3 - T5) % int(self.ecc_table['p'], base=16)
 
            if (T5 % 2) == 1:
                T4 = (T5 + ((T5 + int(self.ecc_table['p'], base=16)) >> 1) - T3) % int(
                    self.ecc_table['p'], base=16)
            else:
                T4 = (T5 + (T5 >> 1) - T3) % int(self.ecc_table['p'], base=16)
 
            T1 = (T1 * T4) % int(self.ecc_table['p'], base=16)
            y3 = (T1 - T2) % int(self.ecc_table['p'], base=16)
 
            form = '%%0%dx' % self.para_len
            form = form * 3
            return form % (x3, y3, z3)
 
    def _add_point(self, P1, P2):
        if P1 == '0' * 128:
            return '%s%s' % (P2, '1')
        if P2 == '0' * 128:
            return '%s%s' % (P1, '1')
        len_2 = 2 * self.para_len
        l1 = len(P1)
        l2 = len(P2)
        if (l1 < len_2) or (l2 < len_2):
            return None
        else:
            X1 = int(P1[0:self.para_len], 16)
            Y1 = int(P1[self.para_len:len_2], 16)
            if l1 == len_2:
                Z1 = 1
            else:
                Z1 = int(P1[len_2:], 16)
            x2 = int(P2[0:self.para_len], 16)
            y2 = int(P2[self.para_len:len_2], 16)
 
            T1 = (Z1 * Z1) % int(self.ecc_table['p'], base=16)
            T2 = (y2 * Z1) % int(self.ecc_table['p'], base=16)
            T3 = (x2 * T1) % int(self.ecc_table['p'], base=16)
            T1 = (T1 * T2) % int(self.ecc_table['p'], base=16)
            T2 = (T3 - X1) % int(self.ecc_table['p'], base=16)
            T3 = (T3 + X1) % int(self.ecc_table['p'], base=16)
            T4 = (T2 * T2) % int(self.ecc_table['p'], base=16)
            T1 = (T1 - Y1) % int(self.ecc_table['p'], base=16)
            Z3 = (Z1 * T2) % int(self.ecc_table['p'], base=16)
            T2 = (T2 * T4) % int(self.ecc_table['p'], base=16)
            T3 = (T3 * T4) % int(self.ecc_table['p'], base=16)
            T5 = (T1 * T1) % int(self.ecc_table['p'], base=16)
            T4 = (X1 * T4) % int(self.ecc_table['p'], base=16)
            X3 = (T5 - T3) % int(self.ecc_table['p'], base=16)
            T2 = (Y1 * T2) % int(self.ecc_table['p'], base=16)
            T3 = (T4 - X3) % int(self.ecc_table['p'], base=16)
            T1 = (T1 * T3) % int(self.ecc_table['p'], base=16)
            Y3 = (T1 - T2) % int(self.ecc_table['p'], base=16)
 
            form = '%%0%dx' % self.para_len
            form = form * 3
            return form % (X3, Y3, Z3)
    
    def _convert_jacb_to_nor(self, Point):
        len_2 = 2 * self.para_len
        x = int(Point[0:self.para_len], 16)
        y = int(Point[self.para_len:len_2], 16)
        z = int(Point[len_2:], 16)
        z_inv = pow(
            z, int(self.ecc_table['p'], base=16) - 2, int(self.ecc_table['p'], base=16))
        z_invSquar = (z_inv * z_inv) % int(self.ecc_table['p'], base=16)
        z_invQube = (z_invSquar * z_inv) % int(self.ecc_table['p'], base=16)
        x_new = (x * z_invSquar) % int(self.ecc_table['p'], base=16)
        y_new = (y * z_invQube) % int(self.ecc_table['p'], base=16)
        z_new = (z * z_inv) % int(self.ecc_table['p'], base=16)
        if z_new == 1:
            form = '%%0%dx' % self.para_len
            form = form * 2
            return form % (x_new, y_new)
        else:
            return None

bits = [40, 140]
ecc_table = [ecc_table1, ecc_table2]

for _ in range(2):
    sk = getrandbits(bits[_])
    key = hashlib.md5(str(sk).encode()).digest()
    aes = AES.new(key, AES.MODE_ECB)
    msg = aes.encrypt(msg)
    sk = sk % int(ecc_table[_]['n'], 16)
    sm2 = TSM2(sk, ecc_table[_])
    pk = sm2.pk
    print(pk)

print(msg)
'''
749180747cd167402449aa2ffdfacd87a3d8ff9a3e00670e9501cefc96c712f16ddca6ba23ed4c9bbbd23bb5d4f16627cd30a82f64a75195c4c3337c17041e89
0000000000000000000000000000002030737f5765b0ca638572035d1143a2210000000000000000000000000000001c69be4b4c104fe2c888be4b659620d933
b'\xb6\xe2\xd0\x00\xdb=3\xf2\xc3\x10\x9b\xb7\x04\xb2\t\x92\xed$Z\xf6.\xd8\xdd/\xad\x03M\xaec\xce\xa9Z\x00}\xfd\xe2T\x859\x0b\x1f\xd3\xd4H^\xfa8\xb7'
'''

特别长的题目，不过其实注意到函数：

def _convert_jacb_to_nor(self, Point)

看名字，这个函数的作用应该是把雅可比坐标转换成普通坐标，那么前面的很长的点加、倍乘的实现，估计也就是雅可比坐标下的ECC点运算。

所以其实很长一段代码都可以忽略，其实这个题就是定义了两个ECC，每张ecc_table中，p、a、b就是椭圆曲线方程中的各参数，n是曲线的阶，g是曲线上的一个生成元。

明白了这一点后再看题目主要任务：

bits = [40, 140]
ecc_table = [ecc_table1, ecc_table2]

for _ in range(2):
    sk = getrandbits(bits[_])
    key = hashlib.md5(str(sk).encode()).digest()
    aes = AES.new(key, AES.MODE_ECB)
    msg = aes.encrypt(msg)
    sk = sk % int(ecc_table[_]['n'], 16)
    sm2 = TSM2(sk, ecc_table[_])
    pk = sm2.pk
    print(pk)

print(msg)

可以看到，他其实给了两个曲线上各自生成元的sk倍点坐标，sk分别为40bit和140bit，然后他将sk分别作为AES密钥，连续对flag进行两次加密，最后给出密文。那么我们想要得到flag，也就要解出两个sk然后AES解密。

也就是说，我们实际上是要完成两个不同曲线上的DLP问题。而略作测试就会发现：

• 对于曲线一，由于私钥sk仅有40bit，因此可以考虑bsgs
• 对于曲线二，order=p，因此smart attack解决问题

而实际上bsgs也需要一定的时间，耐心一点就好了。

而最后还有一个小细节需要注意：每一次加密中，sk模了曲线阶n，而第二个曲线阶略小于140比特，因此求出DLP后还需要爆破几位才能得到正确sk。

exp：

from tqdm import *
from Crypto.Util.number import *
from Crypto.Cipher import AES
import hashlib

ecc_table1 = {
    'n': 'FFFFFFFEFFFFFFFFFFFFFFFFFFFFFFFF7203DF6B21C6052B53BBF40939D54123',
    'p': 'FFFFFFFEFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF00000000FFFFFFFFFFFFFFFF',
    'gx': '32c4ae2c1f1981195f9904466a39c9948fe30bbff2660be1715a4589334c74c7',
    'gy': 'bc3736a2f4f6779c59bdcee36b692153d0a9877cc62a474002df32e52139f0a0',
    'a': 'FFFFFFFEFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF00000000FFFFFFFFFFFFFFFC',
    'b': '28E9FA9E9D9F5E344D5A9E4BCF6509A7F39789F515AB8F92DDBCBD414D940E93',
}

ecc_table2 = {
    'n': '00000000000000000000000000000024c0a1eef669c78f5a60af4b7eece0cec3',
    'p': '00000000000000000000000000000024c0a1eef669c78f5a60af4b7eece0cec3',
    'gx': '00000000000000000000000000000001333d7197c5da65ad36ebbb3589634ad6',
    'gy': '00000000000000000000000000000007b81a5934e5b39c8d36f449527767209b',
    'a': '00000000000000000000000000000024c0a1eef669c78f5a60af4b7e22c0cec3',
    'b': '00000000000000000000000000000024c0a1eef669c78f5a60af04915ce0cec3',
}

ecc1 = EllipticCurve(GF(int(ecc_table1['p'],16)),[int(ecc_table1['a'],16),int(ecc_table1['b'],16)])
ecc2 = EllipticCurve(GF(int(ecc_table2['p'],16)),[int(ecc_table2['a'],16),int(ecc_table2['b'],16)])
g1 = ecc1((int(ecc_table1['gx'],16),int(ecc_table1['gy'],16)))
g2 = ecc2((int(ecc_table2['gx'],16),int(ecc_table2['gy'],16)))

k1g1 = ecc1(int("749180747cd167402449aa2ffdfacd87a3d8ff9a3e00670e9501cefc96c712f1",16),int("6ddca6ba23ed4c9bbbd23bb5d4f16627cd30a82f64a75195c4c3337c17041e89",16))
k2g2 = ecc2(int("2030737f5765b0ca638572035d1143a221",16),int("1c69be4b4c104fe2c888be4b659620d933",16))
enc = b'\xb6\xe2\xd0\x00\xdb=3\xf2\xc3\x10\x9b\xb7\x04\xb2\t\x92\xed$Z\xf6.\xd8\xdd/\xad\x03M\xaec\xce\xa9Z\x00}\xfd\xe2T\x859\x0b\x1f\xd3\xd4H^\xfa8\xb7'

#part1 bsgs
'''
dic = {}
for b in trange(2**20):
    dic[b*g1] = b

cc = dic.keys()
for a in trange(464193,2**20):
    temp = k1g1 - a*2^20*g1
    if(temp in cc):
        print(a)
        print(dic[temp])
'''
#d = a*2^20+b
a = 464193
b = 319430
k1 = a*2^20+b
#print(k1*g1 == k1g1)

#part2 smart_attack
def SmartAttack(P,Q,p):
    E = P.curve()
    Eqp = EllipticCurve(Qp(p, 2), [ ZZ(t) + randint(0,p)*p for t in E.a_invariants() ])

    P_Qps = Eqp.lift_x(ZZ(P.xy()[0]), all=True)
    for P_Qp in P_Qps:
        if GF(p)(P_Qp.xy()[1]) == P.xy()[1]:
            break

    Q_Qps = Eqp.lift_x(ZZ(Q.xy()[0]), all=True)
    for Q_Qp in Q_Qps:
        if GF(p)(Q_Qp.xy()[1]) == Q.xy()[1]:
            break

    p_times_P = p*P_Qp
    p_times_Q = p*Q_Qp

    x_P,y_P = p_times_P.xy()
    x_Q,y_Q = p_times_Q.xy()

    phi_P = -(x_P/y_P)
    phi_Q = -(x_Q/y_Q)
    k = phi_Q/phi_P
    return ZZ(k)

k2 = int(SmartAttack(g2, k2g2, ecc2.order()))
#print(k2*g2 == k2g2)

#part3 AES
for i in range(2**7):
    temp = k2 + i*g2.order()
    #print(temp*g2 == k2g2)
    key2 = hashlib.md5(str(temp).encode()).digest()
    aes2 = AES.new(key2, AES.MODE_ECB)
    dec2 = aes2.decrypt(enc)

    key1 = hashlib.md5(str(k1).encode()).digest()
    aes1 = AES.new(key1, AES.MODE_ECB)
    dec1 = aes1.decrypt(dec2)

    if(b"flag" in dec1):
        print(dec1)
        break

#flag{609a091e-6666-b7c9-edc2-c721958409e2}

Rosita

题目来源：巅峰极客 2023

题目：

# Problem by rec, without any sleep at all.
from Crypto.Util.number import bytes_to_long as b2l
from hashlib import sha256
from os import urandom
from secret import p, a, b, flag

ECC = EllipticCurve(GF(p), [a, b])
R, E, C = [ECC.random_point() for _ in range(3)]

pad = lambda m: urandom(8) + m + b'\x00' * (ZZ(p).nbits() // 8 - len(m) - 8 - 1)
out = list()
for i in range(len(flag)):
    m = pad(chr(flag[i]).encode())
    nonce = urandom(16)
    sh = sha256(nonce + m).digest()
    
    Q = b2l(m)*R + b2l(nonce)*E + b2l(sh)*C
    out.append(Q.xy())

with open('out.tuo', 'w') as f:
    f.write(str(out))

题目内容很简短，流程如下：

• 用自定义参数a、b、p生成一条曲线ECC，并在上面随机取了三点R、E、C。曲线参数与REC三点坐标均未知
• 对flag的每一个字节填充，填充方式为：

pad = lambda m: urandom(8) + m + b'\x00' * (ZZ(p).nbits() // 8 - len(m) - 8 - 1)

• 每一次加密生成一个临时密钥nonce，并依据nonce与填充后的m计算sh
• 计算：

$$Q_i = m_i*R + nonce_i*E + sh_i*C$$

• 给出所有Qi坐标，要求还原m

还原ecc参数

第一步肯定是要先找到原曲线的a、b、p参数的，这里我采用消元办法，思路如下：

我们知道椭圆曲线的点满足方程：

$$y^2 = x^3+ax+b \quad(mod\;p)$$

而Qi均为曲线上的点，所以我们取其中前四个点的坐标，就有四组方程：

$$y_1^2 = x_1^3+ax_1+b \quad(mod\;p)$$ $$y_2^2 = x_2^3+ax_2+b \quad(mod\;p)$$ $$y_3^2 = x_3^3+ax_3+b \quad(mod\;p)$$ $$y_4^2 = x_4^3+ax_4+b \quad(mod\;p)$$

要求从四组方程中还原a、b、p，其实跟无参数的LCG很像，主要目的就是消元，求gcd得到p，然后再回到模方程中求解a、b。

消元过程如下，我们取前两组方程：

$$y_1^2 - x_1^3 =ax_1+b \quad(mod\;p)$$ $$y_2^2 - x_2^3= ax_2+b \quad(mod\;p)$$

作差消掉b：

$$y_1^2 - x_1^3 - y_2^2 + x_2^3=ax_1 - ax_2 \quad(mod\;p)$$

把a单独提出来：

$$(y_1^2 - x_1^3 - y_2^2 + x_2^3)(x_1-x_2)^{-1} = a \quad(mod\;p)$$

同理，我们取第二三组方程也能得到：

$$(y_2^2 - x_2^3 - y_3^2 + x_3^3)(x_2-x_3)^{-1} = a \quad(mod\;p)$$

那么作差就有：

$$(y_1^2 - x_1^3 - y_2^2 + x_2^3)(x_1-x_2)^{-1} - (y_2^2 - x_2^3 - y_3^2 + x_3^3)(x_2-x_3)^{-1} = 0 \quad(mod\;p)$$

那么同样再取一组求gcd，并消去小因子就能得到p，得到p过后还原ab就很轻松。

当然，有更直接的梭哈办法就是直接定义在ZZ域上，用四个点的方程求Groebner就可以轻松还原a、b、p参数。

Smart Attack

有了参数后我们就有了原曲线，经过观察与测试可以知道这个曲线有一个很特殊的性质：

$$E.order() = p$$

有这个性质的曲线可以用Smart attack轻松解决离散对数问题，那么接下来的思路就是如何利用DLP来转化问题。

那么首先，我们期望还原的是明文的每个字符mi，想要解决这一点就要从Qi的生成方式入手：

$$Q_i = m_i*R + nonce_i*E + sh_i*C$$

给了多组Qi相关的线性等式，可以想到应该是要造格，但是造格之前我们需要想办法把ECC相关的计算式转到数域上，而由于DLP在这个曲线上能用Smart Attack轻松求解，因此我们可以考虑如下方式，将所以Qi的计算式转化到数域上：

令G为曲线E的生成元，则G的阶等于曲线的阶，在这条曲线上也就等于p。然后我们就可以把Q、R、E、C等点均用生成元的倍点表示如下：

$$Q_i = d_iG$$ $$R = rG,E=eG,C=cG$$

代回到Qi的计算式中，就有：

$$d_iG = m_i*rG + nonce_i*eG + sh_i*cG$$

即：

$$d_iG = (m_i*r + nonce_i*e + sh_i*c)G$$

而由于G的阶就是曲线的阶，因此有：

$$d_i \equiv m_i*r + nonce_i*e + sh_i*c \quad(mod\;order)$$

又因为曲线阶就等于p，所以有：

$$d_i \equiv m_i*r + nonce_i*e + sh_i*c \quad(mod\;p)$$

到这里，我们就完全脱离了曲线，而转化为了一个数域上的问题，也就是知道多组上述等式，如何求解mi。

格

现在我们有多组如下等式：

$$d_i \equiv m_i*r + nonce_i*e + sh_i*c \quad(mod\;p)$$

而注意到m_i的填充方式很奇怪：

pad = lambda m: urandom(8) + m + b'\x00' * (ZZ(p).nbits() // 8 - len(m) - 8 - 1)

他对flag逐个字节加密，因此len(m)是1，代入我们计算出的p可以得到，这种方式相当于是：

pad(m) = urandom(8) + m + b'\x00'*54

而尾缀这些b’\x00’，其实也就相当于把(urandom(8) + m)左移54x8位，其实也就是乘2^432，所以上述等式可以进一步写成：

$$d_i \equiv m_i*2^{432}r + nonce_i*e + sh_i*c \quad(mod\;p)$$

我们把这些等式写成模p下的矩阵乘法的形式，便于我们发现格构造的方法：

$$\left( \begin{matrix} m_1 & nonce_1 &sh_1\\ m_2 & nonce_2 &sh_2\\ & ...&\\ m_{72} & nonce_{72} &sh_{72}\\ m_{73} & nonce_{73} &sh_{73}\\ \end{matrix} \right)_{73*3} \left( \begin{matrix} 2^{432}r\\ e\\ c\\ \end{matrix} \right)_{3*1} = \left( \begin{matrix} d_1\\ d_1\\ ...\\ d_{72}\\ d_{73}\\ \end{matrix} \right)_{73*1}$$

而观察到mi均为9字节的短向量，因此应该可以利用格规约。

做到这里就没思路了，这能怎么造格呢？想了好一会儿也没想到。

然后就去找别的师傅的wp，学习到了一种基于正交格的很妙的解法：

巅峰极客2023 Crypto Rosita-CSDN博客

这个思路就是，我们首先把上述矩阵方程写作：

$$M_{73*3} \left( \begin{matrix} 2^{432}r\\ e\\ c\\ \end{matrix} \right)_{3*1} = D_{73*1} \quad(mod\;p)$$

然后，我们假设有一个kx73的矩阵R满足下式：

$$R_{k*73}M_{73*3} = \left( \begin{matrix} 0&0&0\\ &...&\\ 0&0&0\\ \end{matrix} \right)_{k*3} \quad(mod\;p)$$

这也就是说R与M正交，然后把这个R矩阵同时乘在矩阵方程的两侧：

$$R_{k*73}M_{73*3} \left( \begin{matrix} 2^{432}r\\ e\\ c\\ \end{matrix} \right)_{3*1} = R_{k*73}D_{73*1} \quad(mod\;p)$$

显然就有：

$$R_{k*73}D_{73*1} = \left( \begin{matrix} 0\\ ...\\ 0\\ \end{matrix} \right)_{k*1} \quad(mod\;p)$$

也就是说，如果R与M正交，那么R一定满足以下几点：

• R与D也正交

• M的三列均在R的右核空间中

• M的第一列是R的右核空间中的一个短向量

所以我们就要完成以下的对应目标：

• 找到与M正交的矩阵R
• 对R的右核空间进行规约得到M的第一列

找到R

由于我们没有M，所以要找到与M正交的矩阵的话，只能通过R与D也正交这一关系来找。为此我们构造如下的格来找到R：

$$\left( \begin{matrix} E&D\\ 0&p\\ \end{matrix} \right)$$

其中，E是单位矩阵，D是DLP得到的列向量，p是曲线的阶，这里也正好是曲线的模数。

构造这个格是基于：

$$R_{k*73}D_{73*1} = \left( \begin{matrix} 0\\ ...\\ 0\\ \end{matrix} \right)_{k*1} \quad(mod\;p)$$

那么我们对R中的任意一行：

$$(r_1,r_2,...,r_{73})$$

都有：

$$(r_1,r_2,...,r_{73},t) \left( \begin{matrix} E&D\\ 0&p\\ \end{matrix} \right) = (r_1,r_2,...,r_{73},0)$$

因此我们就能用LLL规约出若干短向量，只要规约出的最后一个数字为0，就可以作为R的其中一行加入到R中。所以我们可以对格的最后一列配上一个大系数，从而保证规约出0，这样我们就能从中取70行，得到一个70x73的矩阵R。

求解右核

我们通过刚才的步骤，找到了一个符合要求的70x73的矩阵R，他满足：

$$R_{70*73} \left( \begin{matrix} m_1 & nonce_1 &sh_1\\ m_2 & nonce_2 &sh_2\\ & ...&\\ m_{72} & nonce_{72} &sh_{72}\\ m_{73} & nonce_{73} &sh_{73}\\ \end{matrix} \right)_{73*3} = \left( \begin{matrix} 0&0&0\\ &...&\\ 0&0&0\\ \end{matrix} \right)_{70*3} \quad(mod\;p)$$

而R的右核空间，指的就是满足Rv=0的所有向量v张成的空间，那么显然，M的三列均在R的右核空间中，并且M的第一列是R的右核空间中的一个短向量。

那么我们就可以按如下方式找到M的第一列：

• 求解R的右核空间
• 对R的右核空间进行规约

找到M的第一列后，依次取最低字节就可以还原flag。

这里可能会有以下几个问题：

• 为什么不取73行？

这是因为，取73行后的R矩阵是满秩的，右核空间是0维，而我们需要一个至少一维的右核进行规约。

• 为什么不取更少的比如60、50行？

这是因为，少取若干行之后，右核空间维数也会对应增大，高斯启发式期望的短向量长度会显著减小，从而导致我们的目标向量规约不出。经测试取67以下的就规约不出了。

• 为什么不取71、72行？

这一点我也暂时不理解，因为既然规约出最后一列为0，就说明该行与D模p下正交，就应该可以作为R的一行才对。但事实上，只有前70行可以作规约，这也就是说，必须是不配大系数也能规约出的正交行，才能作为R中的一行，这我就暂时没有理解到，如果有明白的师傅欢迎与我交流！

exp：

from Crypto.Util.number import *
from tqdm import *

#part1 recover a,b,p
def calc(points):
    (x1,y1),(x2,y2),(x3,y3),(x4,y4) = points
    t1 = (y2**2-y1**2-x2**3+x1**3)
    t2 = (y3**2-y2**2-x3**3+x2**3)
    t3 = (y4**2-y3**2-x4**3+x3**3)
    k1p = t1*(x3-x2) - t2*(x2-x1)
    k2p = t2*(x4-x3) - t3*(x3-x2)
    k3p = t1*(x4-x3) - t3*(x2-x1)
    p = GCD(k1p,k2p)
    p = GCD(p,k3p)

    for i in range(2,1000):
        while(p % i == 0):
            p //= i
    a = inverse(x2-x1,p)*t1 % p
    b = (y1**2-x1**3-a*x1) % p

    return a,b,p

with open(r'E:\题\历届赛题\巅峰极客 2023\Rosita\out.tuo', 'r') as f:
    points = eval(f.read())
a,b,p = calc(points[:4])


#Smart's attack
def SmartAttack(P,Q,p):
    E = P.curve()
    Eqp = EllipticCurve(Qp(p, 2), [ ZZ(t) + randint(0,p)*p for t in E.a_invariants() ])

    P_Qps = Eqp.lift_x(ZZ(P.xy()[0]), all=True)
    for P_Qp in P_Qps:
        if GF(p)(P_Qp.xy()[1]) == P.xy()[1]:
            break

    Q_Qps = Eqp.lift_x(ZZ(Q.xy()[0]), all=True)
    for Q_Qp in Q_Qps:
        if GF(p)(Q_Qp.xy()[1]) == Q.xy()[1]:
            break

    p_times_P = p*P_Qp
    p_times_Q = p*Q_Qp

    x_P,y_P = p_times_P.xy()
    x_Q,y_Q = p_times_Q.xy()

    phi_P = -(x_P/y_P)
    phi_Q = -(x_Q/y_Q)
    k = phi_Q/phi_P
    return ZZ(k)

E = EllipticCurve(Zmod(p),[a,b])
G = E.gens()[0]
d = []
for i in trange(len(points)):
    d.append(SmartAttack(G,E(points[i]),p))


#part3 LLL
L = Matrix(ZZ,74,74)
K = 2^100
for i in range(73):
    L[i,i] = 1
    L[i,-1] = d[i]*K
L[-1,-1] = p*K
res = L.LLL()

R = Matrix(ZZ,70,73)
for i in range(70):
    R[i] = res[i][:-1]

bas = R.right_kernel().basis()
bas = Matrix(ZZ,bas)
res = bas.LLL()[0]

flag = b""
for i in res:
    flag += long_to_bytes(int(i&0xff))
print(flag)

#Congratulations! Your flag is: flag{893d041e-c0a2-3145-5320-cdee7d3c87fb}

ecdh

题目来源：强网拟态 2022

题目：

#!/usr/bin/env python3.9
# -*- coding: utf-8 -*-
import random
import string
import socketserver
import signal
from hashlib import sha256
from os import urandom
from gmpy2 import invert
from secret import FLAG

p = 0xc483230557557e3d94b7407f3355b8a5b26bda29119babcb8d72b5a19e10e113
a = 0xb5ecbb93d4e49fe3c93ad770343ec5ab70131151151fcc830f6c658223c92e1d
b = 0xb0d9363d1828cfa7c2aba27b0d483fe808637adf6e3a0a5bbc6cb53fdd1e4d85

Px = 0xe9592b5211516c197f0fd31cf0e28201ea0bcc67f7356d8732ded234045259e3
Py = 0xe825e2821cc58e97816ca9877b7604b05a9a4bbc03110bc2124be49eb2718c23
P = (Px, Py)
zero = (0, 0)

           
MENU = rb'''1.sign
2.get flag
'''



def add(p1, p2):
    if p1 == zero:
        return p2
    if p2 == zero:
        return p1
    p1x, p1y = p1
    p2x, p2y = p2
    if p1x == p2x and (p1y != p2y or p1y == 0):
        return zero
    if p1x == p2x:
        tmp = (3 * p1x * p1x + a) * invert(2 * p1y, p) % p
    else:
        tmp = (p2y - p1y) * invert(p2x - p1x, p) % p
    x = (tmp * tmp - p1x - p2x) % p
    y = (tmp * (p1x - x) - p1y) % p
    return x, y


def mul(p, n):
    r = zero
    tmp = p
    while 0 < n:
        if n & 1 == 1:
            r = add(r, tmp)
        n, tmp = n >> 1, add(tmp, tmp)
    return r


class Task(socketserver.BaseRequestHandler):
    def _recvall(self):
        BUFF_SIZE = 1024
        data = b''
        while True:
            part = self.request.recv(BUFF_SIZE)
            data += part
            if len(part) < BUFF_SIZE:
                break
        return data.strip()

    def send(self, msg, newline=True):
        try:
            if newline:
                msg += b'\n'
            self.request.sendall(msg)
        except:
            pass

    def recv(self, prompt=b'> '):
        self.send(prompt, newline=False)
        return self._recvall()

    def timeout_handler(self, signum, frame):
        raise TimeoutError

    def proof_of_work(self):
        random.seed(urandom(32))
        alphabet = string.ascii_letters + string.digits
        proof = ''.join(random.choices(alphabet, k=32))
        hash_value = sha256(proof.encode()).hexdigest()
        self.send(f'sha256(XXXX+{proof[4:]}) == {hash_value}'.encode())
        nonce = self.recv(prompt=b'Give me XXXX > ')
        if len(nonce) != 4 or sha256(nonce + proof[4:].encode()).hexdigest() != hash_value:
            return True
        return True

    def handle(self):
        try:
            signal.signal(signal.SIGALRM, self.timeout_handler)
            signal.alarm(60)

            if not self.proof_of_work():
                self.send(b'\nWrong!')
                return

            self.secret = random.randint(0, p)
            Q = mul(P, self.secret)

            self.send(b'p: ' + str(p).encode())
            self.send(b'a: ' + str(a).encode())
            self.send(b'P: ' + self.point_to_string(P).encode())
            self.send(b'Q: ' + self.point_to_string(Q).encode())

            signal.alarm(300)
		
            for _ in range(45):
                self.send(MENU, newline=False)
                choice = int(self.recv(prompt=b'> '))
                if choice == 1:
                    self.sign()
                elif choice == 2:
                    self.get_flag()
                else:                
                    break

            self.send(b'Bye!\n')

        except TimeoutError:
            self.send(b'\nTimeout!')
        except Exception:
            self.send(b'Something Wrong!')
        finally:
            self.request.close()

    @staticmethod
    def point_to_string(p):
        return f'({p[0]}, {p[1]})'

    def sign(self):
        self.send(b'Give me your key:')
        x = int(self.recv(prompt=b'X > '))
        y = int(self.recv(prompt=b'Y > '))
        point = (x, y)
        result = mul(point, self.secret)
        self.send(b'result: ' + self.point_to_string(result).encode())

    def get_flag(self):
        s = int(self.recv(prompt=b'Give me the secret > '))
        if s == self.secret:
            self.send(b'Here is your flag:')
            self.send(FLAG)
        else:
            self.send(b'wrong.')


class ForkedServer(socketserver.ForkingMixIn, socketserver.TCPServer):
    pass


if __name__ == '__main__':
    HOST, PORT = '0.0.0.0', 10002
    print(HOST, PORT)
    server = ForkedServer((HOST, PORT), Task)
    server.allow_reuse_address = True
    server.serve_forever()

题目代码比较长，但其实主要任务很容易描述：

• 生成一条给定参数的椭圆曲线
• 生成一个(0,p)之间的随机数当作私钥secret
• 通过proof后，有45次交互机会。输入1，可以输入一个点的坐标，并获得他的secret倍点坐标；输入2，可以核验secret值，如果输入的值等于secret则得到flag

因此目标明确：想办法构造一些点，并利用输入1的交互返回的倍点来求解secret的有关信息。

首先检查题目椭圆曲线的几个参数，可以发现这个曲线并不是易受攻击的椭圆曲线(比如阶光滑、阶为p之类)，那么直接求解DLP肯定是不行的。

然后注意到：

self.send(b'p: ' + str(p).encode())
self.send(b'a: ' + str(a).encode())
self.send(b'P: ' + self.point_to_string(P).encode())
self.send(b'Q: ' + self.point_to_string(Q).encode())

这是通过proof后靶机会发送给我们的信息，可以注意到一个隐晦的提示：为什么不发送参数b给我们？

要知道，参数b其实是题目上面有的。结合题目任务：我们可以输入1，得到自定义的一个点的倍点坐标。紧接着就会注意到：这里根本没有对输入的点在不在曲线上进行检查。

也就是说，我们完全可以输入一个不在题目曲线上的点，靶机端仍然会照常进行倍点运算并返回给我们倍点坐标。那么如何利用这一点呢？

熟悉倍点计算过程的话应该会知道，倍点的计算其实和参数b没有一点关系，也就是说，我们可以自选一个合适的b’，然后构造一条新的曲线：

$$E': \quad y^2 = x^3 + ax + b' \quad(mod\;p)$$

也就是说，由于倍点与b无关，我们输入给靶机这条曲线上的点的话，靶机计算出来的其实也是这条曲线上的倍点。因此我们现在就可以想办法改变b’，从而找到一条合适的曲线。

方法一

什么样的参数b’才算合适呢，首先能想到的肯定是找到的b’应该使曲线的DLP易求。就比如说，假如我们可以找到b’使得曲线E’的阶是光滑的，那么就可以输入一个E’的生成元，并很轻松地由倍点求解DLP得到secret；又或者我们能找到b’使得E’是anomalous curve，也就是E’.order()=p，那么就可以用Smart Attack攻击等等。

大致思路如下：

from Crypto.Util.number import * 

p = 0xc483230557557e3d94b7407f3355b8a5b26bda29119babcb8d72b5a19e10e113
a = 0xb5ecbb93d4e49fe3c93ad770343ec5ab70131151151fcc830f6c658223c92e1d

for i in range(1,100000):
    b = i
    temp = factor(EllipticCurve(GF(p), [a, b]).order())
    print(temp)
    maxp = [int(p) for p, e in temp][-1]
    if(maxp < 2 ^ 40):
        print(b)

但是实际这样做的话，会发现由于order比较大，所以factor执行的很慢，因此很难在短时间内找到一个符合条件的b’。所以要想别的更通用一点的办法，毕竟我们其实有多次交互机会可以用。

方法二

可以想到用CRT来解题，我们假设取了某个b’，得到的E’的阶order是：

$$order = p_1^{k_1}p_2^{k_2}...p_n^{k_n}$$

其中，因子分解按从小到大的递增顺序排列，那么假设我们用factordb获取到其中的一个小素因子a1。(较小但是又不那么小，比如3347、3517这种)

然后，我们取E’的一个生成元G，G的阶就等于曲线的阶order，然后我们获取以下点：

$$G_1 = \frac{order}{a_1}*G$$

那么G1是a1阶元素，这是因为：

$$a_1G_1 = a_1*\frac{order}{a_1}*G = order*G = O$$

而a1是我们分解order后得到的一个小素因子，那么我们发送G1给靶机得到Q1，然后就可以在(0,a1)内爆破出对应的倍数d1，d1满足：

$$d_1 = secret \quad(mod\;a_1)$$

然后我们改变b，并如此重复多次，就可以CRT得到secret在更大模数下的值。也就是说，我们不断调整b，并选取order的一些小因子去求secret模这些小因子下的结果，只要能使得最后的模数乘积大于p，我们就能得到secret的原本值了。

然后还有一个问题就是限时300s且只有44次输入1的交互，所以选的素因子最佳的就是大小为几千的数，选太大会加长爆破时间，选太小可能会注意不到模数重复，并且可能凑不够256比特。同时proof过的太慢的话可以掐掉重来(不过选的数字恰当的话不会这么紧)。

exp：

from Crypto.Util.number import *
from Pwn4Sage.pwn import *
from tqdm import *
from hashlib import sha256
import string
from sympy.ntheory.modular import crt

def proof_of_work():
	table = string.digits + string.ascii_letters
	temp = r.recvuntil(b"sha256(XXXX+")
	temp = r.recvline()
	suffix = temp[:28].decode()
	hex1 = temp[33:].strip().decode()
	for i in tqdm(table):
		for j in table:
			for k in table:
				for m in table:
					temp1 = i+j+k+m
					if(sha256((temp1+suffix).encode()).hexdigest() == hex1):
						r.sendline(temp1.encode())
						return

#part1 Get parameters
r = remote("node5.anna.nssctf.cn",28809)
proof_of_work()
r.recvuntil(b"XXXX > ")

p = int(r.recvline().strip().decode()[3:])
a = int(r.recvline().strip().decode()[3:])
P = eval(r.recvline().strip().decode()[3:])
Q = eval(r.recvline().strip().decode()[3:])


#part2 Get secret % num_i
B = [(5,12841),(9,9463),(12,1567,2371,5981),(14,7417),(16,3413),(27,4127),(28,9901),(36,8581),(41,5659),(45, 5477,7057),(54,1109),(61,409),(64,7243),(66,1049),(67,3469),(68,1979),(69,1987,6967),(73,1873),(77,1153)]
modnum = []
secret = []

for i in trange(len(B)):
	b = B[i][0]
	E = EllipticCurve(Zmod(p),[a,b])
	order = E.order()
	G = E.gens()[0]
	factors = B[i][1:]
	for j in factors:
		r.recvuntil(b"> ")
		r.sendline(b"1")

		sendG = (order//j)*G
		modnum.append(j)
		sendGx = int(sendG[0])
		sendGy = int(sendG[1])
		r.sendline(str(sendGx).encode())
		r.sendline(str(sendGy).encode())

		r.recvuntil(b'result: ')
		temp = eval(r.recvline().strip().decode())
		Qx = int(temp[0])
		Qy = int(temp[1])
		Q = E(Qx,Qy)

		for k in range(j):
			if(k*sendG == Q):
				secret.append(k)
				break

#part3 Get flag using CRT
secret1 = crt(modnum,secret)[0]	
r.recvuntil(b"> ")
r.sendline(b"2")
r.recvuntil(b"the secret > ")
r.sendline(str(secret1).encode())
r.recvline()
print(r.recvline())


#NSSCTF{dabf6daf-d95d-4764-8af0-b5fa9b8dc442}

Tiny ECC

题目来源：CryptoCTF 2021

题目：

#!/usr/bin/env python3
# In the name of Allah

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

def tonelli_shanks(n, p):
	if pow(n, int((p-1)//2), p) == 1:
			s = 1
			q = int((p-1)//2)
			while True:
				if q % 2 == 0:
					q = q // 2
					s += 1
				else:
					break
			if s == 1:
				r1 = pow(n, int((p+1)//4), p)
				r2 = p - r1
				return r1, r2
			else:
				z = 2
				while True:
					if pow(z, int((p-1)//2), p) == p - 1:
						c = pow(z, q, p)
						break
					else:
						z += 1
				r = pow(n, int((q+1)//2), p)
				t = pow(n, q, p)
				m = s
				while True:
					if t == 1:
						r1 = r
						r2 = p - r1
						return r1, r2
					else:
						i = 1
						while True:
							if pow(t, 2**i, p) == 1:
								break
							else:
								i += 1
						b = pow(c, 2**(m-i-1), p)
						r = r * b % p
						t = t * b ** 2 % p
						c = b ** 2 % p
						m = i
	else:
		return False

def random_point(p, a, b):
	while True:
		gx = getRandomRange(1, p-1)
		n = (gx**3 + a*gx + b) % p
		gy = tonelli_shanks(n, p)
		if gy == False:
			continue
		else:
			return (gx, gy[0])

def die(*args):
	pr(*args)
	quit()

def pr(*args):
	s = " ".join(map(str, args))
	sys.stdout.write(s + "\n")
	sys.stdout.flush()

def sc():
	return sys.stdin.readline().strip()

def main():
	border = "+"
	pr(border*72)
	pr(border, "  Dual ECC means two elliptic curve with same coefficients over the ", border)
	pr(border, "  different fields or ring! You should calculate the discrete log   ", border)
	pr(border, "  in dual ECCs. So be smart in choosing the first parameters! Enjoy!", border)
	pr(border*72)

	bool_coef, bool_prime, nbit = False, False, 128
	while True:
		pr(f"| Options: \n|\t[C]hoose the {nbit}-bit prime p \n|\t[A]ssign the coefficients \n|\t[S]olve DLP \n|\t[Q]uit")
		ans = sc().lower()
		if ans == 'a':
			pr('| send the coefficients a and b separated by comma: ')
			COEFS = sc()
			try:
				a, b = [int(_) for _ in COEFS.split(',')]
			except:
				die('| your coefficients are not valid, Bye!!')
			if a*b == 0:
				die('| Kidding me?!! a*b should not be zero!!')
			else:
				bool_coef = True
		elif ans == 'c':
			pr('| send your prime: ')
			p = sc()
			try:
				p = int(p)
			except:
				die('| your input is not valid :(')
			if isPrime(p) and p.bit_length() == nbit and isPrime(2*p + 1):
				q = 2*p + 1
				bool_prime = True
			else:
				die(f'| your integer p is not {nbit}-bit prime or 2p + 1 is not prime, bye!!')
		elif ans == 's':
			if bool_coef == False:
				pr('| please assign the coefficients.')
			if bool_prime == False:
				pr('| please choose your prime first.')
			if bool_prime and bool_coef:
				Ep = CurveOverFp(0, a, b, p)
				Eq = CurveOverFp(0, a, b, q)

				xp, yp = random_point(p, a, b)
				P = Point(xp, yp)

				xq, yq = random_point(q, a, b)
				Q = Point(xq, yq)

				k = getRandomRange(1, p >> 1)
				kP = Ep.mult(P, k)

				l = getRandomRange(1, q >> 1)
				lQ = Eq.mult(Q, l)
				pr('| We know that: ')
				pr(f'| P = {P}')
				pr(f'| k*P = {kP}')
				pr(f'| Q = {Q}')
				pr(f'| l*Q = {lQ}')
				pr('| send the k and l separated by comma: ')
				PRIVS = sc()
				try:
					priv, qriv = [int(s) for s in PRIVS.split(',')]
				except:
					die('| your input is not valid, Bye!!')
				if priv == k and qriv == l:
					die(f'| Congrats, you got the flag: {flag}')
				else:
					die('| sorry, your keys are not correct! Bye!!!')
		elif ans == 'q':
			die("Quitting ...")
		else:
			die("Bye ...")

if __name__ == '__main__':
	main()

题目比较长，但是其实提取出来的任务还是很清楚，靶机提供给了我们三个有用选项：

• 输入a，可以提供ECC的a、b参数，这里要求a、b均不为0
• 输入c，可以提供ECC的模数p，这里要求p是一个128比特的素数且q=2p+1也是一个素数
• 在参数均给定后，输入s，题目会用输入的参数生成以下两条曲线：

$$y^2 = x^3 + ax + b \quad(mod\; p)$$ $$y^2 = x^3 + ax + b \quad(mod\; q)$$

• 并且，题目会分别取这两条曲线上的随机点P、Q，并生成两个随机数k、l，计算倍点kP、lQ，并给出四个点的坐标。如果能够计算出正确的k、l并提交，我们就能获取flag

那么其实也就是解一个离散对数问题，而ECC的所有参数都是我们可以自定义的，那么可操作性非常大。这里我阐述两种思路。

爆破

最直接的思路就是：直接爆破p、q、a、b，使得两条曲线均光滑，这样我们就可以pohlig-hellman轻松解决DLP。这一部分实现如下：

def genprime(nbits):
    while(1):
        a = getPrime(16)
        b = getPrime(16)
        p = int(getPrime(128))
        q = int(2*p + 1)
        if(isPrime(q)):
            Ep = EllipticCurve(GF(p),[a,b])
            Eq = EllipticCurve(GF(q),[a,b])
            factorsEp = list(factor(Ep.order()))
            factorsEq = list(factor(Eq.order()))
            print(int(factorsEp[-1][0]).bit_length(),int(factorsEq[-1][0]).bit_length())
            if(factorsEp[-1][0] < 2^35 and factorsEq[-1][0] < 2^35):
                print(factorsEp)
                print(factorsEq)
                return a,b,p,q
        else:
            continue

这个方法确实用一小段时间就能找到一组合理值。

映射

这个方法我是看maple的博客看到的，原wp指路：

Crypto CTF 2021 WriteUps | 廢文集中區 (maple3142.net)

也就是说，在输入a、b时，他只检查a、b均不为0.并且他的椭圆曲线实现中不会检查曲线是singular curve(4a^3+27b^2=0)的情况(而sage中会)，所以我们可以令a，b均等于pq，这样两条曲线就分别是：

$$y^2 = x^3 \quad(mod\;p)$$ $$y^2 = x^3 \quad(mod\;q)$$

而这样的singular curve可以作椭圆曲线加法群中的点到有限域中的数的映射：

$$(x,y) \rightarrow \frac{y}{x}$$

因此可以直接求解逆元求解离散对数问题。而其实取其他singular curve的参数，我们同样可以找到映射，但是肯定没有这个这么简单。

exp：

from Crypto.Util.number import *
from Pwn4Sage.pwn import *

def genprime(nbits):
    while(1):
        a = getPrime(16)
        b = getPrime(16)
        p = int(getPrime(128))
        q = int(2*p + 1)
        if(isPrime(q)):
            Ep = EllipticCurve(GF(p),[a,b])
            Eq = EllipticCurve(GF(q),[a,b])
            factorsEp = list(factor(Ep.order()))
            factorsEq = list(factor(Eq.order()))
            print(int(factorsEp[-1][0]).bit_length(),int(factorsEq[-1][0]).bit_length())
            if(factorsEp[-1][0] < 2^35 and factorsEq[-1][0] < 2^35):
                print(factorsEp)
                print(factorsEq)
                return a,b,p,q
        else:
            continue

r = remote("node4.anna.nssctf.cn",28697)

#part1 send p
r.recvuntil(b"[Q]uit")
r.sendline(b"c")
r.recvuntil(b"prime:")
while(1):
    p = int(getPrime(128))
    q = int(2*p+1)
    if(isPrime(q)):
        break
r.sendline(str(p).encode())

#part2 send a,b
r.recvuntil(b"[Q]uit")
r.sendline(b"a")
r.recvuntil(b"by comma:")
a = p*q
b = p*q
r.sendline(str(a).encode() + b"," + str(b).encode())

#part3 get point
r.recvuntil(b"[Q]uit")
r.sendline(b"s")
r.recvuntil(b" = (")
P = eval("(" + r.recvline().strip().decode())
r.recvuntil(b" = (")
kP = eval("(" + r.recvline().strip().decode())
r.recvuntil(b" = (")
Q = eval("(" + r.recvline().strip().decode())
r.recvuntil(b" = (")
lQ = eval("(" + r.recvline().strip().decode())
r.recvuntil(b"comma: ")

#part4 solveDLP Getflag
def mapping(point,prime):
    x = point[0]
    y = point[1]
    return y*inverse(x,prime) % prime
map_P = mapping(P,p)
map_kP = mapping(kP,p)
map_Q = mapping(Q,q)
map_lQ = mapping(lQ,q)

k = inverse(map_kP,p)*map_P % p
l = inverse(map_lQ,q)*map_Q % q
r.sendline(str(k).encode() + b"," + str(l).encode())
r.recvline()
print(r.recvline())

#NSSCTF{c375f851-5f15-438d-9726-fd20b8903ebe}