上个月排查一段聚类代码，距离矩阵里冒出了 -3.7e-9 这样的负数。距离平方不可能为负，所以问题不在聚类算法本身，而在那句看起来人畜无害的‖x‖² + ‖y‖² - 2x·y。这个式子没错，数学上严丝合缝，它就是范数平方的展开公式在距离计算里的标准用法。错的是我把它用在了两个几乎完全重合的向量上，浮点减法把有效位数吃干净了。

范数平方的展开公式，说穿了就是把‖x±y‖²这种"先做加减、再做范数"的表达式，拆成‖x‖²、‖y‖²和交叉项2⟨x,y⟩的组合。它在机器学习、数值线性代数、最优化、核方法里出现的频率高得离谱：最小二乘的正规方程是它推出来的，RBF 核的定义是它写出来的，PCA 的方差目标靠它变换，K-means 的簇心分配更是直接用它省掉了一半计算量。

真正上手用的时候，纸面上看不出来的细节有一堆。什么条件下才能展开、展开之后数值上会发生什么、复数域要不要取实部、p 范数为什么不能这么玩、加权范数对矩阵有什么硬性要求。这篇就把推导、代码和踩过的坑一次性理清楚，适合写过梯度下降、跑过 SVM、手撸过距离矩阵的人，也适合正在补线性代数的人把它当成一条实打实的应用线索。

1. 展开公式的骨架：范数平方为什么能拆成三项

1.1 能展开的前提是范数由内积诱导

很多人一上来就背‖x+y‖² = ‖x‖² + 2⟨x,y⟩ + ‖y‖²，但没注意这个式子成立是有前提的。范数的定义只要求三条：正定性（‖x‖ ≥ 0，且‖x‖ = 0当且仅当x = 0）、齐次性（‖αx‖ = |α|‖x‖）、三角不等式（‖x+y‖ ≤ ‖x‖ + ‖y‖）。这三条里面没有一条能推出展开公式。

能推出展开公式的是更强的一个条件：这个范数是由内积诱导出来的。也就是说，存在一个内积⟨·,·⟩，使得‖x‖² = ⟨x,x⟩。实空间上的内积要满足双线性（对两个位置都线性）、对称性（⟨x,y⟩ = ⟨y,x⟩）和正定性（⟨x,x⟩ ≥ 0，等号只在零向量处成立）。

反过来问：随便给一个范数，都能找到内积让它变成"内积诱导范数"吗？不能。Jordan-von Neumann 定理给出了判别标准——当且仅当这个范数满足平行四边形恒等式‖x+y‖² + ‖x-y‖² = 2‖x‖² + 2‖y‖²。这条定理把"能不能展开"变成了一道可以验算的题目，第 2 节会拿具体反例走一遍。

我习惯用一句话记住这个层次关系：内积是比较两个向量的"互相投影量"，范数平方是单个向量的"长度平方"。展开公式告诉我们，长度平方不是孤立的数字，它和互相投影量之间有确定的代数关系。这个关系在有限维实空间里最常见的载体就是 2-范数，也就是欧几里得范数。

1.2 三项式展开的逐项推导

推一遍不费什么时间，但推过之后对每一项的来源会清楚很多。把x+y当作内积的第一个参数：

‖x+y‖² = ⟨x+y, x+y⟩

对第一个位置用线性，拆成两项：

= ⟨x, x+y⟩ + ⟨y, x+y⟩

再对第二个位置用线性（实空间下对称，不用管顺序）：

= ⟨x,x⟩ + ⟨x,y⟩ + ⟨y,x⟩ + ⟨y,y⟩

最后用对称性⟨x,y⟩ = ⟨y,x⟩合并交叉项：

‖x+y‖² = ‖x‖² + 2⟨x,y⟩ + ‖y‖²

把y换成-y，利用⟨x,-y⟩ = -⟨x,y⟩，立刻得到差的形式：

‖x-y‖² = ‖x‖² - 2⟨x,y⟩ + ‖y‖²

记忆上有个偷懒的办法：这两个式子和初中的(a+b)² = a² + 2ab + b²、(a-b)² = a² - 2ab + b²结构完全一样，只是把ab换成了⟨x,y⟩。唯一需要小心的是⟨x,y⟩可以是负数，而‖x‖²永远非负。

