SignedWadMath.sol 以 int256 统一18 位小数定点数（wad）格式，通过 Padé 有理逼近 + 全程 unchecked assembly 在 EVM 上实现 exp/ln/pow 高等数学运算，服务于 DeFi 连续复利、指数衰减与 bonding curve 等场景。

## 一、概述

**版本说明**：[solmate]：main分支 [commit: 89365b8]，[forge-std]：v1.15.0

**源码链接**：https://github.com/RevelationOfTuring/solmate/blob/main/src/utils/SignedWadMath.sol

`SignedWadMath` 是 solmate 的有符号 18 位小数定点数（wad）数学库，提供单位转换、定点乘除法、幂运算、指数函数和自然对数共 11 个函数。所有函数均为 **free function**（不属于任何合约或 library），可直接 import 使用。

**什么是 wad？**

wad 是 DeFi 中的标准定点数表示法，源自 MakerDAO 的 ds-math 库。1 wad = 1e18，即用 18 位小数表示实数：

```
实数 1.5   → 存储为 1500000000000000000（1.5 × 1e18）
实数 0.001 → 存储为 1000000000000000（0.001 × 1e18）
```

**解决的核心问题**：

Solidity 没有原生浮点数，所有计算都是整数运算。wad 库通过统一的 18 位小数定点格式，让智能合约可以进行高精度的数学运算（乘、除、幂、指数、对数），同时保持 gas 效率。

**设计亮点**：

- 整个库统一使用 `int256` 作为输入输出类型（包括 `toWadUnsafe` 等转换函数），使所有函数可以直接串联调用，无需在调用链中反复进行 `int256`/`uint256` 转型
- `wadExp` 和 `wadLn` 使用 Remco Bloemen 的算法，将 10^18 定点数转为 2^96 二进制定点数后用 Padé 有理逼近，全程 `unchecked` + `assembly`，极致 gas 优化

## 二、适用场景

| 适合 | 不适合 |
|------|--------|
| 连续复利计算（DeFi 借贷协议） | 只需简单整数加减的场景 |
| 代币排放曲线（指数衰减/增长） | 需要无符号定点数的场景（本库为有符号） |
| 线性释放（vesting）中的时间→天数转换 | 对精度要求超过 18 位小数的场景 |
| 定价公式中的幂运算（如 bonding curve） | 输入范围可能导致 wadExp 溢出的场景（x > 135.3 × 1e18） |
| 任何需要 ln/exp/pow 的链上数学运算 | 不需要有符号数的纯正数运算（可用更简单的库） |

## 三、函数总览

```
SignedWadMath（11 个 free function）
│
├── 单位转换（3 个）
│   ├── toWadUnsafe(uint256) → int256           ← 整数 → wad（× 1e18）
│   ├── toDaysWadUnsafe(uint256) → int256       ← 秒 → wad 天数
│   └── fromDaysWadUnsafe(int256) → uint256     ← wad 天数 → 秒
│
├── 定点乘除法（4 个）
│   ├── unsafeWadMul(int256, int256) → int256   ← wad 乘法（不检查溢出）
│   ├── unsafeWadDiv(int256, int256) → int256   ← wad 除法（不检查溢出/除零）
│   ├── wadMul(int256, int256) → int256         ← wad 乘法（带溢出检查）
│   └── wadDiv(int256, int256) → int256         ← wad 除法（带溢出检查）
│
├── 高等数学（3 个）
│   ├── wadPow(int256, int256) → int256         ← x^y = e^(ln(x) × y)
│   ├── wadExp(int256) → int256                 ← e^x（Padé 逼近）
│   └── wadLn(int256) → int256                  ← ln(x)（Padé 逼近）
│
└── 工具函数（1 个）
    └── unsafeDiv(int256, int256) → int256      ← 有符号整数除法（不检查除零）
```

**命名规则**：`unsafe` 前缀/后缀 = 不做溢出或除零检查，由调用方确保输入安全。

## 四、源码逐行解析

### 4.1 toWadUnsafe —— 整数转 wad

```solidity
function toWadUnsafe(uint256 x) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // x × 1e18，直接 mul 不做溢出检查
        r := mul(x, 1000000000000000000)
    }
}
```

**作用**：将无符号整数转为 wad 格式（乘以 1e18）。不做溢出检查。

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `uint256` | 要转换的无符号整数 |

**返回值**：`int256` —— wad 格式的结果。

**注意**：调用方需确保 `x × 1e18` 不溢出 `int256`（即 `x ≤ type(int256).max / 1e18`）。

### 4.2 toDaysWadUnsafe —— 秒转 wad 天数

```solidity
function toDaysWadUnsafe(uint256 x) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // x × 1e18 / 86400
        r := div(mul(x, 1000000000000000000), 86400)
    }
}
```

**作用**：将秒数转为 wad 格式的天数。是进入 wad 数学世界的入口，把链上时间差转为高精度天数喂给 `wadExp` 等函数。

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `uint256` | 秒数 |

**返回值**：`int256` —— wad 格式的天数（1 天 = 1e18）。

**公式**：`r = x × 1e18 / 86400`（1 天 = 86400 秒）

**为什么需要这个函数？** 凡是要把时间差喂给 `wadExp` / `wadMul` / `wadDiv` 做连续数学运算的场景，都需要先转成 wad 天数。否则整数除法会把不足一天的部分截断为 0。

**典型使用场景**：
1. 线性释放（vesting）
2. 代币排放曲线
3. 连续利率计算（DeFi 借贷中按时间连续复利）

### 4.3 fromDaysWadUnsafe —— wad 天数转秒

```solidity
function fromDaysWadUnsafe(int256 x) pure returns (uint256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // x × 86400 / 1e18
        r := div(mul(x, 86400), 1000000000000000000)
    }
}
```

**作用**：将 wad 格式的天数转回秒数。是 wad 数学世界的出口，把计算结果转回链上可用的秒数。

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `int256` | wad 格式的天数 |

**返回值**：`uint256` —— 秒数。

**典型调用链**：

```
秒数 -> toDaysWadUnsafe -> wad天数 -> wadExp/wadMul等运算 -> wad结果 -> fromDaysWadUnsafe -> 秒数
```

### 4.4 unsafeWadMul —— wad 乘法（不检查溢出）

```solidity
function unsafeWadMul(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // (x × y) / 1e18
        // 注：sdiv 是 EVM 的有符号整数除法操作码（Signed DIVision），对应无符号版本是 div
        r := sdiv(mul(x, y), 1000000000000000000)
    }
}
```

**作用**：wad 定点乘法，不检查溢出。

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `int256` | 第一个 wad 操作数 |
| `y` | `int256` | 第二个 wad 操作数 |

**返回值**：`int256` —— wad 格式的乘积。

**数学原理**：

```
a_wad × b_wad = (a × 1e18) × (b × 1e18) = a × b × 1e36
需要除以 1e18 才能得到正确的 wad 结果：a × b × 1e18
```

### 4.5 unsafeWadDiv —— wad 除法（不检查溢出和除零）

```solidity
function unsafeWadDiv(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // (x × 1e18) / y
        r := sdiv(mul(x, 1000000000000000000), y)
    }
}
```

**作用**：wad 定点除法，不检查溢出和除零。`y` 为 0 时返回 0（`sdiv` 除零行为），不会 revert。

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `int256` | 被除数（wad 格式） |
| `y` | `int256` | 除数（wad 格式） |

**返回值**：`int256` —— wad 格式的商。

**数学原理**：

```
(a_wad / b_wad) = (a × 1e18) / (b × 1e18) = a / b
需要先乘 1e18 才能得到正确的 wad 结果：(a / b) × 1e18
```

### 4.6 wadMul —— wad 乘法（带溢出检查）

```solidity
function wadMul(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // 先计算 x * y 存入 r
        r := mul(x, y)

// 组合溢出检查（合并为 1 个表达式，只产生 1 个 JUMPI）：
        if iszero(
            and(
                // 条件 A: or(iszero(x), eq(sdiv(r, x), y))
                //   → x == 0（0 乘任何数不溢出）或 (x*y)/x == y（验证乘法可逆）
                or(iszero(x), eq(sdiv(r, x), y)),
                // 条件 B: or(lt(x, not(0)), sgt(y, 0x80...0))
                //   → x < -1 或 y > type(int256).min
                //   → 排除 x=-1, y=type(int256).min 的二补码陷阱
                or(lt(x, not(0)), sgt(y, 0x8000000000000000000000000000000000000000000000000000000000000000))
            )
        ) {
            revert(0, 0)
        }

// 溢出检查通过，除以 1e18 得到 wad 结果
        r := sdiv(r, 1000000000000000000)
    }
}
```

**作用**：wad 定点乘法，带溢出检查。溢出时 revert（无错误信息）。

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `int256` | 第一个 wad 操作数 |
| `y` | `int256` | 第二个 wad 操作数 |

**返回值**：`int256` —— wad 格式的乘积，溢出时 revert。

**溢出检查的两个条件**：

| 条件 | 含义 | 作用 |
|------|------|------|
| A：`x == 0 || (x*y)/x == y` | 标准乘法溢出检查 | 验证乘法可逆 |
| B：`x < -1 || y > type(int256).min` | 二补码边缘情况检查 | 排除 x=-1, y=min 的陷阱 |

**为什么需要条件 B？**

当 x = -1, y = type(int256).min（即 -2^255）时：
- `(-1) × (-2^255) = 2^255` 超过 int256.max（2^255-1），产生溢出
- 但 EVM 的 `mul` 不区分有符号/无符号，按无符号相乘后取低 256 位

