"""Tensor + 反向模式自动求导(只依赖 numpy)。 每个运算产出一个新 Tensor,并记录: _parents : 参与运算的输入 Tensor _backward : 函数 g -> (对每个 parent 的梯度, ...) backward() 先拓扑排序,再从输出往回逐个调用 _backward,把梯度累加到 .grad。 """ import numpy as np _grad_enabled = True _default_dtype = np.dtype(np.float32) def set_default_dtype(dtype): """所有浮点数据进入框架时统一转成这个类型(默认 float32)。 只有测试里的数值梯度检查需要切到 float64。整数数组(下标)不受影响。""" global _default_dtype _default_dtype = np.dtype(dtype) def get_default_dtype(): return _default_dtype _flops = None # count_flops() 打开时是一个 [计数],否则 None class count_flops: """with count_flops() as c: ... 之后 c.flops = 期间所有矩阵乘的浮点运算次数(乘加算 2 次)。 只统计矩阵乘(模型里绝大部分计算量),不统计逐元素运算;不受机器负载影响,比计时可靠。""" def __enter__(self): global _flops self._prev, _flops = _flops, [0] return self def __exit__(self, *exc): global _flops self.flops = _flops[0] _flops = self._prev class no_grad: """with no_grad(): 内部的运算不建计算图(推理/评估时用,省内存)。""" def __enter__(self): global _grad_enabled self._prev, _grad_enabled = _grad_enabled, False def __exit__(self, *exc): global _grad_enabled _grad_enabled = self._prev def _unbroadcast(grad, shape): """把广播后的梯度求和还原成原始形状。""" while grad.ndim > len(shape): grad = grad.sum(axis=0) for i, s in enumerate(shape): if s == 1 and grad.shape[i] != 1: grad = grad.sum(axis=i, keepdims=True) return grad def _to_array(x): if isinstance(x, Tensor): return x.data arr = np.asarray(x) # 浮点一律转成默认类型,避免 float32 和 float64 混在一起被 numpy 悄悄升成 float64 if arr.dtype.kind == "f" and arr.dtype != _default_dtype: arr = arr.astype(_default_dtype) return arr class Tensor: __array_priority__ = 100 # 让 ndarray + Tensor 走 Tensor.__radd__ def __init__(self, data, requires_grad=False): self.data = _to_array(data) self.requires_grad = requires_grad self.grad = None self._parents = () self._backward = None # ---------- 基础属性 ---------- @property def shape(self): return self.data.shape @property def ndim(self): return self.data.ndim @property def dtype(self): return self.data.dtype def __len__(self): return len(self.data) def __repr__(self): return f"Tensor({self.data}, requires_grad={self.requires_grad})" def numpy(self): return self.data def item(self): return self.data.item() def detach(self): return Tensor(self.data) # ---------- 建图 ---------- @staticmethod def _op(data, parents, backward): """创建运算结果。backward(g) 返回与 parents 一一对应的梯度(可为 None)。""" out = Tensor(data) if _grad_enabled and any(p.requires_grad for p in parents): out.requires_grad = True out._parents = parents out._backward = backward return out def backward(self, grad=None): if grad is None: if self.data.size != 1: raise RuntimeError("非标量输出调用 backward() 需要显式传入 grad") grad = np.ones_like(self.data) # 迭代式拓扑排序(避免深图递归爆栈) topo, visited, stack = [], set(), [(self, False)] while stack: node, expanded = stack.pop() if expanded: topo.append(node) continue if id(node) in visited: continue visited.add(id(node)) stack.append((node, True)) for p in node._parents: if p.requires_grad and id(p) not in visited: stack.append((p, False)) self.grad = grad if self.grad is None else self.grad + grad for node in reversed(topo): if node._backward is None or node.grad is None: continue for p, g in zip(node._parents, node._backward(node.grad)): if g is not None and p.requires_grad: p.grad = g if p.grad is None else p.grad + g # 释放计算图,中间节点的梯度也不再需要 node._parents, node._backward, node.grad = (), None, None # ---------- 逐元素运算 ---------- def __add__(self, other): other = other if isinstance(other, Tensor) else Tensor(other) a, b = self, other return Tensor._op(a.data + b.data, (a, b), lambda g: (_unbroadcast(g, a.shape), _unbroadcast(g, b.shape))) def __mul__(self, other): other = other if isinstance(other, Tensor) else Tensor(other) a, b = self, other return Tensor._op(a.data * b.data, (a, b), lambda g: (_unbroadcast(g * b.data, a.shape), _unbroadcast(g * a.data, b.shape))) def __truediv__(self, other): other = other if isinstance(other, Tensor) else Tensor(other) a, b = self, other return Tensor._op(a.data / b.data, (a, b), lambda g: (_unbroadcast(g / b.data, a.shape), _unbroadcast(-g * a.data / (b.data ** 2), b.shape))) def __pow__(self, p): assert isinstance(p, (int, float)), "只支持标量指数" a = self return Tensor._op(a.data ** p, (a,), lambda g: (g * p * a.data ** (p - 1),)) def __neg__(self): return Tensor._op(-self.data, (self,), lambda g: (-g,)) def __sub__(self, other): return self + (-other if isinstance(other, Tensor) else Tensor(-_to_array(other))) def __radd__(self, other): return self + other def __rmul__(self, other): return self * other def __rsub__(self, other): return (-self) + other def __rtruediv__(self, other): return Tensor(other) / self def exp(self): out = np.exp(self.data) return Tensor._op(out, (self,), lambda g: (g * out,)) def log(self): a = self return Tensor._op(np.log(a.data), (a,), lambda g: (g / a.data,)) def sqrt(self): return self ** 0.5 def tanh(self): out = np.tanh(self.data) return Tensor._op(out, (self,), lambda g: (g * (1 - out ** 2),)) def relu(self): mask = self.data > 0 return Tensor._op(self.data * mask, (self,), lambda g: (g * mask,)) def masked_fill(self, mask, value): mask = _to_array(mask).astype(bool) return Tensor._op(np.where(mask, value, self.data).astype(self.dtype), (self,), lambda g: (np.where(mask, 0, g),)) # ---------- 归约 ---------- def sum(self, axis=None, keepdims=False): a = self def backward(g): if axis is not None and not keepdims: g = np.expand_dims(g, axis) return (np.broadcast_to(g, a.shape).copy(),) return Tensor._op(a.data.sum(axis=axis, keepdims=keepdims), (a,), backward) def mean(self, axis=None, keepdims=False): n = self.data.size if axis is None else np.prod([self.shape[i] for i in np.atleast_1d(axis)]) return self.sum(axis, keepdims) * (1.0 / n) def max(self, axis=None, keepdims=False): a = self out = a.data.max(axis=axis, keepdims=True) def backward(g): if axis is not None and not keepdims: g = np.expand_dims(g, axis) mask = a.data == out return (mask * g / mask.sum(axis=axis, keepdims=True),) # 并列最大值平分梯度 return Tensor._op(out if keepdims else out.squeeze(axis), (a,), backward) # ---------- 形状 ---------- def reshape(self, *shape): a = self return Tensor._op(a.data.reshape(*shape), (a,), lambda g: (g.reshape(a.shape),)) def transpose(self, *axes): axes = axes or tuple(reversed(range(self.ndim))) inv = np.argsort(axes) return Tensor._op(self.data.transpose(axes), (self,), lambda g: (g.transpose(inv),)) @property def T(self): return self.transpose() def __getitem__(self, idx): a = self if isinstance(idx, Tensor): idx = idx.data parts = idx if isinstance(idx, tuple) else (idx,) fancy = any(isinstance(i, (np.ndarray, list)) for i in parts) def backward(g): grad = np.zeros_like(a.data) if fancy: np.add.at(grad, idx, g) # 数组下标可能重复(Embedding 查表),必须累加 else: grad[idx] = g # 普通切片不会重复,直接写,比 add.at 快得多 return (grad,) return Tensor._op(a.data[idx], (a,), backward) # ---------- 矩阵乘 ---------- def __matmul__(self, other): other = other if isinstance(other, Tensor) else Tensor(other) a, b = self, other assert a.ndim >= 2 and b.ndim >= 2, "matmul 需要至少 2 维" out = a.data @ b.data if _flops is not None: _flops[0] += 2 * out.size * a.shape[-1] def backward(g): ga = g @ np.swapaxes(b.data, -1, -2) gb = np.swapaxes(a.data, -1, -2) @ g return _unbroadcast(ga, a.shape), _unbroadcast(gb, b.shape) return Tensor._op(out, (a, b), backward) def concat(tensors, axis=0): """沿 axis 拼接。反向时把梯度按原来的大小切回去。""" tensors = tuple(t if isinstance(t, Tensor) else Tensor(t) for t in tensors) splits = np.cumsum([t.shape[axis] for t in tensors])[:-1] return Tensor._op(np.concatenate([t.data for t in tensors], axis=axis), tensors, lambda g: tuple(np.split(g, splits, axis=axis))) # ---------- 构造函数 ---------- def tensor(data, requires_grad=False): return Tensor(data, requires_grad) def zeros(*shape, requires_grad=False): return Tensor(np.zeros(shape, dtype=_default_dtype), requires_grad) def ones(*shape, requires_grad=False): return Tensor(np.ones(shape, dtype=_default_dtype), requires_grad) def randn(*shape, requires_grad=False): return Tensor(np.random.randn(*shape).astype(_default_dtype), requires_grad)