密码学氛围代码垃圾

加密技术氛围下的垃圾代码

我一直在阅读这篇论文 [2]：

它是一种相对较新的将三次 Pell 曲线与 RSA 方法结合使用的方法。使用 RSA（Rivest、Shamir 和 Adleman）方法，我们生成两个素数（$p$ 和 $q$），并创建一个模数 $N = p \cdot q$。 因此，该方法的安全性基于将模数分解为其原始素数的难度。在 2018 年，Murru 和 Saettone 提出了一种新的 RSA 方法，该方法基于具有 $N = p.q$ 形式的模数的三次 Pell 曲线 [1]。 通过这种方法，我们基于三次 Pell 曲线执行算术运算，而这只有消息的发送者和接收者知道。在 2024 年的一篇论文中，Seck 和 Nitaj [2] 通过 $N = p^r q^s$ 的模数增强了这一点。

在 RSA 中，我们有一个公钥 (e,N) 和一个私钥 (d,N)。 Pell 曲线增强了基本的加密和解密过程，而不是使用指数运算，例如：

C = M^e (mod N)

并且解密：

D = C^d (mod N)

我们使用三次 Pell 曲线运算。 Pell 曲线的基本方程是：

x³+ay³+a²z³-(3axyz) (mod N)

其中 N=p.q [1]（或 N=p^r q^s [2]）。该方法的核心是找到以下方程的解：

x³+ay³+a²z³-(3axyz) = 1 (mod N)

ChatGPT

在理解了论文中的主要方法之后：

我决定尝试一下 ChatGPT 的代码解决方案：

编写一个基于三次
Pell 方程实现 RSA 变体的 Python 程序

不幸的是，它给了我垃圾代码，其中包含太多错误：

##!/usr/bin/env python3
"""
Cubic Pell RSA (核心算术的参考实现)

实现：
  - 环 R = Z_N[α]/(α^3 - a)
  - 元素表示为三元组 (x, y, z)，表示 x + y*α + z*α^2
  - 由 α^3 = a 引起的乘法
  - “单位群”/范数-1 子群算术是 cubic-Pell RSA 变体使用的
  - 类似 RSA 的幂运算：C = P^e, P = C^d

此文件有意将消息编码/解码保留为可插入函数，
因为不同的论文使用不同的编码和参数约束。
"""

from __future__ import annotations
from dataclasses import dataclass
from typing import Tuple, Optional
import secrets
import math

## -------------------------
## 基本数论实用程序
## -------------------------

def egcd(a: int, b: int) -> Tuple[int, int, int]:
    """扩展 GCD：返回 (g, x, y)，其中 ax + by = g = gcd(a,b)。"""
    x0, y0, x1, y1 = 1, 0, 0, 1
    while b != 0:
        q, a, b = a // b, b, a % b
        x0, x1 = x1, x0 - q * x1
        y0, y1 = y1, y0 - q * y1
    return a, x0, y0

def modinv(a: int, n: int) -> int:
    """a 模 n 的模逆。"""
    a %= n
    g, x, _ = egcd(a, n)
    if g != 1:
        raise ValueError(f"no inverse: gcd({a},{n})={g}")
    return x % n

def is_probable_prime(n: int, rounds: int = 40) -> bool:
    """Miller-Rabin 概率素性测试。"""
    if n < 2:
        return False
    small_primes = [2,3,5,7,11,13,17,19,23,29,31,37]
    for p in small_primes:
        if n == p:
            return True
        if n % p == 0:
            return False

    # 写成 n-1 = d * 2^s，其中 d 为奇数
    d = n - 1
    s = 0
    while d % 2 == 0:
        s += 1
        d //= 2

    def check(a: int) -> bool:
        x = pow(a, d, n)
        if x == 1 or x == n - 1:
            return True
        for _ in range(s - 1):
            x = (x * x) % n
            if x == n - 1:
                return True
        return False

    for _ in range(rounds):
        a = secrets.randbelow(n - 3) + 2
        if not check(a):
            return False
    return True

def random_prime(bits: int, *, mod: Optional[int] = None, rem: int = 0) -> int:
    """生成一个 `bits` 位的随机素数；可以选择强制执行 p % mod == rem。"""
    if bits < 2:
        raise ValueError("bits must be >= 2")
    while True:
        x = secrets.randbits(bits)
        x |= (1 << (bits - 1))  # 确保最高位
        x |= 1                  # 确保为奇数
        if mod is not None and x % mod != rem:
            continue
        if is_probable_prime(x):
            return x

## -------------------------
## 三次 Pell 元素
## -------------------------