```
  EVM mul 实际计算：把两个操作数当作无符号整数相乘

-1 的补码 = 0xFFFF...FFFF（无符号值 = 2^256 - 1）：
  ┌──────────────────────────────────────────────────────────────────────┐
  │ 0xFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF   │
  └──────────────────────────────────────────────────────────────────────┘

-2^255 的补码 = 0x8000...0000（无符号值 = 2^255）：
  ┌──────────────────────────────────────────────────────────────────────┐
  │ 0x8000000000000000000000000000000000000000000000000000000000000000   │
  └──────────────────────────────────────────────────────────────────────┘

＝

无符号乘积：(2^256 - 1) × 2^255 = 2^511 - 2^255

完整 512 位结果（二进制）：
  高 256 位                                低 256 位
  ┌──────────────────────────────────┬──────────────────────────────────┐
  │ 0111 1111 ... 1111 1111 1111 1111│ 1000 0000 ... 0000 0000 0000 0000│
  └──────────────────────────────────┴──────────────────────────────────┘
    					↑ 高 256 位丢弃                   ↑ EVM 只保留低 256 位

截断后的低 256 位：
  ┌──────────────────────────────────────────────────────────────────────┐
  │ 0x8000000000000000000000000000000000000000000000000000000000000000   │
  └──────────────────────────────────────────────────────────────────────┘
  = -2^255 = type(int256).min = y
```

- 此时 `sdiv(-2^255, -1) = -2^255 = y`，条件 A 误判为"未溢出"

注：这不是补码层面的计算逻辑，而是 EVM 规范的硬编码特例（Yellow Paper中有定义：https://ethereum.github.io/yellowpaper/paper.pdf ）：

