零知识证明一
原文：Quadratic Arithmetic Programs: from Zero to Hero
简介：本文是Vitalik写于2016年12月，用于介绍零知识证明的数学实现方式的论文。文章思路清晰，通俗易懂，也因此，该文成为区块链行业技术人员学习这方面知识的首选文章之一。
1 术语介绍
解决一个问题需要花费时间。如果解决问题需要的时间与问题的规模之间是多项式关系，则可以称该问题具有多项式复杂度。
P 问题 指的是在多项式时间内可解的问题。
NP 问题 （Non-Deterministic Polynomial Problem，非确定性多项式问题），指不能在多项式内可解，但是可以在多项式时间内验证的问题。
QSP 问题 （QSP - Quadratic Span Program），实现基于布尔电路的 NP 问题的证明和验证。
QAP 问题 （Quadratic Arithmetic Program），实现基于算术电路的 NP 问题的证明和验证，相对于 QSP，QAP 有更好的普适性。
2 详解
最近人们对 zk-SNARKs（零知识证明）背后的技术有很多兴趣，人们越来越多地试图去揭开一些被许多人称为“月球数学”的东西的神秘面纱，因为人们认为它的复杂性非常难以理解。zk-SNARKs的理解确实相当具有挑战性，尤其是由于整个系统需要组装起来才能工作的移动部件太多，但如果我们把这项技术一件一件地分解，那么理解起来就会变得更简单。
这篇文章的目的不是用于完整的介绍 zk-SNARKs，它假定您具有以下背景知识：
1 你知道 zk-SNARKs 和他的大致原理；（译者注：如果不知道 zk-SNARK，建议您可以参考《一个数独引发的惨案：零知识证明》）
2 你有足够的数学知识，能理解一些基本的多项式知识。(如if P(x)+Q(x)=(P + Q)(x)，P和Q代表多项式，如果你对这类多项式表述方式已经非常熟悉，说明你符合继续阅读的要求)。
如上图，可以将以上零知识证明分为由上至下的多个阶段：计算，代数电路，R1CS（即rank-1 constraint system，秩为1的约束系统，可将其理解为一个方程组），QAP（二次算术程序），Linear PCP。首先，zk-SNARK不能直接应用于任何计算问题；相反，您必须将问题转换为操作的正确“形式”。这种形式被称为“二次算术程序”(QAP)，将函数的代码转换成这些代码本身就非常重要。与将函数代码转换为 QAP的过程一起运行的还有另一个过程，这样，如果对代码有输入，就可以创建相应的解决方案(有时称为QAP的“见证”)。这就是本文需要讲述的内容。
在此之后，还有另一个相当复杂的过程来为这个QAP创建实际的“零知识证明”，还有一个单独的过程来验证别人传给你的证据，但是这些细节超出了本文的范围。
在下面示例中，我们将选择一个非常简单的问题：求一个三次方程的解：x^3 + x + 5 == 35(提示:答案是 3)。这个问题很简单，由此产生的QAP不会大到令人生畏，但是也很重要的，你可以由此案例看到所有的功能是如何发挥作用的。
用编程语言描述以上方程如下：
def qeval(x): y = x^3 return x + y + 5
我们在这里使用的简单编程语言支持基本的算术(+、-、*、/)、常量幂指数(如x^7，但不是x^y)和变量赋值，这足够强大到理论上可以在其中进行任何计算(只要计算步骤的数量是有界的;不允许循环)。注意模(%)和比较运算符(<、>、≤≥)不支持,因为没有有效的方法做模或直接比较有限循环群算法。
您可以通过位分解来将语言扩展到模和比较，（例如：13 = 2^3 + 2^2 + 1 = 8 + 4 + 1)作为辅助输入，证明这些分解的正确性，并在二进制电路中进行数学运算；在有限域算法中，执行等式(==)检查也是可行的，实际上更容易一些，但这两个细节我们现在都不讨论。我们可以扩展语言来支持条件句（例如将语句：if x < 5: y = 7; else: y = 9; 转换为算术形式：y = 7 * ( x < 5) + 9 * (x >= 5);）不过请注意，条件的两个“路径”都需要执行，如果您有许多嵌套的条件，那么这会导致大量开销。
现在让我们一步一步地经历这个过程。如果你想自己做任何代码，我在这里用 Python 实现了一段代码（仅用于教育目的；还没有准备好为现实世界的zk-SNARK制作QAPs!）
2.1 扁平化
第一步是一个“压扁”的过程，我们把原来的代码（这可能包含任意复杂的语句和表达式）分解为最简单的表达式，这种表达式有两种形式：
1 x = y（y 可以是变量或数字)
2 x = y (op) z (op 可以+,-,*,/，y 和 z 可以是变量，数字或子表达式)
你可以把这些表述看成是电路中的逻辑门，上述表达式x^3 + x + 5的扁平化过程结果如下：
sym_1 = x * x y = sym_1 * x //相当于实现了幂函数y = x^3 sym_2 = y + x ~out = sym_2 + 5
你可以认为上述的每一行声明都是一个电路中的逻辑门，4个表达式，为4个约束条件或者4个门电路，与原始代码相比，这里我们引入了中间变量sym_1，y和sym_2，还有一个表示输出的冗余变量~out，不难看出“压扁”后的声明序列和原始代码是等价的。
2.2 R1CS转换
现在，我们把它转换成一个称为R1CS（Rand-1 Constraint System）的东西。R1CS是由三个向量(a, b, c)组成的序列，R1CS的解是一个向量s，其中s必须满足方程
s . a * s . b - s . c = 0
其中.代表内积运算。例如，以下是一个令人满意的 R1CS：
a = (5,0,0,0,0,1) b = (1,0,0,0,0,0) c = (0,0,1,0,0,0) s = (1,3,35,9,27,30)
（第一个35=1*5 + 30*1，第二个35=35*1）
上述例子只是一个约束，接下来我们要将每个逻辑门（即“压扁”后的每一个声明语句）转化成一个约束（即一个(a, b, c)三向量组），转化的方法取决于声明是什么运算 (+,-,*,/) 和声明的参数是变量还是数字。在我们这个例子中，除了“压扁”后的五个变量 ('x', '~out', 'sym_1', 'y', 'sym_2') 外，还需要在第一个分量位置处引入一个冗余变量~one来表示数字 1，就我们这个系统而言，一个向量所对应的 6 个分量是(可以是其他顺序，只要对应起来即可)：
s = (1, x, ~out, sym1, y, sym2)
第一个约束（第一个门电路）
sym_1 = x * x，即 x*x - sym_1 = 0
由于s向量是已知的，我们列出以下表格。
即(1·0+x·1+~out·0+sym_1·0+y·0+sym_2·0)*(1·0+x·1+~out·0+sym_1·0+y·0+sym_2·0)-(1·0+x·0+~out·0+sym_1·1+y·0+sym_2·0)=x^2-sym_1=0，正好是第一个约束，所以可以得到如下向量组：
a1 = [0, 1, 0, 0, 0, 0] b1 = [0, 1, 0, 0, 0, 0] c1 = [0, 0, 0, 1, 0, 0]
如果解向量s的第二个标量是3，第四个标量是9，无论其他标量是多少，则最终的点乘会归结为3*3-9=0，即满足方程。同样，如果s的第二个标量是7，第四个标量是49，也会通过检查，第一次检查仅仅是为了验证第一个门的输入和输出的一致性。
第二个约束（第二个门电路）
y = sym_1 * x，即 sym_1 * x - y = 0
可以得到以下向量组：
a2 = [0, 0, 0, 1, 0, 0] b2 = [0, 1, 0, 0, 0, 0] c2 = [0, 0, 0, 0, 1, 0]
第三个约束（第三个门电路）
sym_2 = y + x，加法门需要转换为：(x + y) * 1 - sym_2 = 0
得到以下向量组：
a3 = [0, 1, 0, 0, 1, 0] 对应(x + y) b3 = [1, 0, 0, 0, 0, 0] 对应常量1 c3 = [0, 0, 0, 0, 0, 1]
第四个约束（第四个门电路）
~out = sym_2 + 5，即 (sym_2 + 5) * 1 - ~out = 0
得到以下向量组：
a4 = [5, 0, 0, 0, 0, 1] b4 = [1, 0, 0, 0, 0, 0] c4 = [0, 0, 1, 0, 0, 0]
现在，我们假设x = 3，根据第一个门，得到sym_1 = 9，根据第二个门得到y = 27，根据第三个门，得到sym_2 = 30，根据第四个门得到~out = 35，因此，可以得到：
s = [1, 3, 35, 9, 27, 30]
如果假设不同的 x，都可以得到不同的 s，但所有 s 都可以用来验证(a, b, c)，现在我们得到了四个约束的 R1CS，完整的 R1CS 如下：
A [0, 1, 0, 0, 0, 0] [0, 0, 0, 1, 0, 0] [0, 1, 0, 0, 1, 0] [5, 0, 0, 0, 0, 1] B [0, 1, 0, 0, 0, 0] [0, 1, 0, 0, 0, 0] [1, 0, 0, 0, 0, 0] [1, 0, 0, 0, 0, 0] C [0, 0, 0, 1, 0, 0] [0, 0, 0, 0, 1, 0] [0, 0, 0, 0, 0, 1] [0, 0, 1, 0, 0, 0]
2.3 从R1CS到QAP
下一步是将这个R1CS转换成QAP形式，它实现了完全相同的逻辑，只是使用多项式而不是内积。我们是这样做的：从4组长度为6的3个向量到6组长度为3的多项式，在每个x坐标处求多项式代表一个约束条件。也就是说，如果我们求出x=1处的多项式，我们就得到了第一组向量，如果我们求出x=2处的多项式，我们就得到第二组向量，以此类推。
我们可以用拉格朗日插值来做这个变换。拉格朗日插值法解决的问题是：如果你有一组点(即(x, y)坐标对），然后对这些点做拉格朗日插值得到一个经过所有这些点的多项式。我们通过分解问题：对于每个 x 坐标，我们创建一个多项式，在x座标处有相应的y坐标，在其他坐标我们不感兴趣的x处y坐标为0，然后我们将所有的多项式相加得到最终的多项式。
让我们做一个例子。假设我们想要一个多项式经过(1,3)，(2,2)和(3,4)。我们首先做一个多项式，经过(1,3)，(2,0)和(3,0)。事实证明,一个多项式，“伸出”x = 1 在其他点纵坐标0，很容易得到以下方程：
(x - 2) * (x - 3)
如图所示：
然后，在 y 轴方向“拉伸”，得如下多项式：
(x - 2) * (x - 3) * 3 / ((1 - 2) * (1 - 3))
即
1.5 * x^2 - 7.5 * x + 9
在其他两点做相同得操作，使得他们分别在x=2和x=3时“伸出”，会得到两个相似的多项式
(x-1)*(x-3)*2/((2-1)*(2-3)) (x-1)*(x-2)*4/((3-1)*(3-2))
将这两个多项式和之前得到的那个多项式相加得：
1.5 * x^2 - 5.5 * x + 7
最终函数方程如下图：
上述算法需要O(n^3)时间，因为有n个点，每个点都需要O(n^2)时间将多项式相乘。稍微思考一下，这就可以减少到O(n^2)的时间，再多思考一下，使用快速的傅里叶变换算法等等，它可以进一步减少——这是一个关键的优化，当在zk- spuks中使用的函数通常有成千上万个门时。
拉格朗日插值的相关理论可以参考这边文章，这里直接给出插值公式，通过n个点(x1,y1),(x2,y2),...,(xn,yn) 的n-1阶多项式为：
例如上例中，通过点(1,3), (2,2), (3,4)的多项式为：
现在我们要将A，B，C的6个列向量转化为六组多项式，每组多项式包括三个三阶多项式，我们在每个x点处来评估不同的约束，在这里，我们共有四个约束，因此我们分别用多项式在 x = 1,2,3,4 处来评估这四个向量组。
现在我们使用拉格朗日差值公式来将 R1CS 转化为 QAP 形式。我们先求出四个约束所对应的每个 a 向量的第一个值的多项式，也就是说使用拉格朗日插值定理求过点(1,0), (2,0), (3,0), (4,5)的多项式，类似的我们可以求出其余的四个约束所对应的每个向量的第 i 个值的多项式。这里，直接给出答案：
A polynomials [-5.0, 9.166, -5.0, 0.833] [8.0, -11.333, 5.0, -0.666] [0.0, 0.0, 0.0, 0.0] [-6.0, 9.5, -4.0, 0.5] [4.0, -7.0, 3.5, -0.5] [-1.0, 1.833, -1.0, 0.166] B polynomials [3.0, -5.166, 2.5, -0.333] [-2.0, 5.166, -2.5, 0.333] [0.0, 0.0, 0.0, 0.0] [0.0, 0.0, 0.0, 0.0] [0.0, 0.0, 0.0, 0.0] [0.0, 0.0, 0.0, 0.0] C polynomials [0.0, 0.0, 0.0, 0.0] [0.0, 0.0, 0.0, 0.0] [-1.0, 1.833, -1.0, 0.166] [4.0, -4.333, 1.5, -0.166] [-6.0, 9.5, -4.0, 0.5] [4.0, -7.0, 3.5, -0.5]
这些系数是升序排序的，例如上述第一个多项式是0.833*x^3 - 5*x^2 + 9.166*x - 5。
（计算过程： 将 (1,0), (2,0), (3,0), (4,5)，带入拉格朗日插值公式计算得出，第二个多项式则是将 (1,1),(2,0),(3,1),(4,0) 带入拉格朗日插值公式计算得出, 由此类推）。
这组多项式（加上我将在后面解释的Z多项式）构成了这个特定QAP实例的参数。请注意，对于尝试使用zk-SNARKs验证的每个功能，到此为止的所有工作只需要完成一次；一旦生成了QAP参数，就可以重复使用它们。
让我们尝试在x=1处评估所有这些多项式。 在x=1处评估多项式仅意味着将所有系数相加（对于所有k，1^k = 1），因此并不困难。 我们得到：
如果我们将 x=1 带入上述十八个多项式，可以得到第一个约束的三个向量：
A results at x=1 0 1 0 0 0 0 B results at x=1 0 1 0 0 0 0 C results at x=1 0 0 0 1 0 0
可以看到，我们在这里所拥有的与我们在上面创建的第一个逻辑门的三个向量的集合完全相同。
（以上计算过程： x=1 代入第一个多项式0.833 * x^3 - 5 * x^2 + 9.166 * x - 5 = 0，x=1 代入第二个多项式-0.666 * x^3 + 5 * x^2 - 11.333 * x + 8 =1, 由此类推）
2.4 检查QAP
通过将R1CS转换成QAP，我们可以通过多项式的内积运算来同时检查所有的约束，而不是像R1CS那样单独的检查每一个约束。如下图所示：
因为在这种情况下，点积检验是一系列多项式的加法和乘法，结果本身就是一个多项式。如果得到的多项式，在我们上面用来表示逻辑门的每一个x坐标处的值，等于0，那就意味着所有的检查都通过了；如果结果多项式至少有一个非零值，那么这就意味着进出逻辑门的值是不一致的。（例如，门为y = x * sym_1，但是提供的值为x = 2，sym_1 = 2，y = 5）。
值得注意的是，得到的多项式本身不一定是零，事实上在大多数情况下不是；它可以在不符合任何逻辑门的点上有任何行为，只要在所有符合某些门的点上结果是零。为了验证正确性，我们不计算多项式t = A . s * B . s - C . s在每一点对应一个门；相反，我们把t除以另一个多项式Z，然后检查Z是否均匀地除t，也就是说，t/Z没有余数。
Z定义为(x - 1) * (x - 2) * (x - 3)…-最简单的多项式，在所有对应逻辑门的点上都等于0。这是代数的一个基本事实任何多项式在所有这些点上等于零都必须是这个最小多项式的倍数，如果一个多项式是Z的倍数那么它在任何这些点上的值都是零；这种对等使我们的工作容易得多。
现在，让我们用上面的多项式做内积检验。
首先，我们得到中间多项式：
A . s = [43.0, -73.333, 38.5, -5.166] B . s = [-3.0, 10.333, -5.0, 0.666] C . s = [-41.0, 71.666, -24.5, 2.833]
以上计算过程：
43.0 = -5 * 1 + 8 * 3 + 0 * 35 - 6 * 9 + 4 * 27 - 1 * 30
-73.333 = 9.166 * 1 - 11.333 * 3 + 0 * 35 + 9.5 * 9 - 7 * 27 + 1.833 * 30
...
-3 = 3 * 1 - 2 * 3 + 0 * 35 + 0 * 9 + 0 * 27 + 0 * 30
...
以上多项式经过：A . s * B . s - C . s 计算后得到：
t = [-88.0, 592.666, -1063.777, 805.833, -294.777, 51.5, -3.444]
计算过程：
A · s = [43.0, -73.333, 38.5, -5.166] = -5.166x^3 + 38.5x^2 - 73.333x + 43
B · s = [-3.0, 10.333, -5.0, 0.666] = 0.666x^3 - 5x^2 + 10.333x - 3.0
C · s = [-41.0, 71.666, -24.5, 2.833] = 2.833x^3 - 24.5x^2 + 71.666x - 41.0
A · s * B · s - C · s 就是上面多项式的计算，计算后，按幂从低到高排列系数，得到： [-88.0, 592.666, -1063.777, 805.833, -294.777, 51.5, -3.444]
最小多项式为：
Z = (x - 1) * (x - 2) * (x - 3) * (x - 4)
即：
Z = [24, -50, 35, -10, 1]
现在计算多项式相除：
h = t / Z = [-3.666, 17.055, -3.444]
没有任何余数。
我们有了QAP的解。如果我们试图伪造R1CS中的变量，而这个R1CS推导出了QAP解决方案——比如，将s的最后一个数字设为31，而不是30，我们检查多项式t时将会不通过(在特定情况下,在 x = 3 处的结果将等于 -1 而不是 0)，而且t不会是Z的倍数; 相反，除以t / Z会得到[-5.0, 8.833, -4.5, 0.666]的余数。
注意，以上只是一个非常简单的示例；在现实世界中，加减乘除运算通常伴随着非常规的数字，所以所有的我们知道和喜爱的代数定律还是有用的，但是所有答案是一些有限大小集合的元素，通常是从 0 到 n - 1 范围内的整数 n。例如，如果 n = 13，则 1 / 2 = 7 (7 * 2 = 1，7为2的逆), 3 * 5 = 2 等等。使用有限域算法消除了对舍入误差的担心，并允许系统与椭圆曲线很好地工作，这最终对使 zk-SNARK 协议变得真正安全。
3 源码
用python代码实现了大致流程一并给出
class Point: def __init__(self, x, y): self.x = x self.y = y def __repr__(self): return f'Point({self.x}, {self.y})' def add_poly(L1,L2): #多项式加法，同次项系数相加 R=[] if len(L1)>len(L2):#默认L2比较长 L1,L2=L2,L1 i=0 while i<len(L1): R.append(L1[i]+L2[i])#从低次项开始对应相加 i+=1 R=R+L2[len(L1):len(L2)]#较长的多项式高次项直接复制 return R def subtract_poly(L1,L2): #多项式减法 L2=L2[:]#为了不改变传入的L2 for i in range(len(L2)): L2[i]=-L2[i] return(add_poly(L1,L2)) def multiply_poly(L1,L2):#多项式乘法 if len(L1)>len(L2): L1,L2=L2,L1 zero=[];R=[] for i in L1: T=zero[:]#存储中间产生的结果多项式，每次更新结果多项式的列表长度 for j in L2:#一个单项式乘以多项式的每一项 T.append(i*j) R=add_poly(R,T) zero=zero+[0]# 每一个新的多形式都要比前一个多项式次数高1，列表长度增加，所以多补一个0 return R def divide_poly(L1,L2): if len(L1)<len(L2): return 0,L1 #商为0的情况 d=len(L1)-len(L2) #循环除的次数为d+1 T=L1[:] #被除数 R=[] #商 for i in range(d+1): n=T[len(T)-1]/L2[len(L2)-1]#被除数的最高次项系数除以除数的最高次项的系数 R=[n]+R T1=[0]*(d-i)+[n]#对得到的商补若干个0 T2=multiply_poly(T1,L2)#商乘以除数结果为T2 T=subtract_poly(T,T2)#被除数减去T2 T=T[:len(T)-1]#切片，被除数减少一次方 return R,T def getPolynomialCoef(p1, p2, p3, p4): coef=p1.y/((p1.x-p2.x)*(p1.x-p3.x)*(p1.x-p4.x)) x3coef=coef x2coef=-(p2.x+p3.x+p4.x)*coef x1coef=(p2.x*p3.x+p2.x*p4.x+p3.x*p4.x)*coef x0coef=-p2.x*p3.x*p4.x*coef coef=p2.y/((p2.x-p1.x)*(p2.x-p3.x)*(p2.x-p4.x)) x3coef+=coef x2coef+=-(p1.x+p3.x+p4.x)*coef x1coef+=(p1.x*p3.x+p1.x*p4.x+p3.x*p4.x)*coef x0coef+=-p1.x*p3.x*p4.x*coef coef=p3.y/((p3.x-p1.x)*(p3.x-p2.x)*(p3.x-p4.x)) x3coef+=coef x2coef+=-(p1.x+p2.x+p4.x)*coef x1coef+=(p1.x*p2.x+p1.x*p4.x+p2.x*p4.x)*coef x0coef+=-p1.x*p2.x*p4.x*coef coef=p4.y/((p4.x-p1.x)*(p4.x-p2.x)*(p4.x-p3.x)) x3coef+=coef x2coef+=-(p1.x+p2.x+p3.x)*coef x1coef+=(p1.x*p2.x+p1.x*p3.x+p2.x*p3.x)*coef x0coef+=-p1.x*p2.x*p3.x*coef return x0coef, x1coef, x2coef, x3coef if __name__ == '__main__': A = [ [0, 1, 0, 0, 0, 0], [0, 0, 0, 1, 0, 0], [0, 1, 0, 0, 1, 0], [5, 0, 0, 0, 0, 1] ] B = [ [0, 1, 0, 0, 0, 0], [0, 1, 0, 0, 0, 0], [1, 0, 0, 0, 0, 0], [1, 0, 0, 0, 0, 0] ] C = [ [0, 0, 0, 1, 0, 0], [0, 0, 0, 0, 1, 0], [0, 0, 0, 0, 0, 1], [0, 0, 1, 0, 0, 0] ] print("A\n{}\nB\n{}\nC\n{}".format(A,B,C)) x = 1 a = [] polyA = [] print("\nA polynomials==================") for i in range(6): p1=Point(1, A[0][i]) p2=Point(2, A[1][i]) p3=Point(3, A[2][i]) p4=Point(4, A[3][i]) poly = getPolynomialCoef(p1, p2, p3, p4) print(poly) y=poly[0]+poly[1]*pow(x,1)+poly[2]*pow(x,2)+poly[3]*pow(x,3) a.append(format(y, ".2f")) polyA.append(poly) b = [] polyB = [] print("\nB polynomials==================") for i in range(6): p1=Point(1, B[0][i]) p2=Point(2, B[1][i]) p3=Point(3, B[2][i]) p4=Point(4, B[3][i]) poly = getPolynomialCoef(p1, p2, p3, p4) print(poly) y=poly[0]+poly[1]*pow(x,1)+poly[2]*pow(x,2)+poly[3]*pow(x,3) b.append(format(y, ".2f")) polyB.append(poly) c = [] polyC = [] print("\nC polynomials==================") for i in range(6): p1=Point(1, C[0][i]) p2=Point(2, C[1][i]) p3=Point(3, C[2][i]) p4=Point(4, C[3][i]) poly = getPolynomialCoef(p1, p2, p3, p4) print(poly) y=poly[0]+poly[1]*pow(x,1)+poly[2]*pow(x,2)+poly[3]*pow(x,3) c.append(format(y, ".2f")) polyC.append(poly) print("\na {}".format(a)) print("\nb {}".format(b)) print("\nc {}".format(c)) #print(polyA) s = [1, 3, 35, 9, 27, 30] As = [] AsF = [] for i in range(4): value = s[0]*polyA[0][i]+s[1]*polyA[1][i]+s[2]*polyA[2][i]+s[3]*polyA[3][i]+s[4]*polyA[4][i]+s[5]*polyA[5][i] As.append(value) AsF.append(format(value, ".3f")) print("\nAs {}".format(AsF)) Bs = [] BsF = [] for i in range(4): value = s[0]*polyB[0][i]+s[1]*polyB[1][i]+s[2]*polyB[2][i]+s[3]*polyB[3][i]+s[4]*polyB[4][i]+s[5]*polyB[5][i] Bs.append(value) BsF.append(format(value, ".3f")) print("\nBs {}".format(BsF)) Cs = [] CsF = [] for i in range(4): value = s[0]*polyC[0][i]+s[1]*polyC[1][i]+s[2]*polyC[2][i]+s[3]*polyC[3][i]+s[4]*polyC[4][i]+s[5]*polyC[5][i] Cs.append(value) CsF.append(format(value, ".3f")) print("\nCs {}".format(CsF)) rs = multiply_poly(As, Bs) rs = subtract_poly(rs, Cs) print("As*Bs-Cs\n{}".format(rs)) div = multiply_poly(multiply_poly([-1,1], [-2,1]), multiply_poly([-3,1], [-4,1])) print("divisor {}".format(div)) Quotient, Remainder = divide_poly(rs, div) print("Quotient\n[", end=' ') for item in Quotient: print(format(item, ".3f"), end=' ') print("]\nRemainder\n[", end=' ') for item in Remainder: print(format(item, ".2f"), end=' ') print("]")
参考：
