Thoughts

这篇文章的话特殊一点，记录的题目主要满足以下几条：

我还没有能解出
搜索不到wp
有一些思路，并且题目本身我觉得比较有趣

所以文章命名为Thoughts，对这些题目提供我的一些思路，希望对各位师傅有一些启发最终解决题目。如果其中的题目有师傅解出来了或者能找到wp，欢迎随时联系我！

DangerousRSA

题目来源：广东强网杯 2021 个人组

题目描述：

1	It is very dangerous

题目：

from Crypto.Util.number import *

n=1130559951242142418553461557920917245278152787036994964813122956739956040733223190585188299825937660853090732755660557900622564346469528239989519875434612493124663925800482373019421482927567375862528994310090666571233446002163478309309649223886467303834842172338454971591276606602163799133000582891965785243482166772793400065793927508012418516719549212969303944155066541587875789636798394581666621711612432864870236407912486149901575077993702841083287354610398446101725184000000000000000000000000000000149

e=65537

m=bytes_to_long(flag)

def zhuanzhuan(n):
    
    if n < 2: return n

    return 9* zhuanzhuan(n - 1) + zhuanzhuan(n - 2)

g = 1766592252951269548814352411327322447034194093921285273353477875204196603230641896039854934719468650093602325707751568

hint_0=zhuanzhuan(g)%50254221567892101225665855


def encrypt(hint_0):
    #n_1,n_2 both are int   
    if   (n_1+hint_0)**2+1712463104290490642*n_1+4730060142800357260761833404838169046987798-66666666*n_2**2=0:

        d=55607604802787651645008458776960381905274761803290637507632188635838921065681450951091689278291418693230528750002533602907652688436593464630306292523277204515196410606067221890882848225521819884919432609641732036959698443050773592482567332594684068998082886109352105201496925283866221090767083891241755981564899917062843679988341936588610298860835342859782896551931353504004447312763256966682142206683209406463516580194846276027818291321875085607298572702721743561*n_1+n_2-69462057634597116805075535702860983356198049

   
    return d  
    

c=pow(m,e,n)=encrypt(hint_0)

#The n is very dangerous

题目加密流程如下：

定义一个递归函数zhuanzhuan()，并用一个大数g作为该函数的形参计算函数值后，模一个数作为hint0
以hint0定义了一个关于n1、n2的方程
用n1、n2生成一个加密值d，保证其与RSA得到的密文相等

那么其实要解决的就是以下几个问题：

计算zhuanzhuan()函数得到hint0
用hint0解方程得到n1、n2，然后得到RSA密文
分解n从而求解RSA明文

计算zhuanzhuan()函数

zhuanzhuan函数如下：

1
2
3

def zhuanzhuan(n):   
    if n < 2: return n
    return 9* zhuanzhuan(n - 1) + zhuanzhuan(n - 2)

其实就是个加了一项系数的斐波那契，而直接用它进行计算的话，由于g过大，肯定会爆栈，因此要优化一下。

而对斐波那契数列有一种矩阵快速幂计算的办法，当然也就可以用在这个题，具体来说如下。

我们直接把zhuanzhuan(n)叫做f(n)，这样看上去更直接，那么开始我们拥有初值：

$f(0) = 0$ $f(1) = 1$

而后面的所有项都是基于这两项进行递推的，写成矩阵形式就是：

$\left( \begin{matrix} f(2)\\ f(1)\\ \end{matrix} \right) = \left( \begin{matrix} 9&1\\ 1&0\\ \end{matrix} \right) \left( \begin{matrix} f(1)\\ f(0)\\ \end{matrix} \right)$

同理：

$\left( \begin{matrix} f(3)\\ f(2)\\ \end{matrix} \right) = \left( \begin{matrix} 9&1\\ 1&0\\ \end{matrix} \right) \left( \begin{matrix} f(2)\\ f(1)\\ \end{matrix} \right) = \left( \begin{matrix} 9&1\\ 1&0\\ \end{matrix} \right)^2 \left( \begin{matrix} f(1)\\ f(0)\\ \end{matrix} \right)$

所以就有通项：