![image-20260410095934058.png](https://img.learnblockchain.cn/attachments/2026/04/8nKrauB269ed944d037a1.png)

- 条件 B 额外捕获这种二补码陷阱

**为什么用 `lt`/`sgt` 而不是 `!= -1` 和 `!= min`？**

Yul 没有 `neq` 操作码，用 `lt(x, not(0))`（即 x < -1）和 `sgt(y, 0x80...0)`（即 y > min）替代 `iszero(eq(...))`，省 2 个操作码。虽然范围更宽，但只会放过更多安全情况，不会漏掉目标。

**gas 优化技巧**：两个条件合并为 1 个 `and` 表达式，编译后只产生 1 个 `JUMPI` 而非 2 个。

### 4.7 wadDiv —— wad 除法（带溢出检查）

```solidity
function wadDiv(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // 先计算 x * 1e18 存入 r
        r := mul(x, 1000000000000000000)

// 检查：y != 0 且 r / 1e18 == x（确保 x * 1e18 没有溢出）
        // 注：iszero(iszero(y)) 是 Yul 中判断 y != 0 的标准写法
        if iszero(and(iszero(iszero(y)), eq(sdiv(r, 1000000000000000000), x))) {
            revert(0, 0)
        }

// r 除以 y 得到最终结果
        r := sdiv(r, y)
    }
}
```

**作用**：wad 定点除法，带溢出检查。除零或溢出时 revert。

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `int256` | 被除数（wad 格式） |
| `y` | `int256` | 除数（wad 格式） |

**返回值**：`int256` —— wad 格式的商，溢出或除零时 revert。

**溢出检查**：`y != 0 && (x * 1e18) / 1e18 == x`，确保除数不为零且 `x * 1e18` 没有溢出。

### 4.8 wadPow —— wad 幂运算

```solidity
function wadPow(int256 x, int256 y) pure returns (int256) {
    // x^y 等价于 e^(ln(x) * y)
    // 其中，wadLn() 要求 x > 0
    return wadExp((wadLn(x) * y) / 1e18);
}
```

**作用**：wad 定点幂运算。

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `int256` | 底数（wad 格式，必须为正） |
| `y` | `int256` | 指数（wad 格式） |

**返回值**：`int256` —— wad 格式的 x^y。

**数学原理**：`x^y = (e^ln(x))^y = e^(ln(x) × y)`

**注意**：只支持正数底数（x > 0），因为 ln对非正数无定义。

### 4.9 wadExp —— 自然指数函数 e^x

这是整个库最复杂的函数，算法来自 [Remco Bloemen](https://xn--2-umb.com/22/exp-ln/)。

#### 函数签名

```solidity
function wadExp(int256 x) pure returns (int256 r) {
    unchecked { ... }
}
```

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `int256` | wad 格式的指数（代表实数 x/1e18） |

**返回值**：`int256` —— wad 格式的 e^(x/1e18)。

#### 算法流程（5 步）

```
wadExp(x)
  │
  ├── 1. 边界检查
  │     ├── x ≤ -42139678854452767551 → return 0（结果太小）
  │     └── x ≥ 135305999368893231589 → revert（结果溢出）
  │
  ├── 2. 基数转换：10^18 定点 → 2^96 定点
  │     └── x = (x << 78) / 5^18
  │
  ├── 3. 范围缩减：exp(x) = exp(x') × 2^k
  │     ├── k = round(x / ln2)
  │     └── x' = x - k × ln2，使 x' ∈ (-½ln2, ½ln2)
  │
  ├── 4. (6,7) 阶 Padé 有理逼近
  │     ├── p(x)：5 次多项式（分子）
  │     ├── q(x)：6 次多项式（分母）
  │     └── r = p / q
  │
  └── 5. 最终组装
        └── r = r × s × 2^k × 10^18 / 2^96（一步完成）
```

#### 第 1 步：边界检查

```solidity
// 1.1 下界检查
//  wad 是整数，所以最小正值就是 1（代表实数 1/1e18 = 10^-18）
//  如果计算结果的实数值 < 0.5 × 10^-18，四舍五入后就是 0，直接返回 0
//  -42139678854452767551 = ⌊ln(0.5e-18) × 1e18⌋
if (x <= -42139678854452767551) return 0;

// 1.2 上界检查
//  如果 e^(x/1e18) × 1e18 超过 type(int256).max = 2^255-1，int256 就装不下了
//  135305999368893231589 = ⌊ln((2^255-1)/1e18) × 1e18⌋
if (x >= 135305999368893231589) revert("EXP_OVERFLOW");
```

x 的有效范围是 (-42.14, 135.31) × 1e18。

#### 第 2 步：基数转换

```solidity
x = (x << 78) / 5 ** 18;
```

**为什么要转换？** x 是 10^18 定点数，每次乘除都涉及 `/ 1e18`，EVM 做除法很贵。转成 2^96 定点后，除法变成移位（`>> 96`），便宜很多。

注：`SHR` 3 gas，`DIV`/`SDIV` 5 gas，单次省 2 gas。

**转换因子推导**：

```
10^18 定点：实数 1.0 存储为 1e18
2^96 定点：实数 1.0 存储为 2^96

x_new = x × (2^96 / 10^18)
      = x × 2^96 / (2×5)^18
      = x × 2^96 / 2^18 / 5^18
      = x × 2^78 / 5^18
      = (x << 78) / 5^18
```

**为什么选 2^96？**

| 原因 | 说明 |
|------|------|
| 精度 | 2^96 ≈ 7.9 × 10^28，约 29 位十进制精度，远超 wad 的 18 位 |
| 不溢出 | 两个 2^96 数相乘得 2^192，不超过 256 位 |
| 省 gas | `>> 96` 替代 `/ 1e18`，Padé 逼近有十几次乘除，每次都省一点 |

#### 第 3 步：范围缩减

```solidity
// k = round(x / ln2) = floor(x / ln2 + 0.5)
// 常数：54916777467707473351141471128 = ln(2) × 2^96
// 常数：2^95 = 0.5 × 2^96
int256 k = ((x << 96) / 54916777467707473351141471128 + 2 ** 95) >> 96;
// x' = x - k × ln2
x = x - k * 54916777467707473351141471128;
```

**为什么要范围缩减？**

Padé 逼近只在很小的范围内精确。exp(100) 很大，直接逼近误差爆炸。利用数学恒等式：

```
exp(x) = exp(x - k×ln2 + k×ln2)
       = exp(x - k×ln2) × exp(k×ln2)
       = exp(x') × 2^k
```

其中 x' = x - k×ln2。选合适的 k 让 x' 很小，最后乘以 2^k（移位操作，低 gas）还原。

**如何确定 k？**

要让 x' = x - k × ln2 足够小，就看 x 中含有多少个 ln2，将这个数定为 k：

```
k = round(x / ln2)
```

令 x / ln2 的真实值 = k + ε，其中 ε 是误差。由于 round 是四舍五入，所以 |ε| < 0.5。

```
x' = x - k × ln2
   = (k + ε) × ln2 - k × ln2
   = ε × ln2

因为 |ε| < 0.5，所以 |x'| = |ε × ln2| < 0.5 × ln2 = ½ln2
即 x' 落在 (-½ln2, ½ln2) 范围内
```

**为什么用 round 而不是 floor？**

```
round: k = round(x / ln2)  → 误差 |ε| < 0.5  → x' ∈ (-½ln2, ½ln2)，关于原点对称
floor: k = floor(x / ln2)  → 误差 ε ∈ [0, 1)  → x' ∈ [0, ln2)，偏离原点
```

Padé 逼近在原点附近最精确，round 产生的 x' 范围关于原点对称，逼近精度更好。

**round 的工程实现**：

```
k = round(x / ln2) = floor(x / ln2 + 0.5)
```

为了保留小数部分，将整体扩大 2^96 倍做加法，再缩小 2^96 倍实现 floor：

```
k = ((x × 2^96 / LN2_96) + 0.5_96) >> 96
```

其中 `LN2_96 = 54916777467707473351141471128`，`0.5_96 = 2^95`。

通过 x 的范围可以计算出 k 的范围：[-61, 195]。

#### 第 4 步：(6,7) 阶 Padé 有理逼近

**什么是 Padé 逼近？**

Taylor 级数用多项式逼近函数，Padé 用有理函数（两个多项式相除）逼近：

```
Taylor: e^x ≈ 1 + x + x²/2 + x³/6 + ...        ← 多项式
Padé:   e^x ≈ P(x) / Q(x)                      ← 有理函数
```

为什么 Padé 更好？因为 exp(x) 在 x → -∞ 时趋近 0。多项式只能趋近 ±∞，无法趋近 0；而 P(x)/Q(x) 可以（只要分母增长比分子快）。

**系数的确定**：一旦确定逼近目标（e^x）、阶数（(6,7)）、展开点（x = 0），系数就唯一确定。

**Horner 法则**：多项式求值的快速算法，通过嵌套形式减少乘法次数：

```
标准形式：ax³ + bx² + cx + d        ← 需要 6 次乘法 + 3 次加法
Horner： ((a*x + b)*x + c)*x + d   ← 需要 3 次乘法 + 3 次加法
```

以下代码就是用 Horner 法则计算分子 P(x) 和分母 Q(x)，其中系数是预先算好的魔数常量（来自 Remco Bloemen 的论文）。

**分子 p(x) 的计算**（非标准 Horner，借用中间变量 y 做交叉展开）：

```solidity
// y = x + c0 → 1 次
int256 y = x + 1346386616545796478920950773328;
// y = y * x + c1 = x²+c0*x+c1 → 2 次
y = ((y * x) >> 96) + 57155421227552351082224309758442;
// p = y + x + c2 → 2 次
int256 p = y + x - 94201549194550492254356042504812;
// p = p*y + c3（2次*2次）→ 4 次
p = ((p * y) >> 96) + 28719021644029726153956944680412240;  
// 最后一步故意不 >> 96，p 留在 2^192 基数，省一次移位
// 这就是为什么常数4385272521454847904659076985693276也要<< 96 —— 把它也提升到 2^192 基数与 p 对齐再做加法
// p = p*x + c4 → 5 次
p = p * x + (4385272521454847904659076985693276 << 96);
```

**分母 q(x) 的计算**（标准 Horner）：

```solidity
// 每步 * x >> 96 就是 2^96 定点下的乘法。6 步展开 = 6 次多项式
int256 q = x - 2855989394907223263936484059900;
q = ((q * x) >> 96) + 50020603652535783019961831881945;
q = ((q * x) >> 96) - 533845033583426703283633433725380;
q = ((q * x) >> 96) + 3604857256930695427073651918091429;
q = ((q * x) >> 96) - 14423608567350463180887372962807573;
q = ((q * x) >> 96) + 26449188498355588339934803723976023;
```

**计算 r = p / q**：

```solidity
assembly {
    // p(2^192基数) / q(2^96 基数) = r(2^96 基数)
    // gas技巧：Solidity的 / 运算符即使在 unchecked {} 块中，编译器仍然会插入一个"除数是否为零"的检查，
    // 由于 q(x) 是一个 6 次多项式，它的 6 个根（令 q = 0 的解）全部是复数，没有实数根
    // 这是Padé逼近的数学性质保证的——论文作者在选系数时就确保了分母没有实数根
    // 由于 x 取值在 int256 范围内（实数域的子集），而 q(x) 的根全是复数，所以可知 q(x) 在实数域上永远不会为0，这里特意使用Yul的sdiv来跳过除数不为0的检查
    r := sdiv(p, q)
}
```

此时 r 的范围是 (0.09, 0.25) × 2^96。

#### 第 5 步：最终组装

```solidity
// CONST = round(s × 10^18 × 2^99) = 3822833074963236453042738258902158003155416615667
r = int256((uint256(r) * 3822833074963236453042738258902158003155416615667) >> uint256(195 - k));
```

**现在的 r 是什么？**

目前的 r = P(x') / Q(x')，是 exp(x') 的近似值。但由于 P 和 Q 都是首一多项式（最高次系数 = 1），P/Q 的结果差一个常数缩放因子 s ≈ 6.031367120。真正的 exp(x') = r × s。

**如何将 exp(x') 还原成 exp(x)？**

几个遗留问题：
1. 目前 r 在 2^96 基数下，但最终结果是 wad（10^18 基数）
2. 范围缩减时拆出了 2^k，还没乘回来：exp(x) = exp(x') × 2^k = r × s × 2^k

所以一步要完成三个操作：

1. **乘以缩放因子 s ≈ 6.031367120**：补偿首一多项式的缩放
2. **乘以 2^k**：恢复范围缩减中分离的指数部分
3. **基数转换 10^18 / 2^96**：从二进制定点转回 wad

**数学推导**：

```
result = r × s × 2^k × 10^18 / 2^96
       = r × (s × 10^18 × 2^99) × 2^k / 2^(96+99)
       = r × (s × 10^18 × 2^99) / 2^(195-k)
       = r × CONST >> (195 - k)

CONST = round(s × 10^18 × 2^99) = 3822833074963236453042738258902158003155416615667
```

**为什么要凑出 2^195？**

因为 k 的范围是 [-61, 195]，右移量 195-k 的范围是 [0, 256]，始终安全。

### 4.10 wadLn —— 自然对数 ln(x)

同样来自 [Remco Bloemen](https://xn--2-umb.com/22/exp-ln/)，是 wadExp 的逆运算。

#### 函数签名

```solidity
function wadLn(int256 x) pure returns (int256 r) {
    unchecked { ... }
}
```

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `int256` | wad 格式的输入（代表实数 x/1e18，必须为正） |

**返回值**：`int256` —— wad 格式的 ln(x/1e18)。

#### 算法流程（5 步）

```
wadLn(x)
  │
  ├── 1. 域检查：require(x > 0)
  │
  ├── 2. 二分搜索找 r = floor(log2(x))（8 步）
  │
  ├── 3. 范围缩减：将x归一化到 [2^96, 2^97)
  │     ├── k = r - 96
  │     ├── x <<= (159 - k)
  │     └── x = int256(uint256(x) >> 159)
  │
  ├── 4. (8,8) 阶 Padé 有理逼近
  │     ├── p(x)：7 次多项式（分子）
  │     ├── q(x)：7 次多项式（分母）
  │     └── r = p / q
  │
  └── 5. 最终组装
        ├── r *= s_const        ← 补偿缩放因子
        ├── r += ln2_const × k  ← 补偿归一化
        ├── r += offset_const   ← 补偿 ln(2^96/10^18)
        └── r >>= 174           ← 转回 wad
```

#### 第 1 步：域检查

```solidity
require(x > 0, "UNDEFINED");
```

ln 对非正数无定义。

#### 第 2 步：二分搜索找 floor(log2(x))

目的是找到 x 的最高有效位（MSB）在第几位。

**思路**：将 256 位 MSB 位置拆解为 8 个二进制位 (b7~b0)，每步将搜索范围减半：

```
1. b7: x > 2^128-1? → 确定 MSB 在上半还是下半 256 位
2. b6: x>>r > 2^64-1? → 继续二分
3. ...依此类推直到 b0

结果：r = 128×b7 + 64×b6 + ... + 1×b0 = floor(log2(x))
```

```solidity
assembly {
    // b7：判断 x > 2^128 - 1，得出 MSB 在高 128 位 (128~255) 还是低 128 位 (0~127)
    // 如果是：r = 128，否则：r = 0
    r := shl(7, lt(0xffffffffffffffffffffffffffffffff, x))
    // b6：将已确定的位移掉后，判断 剩余部分 > 2^64 - 1
    // 如果是：r += 64，否则：r不变
    r := or(r, shl(6, lt(0xffffffffffffffff, shr(r, x))))
    // b5：将已确定的位移掉后，判断 剩余部分 > 2^32 - 1
    // 如果是：r += 32，否则：r不变
    r := or(r, shl(5, lt(0xffffffff, shr(r, x))))
    // b4：将已确定的位移掉后，判断 剩余部分 > 2^16 - 1
    // 如果是：r += 16，否则：r不变
    r := or(r, shl(4, lt(0xffff, shr(r, x))))
    // b3：将已确定的位移掉后，判断 剩余部分 > 2^8 - 1
    // 如果是：r += 8，否则：r不变
    r := or(r, shl(3, lt(0xff, shr(r, x))))
    // b2：将已确定的位移掉后，判断 剩余部分 > 2^4 - 1
    // 如果是：r += 4，否则：r不变
    r := or(r, shl(2, lt(0xf, shr(r, x))))
    // b1：将已确定的位移掉后，判断 剩余部分 > 2^2 - 1
    // 如果是：r += 2，否则：r不变
    r := or(r, shl(1, lt(0x3, shr(r, x))))
    // b0：将已确定的位移掉后，判断 剩余部分 > 2^1 - 1
    // 如果是：r += 1，否则：r不变
    r := or(r, lt(0x1, shr(r, x)))
    // 此时r = floor(log2(x))
}
```

**原理**：以 b7 为例——如果 x > 2^128 - 1，说明 MSB 在第 128~255 位，r 先设为 128；否则 MSB 在第 0~127 位，r = 0。然后将 x 右移 r 位，在剩余范围内继续二分。

#### 第 3 步：范围缩减（归一化）

先提取 m，使得 x = 2^r × m，m ∈ [1, 2)，并将 m 表示为 2^96 定点数。

**为什么要归一化？**

Padé 逼近只在小范围内精确，需要将 x 归一化到固定区间 [2^96, 2^97)。

```
数学依据：ln(x) = ln(2^r × m) = r × ln(2) + ln(m)
归一化后 Padé 只负责算 ln(m)，r × ln(2) 在第 5 步补回
```

**m 为什么用 2^96 定点数表示？**

同 wadExp 中 x 用 2^96 定点数表示的原因（精度、不溢出、移位省 gas）。将 m 表示为 2^96 定点数后：

```
归一化前：x = 2^r × m
归一化后：x_new = 2^96 × m
令 k = r - 96（即 r = k + 96），则：
x_new = 2^96 × m = 2^r × m × 2^(96-r) = x × 2^(-k)
```

**r × ln(2) 怎么补回？**

由前面的分解：ln(x) = r × ln(2) + ln(m)，代入 r = k + 96：

```
ln(x) = (k + 96) × ln(2) + ln(m)
归一化后 Padé 只负责算 ln(m)，剩下的 (k + 96) × ln(2) 需要在第 5 步补回
但 floor(log2(x)) 存在变量 r 中，第 4 步被 p/q 覆盖，只剩 k
所以拆成两部分补：k × ln(2) 显式加回，96 × ln(2) 合并到常数 ln(2^96/10^18) 中补
```

**归一化**

```solidity
// 对 x 做归一化，即把 x 从 2^r × m 变成 2^96 × m，即把 MSB 从第 r 位挪到第 96 位
// 这里先左移到顶再右移，而不是直接 x >> (r-96)。原因是：r 可能 < 96，EVM 没法右移一个负数位

// k 为归一化偏移量
int256 k = r - 96;

// 将 x 的 MSB（第 r 位）左移到第 255 位
// 左移位数为 255-r = 255-(k+96) = 159-k
// 目的：把有效位对齐到最高位，为下一步统一右移做准备
x <<= uint256(159 - k);

// 再从第 255 位统一移到第 96 位（255 - 159 = 96）
// 使用 uint256 cast 确保逻辑右移（填充 0），避免算术右移（填充符号位）
x = int256(uint256(x) >> 159);

// 此时，x 被归一化到 [2^96, 2^97) 范围内，即 2^96 定点数下的 [1, 2)
```

#### 第 4 步：(8,8) 阶 Padé 有理逼近

比 wadExp 的 (6,7) 阶更高，因为 ln 在 [1,2) 上的曲率比 exp 在 (-½ln2, ½ln2) 上更大，需要更高阶才能达到同等精度。

**分子 p(x) 的计算**（7 次多项式）：

```solidity
int256 p = x + 3273285459638523848632254066296;
p = ((p * x) >> 96) + 24828157081833163892658089445524;
p = ((p * x) >> 96) + 43456485725739037958740375743393;
p = ((p * x) >> 96) - 11111509109440967052023855526967;
p = ((p * x) >> 96) - 45023709667254063763336534515857;
p = ((p * x) >> 96) - 14706773417378608786704636184526;
// 最后一步故意不 >> 96，p 留在 2^192 基数
// 目的：最终要算 r = p / q。这样做省了一次移位操作
p = p * x - (795164235651350426258249787498 << 96);
```

**分母 q(x) 的计算**（7 次多项式）：

```solidity
int256 q = x + 5573035233440673466300451813936;
q = ((q * x) >> 96) + 71694874799317883764090561454958;
q = ((q * x) >> 96) + 283447036172924575727196451306956;
q = ((q * x) >> 96) + 401686690394027663651624208769553;
q = ((q * x) >> 96) + 204048457590392012362485061816622;
q = ((q * x) >> 96) + 31853899698501571402653359427138;
// q 在 2^96 基数
q = ((q * x) >> 96) + 909429971244387300277376558375;
```

**计算 r = p / q**：

```solidity
assembly {
    // 同 wadExp 中的 gas 技巧：q(x) 的根全是复数，实数域上不为零，
    // 使用 Yul 的 sdiv 跳过 Solidity 编译器自动插入的除零检查
    r := sdiv(p, q)
}
```

此时 r 的范围是 (0, 0.125) × 2^96。

#### 第 5 步：最终组装

**现在的 r 是什么？**

前面计算的 p(x) 和 q(x) 都是首一多项式，结果差一个缩放因子 s ≈ 5.549。真正的 ln(m) = r × s。

**从分解公式推导最终组装**

目标：ln(x / 10^18)

由 x = 2^floor(log2(x)) × m = 2^(k+96) × m，得：

```
ln(x / 10^18) = ln(2^(k+96) × m / 10^18)
              = ln(m) + (k + 96) × ln(2) - ln(10^18)
              = ln(m) + k × ln(2) + 96 × ln(2) - ln(10^18)
              = ln(m) + k × ln(2) + ln(2^96 / 10^18)
```

三项分别对应：
- **r × s → ln(m)**：补偿缩放因子，得到真正的 ln(m)
- **k × ln(2)** → 补回 r × ln(2) = (k + 96) × ln(2) 中的 k 份（r 的值已被 p/q 覆盖，目前只有 k 可用）
- **ln(2^96 / 10^18)** → 常数补回，包含两项：
  - 96 × ln(2)：(k + 96) × ln(2) 中 k 份已显式补回，这里补剩余的 96 份
  - -ln(10^18)：输入 x 是 wad 定点，目标是 ln(x/10^18)，需减掉多出的 ln(10^18)

**代码实现：在 5^18 × 2^192 基数下统一计算**

数学目标：`result_wad = (r × s + k × ln(2) + ln(2^96 / 10^18)) × 10^18`

但不能直接乘 10^18（r 在 2^96 基数，直接乘会损失精度或溢出）。

技巧：10^18 = 5^18 × 2^18，选 5^18 × 2^192 作为中间基数，算完后只需 >> 174（= 192 - 18）就能转回 wad（5^18 × 2^18 = 10^18）。

```solidity
// 补偿缩放因子 s，并将基数提高到 5^18 × 2^192
// 常数 1677202110996718588342820967067443963516166 = s × 5^18 × 2^96（r 是 2^96 基数，乘完变为 5^18 × 2^192 基数）
r *= 1677202110996718588342820967067443963516166;

// 补回 r × ln(2) = (k + 96) × ln(2) 中的 k 份，并提升到 5^18 × 2^192 基数下
// 常数 16597577552685614221487285958193947469193820559219878177908093499208371 = ln(2) × 5^18 × 2^192
r += 16597577552685614221487285958193947469193820559219878177908093499208371 * k;

// 补回 ln(2^96 / 10^18)（= 96×ln(2) - ln(10^18)，即 (k+96)×ln(2) 中剩余的 96 份和基数转换）
// 常数 600920179829731861736702779321621459595472258049074101567377883020018308 = ln(2^96 / 10^18) × 5^18 × 2^192
r += 600920179829731861736702779321621459595472258049074101567377883020018308;

// 从 5^18 × 2^192 基数转回 wad（10^18 = 5^18 × 2^18）
// 5^18 × 2^192 / 2^174 = 5^18 × 2^18 = 10^18
r >>= 174;
```
### 4.11 unsafeDiv —— 有符号整数除法（不检查除零）

```solidity
function unsafeDiv(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        r := sdiv(x, y)
    }
}
```

**作用**：有符号整数除法，不检查除零。y 为 0 时返回 0（`sdiv` 除零行为），不会 revert。

**参数**：

| 参数 | 类型 | 含义 |
|------|------|------|
| `x` | `int256` | 被除数 |
| `y` | `int256` | 除数 |

**返回值**：`int256` —— 商。

**注意**：这不是 wad 除法，是普通整数除法（没有 × 1e18 的步骤）。

## 五、wadExp / wadLn 算法流程图

### wadExp 流程

```
输入: x (wad 格式)
  │
  ├── x ≤ -42139678854452767551 ?
  │     └── YES → return 0
  │
  ├── x ≥ 135305999368893231589 ?
  │     └── YES → revert("EXP_OVERFLOW")
  │
  ├── 基数转换: x = (x << 78) / 5^18
  │     └── 10^18 定点 → 2^96 定点
  │
  ├── 范围缩减:
  │     ├── k = round(x / ln2)                    ← k ∈ [-61, 195]
  │     └── x' = x - k × ln2                      ← x' ∈ (-½ln2, ½ln2)
  │
  ├── Padé 逼近:
  │     ├── p = P(x')  (5 次多项式，2^192 基数)
  │     ├── q = Q(x')  (6 次多项式，2^96 基数)
  │     └── r = p / q  (2^96 基数)
  │
  └── 最终组装:
        └── r = r × CONST >> (195 - k)
              ↓
        输出: e^(x/1e18) (wad 格式)
```

### wadLn 流程

```
输入: x (wad 格式, x > 0)
  │
  ├── require(x > 0)
  │
  ├── 二分搜索: r = floor(log2(x))
  │     └── 8 步二分，确定 MSB 位置
  │
  ├── 归一化:
  │     ├── k = r - 96
  │     ├── x <<= (159 - k)                        ← MSB 对齐到第 255 位
  │     └── x = uint256(x) >> 159                  ← MSB 移到第 96 位
  │           └── x ∈ [2^96, 2^97)
  │
  ├── Padé 逼近:
  │     ├── p = P(x)  (7 次多项式，2^192 基数)
  │     ├── q = Q(x)  (7 次多项式，2^96 基数)
  │     └── r = p / q  (2^96 基数)
  │
  └── 最终组装 (5^18 × 2^192 基数):
        ├── r *= s_const                           ← 补偿缩放因子
        ├── r += ln2_const × k                     ← 补偿归一化
        ├── r += offset_const                      ← 补偿 ln(2^96/10^18)
        └── r >>= 174                              ← 转回 wad
              ↓
        输出: ln(x/1e18) (wad 格式)
```

## 六、设计思想

### 6.1 极简主义

- 所有函数都是 free function，不依赖任何合约或 library
- 没有自定义 error、event、modifier，没有 storage
- 统一使用 `int256` 类型，消除调用链中的类型转换

### 6.2 gas 极致优化

| 技巧 | 说明 |
|------|------|
| 2^96 定点数 | `>> 96` 替代 `/ 1e18`，移位比除法便宜 |
| unchecked | wadExp/wadLn 全程 unchecked，省去溢出检查 gas |
| Yul sdiv | 跳过 Solidity 编译器自动插入的除零检查（因为数学上证明了分母不为 0） |
| 合并溢出检查 | wadMul 两个条件合并为 1 个 JUMPI |
| 分子不右移 | p 留在 2^192 基数，p/q 自动得到 2^96 结果，省一次移位 |
| 常数预计算 | 所有常数都是编译期已知的魔数，运行时无额外计算 |

### 6.3 Padé 逼近 vs Taylor 级数

| 对比项 | Taylor 级数 | Padé 逼近 |
|--------|------------|-----------|
| 形式 | 多项式 | 有理函数 P(x)/Q(x) |
| 逼近 exp(x) 在 x→-∞ 的行为 | 多项式趋向 ±∞，无法趋近 0 | 分母增长快于分子，可趋近 0 ✅ |
| 同等精度所需项数 | 多 | 少 ✅ |
| EVM 实现复杂度 | 简单 | 稍复杂（需要除法） |

### 6.4 安全 vs 性能的取舍

库提供了 safe/unsafe 两套函数，让调用方自行选择：

| 函数 | 溢出检查 | 除零检查 | gas |
|------|----------|----------|-----|
| unsafeWadMul | ❌ | - | 低 |
| wadMul | ✅ | - | 高 |
| unsafeWadDiv | ❌ | ❌ | 低 |
| wadDiv | ✅ | ✅ | 高 |

## 七、安全注意事项

| 风险 | 建议 |
|------|------|
| unsafe 函数不检查溢出/除零 | 调用方需自行保证输入安全，或优先使用 safe 版本 |
| wadExp 输入范围有限 | x 必须在 (-42.14, 135.31) × 1e18 范围内，低于下界返回 0，高于上界revert |
| wadLn 只接受正数 | x ≤ 0 会 revert("UNDEFINED") |
| wadPow 只支持正数底数 | 内部调用 wadLn，底数必须 > 0 |
| wad 精度限制 | 最小正值为 1（代表 10^-18），小于此值的结果会被截断为 0 |
| 整数截断 | wad 乘除法的中间结果除以 1e18 是整数截断（向零取整），非四舍五入 |
| unsafeWadDiv 除零返回 0 | `sdiv(x, 0)` 在 EVM 中返回 0 而非 revert，可能导致静默错误 |

## 八、与同类方案对比

| 对比项 | Solmate SignedWadMath | OpenZeppelin Math | PRBMath |
|--------|----------------------|-------------------|---------|
| 类型 | int256（有符号） | uint256（无符号） | int256 / uint256 两套 |
| 精度 | 18 位小数 | 无定点 | 18 位小数 |
| exp/ln | ✅ Padé 逼近 | ❌ 无 | ✅ Padé 逼近 |
| 幂运算 | ✅ wadPow | ❌ 无 | ✅ pow |
| 形式 | free function | library | library |
| gas | 极低（全 assembly） | 中 | 中低 |
| 安全检查 | safe + unsafe 两套 | 全部带检查 | 全部带检查 |
| 时间转换 | ✅ toDays/fromDays | ❌ | ❌ |

## 九、测试实战

**Mock 合约**：https://github.com/RevelationOfTuring/foundry-solmate/blob/main/src/utils/MockSignedWadMath.sol

**全部 Foundry 测试合约**：https://github.com/RevelationOfTuring/foundry-solmate/blob/main/test/utils/SignedWadMath.t.sol

一、概述

版本说明：[solmate]：main分支 [commit: 89365b8]，[forge-std]：v1.15.0

源码链接：https://github.com/RevelationOfTuring/solmate/blob/main/src/utils/SignedWadMath.sol

SignedWadMath 是 solmate 的有符号 18 位小数定点数（wad）数学库，提供单位转换、定点乘除法、幂运算、指数函数和自然对数共 11 个函数。所有函数均为 free function（不属于任何合约或 library），可直接 import 使用。

什么是 wad？

wad 是 DeFi 中的标准定点数表示法，源自 MakerDAO 的 ds-math 库。1 wad = 1e18，即用 18 位小数表示实数：

实数 1.5   → 存储为 1500000000000000000（1.5 × 1e18）
实数 0.001 → 存储为 1000000000000000（0.001 × 1e18）

解决的核心问题：

设计亮点：

整个库统一使用 int256 作为输入输出类型（包括 toWadUnsafe 等转换函数），使所有函数可以直接串联调用，无需在调用链中反复进行 int256/uint256 转型
wadExp 和 wadLn 使用 Remco Bloemen 的算法，将 10^18 定点数转为 2^96 二进制定点数后用 Padé 有理逼近，全程 unchecked + assembly，极致 gas 优化

二、适用场景

适合	不适合
连续复利计算（DeFi 借贷协议）	只需简单整数加减的场景
代币排放曲线（指数衰减/增长）	需要无符号定点数的场景（本库为有符号）
线性释放（vesting）中的时间→天数转换	对精度要求超过 18 位小数的场景
定价公式中的幂运算（如 bonding curve）	输入范围可能导致 wadExp 溢出的场景（x > 135.3 × 1e18）
任何需要 ln/exp/pow 的链上数学运算	不需要有符号数的纯正数运算（可用更简单的库）

三、函数总览

SignedWadMath（11 个 free function）
│
├── 单位转换（3 个）
│   ├── toWadUnsafe(uint256) → int256           ← 整数 → wad（× 1e18）
│   ├── toDaysWadUnsafe(uint256) → int256       ← 秒 → wad 天数
│   └── fromDaysWadUnsafe(int256) → uint256     ← wad 天数 → 秒
│
├── 定点乘除法（4 个）
│   ├── unsafeWadMul(int256, int256) → int256   ← wad 乘法（不检查溢出）
│   ├── unsafeWadDiv(int256, int256) → int256   ← wad 除法（不检查溢出/除零）
│   ├── wadMul(int256, int256) → int256         ← wad 乘法（带溢出检查）
│   └── wadDiv(int256, int256) → int256         ← wad 除法（带溢出检查）
│
├── 高等数学（3 个）
│   ├── wadPow(int256, int256) → int256         ← x^y = e^(ln(x) × y)
│   ├── wadExp(int256) → int256                 ← e^x（Padé 逼近）
│   └── wadLn(int256) → int256                  ← ln(x)（Padé 逼近）
│
└── 工具函数（1 个）
    └── unsafeDiv(int256, int256) → int256      ← 有符号整数除法（不检查除零）

命名规则：unsafe 前缀/后缀 = 不做溢出或除零检查，由调用方确保输入安全。

四、源码逐行解析

4.1 toWadUnsafe —— 整数转 wad

function toWadUnsafe(uint256 x) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // x × 1e18，直接 mul 不做溢出检查
        r := mul(x, 1000000000000000000)
    }
}

作用：将无符号整数转为 wad 格式（乘以 1e18）。不做溢出检查。

参数：

参数	类型	含义
`x`	`uint256`	要转换的无符号整数

返回值：int256 —— wad 格式的结果。

注意：调用方需确保 x × 1e18 不溢出 int256（即 x ≤ type(int256).max / 1e18）。

4.2 toDaysWadUnsafe —— 秒转 wad 天数

function toDaysWadUnsafe(uint256 x) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // x × 1e18 / 86400
        r := div(mul(x, 1000000000000000000), 86400)
    }
}

作用：将秒数转为 wad 格式的天数。是进入 wad 数学世界的入口，把链上时间差转为高精度天数喂给 wadExp 等函数。

参数：

参数	类型	含义
`x`	`uint256`	秒数

返回值：int256 —— wad 格式的天数（1 天 = 1e18）。

公式：r = x × 1e18 / 86400（1 天 = 86400 秒）

为什么需要这个函数？ 凡是要把时间差喂给 wadExp / wadMul / wadDiv 做连续数学运算的场景，都需要先转成 wad 天数。否则整数除法会把不足一天的部分截断为 0。

典型使用场景：

线性释放（vesting）
代币排放曲线
连续利率计算（DeFi 借贷中按时间连续复利）

4.3 fromDaysWadUnsafe —— wad 天数转秒

function fromDaysWadUnsafe(int256 x) pure returns (uint256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // x × 86400 / 1e18
        r := div(mul(x, 86400), 1000000000000000000)
    }
}