@dataclass(frozen=True)
class CubicPell:
    """
    表示 Z_N[α]/(α^3 - a) 中的一个元素 u = x + y*α + z*α^2。
    """
    x: int
    y: int
    z: int
    N: int
    a: int

    def __post_init__(self):
        object.__setattr__(self, "x", self.x % self.N)
        object.__setattr__(self, "y", self.y % self.N)
        object.__setattr__(self, "z", self.z % self.N)
        object.__setattr__(self, "a", self.a % self.N)

    @staticmethod
    def one(N: int, a: int) -> "CubicPell":
        return CubicPell(1, 0, 0, N, a)

    def mul(self, other: "CubicPell") -> "CubicPell":
        """
        使用 α^3=a 乘 (x1+y1α+z1α^2)(x2+y2α+z2α^2)。

        使用：
          α^3 = a
          α^4 = aα
        乘积系数：
          x' = x1x2 + a(y1z2 + z1y2)
          y' = x1y2 + y1x2 + a(z1z2)
          z' = x1z2 + z1x2 + (y1y2)
        """
        if self.N != other.N or self.a != other.a:
            raise ValueError("mismatched modulus/parameter")

        N = self.N
        a = self.a

        x1, y1, z1 = self.x, self.y, self.z
        x2, y2, z2 = other.x, other.y, other.z

        x = (x1 * x2 + a * (y1 * z2 + z1 * y2)) % N
        y = (x1 * y2 + y1 * x2 + a * (z1 * z2)) % N
        z = (x1 * z2 + z1 * x2 + (y1 * y2)) % N
        return CubicPell(x, y, z, N, a)

    def pow(self, e: int) -> "CubicPell":
        """通过平方和乘法进行快速求幂。"""
        if e < 0:
            raise ValueError("negative exponent not supported in this reference")
        result = CubicPell.one(self.N, self.a)
        base = self
        k = e
        while k:
            if k & 1:
                result = result.mul(base)
            base = base.mul(base)
            k >>= 1
        return result

    def norm(self) -> int:
        """
        范数 N(x+yα+zα^2) = x^3 + a y^3 + a^2 z^3 - 3 a x y z  (mod N)。
        """
        N = self.N
        a = self.a
        x, y, z = self.x, self.y, self.z
        return (pow(x, 3, N) + a * pow(y, 3, N) + (a * a % N) * pow(z, 3, N) - (3 * a * x * y * z)) % N

    def is_norm_one(self) -> bool:
        return self.norm() % self.N == 1 % self.N

## -------------------------
## 类似 Cubic Pell RSA 的密钥
## -------------------------

@dataclass
class PublicKey:
    N: int
    a: int
    e: int

@dataclass
class PrivateKey:
    N: int
    a: int
    d: int
    # 可选地保留因子（方便测试；在实际使用中省略）
    p: Optional[int] = None
    q: Optional[int] = None

def psi_murru_saettone(p: int, q: int) -> int:
    """ψ(N) = (p^2 + p + 1)(q^2 + q + 1)。"""
    return (p*p + p + 1) * (q*q + q + 1)

def keygen(bits: int = 512, *, a: Optional[int] = None, e: int = 65537) -> Tuple[PublicKey, PrivateKey]:
    """
    使用 N=pq 和
    ed ≡ 1 (mod (p^2+p+1)(q^2+q+1)) 为经典 cubic-Pell RSA 变体生成密钥。

    注意：许多构造还需要对 a 进行条件限制（三次非剩余等）。
    如果未提供，我们在此简单地随机选择一个，并强制执行 gcd(e, ψ)=1。
    """
    if bits < 64:
        raise ValueError("use >= 64 bits even for toy demos")

    while True:
        p = random_prime(bits, mod=3, rem=2)  # 在某些三次剩余设置中很常见
        q = random_prime(bits, mod=3, rem=2)
        if p == q:
            continue
        N = p * q
        psi = psi_murru_saettone(p, q)

        if math.gcd(e, psi) != 1:
            continue

        d = modinv(e, psi)

        if a is None:
            # 在 Z_N^* 中均匀随机选择 a（避免 0 和小 gcd）
            while True:
                aa = secrets.randbelow(N - 2) + 2
                if math.gcd(aa, N) == 1:
                    a_val = aa
                    break
        else:
            a_val = a % N
            if a_val == 0 or math.gcd(a_val, N) != 1:
                raise ValueError("a must be nonzero and coprime to N")

        pk = PublicKey(N=N, a=a_val, e=e)
        sk = PrivateKey(N=N, a=a_val, d=d, p=p, q=q)
        return pk, sk

## -------------------------
## 加密/解密
## -------------------------

