本文介绍了平方乘算法,它是一种在O(log n)时间内计算整数幂的算法,相比于朴素算法的O(n)时间复杂度更高效。文章解释了该算法的原理,包括平方指数序列和如何利用指数的二进制表示来选择正确的平方指数序列元素。此外,还探讨了如何将该算法应用于具有分数指数和固定点数的场景,并提供了Python和Solidity的实现示例,最后提到了Uniswap V3中使用平方乘算法的例子。
平方再乘算法在 $O(logn)$ (对数时间)内计算整数指数。 计算指数 $x^n$ 的简单方法是将 x 与自身相乘 n 次,这将需要 $O(n)$ 时间来计算。
假设我们要计算 $x^{20}$。 不是将 x 与自身相乘 20 次,而是从 x 开始并重复平方,直到计算出 $x^{16}$:
$$ \begin{align} x^2&=x\cdot x\ x^4&=x^2\cdot x^2\ x^8&=x^4\cdot x^4\ x^{16}&=x^8\cdot x^8 \end{align} $$
观察到我们可以将 $x^{20}$ 计算为:
$$x^{20}=x^4x^{16}$$
指数以这种方式“相加”的事实在代数中称为幂的乘积规则。
本文是 关于 Uniswap V3 代码库的系列文章 的一部分。 Uniswap V3 使用平方再乘算法将 tick 索引转换为平方根价格。 但是,本文也适用于任何想学习平方再乘算法工作原理的人。
我们将在本文中多次使用术语“平方指数序列”,即 $x^1,x^2,x^4,x^8,\dots$。 每一项都是前一项的平方。 基数 x 可以有任何非零值。 例如,$0.7^1,0.7^2,0.7^4,…$ 是一个有效的平方指数序列。 类似地,$-1.56,-1.56^2,…,-1.56^{32},-1.56^{64}$ 也是一个有效的平方指数序列。
如果我们预先知道基数将是一个固定值,并且我们可能计算的指数具有已知的上限,那么我们可以预先计算平方指数序列,然后进行查找而不是乘法。 例如,如果我们想要计算 $0.7^x$,并且我们预先知道 x 不会大于 63,我们可以预先计算以下值:
$$ \begin{align*} &0.7\ &0.7^2=0.49\ &0.7^4=0.2401\ &0.7^8=0.05764801\ &0.7^{16}=0.00332329\ &0.7^{32}=0.00001104
\end{align*} $$
然后,例如,如果我们想计算 $0.7^{44}$,我们可以将平方指数序列中的适当的预计算值相乘,如下所示:
$$ \begin{align} 0.7^{44}&=0.7^{32}\times0.7^{8}\times0.7^{4}\ 0.00000015&=0.00001104\times0.05764801\times0.2401 \end{align} $$
鼓励读者在他们的计算器上重做 $0.7^{44}$ 以仔细检查这一点。
当我们以这种方式缓存值时,我们必须提前知道我们的应用程序必须处理的指数的上限。
假设我们不是计算 $0.7^x$,而是想要计算 $\sqrt[3]{0.7^x}$? 这等效于 $0.7^{x/3}$。 我们仍然可以使用平方再乘算法,因为指数除以一个常数,所以幂的乘积规则仍然适用。 换句话说,
$$ 0.7^{44/3}=0.7^{32/3}\times0.7^{8/3}\times0.7^{4/3} $$
要计算 $0.7^{44/3}$,我们将更改 0.7 的幂,我们预先计算如下:
然后,要计算 $0.7^{44/3}$,我们将预先计算的指数相乘,如下所示:
$0.7^{44/3} = 0.7^{32/3} \times 0.7^{8/3} \times 0.7^{4/3}$ $0.005346931087848946 = 0.02226960053248208 \times 0.38630302299215685 \times 0.6215328012198205$
同样,鼓励读者自己重做这些计算。
如果想计算 b44,我们如何确定我们需要平方指数序列中的哪些 b 的幂? 换句话说,我们如何快速确定我们需要 $ b^{32}, b^{8}, b^{4} $而不是 $b^{32}, b^{4}, b^{1}$?
给定一个目标和,例如 44,我们如何快速确定 32、8 和 4 是和? 或者,假设我们正在尝试计算 b37。 我们需要的平方指数序列的元素是 $b^{32},b^{4},b^{1}$ - 我们如何快速确定我们需要 $b^{32},b^{4},b^{1}$,而不是 $b^{2},b^{8},b^{16}$?
一种简单的方法是从最大的预计算指数向下进行线性搜索,并减去不超过目标的最大的一个。 例如,如果我们正在计算 544,这意味着我们已经计算了 5 的平方指数序列。我们向下扫描并发现 32 是最接近 44 的值。然后我们再次向下扫描以找到 16,但请注意,添加 32+16 超过了 44,因此我们跳过 16 并继续执行 8,依此类推。
我们可以通过观察到目标指数的二进制表示准确地告诉我们应该使用平方指数序列中的哪些元素,而不是使用上述线性搜索,我们可以更有效率。
通过示例可以最好地说明这一点。
要将二进制数转换为十进制数,我们使用以下公式,其中 $b_i$ 是二进制数中的第 $i$ 位:
$$ \text{decimal_value}=2^{n-1}b_{n-1}\dots+8b_3+4b_2+2b_1+b0=\sum{i=0}^{n-1}2^ib_i $$
13 的二进制值是 1101,因为
$$ \begin{align} 44=8&\cdot1 \ +4&\cdot1 \ +2&\cdot0 \ +1&\cdot1 \end{align} $$
52 的二进制值是 110100,因为
$$ \begin{align} 52 = 32&\cdot1\ +16&\cdot1\ +8&\cdot0\ +4&\cdot1\ +2&\cdot0\ +1&\cdot0 \end{align} $$
因此,如果我们想计算 $b^{13}$,并且我们的平方指数序列是 $b^8,b^4,b^2,b^1$,那么我们可以查看 1011 的二进制表示,并确定我们应该选择 8 指数、4 指数和 1 指数,然后计算
$$ b^{13}=b^8\cdot b^4\cdot b^1 $$
因为
$$ 13=8+4+1 $$
因此,我们可以快速确定 13 可以写成 8 + 4 + 1,因为第 3 位、第 2 位和第 0 位是 1。
我们可以使用以下逻辑检测二进制表示中是否设置了第 n 位:
function nthBitSet(uint256 x, uint8 n) public pure returns (bool isSet) {
isSet = x & uint256(2)**n != 0;
}
考虑 2 的幂的二进制表示:
2 的幂 | 十进制值 | 二进制表示 |
---|---|---|
2^0 | 1 | 00001 |
2^1 | 2 | 00010 |
2^2 | 4 | 00100 |
2^3 | 8 | 01000 |
2^4 | 16 | 10000 |
uint256(2)**n
创建一个只有第 n
位设置为 1
的数字。 例如,如果 n = 3
,这将生成二进制值 1000
。 23=8,1000
是 8 的二进制表示,如上表所示。 如果 x
和 2^n
的按位与不为零,则 isSet
为 true
。
只有当两个对应的位都是 1 时,按位与才返回 1; 否则,它返回 0。 例如,8 & 13 = 8
,因为 8
的二进制表示是 1000
,而 13 的二进制表示是:
$$ \begin{aligned}13 &= 1 \cdot 2^3 \ &+ 1 \cdot 2^2 \ &+ 0 \cdot 2^1 \ &+ 1 \cdot 2^0\end{aligned} $$
或 1101
。 当我们将 1000
与 1101
进行按位与时,我们得到 1000
,因为 8 和 13 在第 3 位都有一个共同的 1 位。
如下图所示:
假设 x = 13 (二进制: 1101) 并且 n = 3:
x = 1101 (十进制 13)
2**n = 1000 (十进制 8)
------------------------------
按位与 = 1000 (十进制 8) != 0, 因此为真。
因此, nthBitSet(13, 3) 返回 true。
由于最终结果 1000
不等于零,那么我们知道 x
和 2^n
之间必须有 1 位的重叠 - 这表明该位必须等于 1。
示例 2: 13 中是否设置了第 1 位? 我们可以为第一位创建一个“掩码”,如下所示:2**1
假设 x = 13 (二进制: 1101) 并且 n = 1:
x = 1101 (十进制 13)
2**n = 0010 (十进制 2)
------------------------------
按位与 = 0000 (十进制 0) == 0, 因此为假。
因此, nthBitSet(13, 1) 返回 false。
由于最终结果是 0000,我们知道第 1 位未设置。
示例 3: 13 中是否设置了第 0 位?
假设 x = 13 (二进制: 1101) 并且 n = 0:
x = 1101 (十进制 13)
2**n = 0001 (十进制 1)
------------------------------
按位与 = 0001 (十进制 1) == 1, 因此为真。
因此, nthBitSet(13, 0) 返回 true。
我们现在将我们所学的一切结合起来,生成以下函数,该函数将 0.7 提高到最高 32 的幂:
pow1 = 0.7
pow2 = 0.48999999999999994 # 0.7^2
pow4 = 0.24009999999999995 # 0.7^4
pow8 = 0.05764800999999997 # 0.7^8
pow16 = 0.0033232930569600965 # 0.7^16
def raise_07_to_p(p):
# 计算 0.7^p
assert p < 32
# 如果 p = 0, 我们返回 1, 这是正确的
# 因为 0.7^0 = 1
accumulator = 1
if p & 1 != 0:
accumulator = accumulator * pow1
if p & 2 != 0:
accumulator = accumulator * pow2
if p & 4 != 0:
accumulator = accumulator * pow4
if p & 8 != 0:
accumulator = accumulator * pow8
if p & 16 != 0:
accumulator = accumulator * pow16
return accumulator
## check correctness
print(0.7**14, raise_07_to_p(14))
print(0.7**18, raise_07_to_p(18))
print(0.7**30, raise_07_to_p(30))
上面的代码通过仅使用 0.7 的幂来避免重复乘法,这些幂对应于指数的二进制表示中的 1。
到目前为止,我们使用了常规/浮点数。 但是,在像智能合约这样的真实系统中,我们经常使用定点数(或 Q 数)。
将定点数(或 Q 数)相乘后,必须将结果“标准化”,即需要将结果除以比例因子。 例如,当我们将两个 18 位小数的定点数相乘时,我们需要将结果除以 1018,否则最终结果将写成 36 位小数。
我们将首先在 Python 中实现平方再乘算法,然后将该代码转换为 Solidity。 Solidity 版本将使用定点数,因为 Solidity 没有浮点数。 下面,我们使用 Python 将 0.7 的幂计算为 128 位定点 Q 数作为参考实现:
from decimal import *
import math
getcontext().prec = 100
for i in [1,2,4,8,16,32]:
print(f"0.7^{i} * 2**128 =", math.floor(Decimal(0.7) ** Decimal(i) * Decimal(2**128)))
当我们运行上面的代码时,我们得到以下常量
(base) ➜ python sqmul.py
0.7^1 * 2**128 = 238197656844656909312789480019373064192
0.7^2 * 2**128 = 166738359791259825940851714385556537344
0.7^4 * 2**128 = 81701796297717304344478436853479174360
0.7^8 * 2**128 = 19616601291081919795097291374069314602
0.7^16 * 2**128 = 1130858027394302849221997646234907220
0.7^32 * 2**128 = 3758172630847077410107089115980415
下面的 Solidity 实现使用 128 位定点数,这意味着数字后有 128 位,或者等效地,小数乘以 2128 放大。>> 128
操作等效于除以 2128。
contract SquareAndMultiply {
// z7_x 表示 0.7^x
uint256 constant z7_1 = 238197656844656909312789480019373064192;
uint256 constant z7_2 = 166738359791259825940851714385556537344;
uint256 constant z7_4 = 81701796297717304344478436853479174360;
uint256 constant z7_8 = 19616601291081919795097291374069314602;
uint256 constant z7_16 = 1130858027394302849221997646234907220;
uint256 constant z7_32 = 3758172630847077410107089115980415;
function powZ7(uint256 e) public pure returns (uint256) {
require(e < 64, "e 太大");
// 如果 e 为零,则将 1 作为 128 位定点数返回
uint256 acc = 1 << 128;
// 2 的幂以十六进制编写,0x10 = 16 并且 0x20 = 32
if (e & 0x1 != 0) acc = (acc * z7_1) >> 128;
if (e & 0x2 != 0) acc = (acc * z7_2) >> 128;
if (e & 0x4 != 0) acc = (acc * z7_4) >> 128;
if (e & 0x8 != 0) acc = (acc * z7_8) >> 128;
if (e & 0x10 != 0) acc = (acc * z7_16) >> 128;
if (e & 0x20 != 0) acc = (acc * z7_32) >> 128;
// 通过在具有浮点类型
// 的语言中执行 acc / 2**128
// 来检查答案,然后检查其是否等于 0.7^e
return acc;
}
}
当将整数提高到整数幂时,通常使用虚拟机的内置运算码会更有效。
但是,此运算码不包括上述的标准化步骤。 因此,如果 bx 中的 b 或 x 中至少有一个是定点数,那么我们无法使用虚拟机的 exp
运算码。
要将 tick 转换为平方根价格,Uniswap V3 必须计算
$$ \texttt{sqrtPrice}=\sqrt{1.0001^\text{i}} $$
(对 sqrtPrice 作为 Q64.96 数进行了一些更正)。 因为
基数是固定的,唯一的变量是指数,并且
tick 的范围是提前知道的。 Uniswap V3 可以(并且确实)使用平方再乘算法。 这将在下一章中解释。
作为预告,这是来自 Uniswap V3 Tickmath 库 的 getSqrtRatioAtTick()
的屏幕截图。 这应该看起来非常类似于我们上面编写的 Solidity 代码。 在较高层次上,该函数检查 absTick
中的位是否已设置,如果已设置,则将 ratio
乘以 1.0001 的预计算幂。
- 原文链接: rareskills.io/post/squar...
- 登链社区 AI 助手,为大家转译优秀英文文章,如有翻译不通的地方,还请包涵~
如果觉得我的文章对您有用,请随意打赏。你的支持将鼓励我继续创作!