作用：将 wad 格式的天数转回秒数。是 wad 数学世界的出口，把计算结果转回链上可用的秒数。

参数：

参数	类型	含义
`x`	`int256`	wad 格式的天数

返回值：uint256 —— 秒数。

典型调用链：

秒数 -> toDaysWadUnsafe -> wad天数 -> wadExp/wadMul等运算 -> wad结果 -> fromDaysWadUnsafe -> 秒数

4.4 unsafeWadMul —— wad 乘法（不检查溢出）

function unsafeWadMul(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // (x × y) / 1e18
        // 注：sdiv 是 EVM 的有符号整数除法操作码（Signed DIVision），对应无符号版本是 div
        r := sdiv(mul(x, y), 1000000000000000000)
    }
}

作用：wad 定点乘法，不检查溢出。

参数：

参数	类型	含义
`x`	`int256`	第一个 wad 操作数
`y`	`int256`	第二个 wad 操作数

返回值：int256 —— wad 格式的乘积。

数学原理：

a_wad × b_wad = (a × 1e18) × (b × 1e18) = a × b × 1e36
需要除以 1e18 才能得到正确的 wad 结果：a × b × 1e18

4.5 unsafeWadDiv —— wad 除法（不检查溢出和除零）

function unsafeWadDiv(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // (x × 1e18) / y
        r := sdiv(mul(x, 1000000000000000000), y)
    }
}