def encrypt_point(P: CubicPell, pk: PublicKey) -> CubicPell:
    """通过幂运算加密*曲线/环元素* P：C = P^e。"""
    if P.N != pk.N or P.a != pk.a:
        raise ValueError("point parameters mismatch")
    return P.pow(pk.e)

def decrypt_point(C: CubicPell, sk: PrivateKey) -> CubicPell:
    """通过幂运算解密：P = C^d。"""
    if C.N != sk.N or C.a != sk.a:
        raise ValueError("cipher parameters mismatch")
    return C.pow(sk.d)

## -------------------------
## 示例/玩具编码
## -------------------------
## 实际方案使用复杂的 Encode2Curve/DecodeFromCurve 风格的映射。
## 为了实现最小的可运行演示，我们将整数 m 编码为范数 1 元素
## 通过暴力破解一个小调整 t，直到 norm(x,y,z)=1（仅适用于玩具大小）。
## 请勿对实际大小使用暴力编码。

def encode_int_toy(m: int, N: int, a: int, max_tries: int = 100000) -> CubicPell:
    """
    玩具编码：搜索 t，使得元素 (m, t, 1) 具有范数 1 mod N。
    这纯粹是为了演示小 N（否则不可行）。
    """
    m %= N
    for _ in range(max_tries):
        t = secrets.randbelow(N)
        P = CubicPell(m, t, 1, N, a)
        if P.is_norm_one():
            return P
    raise RuntimeError("toy encoding failed (try increasing max_tries or use smaller N)")

def decode_int_toy(P: CubicPell) -> int:
    """玩具解码：仅返回 x 坐标（仅适用于玩具编码器）。"""
    return P.x

def demo():
    # 使用微小的大小，以便玩具演示编码能够快速成功。
    pk, sk = keygen(bits=80, e=65537)
    print("N =", pk.N)
    print("a =", pk.a)
    print("e =", pk.e)
    print("d =", sk.d)

    message = 123456789
    P = encode_int_toy(message, pk.N, pk.a)
    assert P.is_norm_one()
    C = encrypt_point(P, pk)
    P2 = decrypt_point(C, sk)
    assert P2.is_norm_one()

    m2 = decode_int_toy(P2)
    print("message =", message)
    print("decoded =", m2)
    print("OK =", (m2 == message))

if __name__ == "__main__":
    demo()

我真的无法开始概述此代码的问题。 即使代码看起来很权威，但方法论是错误的。 即使你修复了所有基本的代码问题（例如缺少类的初始化），你也不会让它工作。 我花了几个小时修复数学和语法上的众多问题，但最终放弃了，因为它完全是垃圾！

正确的代码

因此，我从论文中找到了一些 Sage 代码 [ here ]，并将其转换为纯 Python [ here ]：

## 基于此处的 Sage 代码：https://github.com/mseckept/schemesecknitaj/tree/main

from Crypto.Util.number import *
from Crypto import Random
import Crypto
from random import randint
from math import gcd
import libnum

