完全指南:从逐元素运算到广播、类型转换与错误处理)
NumPy Universal Functionsufunc完全指南从逐元素运算到广播、类型转换与错误处理【免费下载链接】numpyThe fundamental package for scientific computing with Python.项目地址: https://gitcode.com/gh_mirrors/nu/numpyUniversal functions简称 ufunc是 NumPy 高性能数值计算的核心抽象。本文以官方文档 doc/source/user/basics.ufuncs.rst 为主线结合仓库中的 ufunc 参考手册 与numpy/_core下的源码实现系统讲解 ufunc 的定义、五大实例方法reduce、accumulate、reduceat、outer、at、输出类型确定规则、广播机制、类型转换casting规则、内部缓冲区以及浮点错误处理。读完本文你将掌握 ufunc 的完整行为模型并能利用frompyfunc创建自定义 ufunc、用seterr/errstate精确控制异常处理、用can_cast预判类型提升结果。什么是 ufunc逐元素运算的向量化封装Universal functionufunc是作用于ndarray的逐元素element-by-element函数支持数组广播、类型转换等标准特性。本质上ufunc 是一个针对固定数量输入、固定数量输出函数的向量化vectorized包装器。在 NumPy 中ufunc 是numpy.ufunc类的实例大量内建函数由编译后的 C 代码实现。基本 ufunc 操作标量元素而广义 ufuncgeneralized ufunc简称 gufunc的基本元素是子数组向量、矩阵等并在其他维度上进行广播。最典型的对比是逐元素的np.add与作用于向量/矩阵的np.matmul import numpy as np a np.arange(6).reshape(3, 2) a array([[0, 1], [2, 3], [4, 5]]) np.add(a, a) # 逐元素加法 array([[ 0, 2], [ 4, 6], [ 8, 10]]) np.matmul(a, a.T) # 矩阵乘法 (3x2) (2x3) - (3x3) array([[ 1, 3, 5], [ 3, 13, 23], [ 5, 23, 41]])最简单的 ufunc 用法是算术运算符运算符底层即调用对应的加法 ufunc np.array([0,2,3,4]) np.array([1,1,-1,2]) array([1, 3, 2, 6])对于标量 ufunc输入输出都是标量对于广义 ufunc输入输出则是子数组。想要深入了解 gufunc 的签名signature机制可参阅 doc/source/reference/ufuncs.rst 中的signature/axes/axis参数说明。用 frompyfunc 创建自定义 ufunc除了内建 ufunc你还可以通过工厂函数numpy.frompyfunc将任意 Python 函数包装成 ufunc 实例。该函数在底层由 C 实现入口定义于 numpy/_core/src/umath/umathmodule.c 的ufunc_frompyfunc并注册在 multiarray 模块中。其签名见 numpy/_core/multiarray.pyi支持指定输入个数、输出个数以及可选的恒等元素identitydef frompyfunc(func, nin, nout, *, identityNone)例如把一个接受两个 Python 标量的函数包装成 ufunc def my_add(a, b): ... return a b uf np.frompyfunc(my_add, 2, 1) uf([1, 2, 3], [10, 20, 30]) array([11, 22, 33], dtypeobject)注意frompyfunc生成的对象 ufunc 输出 dtype 为object如需更高效的数值型自定义 ufunc可参考 c-info.ufunc-tutorial 编写 C 扩展。Ufunc 的五大方法reduce 系列与原地操作所有 ufunc 都有5 个方法4 个 reduce 类方法reduce、accumulate、reduceat、outer和 1 个原地操作方法at。这些方法只对接收两个输入、返回一个输出的标量 ufunc 有意义其内层循环作用于单个标量值。对不满足条件的 ufunc 调用这些方法会抛出ValueError对输出多于一个的 ufunc 调用reduce会抛出TypeError除非该 ufunc 的循环实现注册了专用的 reduction 循环此时reduce也可用但返回元组。以np.add为例其方法正常工作 np.add.reduce([1, 2, 3]) 6而np.divmod返回两个输出商与余数且未注册 reduction 循环调用其方法即报错 np.divmod.reduce([1, 2, 3]) Traceback (most recent call last): ... TypeError: divmod.reduce is not supported: the resolved loop does not register a reduction loop若某个多输出 ufunc 的循环确实注册了 reduction 循环则其reduce返回一个每个输出一个数组的元组initial参数既可以传单个值广播到每个输出也可以传一个每个输出一个值的元组。如何为自定义 ufunc 添加 reduction 循环参见 c-info.reduction-loop-tutorial。reduce 类的通用关键字axis、dtype、out所有 reduce 类方法都接受axis、dtype、out关键字且输入数组维度必须 ≥ 1。axis指定归约沿哪个轴进行负值从后往前计数。对reduce而言它还可以是int元组同时归约多个轴或None归约所有轴 x np.arange(9).reshape(3,3) x array([[0, 1, 2], [3, 4, 5], [6, 7, 8]]) np.add.reduce(x, 1) # 沿轴 1 归约 array([ 3, 12, 21]) np.add.reduce(x, (0, 1)) # 同时归约两个轴 36dtype解决结果放不进原数组 dtype这一常见问题。例如对单字节整数数组求和结果可能溢出int8。dtype允许你指定归约运算所用的数据类型即输出类型从而保证精度足够 x.dtype dtype(int64) np.multiply.reduce(x, dtypenp.float64) array([ 0., 28., 80.])有一个自动提升例外若对add/multiply做归约而未指定dtype且输入为比numpy.int_更小的整数或布尔类型则会内部提升到int_或uint。其余情况调整归约类型的责任基本在你身上。out提供输出数组多输出 ufunc 提供输出数组元组。若给定了outdtype只影响内部计算最终结果会写入out并按其 dtype 存储 y np.zeros(3, dtypenp.int_) y array([0, 0, 0]) np.multiply.reduce(x, dtypenp.float64, outy) array([ 0, 28, 80])注在 doc/source/reference/ufuncs.rst 的Optional keyword arguments一节还说明了out的更多细节out可以是元组每输出一项允许None表示由 ufunc 分配单输出 ufunc 也可直接传单个数组默认outNone时创建未初始化数组若结果是零维则转换为标量传out...即outEllipsis可避免该转换。ufunc.at基于高级索引的原地操作第五个方法numpy.ufunc.at允许使用高级索引执行原地操作。在用到高级索引的维度上不使用内部缓冲buffering因此高级索引可以多次列出同一元素操作将基于该元素上一次操作的结果继续执行。例如 a np.zeros(5) np.add.at(a, [1, 1, 1], 1) # 对索引 1 累加 3 次 a array([0., 3., 0., 0., 0.])仓库中 numpy/_core/tests/test_mem_overlap.py 的test_ufunc_at_manual对ufunc.at在各种索引下的行为含副本与重叠做了系统性验证test_custom_dtypes.py 则验证了自定义 dtype 上multiply.at的语义。输出类型确定ndarray、array scalar 与子类包装如果 ufunc或其方法的输入参数是ndarray则输出也是ndarray。唯一例外是结果为零维时会被转换为array scalar传入out...或outEllipsis可以避免这种转换。如果部分或全部输入不是ndarray则输出也不一定是ndarray若任一输入定义了__array_ufunc__方法控制权将完全交给该方法即发生ufunc 覆盖override详见 doc/source/reference/arrays.classes.rst。若没有任何输入覆盖 ufunc则所有输出数组会交给那些定义了__array_wrap__方法、且在所有输入除ndarray与标量外中__array_priority__最高的输入对象处理。ndarray的默认__array_priority__为 0.0子类型默认为 0.0而matrix为 10.0。仓库中的实际例子numpy/_core/memmap.py的memmap将__array_priority__设为 -100.0并自定义__array_wrap__numpy/_core/defchararray.py的chararray也重写了__array_wrap__。这三个协议的完整文档位于 numpy/_core/_add_newdocs.py 的ndarray条目__array_ufunc__、__array_wrap__、__array_priority__。所有 ufunc 都可以接收输出参数输出必须是数组或其子类必要时结果会被转换为所提供输出数组的 dtype。如果输出对象自身定义了__array_wrap__则优先调用它而不是输入上找到的那个。相关测试可参考 numpy/_core/tests/test_multiarray.py 中的test_ufunc_binop_bad_array_priority等用例它们验证了不同__array_priority__下二元运算的返回类型行为。广播Broadcasting维度为 1 时步长为 0每个 ufunc 都通过对输入执行核心函数的逐元素运算来产生数组输出元素通常是标量对 gufunc 则可能是向量或更高阶子数组。标准广播规则被应用使得形状不完全相同的输入仍可参与运算。按照这些规则如果输入形状的某个维度大小为 1则沿该维度的所有计算都使用该维度的第一个数据元素。换句话说ufunc 的步进机制stepping machinery在该维度上不步进——即该维度的stride 为 0。这种stride 置 0的实现方式使得广播不需要复制数据而是直接复用同一内存位置这正是广播高效的原因。详细理论可参见 theory.broadcasting。类型转换Casting规则内层循环与 can_cast 表内层循环与 .types 属性每个 ufunc 的核心是一个一维 strided 循环inner loop为特定类型组合实现实际函数。创建 ufunc 时它携带一个静态的内层循环列表及对应的类型签名列表。ufunc 机制根据输入类型从该列表中选择合适的内层循环。你可以通过 ufunc 的.types属性查看哪些类型组合定义了内层循环及其输出类型输出类型使用字符代码简写如dd-d表示双精度输入输出双精度。安全转换与选择算法当 ufunc 没有针对所给输入类型的内层循环实现时就必须对部分或全部输入做类型转换。算法会搜索一个所有输入都能**安全地safely**转换过去的类型签名选中内部列表里第一个匹配项在完成所有必要转换后执行。注意ufunc 过程中产生的内部拷贝即使为了转换被限制在内部缓冲区大小之内缓冲区大小可由用户设置。NumPy 的 ufunc 支持混合类型签名——同一个 ufunc 可以同时支持浮点与整数输入例如numpy.ldexp。用 can_cast 查看安全转换表上述转换规则本质上等价于何时一种 dtype 可以安全地转换为另一种 dtype。该问题可用numpy.can_cast(fromtype, totype)在 Python 中直接判断。can_cast的实现位于 numpy/_core/multiarray.py默认castingsafe与之配套的类型提升 API 还有numpy.result_type、numpy.promote_types与numpy.min_scalar_type自 NumPy 1.6.0 起用于封装输出类型确定机制同样定义在该文件中。以下代码输出 64 位系统上可安全转换表结果依赖平台32 位系统上整数类型尺寸不同表格会略有差异 mark {False: -, True: Y} def print_table(ntypes): ... print(X .join(ntypes)) ... for row in ntypes: ... print(row, end) ... for col in ntypes: ... print(mark[np.can_cast(row, col)], end) ... print() ... print_table(np.typecodes[All]) X ? b h i l q n p B H I L Q N P e f d g F D G S U V O M m ? Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y - Y b - Y Y Y Y Y Y Y - - - - - - - Y Y Y Y Y Y Y Y Y Y Y - Y h - - Y Y Y Y Y Y - - - - - - - - Y Y Y Y Y Y Y Y Y Y - Y i - - - Y Y Y Y Y - - - - - - - - - Y Y - Y Y Y Y Y Y - Y l - - - - Y Y Y Y - - - - - - - - - Y Y - Y Y Y Y Y Y - Y q - - - - Y Y Y Y - - - - - - - - - Y Y - Y Y Y Y Y Y - Y n - - - - Y Y Y Y - - - - - - - - - Y Y - Y Y Y Y Y Y - Y p - - - - Y Y Y Y - - - - - - - - - Y Y - Y Y Y Y Y Y - Y B - - Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y - Y H - - - Y Y Y Y Y - Y Y Y Y Y Y - Y Y Y Y Y Y Y Y Y Y - Y I - - - - Y Y Y Y - - Y Y Y Y Y - - Y Y - Y Y Y Y Y Y - Y L - - - - - - - - - - - Y Y Y Y - - Y Y - Y Y Y Y Y Y - - Q - - - - - - - - - - - Y Y Y Y - - Y Y - Y Y Y Y Y Y - - N - - - - - - - - - - - Y Y Y Y - - Y Y - Y Y Y Y Y Y - - P - - - - - - - - - - - Y Y Y Y - - Y Y - Y Y Y Y Y Y - - e - - - - - - - - - - - - - - - Y Y Y Y Y Y Y Y Y Y Y - - f - - - - - - - - - - - - - - - - Y Y Y Y Y Y Y Y Y Y - - d - - - - - - - - - - - - - - - - - Y Y - Y Y Y Y Y Y - - g - - - - - - - - - - - - - - - - - - Y - - Y Y Y Y Y - - F - - - - - - - - - - - - - - - - - - - Y Y Y Y Y Y Y - - D - - - - - - - - - - - - - - - - - - - - Y Y Y Y Y Y - - G - - - - - - - - - - - - - - - - - - - - - Y Y Y Y Y - - S - - - - - - - - - - - - - - - - - - - - - - Y Y Y Y - - U - - - - - - - - - - - - - - - - - - - - - - - Y Y Y - - V - - - - - - - - - - - - - - - - - - - - - - - - Y Y - - O - - - - - - - - - - - - - - - - - - - - - - - - - Y - - M - - - - - - - - - - - - - - - - - - - - - - - - Y Y Y - m - - - - - - - - - - - - - - - - - - - - - - - - Y Y - Y注意两点表格中的S、U、V字节串/Unicode/void 等虽然列在表中但不能被 ufunc 直接操作。标量-数组混合的特殊规则混合标量-数组运算使用另一套转换规则只有当标量属于与数组根本不同种类的数据即处于 dtype 层级结构的不同分支时标量才会提升数组。这条规则让你可以在代码中放心使用标量常量它们在 ufunc 中按 Python 类型解释而不必担心标量常量的精度会强制你的大数组低精度被上转。内部缓冲区setbufsize 与逐线程配置ufunc 内部使用缓冲区来处理三类数据未对齐misaligned数据、字节交换swapped数据、以及需要从一种 dtype 转换为另一种 dtype 的数据。内部缓冲区的大小可按线程per-thread设置。最多会创建2 * (n_inputs n_outputs)个指定大小的缓冲区来处理所有输入与输出。缓冲区默认大小为 10,000 个元素。只要所有输入数组小于缓冲区大小这些表现异常或类型错误的数组会在计算前被整体拷贝。调整缓冲区大小可能显著改变各类 ufunc 计算的完成速度。设置缓冲区的简单接口是numpy.setbufsize(size)对应读取接口为numpy.getbufsize()。其实现位于 numpy/_core/_ufunc_config.pydef setbufsize(size): if size 0: raise ValueError(buffer size must be non-negative) old _get_extobj_dict()[bufsize] extobj _make_extobj(bufsizesize) _extobj_contextvar.set(extobj) return old从源码可以看到缓冲区大小存储在线程局部_extobj_contextvar一个 contextvar的扩展对象字典中setbufsize返回旧值便于事后恢复。自 NumPy 2.0 起缓冲区大小的作用域与numpy.errstate上下文绑定退出with errstate():块时缓冲区大小也会一并恢复numpy/_core/_ufunc_config.py 的setbufsizedocstring 中有明确示例。错误处理浮点状态寄存器与 seterr 系列ufunc 会触发硬件中的特殊浮点状态寄存器例如除零。只要平台支持这些寄存器在计算过程中会被定期检查。错误处理同样是按线程控制的可通过numpy.seterr和numpy.seterrcall配置。numpy.seterr支持五类浮点异常与六种处理方式详见 numpy/_core/_ufunc_config.py 的实现与文档异常类型触发场景divide除零有限数相除得无穷over溢出结果大到无法表示under下溢结果太接近 0 而损失精度invalid非法操作结果不可表示通常产生 NaN处理方式行为ignore不采取任何行动warn通过 Pythonwarnings模块发出RuntimeWarning默认raise抛出FloatingPointErrorcall调用seterrcall指定的回调函数print直接向stdout打印警告log将错误记录到seterrcall指定的 Log 对象示例 orig_settings np.seterr(allignore) # 先设置到已知状态 np.seterr(overraise) {divide: ignore, over: ignore, under: ignore, invalid: ignore} old_settings np.seterr(allwarn, overraise) np.int16(32000) * np.int16(3) Traceback (most recent call last): File stdin, line 1, in module FloatingPointError: overflow encountered in scalar multiply np.seterr(**old_settings) # 恢复原设置numpy.seterrcall(func)用于设置call模式下的回调函数f(err, flag)其中flag的低 4 位分别表示 divide/over/under/invalid或log模式下带write(msg)方法的日志对象配套的查询函数是geterr与geterrcall。回调与日志示例可直接在该模块的 docstring 中运行验证。更推荐的做法是使用上下文管理器numpy.errstate自 1.17 起也可作为装饰器使用NumPy 2.0 起完全线程与 asyncio 安全但同一实例不可重复进入也不应用于装饰异步函数 with np.errstate(divideignore): ... np.arange(3) / 0. array([nan, inf, inf]) with np.errstate(invalidraise): ... np.sqrt(-1) Traceback (most recent call last): File stdin, line 2, in module FloatingPointError: invalid value encountered in sqrtseterr/seterrcall/errstate都通过_make_extobj构造扩展对象并写入_extobj_contextvarnumpy/_core/_ufunc_config.py 第 100-108 行因此它们是逐线程生效的。相关行为测试见 numpy/_core/tests/test_errstate.py其中覆盖了非法模式值抛错、自引用数组等边界情况。覆盖 ufunc 行为array_ufunc协议类包括 ndarray 子类可以通过定义特殊方法__array_ufunc__、__array_wrap__和__array_priority__来覆盖 ufunc 作用于它们时的行为。这是 NumPy 1.13 引入的正式分派协议__array_ufunc__(self, ufunc, method, /, *inputs, **kwargs)当输入包含该类的实例时ufunc 的调用包括其所有方法如reduce、at会被完整转发到这里实现完全自定义的分派逻辑__array_wrap__决定结果的包装/返回类型__array_priority__在多个输入都实现包装方法时决定优先级。关于该协议的完整细节含 ndarray 子类示例参见仓库中的 doc/source/reference/arrays.classes.rst 与 basics.dispatch以及ndarray各协议的文档numpy/_core/_add_newdocs.py 第 3139-3156 行附近。总结与实践建议性能优先使用内建 ufunc 而非 Python 循环ufunc 的广播通过 stride 置 0 实现零拷贝类型转换则通过受控的内部缓冲区完成。归约精度对add/multiply之外的自定义归约务必显式指定dtype避免精度溢出利用out复用预分配内存。错误控制用errstate上下文而非全局seterr临时切换浮点异常行为用seterrcall捕获错误详情。类型预判编写通用代码前可用np.can_cast、np.result_type、np.promote_typesnumpy/_core/multiarray.py推演类型提升结果避免意外上转。扩展需要自定义逐元素函数时np.frompyfunc是最快上手的方式追求性能则参考 c-info.ufunc-tutorial 编写 C 扩展并可通过注册专用 reduction 循环让自定义 ufunc 支持reducec-info.reduction-loop-tutorial。进一步阅读Ufunc 参考手册含where、axes、axis、keepdims、casting、order等可选关键字参数的完整说明、广播基础、dtype 层级与字符代码。【免费下载链接】numpyThe fundamental package for scientific computing with Python.项目地址: https://gitcode.com/gh_mirrors/nu/numpy创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考