顺手补一个几何解释。在欧几里得空间里，⟨x,y⟩ = ‖x‖‖y‖cosθ，θ 是两个向量的夹角。代进差的形式就得到‖x-y‖² = ‖x‖² + ‖y‖² - 2‖x‖‖y‖cosθ，这就是余弦定理。所以展开公式不是什么新东西，它和余弦定理是同一件事的两种写法，只是一个用内积表达，一个用角度和长度表达。

1.3 两个必背推论：相加得平行四边形，相减得极化

把和形式与差形式分别相加、相减，能直接挖出两条极其有用的推论。

先相加：

‖x+y‖² + ‖x-y‖² = 2‖x‖² + 2‖y‖²

交叉项+2⟨x,y⟩和-2⟨x,y⟩正好抵消。这就是平行四边形恒等式。几何直觉很直白：平行四边形的两条对角线长度的平方和，等于四条边长度的平方和。它的实际用途是判别——只靠范数值就能判断一个范数是不是由内积诱导的，不需要去构造内积。

再相减：

‖x+y‖² - ‖x-y‖² = 4⟨x,y⟩ ⟨x,y⟩ = (1/4)(‖x+y‖² - ‖x-y‖²)

这是极化恒等式在实空间的形式。它的意义在于：如果你手头只有距离信息（比如从数据库里拿到的只有距离矩阵，原始向量拿不到），你依然可以通过它把内积恢复出来，进而构造 Gram 矩阵，做 MDS、核 PCA 这类降维。经典 MDS 里那个G = -0.5 · J · D² · J的公式，根就在这里，J = I - (1/n)11ᵀ是中心化矩阵。

2. 差形式的工程红利：K-means 与近邻检索都在吃这个红利

2.1 为什么差形式是距离计算的默认工具

差形式‖x-y‖² = ‖x‖² - 2⟨x,y⟩ + ‖y‖²的价值不在于好看，而在于它把"两个向量的差"这个耦合量，拆成了三个各自独立的量。‖x‖²只跟 x 有关，‖y‖²只跟 y 有关，⟨x,y⟩是一次点积。这意味着所有‖y‖²可以在预处理阶段全部算好存起来，查询的时候只需要算一次点积。

K-means 是把这个红利吃到极致的例子。分配步要做的是：

argmin_k ‖x - μ_k‖²

直接用展开：

argmin_k (‖x‖² - 2⟨x, μ_k⟩ + ‖μ_k‖²)

‖x‖²对所有候选簇心 k 都是同一个数，跟 argmin 完全无关，直接扔掉：

argmin_k (‖μ_k‖² - 2⟨x, μ_k⟩)

簇心的范数平方在更新簇心之后算一次就行，整个分配阶段只剩一次矩阵乘法X @ Mᵀ。n 个样本、d 维特征、K 个簇心，复杂度从逐对计算变成了一个大 GEMM。这个技巧在 scikit-learn 的实现里到处都是，手写 K-means 的时候一定要用上，不然速度差一个数量级。

同样的套路在 kNN、推荐系统的召回层、注意力机制里的相似度计算中都能套。只要出现"和大量固定向量的距离比较"，先想想能不能把其中一部分项提到循环外面。

2.2 平行四边形恒等式：拿 l1 范数验一次

判断一个范数能不能展开，最省事的办法就是拿两个简单向量试平行四边形恒等式。取x = (1, 0)、y = (0, 1)，两个都是二维单位向量。

先看 1-范数（各元素绝对值之和）。‖x+y‖₁ = |1| + |1| = 2，平方得 4；‖x-y‖₁ = |1| + |-1| = 2，平方得 4；左边加起来是 8。右边2‖x‖₁² + 2‖y‖₁² = 2·1 + 2·1 = 4。8 ≠ 4，l1 范数不满足平行四边形恒等式，所以它不由任何内积诱导。

再看 ∞-范数（最大绝对值）。‖x+y‖∞ = 1，‖x-y‖∞ = 1，左边是 1 + 1 = 2。右边还是 4。仍然不等。

最后看 2-范数。‖x+y‖₂ = √2，平方得 2；‖x-y‖₂ = √2，平方得 2；左边是 4。右边2·1 + 2·1 = 4。相等。