def is_cube(a, p):
## 测试 a 是否为立方剩余 mod p，p 为素数
    return power_mod(a, (p-1)//gcd(3, p-1), p) == 1

def is_square(x, p):
## 测试 x 是否为二次剩余 mod p，p 为素数
    if (libnum.has_sqrtmod(x,{p: 1})):
       return True
    else:
       return False

def rand_prime_congru_1mod3_and_3mod4(size):
## 随机生成一个素数 p，2 <= p < 2^size，其中 p = 1 mod 3 且 p = 3 mod 4

    p = Crypto.Util.number.getPrime(size, randfunc=Crypto.Random.get_random_bytes)
    while (p%3 != 1) or (p%4 != 3):
        p = Crypto.Util.number.getPrime(size, randfunc=Crypto.Random.get_random_bytes)
    return p

def rand_pq(size):
## 随机生成两个不同的素数 p, q = 1 mod 3 且 p, q = 3 mod 4

    p = rand_prime_congru_1mod3_and_3mod4(size)
    q = p
    while q == p:
        q = rand_prime_congru_1mod3_and_3mod4(size)
    return (p, q)

## 使用 Hensel 引理解二次方程 ax^2+bx+c = 0 mod p**r 的解

def extend_solution_quadratic_equation_mod_pr(coeffs, sol, p, r):
## 使用 Hensel 引理将二次方程 ax^2+bx+c = 0 的解 *sol* 从 (mod p) 扩展到 (mod p^r)

    a, b, c = coeffs; x0 = int(sol)
    f = lambda x : a *x**2 + b*x + c
    fprime = lambda x : 2*a*x + b
    # 我们应用 Hensel 引理
    x0 = int(sol)
    for _ in range(2, r+1):
        x0 = (x0 - f(x0)*inverse_mod(fprime(x0), p**r )) % p**r
    return x0

def solver_quadratic_equation_mod_pr(coeffs, p, r):
## 使用 Hensel 引理解二次方程 ax^2+bx+c = 0 (mod p^r) 的解
    a, b, c = coeffs;

    if  p%4 == 3:
        x1 = Mod((-b- power_mod(Mod(b**2-4*a*c, p), (p+1)//4, p)) * inverse_mod(2*a, p), p)
        x2 = Mod((-b + power_mod(Mod(b**2-4*a*c, p), (p+1)//4, p)) * inverse_mod(2*a, p), p)
    else:
        x1 = Mod((-b - Mod(b**2-4*a*c, p).sqrt()) * inverse_mod(2*a, p), p)
        x2 = Mod((-b + Mod(b**2-4*a*c, p).sqrt()) * inverse_mod(2*a, p), p)

    x1r = extend_solution_quadratic_equation_mod_pr(coeffs, x1, p, r)
    x2r = extend_solution_quadratic_equation_mod_pr(coeffs, x2, p, r)
    return (x1r, x2r)

def get_sol_with_cubic_pell_point_mod_pr_qs(C, p, q, r, s):
## 给定三次 Pell 方程上的一个点，确定 mod p^r 和 mod q^s 的解
    xc, yc, zc = C; (a, b, c) = zc**3, yc**3-3*xc*yc*zc, xc**3 - 1
    sol_mod_pr = solver_quadratic_equation_mod_pr((a, b, c), p, r)
    sol_mod_qs = solver_quadratic_equation_mod_pr((a, b, c), q, s)
    return (sol_mod_pr, sol_mod_qs)

def is_solution_mod_pr(C, sol, p, r):
## 测试 *sol* 是否为 xc^3 + x*yc^3 + x^2 * zc^3 - 3x*xc*yc*zc = 1 mod p^r 的解
    xc, yc, zc = C; (a, b, c) = zc**3, yc**3-3*xc*yc*zc, xc**3 - 1
    f = lambda x : a *x**2 + b*x + c
    return (f(sol) % p**r ) == 0

def get_sols_with_ctheo(xp, yq, p, q, r, s):
##  确定 N=p^rq^s 的解 Y mod N，其中 Y = xp mod p^r 且 Y = yq mod q^s

    N = p**r  * q**s
    X = (xp*q**s *inverse_mod(q**s , p**r ) + yq*p**r *inverse_mod(p**r , q**s ))
    return X % N

def get_possible_a_values(C, p, q, r, s):
## 确定 a 的所有值，使得 xc^3 + a*yc^3 + a^2 * zc^3 - 3a*xc*yc*zc = 1 mod N，其中 N = p^rq^s

    N = p**r * q**s; (xc, yc, zc) = C
    sol_mod_pr, sol_mod_qs = get_sol_with_cubic_pell_point_mod_pr_qs(C, p, q, r, s)
    a1 =  get_sols_with_ctheo(sol_mod_pr[0], sol_mod_qs[0], p, q, r, s)
    a2 =  get_sols_with_ctheo(sol_mod_pr[1], sol_mod_qs[1], p, q, r, s)
    a3 =  get_sols_with_ctheo(sol_mod_pr[0], sol_mod_qs[1], p, q, r, s)
    a4 =  get_sols_with_ctheo(sol_mod_pr[1], sol_mod_qs[0], p, q, r, s)
    possible_a = (a1, a2, a3, a4)

    return possible_a

def inverse_mod(a,b):
   return(pow(a,-1,b))

def Mod(a,b):
    return(a%b)

def power_mod(a,b,N):
    return(pow(a,b,N))

## 三次 Pell 曲线的算术运算

def add(pt1, pt2, a, N):
## 在 Z/NZ 中，将三次 Pell 曲线 C_a(N) = x^3+ay^3+a^2z^3-3axyz =  1 上的两个点 pt1 和 pt2 相加

    (x1, y1, z1) =  pt1; (x2, y2, z2) = pt2
    x3 = (x1*x2 + a*(y2*z1 + y1*z2)) % N
    y3 = (x2*y1 + x1*y2 + a*z1*z2) % N
    z3 = (y1*y2 + x2*z1 + x1*z2) % N
    return (x3, y3, z3)

def exponent(pt, e, a, N):
## 在三次 Pell 曲线 C_a(N) = x^3+ay^3+a^2z^3-3axyz =  1 中，从左到右进行标量乘法
    (x1, y1, z1) = pt; (x2, y2, z2) = (1, 0, 0)
    for bit in bin(e)[2:]:
        (x2, y2, z2) = add((x2, y2, z2), (x2, y2, z2), a, N)
        if bit == '1':
            (x2, y2, z2) = add((x2, y2, z2), (x1, y1, z1), a, N)

    return (x2, y2, z2)

def keygen(size, r, s):

    p, q = rand_pq(size)
    N = p**r  * q**s
    ordn = p*q*(p**2+p+1)*(q**2+q+1)*(p - 1)**2*(q - 1)**2
    e = randint(1, N)
    while gcd(e, ordn) != 1:
        e = randint(1, N)

    d1 = inverse_mod(e, p**(2*(r-1))*q**(2*(s-1))*(p**2+p+1)*(q**2+q+1))
    d2 = inverse_mod(e, p**(2*(r-1))*q**(2*(s-1))*(p - 1)**2*(q - 1)**2)
    d3 = inverse_mod(e, p**(2*(r-1))*q**(2*(s-1))*(p**2+p+1)*(q - 1)**2)
    d4 = inverse_mod(e, p**(2*(r-1))*q**(2*(s-1))*(p - 1)**2*(q**1+q+1))
    return (N, e), (p, q, (r, s), (d1, d2, d3, d4))

def encipher(m, pk):

    (N, e) = pk; (xm, ym) = m
    a = ((1 - xm**3) * inverse_mod(ym**3, N)) % N
    C = exponent((xm, ym, 0), e, a, N)
    return C

def decipher(ctx, sk):

    (p, q, (r, s), (d1, d2, d3, d4)) = sk
    N =  p**r *q**s
    possible_a = get_possible_a_values(ctx, p, q, r, s)
    for a in possible_a:
        if (not is_cube(a, p)) and (not is_cube(a, q)):
            d = d1
        elif is_cube(a, p) and is_cube(a, q):
            d = d2
        elif (not is_cube(a, p)) and is_cube(a, q):
            d = d3
        else:
            d = d4
        (xm, ym, zm) = exponent(ctx, d, a, N)
        if zm == 0:
            return (xm, ym)

r, s = 2, 3

msg="Hello"
bits=60

if (len(sys.argv)>1):
        msg=str(sys.argv[1])
if (len(sys.argv)>2):
        bits=int(sys.argv[2])

pk, sk = keygen(bits, r, s)
N = pk[0]
print(f"pk = {pk}\nsk = {sk}\n")
print(f"Bits = {bits}\nMessage = {msg}\n")

m=  bytes_to_long(msg.encode('utf-8'))
m = (m, m);
print(f"\nMessage: {m}\n")
ctx = encipher(m, pk);
print(f"Cipher = {ctx}\n")
ptx = decipher(ctx, sk=sk);
print(f"Decipher= {ptx}\n")

print ("\n==解密后的消息===\nmsg=%s" %long_to_bytes(ptx[0]) )

当我使用 “Hello” 和 60 位消息运行时，我们得到了正确的结果 [ here ]：

pk = (1080928831331135013779446457169419532527431370087255033221341860232168568083200088298968271, 839318919357785876635889999468890161415510066915672848798116175735363773849202484382355727)
sk = (848652422799117451, 1144930513775975191, (2, 3), (300195589744604718162001536858972385437891337444414027571033863024674871239379084432567260146654173021529193784230753898560163461195579737766556175825451362007554735939673785303383, 984829125046219908975956622365077885396089251496387714601240798419375876355803485567553293047859329797113199413345859391395499614544656634154761849589025472160237605937157390439663, 106890682378044000985457050811691152002571563022688174143901212098399223677253105864658595329691625764223465951446088797750147128019601693436000919243413417127766262895098848537563, 1241964101001002944870754602380875718664425735833398559540944784530123751044201115593823554376019993620191152959498221676597950534828743348260213445237809115267163))

Bits = 60
Message = Hello

Message: (310939249775, 310939249775)

Cipher = (469485677348104705702147635846339994152100534842645981803425124585832145208328242075114880, 881871471207622001832864371562012465373436593613927667095664977168039698476294479248910587, 628516305167444377466415415147842177803171283952167154498329097616649062567788677885702272)

Decipher= (310939249775, 310939249775)

==解密后的消息===
msg=b'Hello'

你可以在这里尝试代码：

https://asecuritysite.com/rsa/rsa_pell

参考文献

[1] Murru, N., & Saettone, F. M. (2017, September). A novel

• 原文链接： medium.com/asecuritysite...
• 登链社区 AI 助手，为大家转译优秀英文文章，如有翻译不通的地方，还请包涵～