$\left( \begin{matrix} f(n)\\ f(n-1)\\ \end{matrix} \right) = \left( \begin{matrix} 9&1\\ 1&0\\ \end{matrix} \right)^{n-1} \left( \begin{matrix} f(1)\\ f(0)\\ \end{matrix} \right)$

而由于hint0最后也要模一个素数，因此完全可以将整个矩阵都定义在模数的有限域下进行计算，所以能用矩阵快速幂计算出f(g)，从而得到hint0。

解方程得到n1、n2

n1、n2的方程如下：

1	(n_1+hint_0)*2+1712463104290490642n_1+4730060142800357260761833404838169046987798-66666666n_2*2=0

槽点不少，比如这里都不是”==”而是”=”。。。而现在由于我们已经有hint0，所以把他展开再合并同类项能有以下形式：

$n_1^2 + tn_1 -Dn_2^2 + c = 0$

然后其实能看出有佩尔方程的雏形，试一下果然发现这几个数都是凑好的，也就是说我们可以把方程写成如下佩尔方程的标准形式：

$(n_1+k)^2 -Dn_2^2 = 1$

那么我们可以求解这个佩尔方程得到多组可能的n1、n2，从而得到多个可能的RSA密文。

但是，没能解出来这个题的问题也就在这里，我会在最后一部分展开详述。

分解n

假设解佩尔方程我们得到了正确的RSA密文，那么现在我们显然是需要获得n的分解才能解出明文的。但是题目里关于n一共就两点信息：

一个是n本身：

n=1130559951242142418553461557920917245278152787036994964813122956739956040733223190585188299825937660853090732755660557900622564346469528239989519875434612493124663925800482373019421482927567375862528994310090666571233446002163478309309649223886467303834842172338454971591276606602163799133000582891965785243482166772793400065793927508012418516719549212969303944155066541587875789636798394581666621711612432864870236407912486149901575077993702841083287354610398446101725184000000000000000000000000000000149

另一个是对n的描述：

1	#The n is very dangerous

因此n肯定是由于质数选取不当等原因，导致其可用某种方法简单分解。而观察到n的低位为：

1	......000000000000000000000000000000149

低位由很多0再加上一个149组成，而n的低位也就是pq乘积的低位，并且149是一个质数，那么这也就是说，p、q之中有一个的十进制低位应该是如下形式：

1	......000000000000001

那么就可以猜测他是p-1光滑数，事实上也确实可以用pohlig-hellman获得n的分解。

看上去题目好像已经顺利解决了，但是接下来就是问题：

问题

问题就出在佩尔方程那一步。我们已知RSA密文是如下方式产生的：

d=55607604802787651645008458776960381905274761803290637507632188635838921065681450951091689278291418693230528750002533602907652688436593464630306292523277204515196410606067221890882848225521819884919432609641732036959698443050773592482567332594684068998082886109352105201496925283866221090767083891241755981564899917062843679988341936588610298860835342859782896551931353504004447312763256966682142206683209406463516580194846276027818291321875085607298572702721743561*n_1+n_2-69462057634597116805075535702860983356198049

写直观一点如下：

$d = kn_1 + n_2 - t$

这一步本身并没有什么问题，因为k、t肯定是用RSA得到的密文以及佩尔方程的某一组解凑出来的两个正数。那么显然，随着n1、n2增大，得到的密文d也会增大。

而已知：

1	c=pow(m,e,n)=encrypt(hint_0)

这里也用的是”=”而不是”==”，无力吐槽。。。总之就是说，我们得到的密文应该是在(0,n)范围内的。所以我们需要先找到佩尔方程的最小正整数解，然后逐个解搜索，一直到得到的d不落在(0,n)范围内就退出。

但是问题就在于，求出的一组最小正整数解(n1,n2)得到的d也大于了n？什么情况？

这里简单阐述一下佩尔方程的求解思路，用连分数解决的思路很易懂，如下：

方程基本形式是：

$x^2-Dy^2 = 1$

就有：

$x^2 = Dy^2 + 1 \approx Dy^2$

所以：