结论很清楚：只有 2-范数这一类"由内积导出"的范数才能享受展开公式带来的全部便利。这个结论的实际影响非常大，第 6 节会展开讲。

2.3 极化恒等式在只有距离数据时的救场

极化恒等式解决的是一个很实际的困境：你拿到的数据不是坐标，而是距离。比如用户之间的相似度表、分子结构的构象距离矩阵、城市之间的路网距离表。这些数据里没有"向量"的概念，但你又想对它们做降维、聚类、可视化。

思路是先把距离转成 Gram 矩阵，再对 Gram 矩阵做特征分解。设距离矩阵D，元素D_ij = ‖x_i - x_j‖²。如果这些点来自某个欧几里得空间，那么可以证明中心化之后的 Gram 矩阵满足：

G = -0.5 · J · D · J

其中J = I - (1/n) · 1 · 1ᵀ。得到的G_ij就等于中心化后的内积⟨x_i - x̄, x_j - x̄⟩。接下来对 G 做特征分解，取前几个最大特征值对应的特征向量乘上特征值的平方根，就是经典 MDS 的坐标输出。

这里有一个坑必须提醒：G如果是数值稳定的距离矩阵算出来的，特征值可能有很小的负数。理论上G是半正定的，实际算出来会有-1e-12这样的值。做特征值开方的时候不能直接np.sqrt(eigvals)，会得到nan，要先np.clip(eigvals, 0, None)。我在做相似度可视化的时候被这个坑绊过两次，第二次才想起来加 clip。

3. 矩阵版本：从 ‖Xw - y‖² 到正规方程的完整链路

3.1 手工展开的每一步

最小二乘的目标函数f(w) = ‖Xw - y‖²，其中X ∈ ℝ^{n×d}、w ∈ ℝ^d、y ∈ ℝ^n。Xw - y是一个 n 维向量，直接对向量做内积展开：

f(w) = (Xw - y)ᵀ(Xw - y)

把两个括号按分配律摊开，得到四项：

= (Xw)ᵀ(Xw) - (Xw)ᵀy - yᵀ(Xw) + yᵀy

第一项把(Xw)ᵀ = wᵀXᵀ代进去，就是wᵀXᵀXw。中间两项需要停一下想清楚：(Xw)ᵀy是一个标量，标量的转置等于它本身，所以(Xw)ᵀy = yᵀ(Xw)，两项完全相等，合起来是-2yᵀXw。最后：

f(w) = wᵀXᵀXw - 2yᵀXw + yᵀy

对照向量版本的展开‖x+y‖² = ‖x‖² + 2⟨x,y⟩ + ‖y‖²，可以发现结构完全一致。把x看作Xw、y看作-y，交叉项就是2⟨Xw, -y⟩ = -2yᵀXw。矩阵情形没有引入任何新东西，只是把内积写成了转置乘法。

维度核对一下，防止抄错公式：XᵀX是d×d，wᵀXᵀXw是1×1标量；Xᵀy是d×1，yᵀXw也是标量，没问题。

3.2 对 w 求导，得到正规方程

对w求梯度。这里用到两条矩阵求导规则，都在标准表格里查得到：

• ∂(wᵀAw)/∂w = (A + Aᵀ)w，A对称时简化为2Aw
• ∂(cᵀw)/∂w = c

令A = XᵀX。它一定是对称的，因为(XᵀX)ᵀ = XᵀX。第一项求导得2XᵀXw。第二项里c = -2Xᵀy，求导得-2Xᵀy。第三项yᵀy与w无关，导数为零。合起来：

∂f/∂w = 2XᵀXw - 2Xᵀy

令它等于零，两边除以 2：

XᵀXw = Xᵀy

这就是正规方程。当X列满秩时，XᵀX可逆，解析解是w = (XᵀX)⁻¹Xᵀy。

这里有第一个真正的工程坑：理论上可以这么写，但实际代码里千万别显式求逆。原因是条件数会在这一步被平方，cond(XᵀX) = cond(X)²。如果X的条件数是 1e6，XᵀX就是 1e12，靠近 float64 有效精度的边界，解出来的w可能完全没有意义。稳妥的做法是走 QR 分解（np.linalg.lstsq内部就是干这个的）或者 SVD，直接对X操作，条件数不平方。我见过有人手写np.linalg.inv(X.T @ X) @ X.T @ y，在特征共线时解出的系数能偏到天上去，换成lstsq立刻就稳了。