作用：wad 定点除法，不检查溢出和除零。y 为 0 时返回 0（sdiv 除零行为），不会 revert。

参数：

参数	类型	含义
`x`	`int256`	被除数（wad 格式）
`y`	`int256`	除数（wad 格式）

返回值：int256 —— wad 格式的商。

数学原理：

(a_wad / b_wad) = (a × 1e18) / (b × 1e18) = a / b
需要先乘 1e18 才能得到正确的 wad 结果：(a / b) × 1e18

4.6 wadMul —— wad 乘法（带溢出检查）

function wadMul(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // 先计算 x * y 存入 r
        r := mul(x, y)

        // 组合溢出检查（合并为 1 个表达式，只产生 1 个 JUMPI）：
        if iszero(
            and(
                // 条件 A: or(iszero(x), eq(sdiv(r, x), y))
                //   → x == 0（0 乘任何数不溢出）或 (x*y)/x == y（验证乘法可逆）
                or(iszero(x), eq(sdiv(r, x), y)),
                // 条件 B: or(lt(x, not(0)), sgt(y, 0x80...0))
                //   → x &lt; -1 或 y > type(int256).min
                //   → 排除 x=-1, y=type(int256).min 的二补码陷阱
                or(lt(x, not(0)), sgt(y, 0x8000000000000000000000000000000000000000000000000000000000000000))
            )
        ) {
            revert(0, 0)
        }

        // 溢出检查通过，除以 1e18 得到 wad 结果
        r := sdiv(r, 1000000000000000000)
    }
}

作用：wad 定点乘法，带溢出检查。溢出时 revert（无错误信息）。

参数：

参数	类型	含义
`x`	`int256`	第一个 wad 操作数
`y`	`int256`	第二个 wad 操作数

返回值：int256 —— wad 格式的乘积，溢出时 revert。

溢出检查的两个条件：

条件	含义	作用
A：`x == 0		(x*y)/x == y`
B：`x < -1		y > type(int256).min`

为什么需要条件 B？

当 x = -1, y = type(int256).min（即 -2^255）时：