$\sqrt{D} \approx \frac{x}{y}$

因此对根号D连分数展开，检验是否有x^2-Dy^2 = 1，第一次得到的就是最小正整数解。正常来说这样解是没什么问题的，但是这题里D=66666666，较大，这样解是不是会因为根号D的精度问题漏解？

然后我找到了这一篇参考文章：

连分数法解佩尔方程特解_连分数求解pell方程-CSDN博客

看了这个文章会发现，一个数的根号如果是无理数，那么这个无理数的连分数其实是可以用一个有循环节的连分数表示的，并且有算法能精确的计算出这个连分数，因为它不需要用到无理数的小数部分。

我照着这个方式实现了一遍根号表示的无理数连分数展开，这里打印连分数的前一百项示例：

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

def sol_Pell(D):
    a = []
    b = []
    c = []

    #init
    D_ = int(iroot(D,2)[0])
    b.append(D_)
    c.append(D-D_**2)

    for i in range(100):
        a.append((D_ + b[i])//c[i])
        b.append(a[i]*c[i] - b[i])
        c.append((D-b[i+1]**2) // c[i])
    a.insert(0,D_)
    print(a)

D = 66666666
sol_Pell(D)

然后与我们常用的连分数求解佩尔方程对比一下：

def solve_pell(N, numTry = 10000):
    ring = RealField(1000)
    cf = continued_fraction(sqrt(ring(N)))
    print(cf[:100])
    for i in range(numTry):
        denom = cf.denominator(i)
        numer = cf.numerator(i)
        if numer^2 - N * denom^2 == 1:
            return numer, denom
    return None, None

打印结果可以发现，这两种方法对于根号66666666的连分数展开是完全一样的，这也就是说精度已经足够展开得到精确的连分数循环表示。而由于参考文章中的一个定理：

那么也就是说，我们得到的n1、n2的确就是最小正整数解。而既然最小正整数解得到的d都大于n，那么由于d单增，之后迭代得到的其他所有正整数解都肯定大于n，因此肯定都不是正确解。

然后我的思路是：一组正整数解其实对应四组解，比如如果(x,y)是佩尔方程的解的话，那么以下三组肯定也都是佩尔方程的解：

$(x,-y),(-x,y),(-x,-y)$

但是显而易见，n1乘了一个大系数，所以要想d落在(0,n)，就只能改变n2为负数试试。但是肯定也能想到杯水车薪，因为n2完全影响不了前面那个大数的比特位数。

那到底哪里出了问题？我猜测了一下是不是n2、n1打反了，因为如果调整位置的话确实就可以刚好落在(0,n)内，但是事实上也出不了结果。

那难道数字打错了或者我还有别的疏忽？这就不得而知了。

总之还是放个exp，毕竟也许代码有点什么错误导致出不了结果，有兴趣的师傅可以看看：

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

n = 1130559951242142418553461557920917245278152787036994964813122956739956040733223190585188299825937660853090732755660557900622564346469528239989519875434612493124663925800482373019421482927567375862528994310090666571233446002163478309309649223886467303834842172338454971591276606602163799133000582891965785243482166772793400065793927508012418516719549212969303944155066541587875789636798394581666621711612432864870236407912486149901575077993702841083287354610398446101725184000000000000000000000000000000149
e = 65537
g = 1766592252951269548814352411327322447034194093921285273353477875204196603230641896039854934719468650093602325707751568
p = 50254221567892101225665855

#part1 factor n
a = 2
m = 2
while True:
    a = pow(a, m, n)
    pp = GCD(a-1, n)
    if pp != 1 and pp != n:
        break
    m += 1
qq = int(n // pp)
phi = (pp-1)*(qq-1)
d = inverse(e, phi)


#part2 Get hint0
L = Matrix(Zmod(p),[[9,1],
                    [1,0]]) 
init = vector(Zmod(p),[1,0])
hint0 = int((L^(g-1)*init)[0])


#part3 pell
#(n_1+hint_0)**2+1712463104290490642*n_1+4730060142800357260761833404838169046987798-66666666*n_2**2=0:
temp = hint0**2+4730060142800357260761833404838169046987798+1
k1 = (2*hint0+1712463104290490642)//2
#print(k1**2 == temp)
#(n1+k1)^2 -D*n2^2 = 1

def solve_pell(N, numTry = 10000):
    ring = RealField(10000)
    cf = continued_fraction(sqrt(ring(N)))
    for i in range(numTry):
        denom = cf.denominator(i)
        numer = cf.numerator(i)
        if numer^2 - N * denom^2 == 1:
            return numer, denom
    return None, None

D = 66666666
n1_,n2 = solve_pell(D)
n1 = n1_ - k1
#print((n1+hint0)**2+1712463104290490642*n1+4730060142800357260761833404838169046987798-66666666*n2**2==0)
k = 55607604802787651645008458776960381905274761803290637507632188635838921065681450951091689278291418693230528750002533602907652688436593464630306292523277204515196410606067221890882848225521819884919432609641732036959698443050773592482567332594684068998082886109352105201496925283866221090767083891241755981564899917062843679988341936588610298860835342859782896551931353504004447312763256966682142206683209406463516580194846276027818291321875085607298572702721743561
t = 69462057634597116805075535702860983356198049


#part4 Get flag
c = k*n1 + n2 - t
flag = long_to_bytes(int(pow(c,d,n)))
print(flag)

111

已解出，可见：

Crypto趣题-DLP | 糖醋小鸡块的blog (tangcuxiaojikuai.xyz)

Crand_all

题目来源：ASISCTF 2023

题目：

#!/usr/bin/env python3

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

def gen_prime(nbit):
	while True:
		p = 0
		for i in range(nbit, nbit>>1, -1):
			p += getRandomRange(0, 2) * 2 ** i
		p += getRandomNBitInteger(nbit>>3)
		if isPrime(p):
			return p

nbit = 512
p, q, r = [gen_prime(nbit) for _ in '012']
e, n = 65537, p * q * r

m = bytes_to_long(flag)
c = pow(m, e, n)

print(f'n = {n}')
print(f'c = {c}')

n = 306274031191812861390142322943384171567281235901023281150152366709090671749498337889684025397615153433878170972936066715221203963189667691979936604050835848397126937722483985215998814393176118687506861486835225543585451024226910726976901506973151291379128731667011697266671380831373270759872483773774444157994113393807135699243427889953159841567700099682080305355884531709396908787777504205775972655070971179462480070831508541674970198913788578939945648909705829
c = 223849092404967389418032074311452301651293919558152640065438289744845837981617978133447578227300156302182335073619698990940597295177229704048168320777498240448811390343120956888620737076325442001666883042713739861144512318389043548008886406065462992332747911568586879243279225955299270912019379515768839390586943094509691563525269138915388095236050068130199760650753758898265542401663546938899693685809537450793106770047460151231224343501919710408682427560812073

题目非常直白，给了一个用特殊方式生成素数的函数gen_prime，并用这个函数生成了三个素数p、q、r作为RSA私钥。要求还原明文。

所以重点在于发现gen_prime有什么漏洞，梳理一下他的加密方式其可以用如下方式表述

初始时，素数为512比特0
然后，对于高256位的每一位，随机填1或者0
最后，再加上一个64比特的数，核验为素数后返回

那么p、q、r也就都可以写成如下表达式：

$prime = 2^{256}x + e$

其中，x是一个大约为256比特的数，e是一个64比特的数。那么对于p、q、r的乘积就有：

$n = pqr = (2^{256}x_1 + e_1)(2^{256}x_2 + e_2)(2^{256}x_3 + e_3)$

而e1e2e3的乘积仅有192比特左右，因此n的低位其实就是e1e2e3的乘积，用factordb分解后进行合理组合，就能得到e1、e2、e3分别的值。

接下来就是如何由e1、e2、e3，也就是p、q、r的各自低位来还原p、q、r。这里我考虑的是copper。但是如果仅考虑单元copper的话，可以想到由于是三素数RSA，beta最多也只能取0.33比较合适，但是这样小根的理论上界只能有150-160比特，远远达不到我们需要的256比特。

然后我也尝试了二元和三元copper，会发现即使格的维数调的很大、耗时很久的情况下依然也出不了结果，并且多元copper的根的范围我并不是很会算，所以这个做法的可行性也打个问号。

这一部分我的代码是：

from Crypto.Util.number import *
from factordb.factordb import FactorDB
from tqdm import *
import itertools
from gmpy2 import iroot
 
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 []


n = 306274031191812861390142322943384171567281235901023281150152366709090671749498337889684025397615153433878170972936066715221203963189667691979936604050835848397126937722483985215998814393176118687506861486835225543585451024226910726976901506973151291379128731667011697266671380831373270759872483773774444157994113393807135699243427889953159841567700099682080305355884531709396908787777504205775972655070971179462480070831508541674970198913788578939945648909705829
c = 223849092404967389418032074311452301651293919558152640065438289744845837981617978133447578227300156302182335073619698990940597295177229704048168320777498240448811390343120956888620737076325442001666883042713739861144512318389043548008886406065462992332747911568586879243279225955299270912019379515768839390586943094509691563525269138915388095236050068130199760650753758898265542401663546938899693685809537450793106770047460151231224343501919710408682427560812073

n_low = 1085715048967907895301761852456620799927891682857410093669
f = FactorDB(n_low)
f.connect()
fs = f.get_factor_list()
#print(fs)

p_low = 9525340967214022543
q_low = 540943066242511*71*293
r_low = 13989332051*724037101

PR.<p_high,q_high,r_high> = PolynomialRing(Zmod(n))
f = (p_high*2^256+p_low)*(q_high*2^256+q_low)*(r_high*2^256+r_low)
bound = (2^260, 2^260,2^260)
res = small_roots(f, bound, m=7,d=3)
print(res)

然后我就在想既然有一个等式：

$n = (2^{256}x_1 + e_1)(2^{256}x_2 + e_2)(2^{256}x_3 + e_3)$

以及每个xi的范围，有没有可能说要自己构造格，我尝试了一下构造下面的格：

$L = \left(\begin{matrix} K^3 &1&&&&&&&& \\ K^2e_1 &&1&&&&&&& \\ K^2e_2 &&&1&&&&&&\\ K^2e_3 &&&&1&&&&& \\ Ke_1e_2 &&&&&1&&&&\\ Ke_2e_3 &&&&&&1&&&\\ Ke_1e_3 &&&&&&&1&&\\ e_1e_2e_3 &&&&&&&&1&\\ n &&&&&&&&&1\\ \end{matrix} \right)$

其中K=2^256，这个格具有如下线性关系：

$(x_1x_2x_3,x_2x_3,x_1x_3,x_1x_2,x_3,x_1,x_2,1,-1) \left(\begin{matrix} K^3 &1&&&&&&&& \\ K^2e_1 &&1&&&&&&& \\ K^2e_2 &&&1&&&&&&\\ K^2e_3 &&&&1&&&&& \\ Ke_1e_2 &&&&&1&&&&\\ Ke_2e_3 &&&&&&1&&&\\ Ke_1e_3 &&&&&&&1&&\\ e_1e_2e_3 &&&&&&&&1&\\ n &&&&&&&&&1\\ \end{matrix} \right) = (0,x_1x_2x_3,x_2x_3,x_1x_3,x_1x_2,x_3,x_1,x_2,1,-1)$

然后将目标向量配成相等数量级，并给格的第一列乘上大系数T以保证规约出0。但是结果表明我们需要的目标向量无法规约出。这一部分的代码如下：

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

n = 306274031191812861390142322943384171567281235901023281150152366709090671749498337889684025397615153433878170972936066715221203963189667691979936604050835848397126937722483985215998814393176118687506861486835225543585451024226910726976901506973151291379128731667011697266671380831373270759872483773774444157994113393807135699243427889953159841567700099682080305355884531709396908787777504205775972655070971179462480070831508541674970198913788578939945648909705829
c = 223849092404967389418032074311452301651293919558152640065438289744845837981617978133447578227300156302182335073619698990940597295177229704048168320777498240448811390343120956888620737076325442001666883042713739861144512318389043548008886406065462992332747911568586879243279225955299270912019379515768839390586943094509691563525269138915388095236050068130199760650753758898265542401663546938899693685809537450793106770047460151231224343501919710408682427560812073

n_low = 1085715048967907895301761852456620799927891682857410093669
f = FactorDB(n_low)
f.connect()
fs = f.get_factor_list()
#print(fs)

e1 = 9525340967214022543
e2 = 540943066242511*71*293
e3 = 13989332051*724037101

K = 2^256
T = 2^1000
L = Matrix(ZZ,9,10)
L[0,0] = K^3
L[1,0] = K^2*e1
L[2,0] = K^2*e2
L[3,0] = K^2*e3
L[4,0] = K*e1*e2
L[5,0] = K*e2*e3
L[6,0] = K*e1*e3
L[7,0] = e1*e2*e3
L[8,0] = -n
for i in range(9):
    L[i,0] *= T
L[0,1] = 1
L[1,2] = K
L[2,3] = K
L[3,4] = K
L[4,5] = K^2
L[5,6] = K^2
L[6,7] = K^2
L[7,8] = K^3
L[8,9] = K^3

res = L.LLL()

之后我在想，其实我们还有一些已知信息没有用上，这些已知信息仍然是基于p、q、r的特殊组成方式而言的，由于他们的等式：

$n = pqr = (2^{256}x_1 + e_1)(2^{256}x_2 + e_2)(2^{256}x_3 + e_3)$

把2^256视作是移位，那么其实由n的各个比特，我们可以知道以下几个额外信息：

x1x2x3的高444位
x1e2e3+x2e1e3+x3e1e2的低256位
x1x3e2+x2x3e1+x1x2e3的中间128位

而这些信息应该对我们造copper或者造格要有一些优化。我们现在知道，如果我们想要求单元copper，其能接受的最大上界最多到160比特左右，那么由这些信息，如果我们能确定x1、x2、x3的100个左右的比特位，我们就可以单元copper求解出n的分解，因此这几个信息目前看来我觉得有用的是：

x1x2x3的高444位
x1e2e3+x2e1e3+x3e1e2的低256位

这个能否用来近似确定x1、x2、x3的分别的高100位或者低100位左右？但目前我还没想出具体怎么操作。

然后还有一些比较离奇的想法：

第一个就是直接z3，试了过后发现z3虽然一般来说很强，但确实也还没有这么强。代码如下，这几个约束想直接解是解不了的(数据是我生成的测试数据)：

from random import randint
from z3 import *

e = [10333799336736048411, 10048938500043630937, 13770567448032700553]
#x1 = 134530121267138061692122246190685342329406941720572324580169098024935310581168
#x2 = 123427567449792804014455419800830894394580447784489056793113376627868537160216
#x3 = 54022871114718186939252665659796781248148877214115427587582690293386954260174
n = 1392662992231878453618153650352932233685479923380378591186749945448215216646971207443890669678951977168432541888832218967441146964183473856940088594238372985102279239596228498511057175872182183657147834996556681844161888332641986487531352611082243721457395282059148273836559101486130523166444160903027354883310771050656525852022898214928288842219786550934402364392466948737121570741902353211576548652037840953521222373985396759933500395549254979602390182860639739

#z3solve
def z3_solve():
    x1 = Int("x1")
    x2 = Int("x2")
    x3 = Int("x3")
    solver = Solver()

    t1 = n >> (1092)
    solver.add((x1*x2*x3 / 2**(768 - (1536-1092))) == t1)

    t2 = ((n >> 256) & ((1 << 256) - 1))
    solver.add(((e[1]*e[2]*x1)+(e[0]*e[2]*x2)+(e[0]*e[1]*x3)) % 2**256 == t2)

    K = 2**256
    solver.add((K*x1+e[0])*(K*x2+e[1])*(K*x3+e[2])  == n)

    solver.add(x1 < K)
    solver.add(x2 < K)
    solver.add(x3 < K)

    if solver.check()==sat:
        print(solver.model())

z3_solve()