3.3 加上正则项之后，可解性从"看运气"变成"一定成立"

岭回归的目标函数是：

f(w) = ‖Xw - y‖² + λ‖w‖²

对两项分别展开。第一项刚才算过，第二项λ‖w‖² = λwᵀw。合起来：

f(w) = wᵀXᵀXw - 2yᵀXw + yᵀy + λwᵀIw = wᵀ(XᵀX + λI)w - 2yᵀXw + yᵀy

求导令零：

2(XᵀX + λI)w - 2Xᵀy = 0 w = (XᵀX + λI)⁻¹ Xᵀy

关键性质：XᵀX是半正定的，特征值都 ≥ 0；加上λI之后（λ > 0），所有特征值都 ≥ λ > 0，矩阵严格正定，必然可逆。这意味着即使X列不满秩（特征数比样本还多，或者特征之间有完全共线），岭回归也有唯一解。这一条是岭回归能在高维场景下稳定工作的根本原因，不是"效果更好"这种模糊说法。

顺带解释一个常见疑问：为什么惩罚项用‖w‖²而不是‖w‖？因为λ‖w‖是 L1 正则（LASSO），它在w = 0处不可导，没有闭式解，得用近端梯度或者坐标下降迭代。而λ‖w‖²展开后是光滑的二次型，一次求导就落到线性方程组上，d×d的方程组一次解完。这是纯粹的代数取舍，跟"哪个正则化效果更好"是两个层面的问题。

4. 核方法里的那次展开：从 ‖φ(x) - φ(y)‖² 到 Gram 矩阵

4.1 核技巧为什么绕不开这个展开

核方法的核心操作是：把数据映射到高维（甚至无穷维）特征空间φ(x)，然后只通过内积K(x,y) = ⟨φ(x), φ(y)⟩来做运算，从不显式构造φ(x)。这省下了巨大的存储和计算，但也带来一个问题——很多算法的目标函数里出现的是距离，不是内积。距离能不能只用核函数值表示？

能，靠的就是差形式的展开：

‖φ(x) - φ(y)‖² = ‖φ(x)‖² - 2⟨φ(x), φ(y)⟩ + ‖φ(y)‖² = K(x,x) - 2K(x,y) + K(y,y)

这一步把"高维特征空间里的距离"完整地翻译成了三个核函数值。如果核是 RBF，K(x,x) = 1恒成立，公式进一步简化成2 - 2K(x,y)。这个简化形式在实现里到处都是。

用得到它的地方不少：核 K-means 的分配步、RBF 核 SVM 的对偶求解、谱聚类里的相似度矩阵、MMD（最大均值差异）的计算。以 MMD 为例，两个分布的核均值嵌入分别是μ_P和μ_Q，MMD 的平方就是‖μ_P - μ_Q‖²，展开后变成E[k(x,x')] - 2E[k(x,y)] + E[k(y,y')]，三个期望都可以用样本均值估计，完全绕开了无穷维嵌入的构造。生成模型评估里常用的那个 MMD 指标，代码核心就是这一行展开。

4.2 RBF 核的定义本身就是一次展开

RBF 核写成K(x,y) = exp(-γ‖x-y‖²)，把差形式代进去：

K(x,y) = exp(-γ(‖x‖² - 2⟨x,y⟩ + ‖y‖²)) = exp(-γ‖x‖²) · exp(2γ⟨x,y⟩) · exp(-γ‖y‖²)

这个分解形式说明了一件很有意思的事：RBF 核可以理解成"先把每个点按exp(-γ‖x‖²)缩放到单位球面上，再做指数变换的内积"。把‖x‖归一化之后，exp(-γ‖x‖²)对每个点都变成相同的常数，两个点的核值就只取决于它们的夹角余弦，也就是exp(2γ⟨x,y⟩)。