(-1) × (-2^255) = 2^255 超过 int256.max（2^255-1），产生溢出
但 EVM 的 mul 不区分有符号/无符号，按无符号相乘后取低 256 位

  EVM mul 实际计算：把两个操作数当作无符号整数相乘

  -1 的补码 = 0xFFFF...FFFF（无符号值 = 2^256 - 1）：
  ┌──────────────────────────────────────────────────────────────────────┐
  │ 0xFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF   │
  └──────────────────────────────────────────────────────────────────────┘

  ×

  -2^255 的补码 = 0x8000...0000（无符号值 = 2^255）：
  ┌──────────────────────────────────────────────────────────────────────┐
  │ 0x8000000000000000000000000000000000000000000000000000000000000000   │
  └──────────────────────────────────────────────────────────────────────┘

  ＝

  无符号乘积：(2^256 - 1) × 2^255 = 2^511 - 2^255

  完整 512 位结果（二进制）：
  高 256 位                                低 256 位
  ┌──────────────────────────────────┬──────────────────────────────────┐
  │ 0111 1111 ... 1111 1111 1111 1111│ 1000 0000 ... 0000 0000 0000 0000│
  └──────────────────────────────────┴──────────────────────────────────┘
    					↑ 高 256 位丢弃                   ↑ EVM 只保留低 256 位

  截断后的低 256 位：
  ┌──────────────────────────────────────────────────────────────────────┐
  │ 0x8000000000000000000000000000000000000000000000000000000000000000   │
  └──────────────────────────────────────────────────────────────────────┘
  = -2^255 = type(int256).min = y

此时 sdiv(-2^255, -1) = -2^255 = y，条件 A 误判为"未溢出"

注：这不是补码层面的计算逻辑，而是 EVM 规范的硬编码特例（Yellow Paper中有定义：https://ethereum.github.io/yellowpaper/paper.pdf ）：

条件 B 额外捕获这种二补码陷阱

为什么用 lt/sgt 而不是 != -1 和 != min？

Yul 没有 neq 操作码，用 lt(x, not(0))（即 x < -1）和 sgt(y, 0x80...0)（即 y > min）替代 iszero(eq(...))，省 2 个操作码。虽然范围更宽，但只会放过更多安全情况，不会漏掉目标。

gas 优化技巧：两个条件合并为 1 个 and 表达式，编译后只产生 1 个 JUMPI 而非 2 个。

4.7 wadDiv —— wad 除法（带溢出检查）

function wadDiv(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        // 先计算 x * 1e18 存入 r
        r := mul(x, 1000000000000000000)

        // 检查：y != 0 且 r / 1e18 == x（确保 x * 1e18 没有溢出）
        // 注：iszero(iszero(y)) 是 Yul 中判断 y != 0 的标准写法
        if iszero(and(iszero(iszero(y)), eq(sdiv(r, 1000000000000000000), x))) {
            revert(0, 0)
        }

        // r 除以 y 得到最终结果
        r := sdiv(r, y)
    }
}

作用：wad 定点除法，带溢出检查。除零或溢出时 revert。

参数：

参数	类型	含义
`x`	`int256`	被除数（wad 格式）
`y`	`int256`	除数（wad 格式）

返回值：int256 —— wad 格式的商，溢出或除零时 revert。

溢出检查：y != 0 && (x * 1e18) / 1e18 == x，确保除数不为零且 x * 1e18 没有溢出。

4.8 wadPow —— wad 幂运算

function wadPow(int256 x, int256 y) pure returns (int256) {
    // x^y 等价于 e^(ln(x) * y)
    // 其中，wadLn() 要求 x > 0
    return wadExp((wadLn(x) * y) / 1e18);
}

作用：wad 定点幂运算。

参数：

参数	类型	含义
`x`	`int256`	底数（wad 格式，必须为正）
`y`	`int256`	指数（wad 格式）

返回值：int256 —— wad 格式的 x^y。

数学原理：x^y = (e^ln(x))^y = e^(ln(x) × y)

注意：只支持正数底数（x > 0），因为 ln对非正数无定义。

4.9 wadExp —— 自然指数函数 e^x

这是整个库最复杂的函数，算法来自 Remco Bloemen。

函数签名

function wadExp(int256 x) pure returns (int256 r) {
    unchecked { ... }
}

参数：

参数	类型	含义
`x`	`int256`	wad 格式的指数（代表实数 x/1e18）

返回值：int256 —— wad 格式的 e^(x/1e18)。

算法流程（5 步）

wadExp(x)
  │
  ├── 1. 边界检查
  │     ├── x ≤ -42139678854452767551 → return 0（结果太小）
  │     └── x ≥ 135305999368893231589 → revert（结果溢出）
  │
  ├── 2. 基数转换：10^18 定点 → 2^96 定点
  │     └── x = (x &lt;&lt; 78) / 5^18
  │
  ├── 3. 范围缩减：exp(x) = exp(x') × 2^k
  │     ├── k = round(x / ln2)
  │     └── x' = x - k × ln2，使 x' ∈ (-½ln2, ½ln2)
  │
  ├── 4. (6,7) 阶 Padé 有理逼近
  │     ├── p(x)：5 次多项式（分子）
  │     ├── q(x)：6 次多项式（分母）
  │     └── r = p / q
  │
  └── 5. 最终组装
        └── r = r × s × 2^k × 10^18 / 2^96（一步完成）

第 1 步：边界检查

// 1.1 下界检查
//  wad 是整数，所以最小正值就是 1（代表实数 1/1e18 = 10^-18）
//  如果计算结果的实数值 &lt; 0.5 × 10^-18，四舍五入后就是 0，直接返回 0
//  -42139678854452767551 = ⌊ln(0.5e-18) × 1e18⌋
if (x &lt;= -42139678854452767551) return 0;

// 1.2 上界检查
//  如果 e^(x/1e18) × 1e18 超过 type(int256).max = 2^255-1，int256 就装不下了
//  135305999368893231589 = ⌊ln((2^255-1)/1e18) × 1e18⌋
if (x >= 135305999368893231589) revert("EXP_OVERFLOW");

x 的有效范围是 (-42.14, 135.31) × 1e18。

第 2 步：基数转换

x = (x &lt;&lt; 78) / 5 ** 18;

为什么要转换？ x 是 10^18 定点数，每次乘除都涉及 / 1e18，EVM 做除法很贵。转成 2^96 定点后，除法变成移位（>> 96），便宜很多。

注：SHR 3 gas，DIV/SDIV 5 gas，单次省 2 gas。

转换因子推导：

10^18 定点：实数 1.0 存储为 1e18
2^96 定点：实数 1.0 存储为 2^96

x_new = x × (2^96 / 10^18)
      = x × 2^96 / (2×5)^18
      = x × 2^96 / 2^18 / 5^18
      = x × 2^78 / 5^18
      = (x &lt;&lt; 78) / 5^18

为什么选 2^96？

原因	说明
精度	2^96 ≈ 7.9 × 10^28，约 29 位十进制精度，远超 wad 的 18 位
不溢出	两个 2^96 数相乘得 2^192，不超过 256 位
省 gas	`>> 96` 替代 `/ 1e18`，Padé 逼近有十几次乘除，每次都省一点

第 3 步：范围缩减

// k = round(x / ln2) = floor(x / ln2 + 0.5)
// 常数：54916777467707473351141471128 = ln(2) × 2^96
// 常数：2^95 = 0.5 × 2^96
int256 k = ((x &lt;&lt; 96) / 54916777467707473351141471128 + 2 ** 95) >> 96;
// x' = x - k × ln2
x = x - k * 54916777467707473351141471128;

为什么要范围缩减？

Padé 逼近只在很小的范围内精确。exp(100) 很大，直接逼近误差爆炸。利用数学恒等式：

exp(x) = exp(x - k×ln2 + k×ln2)
       = exp(x - k×ln2) × exp(k×ln2)
       = exp(x') × 2^k

其中 x' = x - k×ln2。选合适的 k 让 x' 很小，最后乘以 2^k（移位操作，低 gas）还原。

如何确定 k？

要让 x' = x - k × ln2 足够小，就看 x 中含有多少个 ln2，将这个数定为 k：

k = round(x / ln2)

令 x / ln2 的真实值 = k + ε，其中 ε 是误差。由于 round 是四舍五入，所以 |ε| < 0.5。

x' = x - k × ln2
   = (k + ε) × ln2 - k × ln2
   = ε × ln2

因为 |ε| &lt; 0.5，所以 |x'| = |ε × ln2| &lt; 0.5 × ln2 = ½ln2
即 x' 落在 (-½ln2, ½ln2) 范围内

为什么用 round 而不是 floor？

round: k = round(x / ln2)  → 误差 |ε| &lt; 0.5  → x' ∈ (-½ln2, ½ln2)，关于原点对称
floor: k = floor(x / ln2)  → 误差 ε ∈ [0, 1)  → x' ∈ [0, ln2)，偏离原点

Padé 逼近在原点附近最精确，round 产生的 x' 范围关于原点对称，逼近精度更好。

round 的工程实现：

k = round(x / ln2) = floor(x / ln2 + 0.5)

为了保留小数部分，将整体扩大 2^96 倍做加法，再缩小 2^96 倍实现 floor：

k = ((x × 2^96 / LN2_96) + 0.5_96) >> 96

其中 LN2_96 = 54916777467707473351141471128，0.5_96 = 2^95。

通过 x 的范围可以计算出 k 的范围：[-61, 195]。

第 4 步：(6,7) 阶 Padé 有理逼近

什么是 Padé 逼近？

Taylor 级数用多项式逼近函数，Padé 用有理函数（两个多项式相除）逼近：

Taylor: e^x ≈ 1 + x + x²/2 + x³/6 + ...        ← 多项式
Padé:   e^x ≈ P(x) / Q(x)                      ← 有理函数

为什么 Padé 更好？因为 exp(x) 在 x → -∞ 时趋近 0。多项式只能趋近 ±∞，无法趋近 0；而 P(x)/Q(x) 可以（只要分母增长比分子快）。

系数的确定：一旦确定逼近目标（e^x）、阶数（(6,7)）、展开点（x = 0），系数就唯一确定。

Horner 法则：多项式求值的快速算法，通过嵌套形式减少乘法次数：

标准形式：ax³ + bx² + cx + d        ← 需要 6 次乘法 + 3 次加法
Horner： ((a*x + b)*x + c)*x + d   ← 需要 3 次乘法 + 3 次加法

以下代码就是用 Horner 法则计算分子 P(x) 和分母 Q(x)，其中系数是预先算好的魔数常量（来自 Remco Bloemen 的论文）。

分子 p(x) 的计算（非标准 Horner，借用中间变量 y 做交叉展开）：

// y = x + c0 → 1 次
int256 y = x + 1346386616545796478920950773328;
// y = y * x + c1 = x²+c0*x+c1 → 2 次
y = ((y * x) >> 96) + 57155421227552351082224309758442;
// p = y + x + c2 → 2 次
int256 p = y + x - 94201549194550492254356042504812;
// p = p*y + c3（2次*2次）→ 4 次
p = ((p * y) >> 96) + 28719021644029726153956944680412240;  
// 最后一步故意不 >> 96，p 留在 2^192 基数，省一次移位
// 这就是为什么常数4385272521454847904659076985693276也要&lt;&lt; 96 —— 把它也提升到 2^192 基数与 p 对齐再做加法
// p = p*x + c4 → 5 次
p = p * x + (4385272521454847904659076985693276 &lt;&lt; 96);

分母 q(x) 的计算（标准 Horner）：

// 每步 * x >> 96 就是 2^96 定点下的乘法。6 步展开 = 6 次多项式
int256 q = x - 2855989394907223263936484059900;
q = ((q * x) >> 96) + 50020603652535783019961831881945;
q = ((q * x) >> 96) - 533845033583426703283633433725380;
q = ((q * x) >> 96) + 3604857256930695427073651918091429;
q = ((q * x) >> 96) - 14423608567350463180887372962807573;
q = ((q * x) >> 96) + 26449188498355588339934803723976023;

计算 r = p / q：

assembly {
    // p(2^192基数) / q(2^96 基数) = r(2^96 基数)
    // gas技巧：Solidity的 / 运算符即使在 unchecked {} 块中，编译器仍然会插入一个"除数是否为零"的检查，
    // 由于 q(x) 是一个 6 次多项式，它的 6 个根（令 q = 0 的解）全部是复数，没有实数根
    // 这是Padé逼近的数学性质保证的——论文作者在选系数时就确保了分母没有实数根
    // 由于 x 取值在 int256 范围内（实数域的子集），而 q(x) 的根全是复数，所以可知 q(x) 在实数域上永远不会为0，这里特意使用Yul的sdiv来跳过除数不为0的检查
    r := sdiv(p, q)
}

此时 r 的范围是 (0.09, 0.25) × 2^96。

第 5 步：最终组装

// CONST = round(s × 10^18 × 2^99) = 3822833074963236453042738258902158003155416615667
r = int256((uint256(r) * 3822833074963236453042738258902158003155416615667) >> uint256(195 - k));

现在的 r 是什么？

如何将 exp(x') 还原成 exp(x)？

几个遗留问题：

目前 r 在 2^96 基数下，但最终结果是 wad（10^18 基数）
范围缩减时拆出了 2^k，还没乘回来：exp(x) = exp(x') × 2^k = r × s × 2^k

所以一步要完成三个操作：

乘以缩放因子 s ≈ 6.031367120：补偿首一多项式的缩放
乘以 2^k：恢复范围缩减中分离的指数部分
基数转换 10^18 / 2^96：从二进制定点转回 wad

数学推导：

result = r × s × 2^k × 10^18 / 2^96
       = r × (s × 10^18 × 2^99) × 2^k / 2^(96+99)
       = r × (s × 10^18 × 2^99) / 2^(195-k)
       = r × CONST >> (195 - k)

CONST = round(s × 10^18 × 2^99) = 3822833074963236453042738258902158003155416615667

为什么要凑出 2^195？

因为 k 的范围是 [-61, 195]，右移量 195-k 的范围是 [0, 256]，始终安全。

4.10 wadLn —— 自然对数 ln(x)

同样来自 Remco Bloemen，是 wadExp 的逆运算。

函数签名

function wadLn(int256 x) pure returns (int256 r) {
    unchecked { ... }
}

参数：

参数	类型	含义
`x`	`int256`	wad 格式的输入（代表实数 x/1e18，必须为正）

返回值：int256 —— wad 格式的 ln(x/1e18)。

算法流程（5 步）

wadLn(x)
  │
  ├── 1. 域检查：require(x > 0)
  │
  ├── 2. 二分搜索找 r = floor(log2(x))（8 步）
  │
  ├── 3. 范围缩减：将x归一化到 [2^96, 2^97)
  │     ├── k = r - 96
  │     ├── x &lt;&lt;= (159 - k)
  │     └── x = int256(uint256(x) >> 159)
  │
  ├── 4. (8,8) 阶 Padé 有理逼近
  │     ├── p(x)：7 次多项式（分子）
  │     ├── q(x)：7 次多项式（分母）
  │     └── r = p / q
  │
  └── 5. 最终组装
        ├── r *= s_const        ← 补偿缩放因子
        ├── r += ln2_const × k  ← 补偿归一化
        ├── r += offset_const   ← 补偿 ln(2^96/10^18)
        └── r >>= 174           ← 转回 wad

第 1 步：域检查

require(x > 0, "UNDEFINED");

ln 对非正数无定义。

第 2 步：二分搜索找 floor(log2(x))

目的是找到 x 的最高有效位（MSB）在第几位。

思路：将 256 位 MSB 位置拆解为 8 个二进制位 (b7~b0)，每步将搜索范围减半：

1. b7: x > 2^128-1? → 确定 MSB 在上半还是下半 256 位
2. b6: x>>r > 2^64-1? → 继续二分
3. ...依此类推直到 b0

结果：r = 128×b7 + 64×b6 + ... + 1×b0 = floor(log2(x))

assembly {
    // b7：判断 x > 2^128 - 1，得出 MSB 在高 128 位 (128~255) 还是低 128 位 (0~127)
    // 如果是：r = 128，否则：r = 0
    r := shl(7, lt(0xffffffffffffffffffffffffffffffff, x))
    // b6：将已确定的位移掉后，判断 剩余部分 > 2^64 - 1
    // 如果是：r += 64，否则：r不变
    r := or(r, shl(6, lt(0xffffffffffffffff, shr(r, x))))
    // b5：将已确定的位移掉后，判断 剩余部分 > 2^32 - 1
    // 如果是：r += 32，否则：r不变
    r := or(r, shl(5, lt(0xffffffff, shr(r, x))))
    // b4：将已确定的位移掉后，判断 剩余部分 > 2^16 - 1
    // 如果是：r += 16，否则：r不变
    r := or(r, shl(4, lt(0xffff, shr(r, x))))
    // b3：将已确定的位移掉后，判断 剩余部分 > 2^8 - 1
    // 如果是：r += 8，否则：r不变
    r := or(r, shl(3, lt(0xff, shr(r, x))))
    // b2：将已确定的位移掉后，判断 剩余部分 > 2^4 - 1
    // 如果是：r += 4，否则：r不变
    r := or(r, shl(2, lt(0xf, shr(r, x))))
    // b1：将已确定的位移掉后，判断 剩余部分 > 2^2 - 1
    // 如果是：r += 2，否则：r不变
    r := or(r, shl(1, lt(0x3, shr(r, x))))
    // b0：将已确定的位移掉后，判断 剩余部分 > 2^1 - 1
    // 如果是：r += 1，否则：r不变
    r := or(r, lt(0x1, shr(r, x)))
    // 此时r = floor(log2(x))
}

原理：以 b7 为例——如果 x > 2^128 - 1，说明 MSB 在第 128~~255 位，r 先设为 128；否则 MSB 在第 0~~127 位，r = 0。然后将 x 右移 r 位，在剩余范围内继续二分。

第 3 步：范围缩减（归一化）

先提取 m，使得 x = 2^r × m，m ∈ [1, 2)，并将 m 表示为 2^96 定点数。

为什么要归一化？

Padé 逼近只在小范围内精确，需要将 x 归一化到固定区间 [2^96, 2^97)。

数学依据：ln(x) = ln(2^r × m) = r × ln(2) + ln(m)
归一化后 Padé 只负责算 ln(m)，r × ln(2) 在第 5 步补回

m 为什么用 2^96 定点数表示？

同 wadExp 中 x 用 2^96 定点数表示的原因（精度、不溢出、移位省 gas）。将 m 表示为 2^96 定点数后：

归一化前：x = 2^r × m
归一化后：x_new = 2^96 × m
令 k = r - 96（即 r = k + 96），则：
x_new = 2^96 × m = 2^r × m × 2^(96-r) = x × 2^(-k)

r × ln(2) 怎么补回？

由前面的分解：ln(x) = r × ln(2) + ln(m)，代入 r = k + 96：

ln(x) = (k + 96) × ln(2) + ln(m)
归一化后 Padé 只负责算 ln(m)，剩下的 (k + 96) × ln(2) 需要在第 5 步补回
但 floor(log2(x)) 存在变量 r 中，第 4 步被 p/q 覆盖，只剩 k
所以拆成两部分补：k × ln(2) 显式加回，96 × ln(2) 合并到常数 ln(2^96/10^18) 中补

归一化

// 对 x 做归一化，即把 x 从 2^r × m 变成 2^96 × m，即把 MSB 从第 r 位挪到第 96 位
// 这里先左移到顶再右移，而不是直接 x >> (r-96)。原因是：r 可能 &lt; 96，EVM 没法右移一个负数位

// k 为归一化偏移量
int256 k = r - 96;

// 将 x 的 MSB（第 r 位）左移到第 255 位
// 左移位数为 255-r = 255-(k+96) = 159-k
// 目的：把有效位对齐到最高位，为下一步统一右移做准备
x &lt;&lt;= uint256(159 - k);

// 再从第 255 位统一移到第 96 位（255 - 159 = 96）
// 使用 uint256 cast 确保逻辑右移（填充 0），避免算术右移（填充符号位）
x = int256(uint256(x) >> 159);

// 此时，x 被归一化到 [2^96, 2^97) 范围内，即 2^96 定点数下的 [1, 2)

第 4 步：(8,8) 阶 Padé 有理逼近

比 wadExp 的 (6,7) 阶更高，因为 ln 在 [1,2) 上的曲率比 exp 在 (-½ln2, ½ln2) 上更大，需要更高阶才能达到同等精度。

分子 p(x) 的计算（7 次多项式）：

int256 p = x + 3273285459638523848632254066296;
p = ((p * x) >> 96) + 24828157081833163892658089445524;
p = ((p * x) >> 96) + 43456485725739037958740375743393;
p = ((p * x) >> 96) - 11111509109440967052023855526967;
p = ((p * x) >> 96) - 45023709667254063763336534515857;
p = ((p * x) >> 96) - 14706773417378608786704636184526;
// 最后一步故意不 >> 96，p 留在 2^192 基数
// 目的：最终要算 r = p / q。这样做省了一次移位操作
p = p * x - (795164235651350426258249787498 &lt;&lt; 96);

分母 q(x) 的计算（7 次多项式）：

int256 q = x + 5573035233440673466300451813936;
q = ((q * x) >> 96) + 71694874799317883764090561454958;
q = ((q * x) >> 96) + 283447036172924575727196451306956;
q = ((q * x) >> 96) + 401686690394027663651624208769553;
q = ((q * x) >> 96) + 204048457590392012362485061816622;
q = ((q * x) >> 96) + 31853899698501571402653359427138;
// q 在 2^96 基数
q = ((q * x) >> 96) + 909429971244387300277376558375;

计算 r = p / q：

assembly {
    // 同 wadExp 中的 gas 技巧：q(x) 的根全是复数，实数域上不为零，
    // 使用 Yul 的 sdiv 跳过 Solidity 编译器自动插入的除零检查
    r := sdiv(p, q)
}

此时 r 的范围是 (0, 0.125) × 2^96。

第 5 步：最终组装

现在的 r 是什么？

前面计算的 p(x) 和 q(x) 都是首一多项式，结果差一个缩放因子 s ≈ 5.549。真正的 ln(m) = r × s。

从分解公式推导最终组装

目标：ln(x / 10^18)

由 x = 2^floor(log2(x)) × m = 2^(k+96) × m，得：

ln(x / 10^18) = ln(2^(k+96) × m / 10^18)
              = ln(m) + (k + 96) × ln(2) - ln(10^18)
              = ln(m) + k × ln(2) + 96 × ln(2) - ln(10^18)
              = ln(m) + k × ln(2) + ln(2^96 / 10^18)

三项分别对应：

r × s → ln(m)：补偿缩放因子，得到真正的 ln(m)
k × ln(2) → 补回 r × ln(2) = (k + 96) × ln(2) 中的 k 份（r 的值已被 p/q 覆盖，目前只有 k 可用）
ln(2^96 / 10^18) → 常数补回，包含两项：
- 96 × ln(2)：(k + 96) × ln(2) 中 k 份已显式补回，这里补剩余的 96 份
- -ln(10^18)：输入 x 是 wad 定点，目标是 ln(x/10^18)，需减掉多出的 ln(10^18)

代码实现：在 5^18 × 2^192 基数下统一计算

数学目标：result_wad = (r × s + k × ln(2) + ln(2^96 / 10^18)) × 10^18

但不能直接乘 10^18（r 在 2^96 基数，直接乘会损失精度或溢出）。

技巧：10^18 = 5^18 × 2^18，选 5^18 × 2^192 作为中间基数，算完后只需 >> 174（= 192 - 18）就能转回 wad（5^18 × 2^18 = 10^18）。

// 补偿缩放因子 s，并将基数提高到 5^18 × 2^192
// 常数 1677202110996718588342820967067443963516166 = s × 5^18 × 2^96（r 是 2^96 基数，乘完变为 5^18 × 2^192 基数）
r *= 1677202110996718588342820967067443963516166;

// 补回 r × ln(2) = (k + 96) × ln(2) 中的 k 份，并提升到 5^18 × 2^192 基数下
// 常数 16597577552685614221487285958193947469193820559219878177908093499208371 = ln(2) × 5^18 × 2^192
r += 16597577552685614221487285958193947469193820559219878177908093499208371 * k;

// 补回 ln(2^96 / 10^18)（= 96×ln(2) - ln(10^18)，即 (k+96)×ln(2) 中剩余的 96 份和基数转换）
// 常数 600920179829731861736702779321621459595472258049074101567377883020018308 = ln(2^96 / 10^18) × 5^18 × 2^192
r += 600920179829731861736702779321621459595472258049074101567377883020018308;

// 从 5^18 × 2^192 基数转回 wad（10^18 = 5^18 × 2^18）
// 5^18 × 2^192 / 2^174 = 5^18 × 2^18 = 10^18
r >>= 174;

4.11 unsafeDiv —— 有符号整数除法（不检查除零）

function unsafeDiv(int256 x, int256 y) pure returns (int256 r) {
    /// @solidity memory-safe-assembly
    assembly {
        r := sdiv(x, y)
    }
}

作用：有符号整数除法，不检查除零。y 为 0 时返回 0（sdiv 除零行为），不会 revert。

参数：

参数	类型	含义
`x`	`int256`	被除数
`y`	`int256`	除数

返回值：int256 —— 商。

注意：这不是 wad 除法，是普通整数除法（没有 × 1e18 的步骤）。

五、wadExp / wadLn 算法流程图

wadExp 流程

输入: x (wad 格式)
  │
  ├── x ≤ -42139678854452767551 ?
  │     └── YES → return 0
  │
  ├── x ≥ 135305999368893231589 ?
  │     └── YES → revert("EXP_OVERFLOW")
  │
  ├── 基数转换: x = (x &lt;&lt; 78) / 5^18
  │     └── 10^18 定点 → 2^96 定点
  │
  ├── 范围缩减:
  │     ├── k = round(x / ln2)                    ← k ∈ [-61, 195]
  │     └── x' = x - k × ln2                      ← x' ∈ (-½ln2, ½ln2)
  │
  ├── Padé 逼近:
  │     ├── p = P(x')  (5 次多项式，2^192 基数)
  │     ├── q = Q(x')  (6 次多项式，2^96 基数)
  │     └── r = p / q  (2^96 基数)
  │
  └── 最终组装:
        └── r = r × CONST >> (195 - k)
              ↓
        输出: e^(x/1e18) (wad 格式)

wadLn 流程

输入: x (wad 格式, x > 0)
  │
  ├── require(x > 0)
  │
  ├── 二分搜索: r = floor(log2(x))
  │     └── 8 步二分，确定 MSB 位置
  │
  ├── 归一化:
  │     ├── k = r - 96
  │     ├── x &lt;&lt;= (159 - k)                        ← MSB 对齐到第 255 位
  │     └── x = uint256(x) >> 159                  ← MSB 移到第 96 位
  │           └── x ∈ [2^96, 2^97)
  │
  ├── Padé 逼近:
  │     ├── p = P(x)  (7 次多项式，2^192 基数)
  │     ├── q = Q(x)  (7 次多项式，2^96 基数)
  │     └── r = p / q  (2^96 基数)
  │
  └── 最终组装 (5^18 × 2^192 基数):
        ├── r *= s_const                           ← 补偿缩放因子
        ├── r += ln2_const × k                     ← 补偿归一化
        ├── r += offset_const                      ← 补偿 ln(2^96/10^18)
        └── r >>= 174                              ← 转回 wad
              ↓
        输出: ln(x/1e18) (wad 格式)

六、设计思想

6.1 极简主义

所有函数都是 free function，不依赖任何合约或 library
没有自定义 error、event、modifier，没有 storage
统一使用 int256 类型，消除调用链中的类型转换

6.2 gas 极致优化

技巧	说明
2^96 定点数	`>> 96` 替代 `/ 1e18`，移位比除法便宜
unchecked	wadExp/wadLn 全程 unchecked，省去溢出检查 gas
Yul sdiv	跳过 Solidity 编译器自动插入的除零检查（因为数学上证明了分母不为 0）
合并溢出检查	wadMul 两个条件合并为 1 个 JUMPI
分子不右移	p 留在 2^192 基数，p/q 自动得到 2^96 结果，省一次移位
常数预计算	所有常数都是编译期已知的魔数，运行时无额外计算

6.3 Padé 逼近 vs Taylor 级数

对比项	Taylor 级数	Padé 逼近
形式	多项式	有理函数 P(x)/Q(x)
逼近 exp(x) 在 x→-∞ 的行为	多项式趋向 ±∞，无法趋近 0	分母增长快于分子，可趋近 0 ✅
同等精度所需项数	多	少 ✅
EVM 实现复杂度	简单	稍复杂（需要除法）

6.4 安全 vs 性能的取舍

库提供了 safe/unsafe 两套函数，让调用方自行选择：

函数	溢出检查	除零检查	gas
unsafeWadMul	❌	-	低
wadMul	✅	-	高
unsafeWadDiv	❌	❌	低
wadDiv	✅	✅	高

七、安全注意事项

风险	建议
unsafe 函数不检查溢出/除零	调用方需自行保证输入安全，或优先使用 safe 版本
wadExp 输入范围有限	x 必须在 (-42.14, 135.31) × 1e18 范围内，低于下界返回 0，高于上界revert
wadLn 只接受正数	x ≤ 0 会 revert("UNDEFINED")
wadPow 只支持正数底数	内部调用 wadLn，底数必须 > 0
wad 精度限制	最小正值为 1（代表 10^-18），小于此值的结果会被截断为 0
整数截断	wad 乘除法的中间结果除以 1e18 是整数截断（向零取整），非四舍五入
unsafeWadDiv 除零返回 0	`sdiv(x, 0)` 在 EVM 中返回 0 而非 revert，可能导致静默错误

八、与同类方案对比

对比项	Solmate SignedWadMath	OpenZeppelin Math	PRBMath
类型	int256（有符号）	uint256（无符号）	int256 / uint256 两套
精度	18 位小数	无定点	18 位小数
exp/ln	✅ Padé 逼近	❌ 无	✅ Padé 逼近
幂运算	✅ wadPow	❌ 无	✅ pow
形式	free function	library	library
gas	极低（全 assembly）	中	中低
安全检查	safe + unsafe 两套	全部带检查	全部带检查
时间转换	✅ toDays/fromDays	❌	❌

九、测试实战

Mock 合约：https://github.com/RevelationOfTuring/foundry-solmate/blob/main/src/utils/MockSignedWadMath.sol

全部 Foundry 测试合约：https://github.com/RevelationOfTuring/foundry-solmate/blob/main/test/utils/SignedWadMath.t.sol

学分: 0
分类: 框架/库
标签: solmate Solidity Foundry SignedWadMath wad exp