这个视角解释了两个常被提到的性质。第一，RBF 核下所有样本点都落在特征空间的同一个球面上（K(x,x) = 1），所以K(x,x) - 2K(x,y) + K(y,y) = 2 - 2K(x,y)，距离完全由核值决定。第二，γ控制的是"这个球面上的距离尺度"，γ越大，不同方向上的点被推得越开，模型越容易过拟合；γ越小，所有点挤在一起，模型欠拟合。调γ的时候心里有这张图，比盲搜参数靠谱得多。

4.3 距离矩阵的 GEMM 加速写法

把差形式用在整批数据上，就得到了距离矩阵的一次性算法。设X ∈ ℝ^{n×d}，第 i 行是x_i，则：

D_ij = ‖x_i - x_j‖² = ‖x_i‖² + ‖x_j‖² - 2 x_iᵀx_j

令a = [‖x_1‖², ..., ‖x_n‖²]ᵀ，用外积和矩阵乘法写出来就是：

D = a · 1ᵀ + 1 · aᵀ - 2 X Xᵀ

三行代码：

import numpy as np def pairwise_sqdist(X): # X: (n, d)，按行存放样本 sq = np.einsum('ij,ij->i', X, X) # 每行的范数平方，避免生成 n×d 临时数组 D = sq[:, None] + sq[None, :] - 2.0 * (X @ X.T) np.maximum(D, 0.0, out=D) # 兜底，原因见第 5 节 return D

为什么必须这么写？逐对循环的复杂度是O(n²d)，GEMM 版本也是O(n²d)，复杂度没变，但常数差得非常远。X @ X.T走的是 BLAS 的三级运算，能吃到多线程、分块、寄存器复用、SIMD 全套优化；Python 层的双重循环每对都要重新读一遍内存、走一次解释器开销。n = 10000、d = 128 的时候，循环版本要跑几分钟，GEMM 版本几秒就结束了。差的那个数量级全在工程实现上。

np.einsum('ij,ij->i', X, X)这个写法比(X**2).sum(axis=1)略快一点，原因是它不需要先构造一个完整的n×d平方数组再求和，省掉一次内存来回。数据量大的时候这点差别会变得可观。

内存不够的话，把X按行分块，一块一块地和全量算X_block @ X.T，结果存进D的对应行段。分块大小控制在让X_block @ X.T这个中间结果不超过几百 MB 就行，太大反而会因为换页变慢。

5. 数值上的三个坑：溢出、灾难性抵消、负距离

5.1 ‖x‖² 自己就可能溢出或下溢

算范数平方最直接的写法是np.sum(x * x)，但这个写法在数据量级极端的时候会出事。假设x里最大的元素是 1e200，平方就是 1e400，而 float64 的最大值大约是 1.8e308，直接变成inf。反过来，如果元素量级是 1e-200，平方是 1e-400，低于 float64 的最小正规数，直接变成 0，整个范数失真。

BLAS 的dnrm2函数用的是缩放法：先把所有元素除以最大绝对值，算完平方和再乘回来。这样中间结果始终在[0, 1]区间内，彻底避免溢出和下溢。

def safe_norm_sq(x): """返回 ‖x‖²，量级极端时也不会溢出""" scale = np.max(np.abs(x)) if scale == 0.0: return 0.0 return (scale ** 2) * np.sum((x / scale) ** 2)

这里有个容易忽略的差别：np.linalg.norm(x)对内存连续的数组会调用 BLAS 的nrm2，是安全的；而np.sqrt(x @ x)或者np.sum(x*x)是不安全的。深度学习里做 embedding 归一化的时候，很多人的写法是x / np.linalg.norm(x, axis=1, keepdims=True)，这个没问题；但如果改成x / np.sqrt((x**2).sum(axis=1, keepdims=True))，在某些量级极端的样本上就会得到inf或者nan。这个坑我在处理用户行为序列 embedding 的时候踩过一次，单条样本的某个维度数值特别大，导致整行归一化全废。

另外提一句，np.linalg.norm(x)**2是先开方再平方，中间多了一次舍入，精度比直接算平方和略差一点点（大概是最后一位的量级差），但安全性足够，日常用没问题。真正要抠精度的时候再换成np.dot(x, x)，前提是你确认过量级不极端。

5.2 灾难性抵消：负距离是怎么冒出来的

回到开头那个负数的距离。展开式D_ij = ‖x_i‖² + ‖x_j‖² - 2x_iᵀx_j在数学上恒等于‖x_i - x_j‖²，但在浮点运算里，当x_i和x_j几乎重合时，等式右边三项都接近同一个数，结果接近零。浮点减法中两个相近的数相减，前面的有效位会被完全抵消掉。

具体估一下误差。假设‖x_i‖ ≈ ‖x_j‖ ≈ R，两点之间的实际距离是 δ，那么三项的量级都是 R²。浮点乘法和加法的相对误差大约是机器精度eps ≈ 2.2e-16，绝对误差约eps · R²。折算到结果上的相对误差是：

相对误差 ≈ eps · R² / δ²

代入一组数字感受一下。如果δ/R = 1e-8，也就是两点距离只有自身长度的亿分之一，那么R²/δ² = 1e16，相对误差约2.2e-16 × 1e16 = 2.2。相对误差超过 100%，说明算出来的结果连符号都没保证，出现-3.7e-9完全在意料之中。

这个现象在哪些场景会真实发生？去重任务里，同一个物品的两个近似 embedding 之间的距离极小；K-means 迭代后期，样本离它归属的簇心可能已经非常近；时间序列里相邻两个采样点相差无几。这些都是高发区。

还有一类更隐蔽的情况：float32下机器精度是1.2e-7，比float64大了 9 个数量级。也就是说δ/R = 1e-3的时候，float32算出的距离平方就已经失真了。而δ/R = 1e-3在归一化之后的向量里是很常见的数值，很多人用float32存 embedding 同时用展开式算距离，结果在一开始就偏了，只是偏差小到没有引起注意。

5.3 补救措施与阈值判定

抑制灾难性抵消最有效的一招是中心化。把X的每一列减去该列的均值，让每个样本向量的范数降到和典型距离同一量级。这样R²/δ²会大幅下降，误差随之减小。中心化不改变任何成对距离（所有点整体平移，距离不变），所以对 K-means、kNN、聚类这些任务毫无副作用。

其他几招按代价从低到高排列：

• 升精度：float32换成float64，eps直接降 9 个数量级，代价是内存翻倍、GEMM 速度减半左右。
• 缩放：除以np.abs(X).max()，把量级压到 1 附近。
• 候选重算：先用 GEMM 快手算一遍全量距离，筛出距离小于某个阈值的候选对，再用x_i - x_j直接算这些对的距离平方。候选对通常只占少数，直接算法没有抵消问题。
• 分块混合：大块之间走 GEMM，块内元素少的时候走直接法。

不管用哪种，都要加一道检查：

def audit_and_clip(D, tol=1e-9): """检查距离矩阵是否出现明显负数，报警后再裁剪""" scale = max(1.0, float(np.max(np.abs(D)))) neg_mask = D < -tol * scale if neg_mask.any(): print(f"警告：{int(neg_mask.sum())} 个负距离平方，" f"最小值 {D.min():.3e}，矩阵量级 {scale:.3e}") return np.maximum(D, 0.0)

np.maximum(D, 0.0)只是把负数拍成 0，属于遮羞布，会掩盖真实的精度丢失，所以一定要先报警再裁剪。如果不报警直接裁剪，后面聚类结果出现诡异行为的时候会完全查不到原因。

一个经验性判据：如果max|D| / min(|D_ij|, D_ij > 0)超过 1e10，这个距离矩阵的精度就值得怀疑了，至少要抽查几个最小值对应的向量对，用直接法验证一遍。

6. 什么时候不能展开：p 范数、复数域与加权范数

6.1 p ≠ 2 时，展开公式直接失效

p 范数的定义是‖x‖_p = (Σ|x_i|^p)^{1/p}。平方之后是(Σ|x_i|^p)^{2/p}，这个表达式不是二次型，根本没有"交叉项等于 2 倍内积"的结构。第 2.2 节已经用x = (1,0)、y = (0,1)在 l1 和 l∞ 上验过平行四边形恒等式不成立，所以这两个范数不由内积诱导，展开公式也就无从谈起。

这个事实带来一连串的后果，值得单独列出来：

• 距离矩阵不能用‖x‖² + ‖y‖² - 2⟨x,y⟩加速，l1 距离只能老老实实逐对算，没有 GEMM 版本。
• 核技巧失效。核方法依赖内积结构，l1 距离下没有对应的特征映射和核函数。
• 最小二乘没有闭式解。‖Xw - y‖₁的目标函数不可微且非二次，只能用线性规划或者迭代重加权最小二乘，慢很多。
• 极化恒等式没有对应形式，无法从 l1 距离恢复任何内积矩阵。

所以工程上默认用 2 范数，并不是因为"欧几里得距离在语义上最正确"，而是因为它带来的代数便利性太高了。如果你确实需要 L1 距离（高维稀疏文本、需要抗离群点的场景），就得接受这些代价。我不会说 L1 不值得用，只是要在选型的时候清楚自己放弃了什么。

6.2 复数域：交叉项必须取实部

复数内积的定义和实空间不一样。常见约定是⟨x,y⟩ = Σ x̄_i y_i（第一个位置共轭），它满足共轭对称⟨x,y⟩ = conj(⟨y,x⟩)，而不是实空间的对称。把和形式展开：

‖x+y‖² = ‖x‖² + ⟨x,y⟩ + ⟨y,x⟩ + ‖y‖² = ‖x‖² + ⟨x,y⟩ + conj(⟨x,y⟩) + ‖y‖² = ‖x‖² + 2Re⟨x,y⟩ + ‖y‖²

最后一步用到了z + conj(z) = 2Re(z)。

所以复数域的交叉项是2Re⟨x,y⟩，是一个实数，不是复数。这一步有太多人写错：直接写2⟨x,y⟩得到复数结果，然后往np.maximum里面塞，numpy 会直接抛类型错误，或者在np.abs之后静默地给出错误答案。处理频谱特征、复数神经网络、量子模拟这类数据时一定要留意。

极化恒等式在复数域也复杂得多，需要四项：

⟨x,y⟩ = (1/4)(‖x+y‖² - ‖x-y‖² + i‖x+iy‖² - i‖x-iy‖²)

你会发现多了两个带虚数单位系数的项。这不是公式凑出来的，而是因为复数内积有两个实自由度（实部和虚部），单靠‖x+y‖²和‖x-y‖²两个实数只能确定一个。信号处理里那些看起来"凭空多出来"的项，基本都来自这里。

6.3 加权范数：矩阵对称正定不是可选项

加权范数写作‖x‖_W² = xᵀWx，展开是：

‖x+y‖²_W = ‖x‖²_W + 2xᵀWy + ‖y‖²_W

结构上和标准形式一样，只是内积换成了⟨x,y⟩_W = xᵀWy。

这里的关键要求是W必须对称正定。对称保证了交叉项的两个位置相等（xᵀWy = yᵀWx）；正定保证了xᵀWx > 0对所有非零x成立，也就是真的构成范数。

如果W只是半正定（有零特征值），得到的是伪范数：存在非零向量x使得xᵀWx = 0，三角不等式也可能失效。这在实践中非常常见——样本量小于特征数时，样本协方差矩阵Σ一定奇异，秩最多是 n-1。

Mahalanobis 距离就是这个形式：

d(x,y)² = (x-y)ᵀ Σ⁻¹ (x-y)

展开之后需要Σ⁻¹而不是Σ，而Σ必须可逆。工程上的处理办法有两条：加正则Σ + λI（这就是 Ledoit-Wolf 收缩估计的直觉来源），或者用伪逆Σ⁺（对零特征值方向不做惩罚）。选哪条取决于业务背景——如果那些零特征值方向本身就是噪声，加正则更合适；如果它们是稀疏但真实的信号，伪逆更保守。

用这个公式之前，我还习惯做两件事：一是把W对称化（W = (W + Wᵀ) / 2），防止矩阵因为数值误差产生极小的非对称分量；二是用np.linalg.cholesky试一下能否分解成功，能分解就说明是正定的，否则会抛LinAlgError，比用特征值判断省事得多。

最后补一个我在实际代码里养成的习惯：只要项目里出现"距离"两个字，先想清楚三件事——用的是不是 2 范数、向量有没有中心化、精度够不够。这三件事的答案决定了后面所有优化能不能上。很多看起来像是"算法不收敛"的问题，追到底都是其中一个没答对。