NumPy 核心数学库#

numpy 核心数学库 (npymath) 是朝这个方向迈出的第一步.这个库包含了大多数与数学相关的 C99 功能,可以在 C99 支持不佳的平台上使用.核心数学函数的 API 与 C99 的相同,除了 npy_* 前缀.

可用函数在 <numpy/npy_math.h> 中定义 - 如有疑问,请参考此头文件.

备注

正在努力使 npymath 更小（因为编译器的 C99 兼容性随着时间的推移有所改善）并且更容易作为供应商或仅头依赖使用.这将避免使用可能与下游包或终端用户使用的编译器不匹配的编译器构建静态库的问题.详情请参见 gh-20880.

浮点分类#

NPY_NAN#: 这个宏定义为 NaN（非数字）,并保证符号位未设置（’正’ NaN）.相应的单精度和扩展精度宏可以通过后缀 F 和 L 获得.

NPY_INFINITY#: 这个宏定义为正无穷大.对应的单精度和扩展精度宏可以通过后缀 F 和 L 获得.

NPY_PZERO#: 此宏定义为正零.对应的单精度和扩展精度宏可以通过后缀 F 和 L 获得.

NPY_NZERO#: 这个宏定义为负零（即符号位设置）.相应的单精度和扩展精度宏可以通过后缀 F 和 L 获得.

npy_isnan(x)#: 这是C99 isnan 的一个别名:适用于单精度、双精度和扩展精度,如果 x 是 NaN,则返回非 0 值.

npy_isfinite(x)#: 这是 C99 isfinite 的一个别名:适用于单精度、双精度和扩展精度,如果 x 既不是 NaN 也不是无穷大,则返回非 0 值.

npy_isinf(x)#: 这是C99 isinf 的一个别名:适用于单精度、双精度和扩展精度,如果 x 是无穷大（正和负）,则返回非 0 值.

npy_signbit(x)#: 这是 C99 signbit 的一个别名:适用于单精度、双精度和扩展精度,如果 x 设置了符号位（即数字为负）,则返回非 0 值.

npy_copysign(x, y)#: 这是 C99 copysign 的别名:返回 x 并带有与 y 相同的符号.适用于任何值,包括 inf 和 nan.单精度和扩展精度可通过后缀 f 和 l 获得.

有用的数学常数#

以下数学常数在 npy_math.h 中可用.单精度和扩展精度也可以分别通过添加 f 和 l 后缀来使用.

NPY_E#: 自然对数的底数 (\(e\))

NPY_LOG2E#: 欧拉常数的以2为底的对数 (\(\frac{\ln(e)}{\ln(2)}\))

NPY_LOG10E#: 以10为底的对数欧拉常数 (\(\frac{\ln(e)}{\ln(10)}\)).

NPY_LOGE2#: 自然对数2 (\(\ln(2)\))

NPY_LOGE10#: 10的自然对数 (\(\ln(10)\))

NPY_PI#: Pi (\(\pi\))

NPY_PI_2#: Pi 除以 2 (\(\frac{\pi}{2}\)).

NPY_PI_4#: Pi 除以 4 (\(\frac{\pi}{4}\)).

NPY_1_PI#: pi的倒数（\(\frac{1}{\pi}\)）

NPY_2_PI#: 两倍于 pi 的倒数 (\(\frac{2}{\pi}\)).

NPY_EULER#

欧拉常数: \(\lim_{n\rightarrow\infty}({\sum_{k=1}^n{\frac{1}{k}}-\ln n})\)

低级浮点数操作#

这些可以用于精确的浮点数比较.

double npy_nextafter(double x, double y)#: 这是对 C99 nextafter 的别名:从 x 向 y 方向返回下一个可表示的浮点值.单精度和扩展精度可通过后缀 f 和 l 获得.

double npy_spacing(double x)#: 这是一个等同于 Fortran 内置函数的函数.返回 x 和从 x 开始的下一个可表示浮点值之间的距离,例如 spacing(1) == eps.nan 和 +/- inf 的 spacing 返回 nan.单精度和扩展精度可通过后缀 f 和 l 获得.

void npy_set_floatstatus_divbyzero()#: 设置除以零的浮点异常

void npy_set_floatstatus_overflow()#: 设置溢出浮点异常

void npy_set_floatstatus_underflow()#: 设置下溢浮点异常

void npy_set_floatstatus_invalid()#: 设置无效的浮点异常

int npy_get_floatstatus()#

获取浮点状态.返回一个包含以下可能标志的位掩码:

NPY_FPE_DIVIDEBYZERO
NPY_FPE_OVERFLOW
NPY_FPE_UNDERFLOW
NPY_FPE_INVALID

注意,:c:func:npy_get_floatstatus_barrier 是首选的,因为它防止了激进的编译器优化相对于设置状态的代码重新排序调用,这可能导致不正确的结果.

int npy_get_floatstatus_barrier(char*)#

获取浮点状态.传递一个指向局部变量的指针,以防止激进的编译器优化重新排序此函数调用相对于设置状态的代码,这可能导致不正确的结果.

返回一个包含以下可能标志的位掩码:

NPY_FPE_DIVIDEBYZERO
NPY_FPE_OVERFLOW
NPY_FPE_UNDERFLOW
NPY_FPE_INVALID

在 1.15.0 版本加入.

int npy_clear_floatstatus()#

清除浮点状态.返回之前的状态掩码.

注意,:c:func:npy_clear_floatstatus_barrier 是首选的,因为它防止了激进的编译器优化重新排序相对于设置状态的代码的调用,这可能导致不正确的结果.

int npy_clear_floatstatus_barrier(char*)#: 清除浮点状态.传递一个指向局部变量的指针,以防止编译器激进的优化重排此函数调用.返回之前的状态掩码.

在 1.15.0 版本加入.

支持复数#

已添加类似 C99 的复数函数.如果您希望实现可移植的 C 扩展,可以使用这些函数.自 NumPy 2.0 以来,我们使用 C99 复数类型作为底层类型:

typedef double _Complex npy_cdouble;
typedef float _Complex npy_cfloat;
typedef long double _Complex npy_clongdouble;

MSVC 本身不支持 _Complex 类型,但通过提供自己的实现,增加了对 C99 complex.h 头文件的支持.因此,在 MSVC 下,将使用等效的 MSVC 类型:

typedef _Dcomplex npy_cdouble;
typedef _Fcomplex npy_cfloat;
typedef _Lcomplex npy_clongdouble;

因为 MSVC 仍然不支持用于初始化复数的 C99 语法,您需要限制使用与 C90 兼容的语法,例如:

/* a = 1 + 2i \*/
npy_complex a = npy_cpack(1, 2);
npy_complex b;

b = npy_log(a);

一些实用工具也被添加到了 numpy/npy_math.h 中,以便检索或设置复数的实部或虚部:

npy_cdouble c;
npy_csetreal(&c, 1.0);
npy_csetimag(&c, 0.0);
printf("%d + %di\n", npy_creal(c), npy_cimag(c));

在 2.0.0 版本发生变更: numpy 所有复数类型的底层 C 类型已更改为使用 C99 复数类型.直到现在,以下内容一直用于表示复数类型:

typedef struct { double real, imag; } npy_cdouble;
typedef struct { float real, imag; } npy_cfloat;
typedef struct {npy_longdouble real, imag;} npy_clongdouble;

使用 struct 表示确保了复数可以在所有平台上使用,即使是没有内置复数类型支持的平台.这也意味着必须与 NumPy 一起发布一个静态库,为下游包提供一个 C99 兼容层.然而,近年来对原生复数类型的支持有了极大的改进,MSVC 在2019年增加了对 complex.h 头文件的内置支持.

为了简化跨版本兼容性,已添加使用新集合API的宏.

#define NPY_CSETREAL(z, r) npy_csetreal(z, r)
#define NPY_CSETIMAG(z, i) npy_csetimag(z, i)

兼容层也在 numpy/npy_2_complexcompat.h 中提供.它会检查宏是否存在,如果不存在则回退到 1.x 语法.

#include <numpy/npy_math.h>

#ifndef NPY_CSETREALF
#define NPY_CSETREALF(c, r) (c)->real = (r)
#endif
#ifndef NPY_CSETIMAGF
#define NPY_CSETIMAGF(c, i) (c)->imag = (i)
#endif

我们建议所有需要此功能的下游包将兼容层代码复制粘贴到他们自己的源代码中并使用它,以便他们可以继续支持 NumPy 1.x 和 2.x 而不会出现问题.还请注意,``complex.h`` 头文件包含在 numpy/npy_common.h 中,这使得 complex 成为一个保留关键字.

在扩展中链接核心数学库#

要在你自己的 Python 扩展中使用 NumPy 作为静态库的核心数学库,你需要将 npymath 编译和链接选项添加到你的扩展中.具体步骤取决于你使用的构建系统.一般的步骤是:

将 numpy 包含目录（= np.get_include() 的值）添加到您的包含目录中,
npymath 静态库位于 lib 目录中,紧邻 numpy 的包含目录（即 pathlib.Path(np.get_include()) / '..' / 'lib'）.将其添加到您的库搜索目录中.
链接 libnpymath 和 libm.

备注

请记住,当你在进行交叉编译时,你必须使用为目标平台构建的 numpy,而不是构建机器的本地版本.否则你会选择一个为错误架构构建的静态库.

当你使用 numpy.distutils （已弃用）构建时,然后在你的 setup.py 中使用这个:

>>> from numpy.distutils.misc_util import get_info
>>> info = get_info('npymath')
>>> _ = config.add_extension('foo', sources=['foo.c'], extra_info=info)

换句话说,使用 info 的方式与使用 blas_info 等完全相同.

当你使用 Meson 构建时,请使用:

# Note that this will get easier in the future, when Meson has
# support for numpy built in; most of this can then be replaced
# by `dependency('numpy')`.
incdir_numpy = run_command(py3,
  [
    '-c',
    'import os; os.chdir(".."); import numpy; print(numpy.get_include())'
  ],
  check: true
).stdout().strip()

inc_np = include_directories(incdir_numpy)

cc = meson.get_compiler('c')
npymath_path = incdir_numpy / '..' / 'lib'
npymath_lib = cc.find_library('npymath', dirs: npymath_path)

py3.extension_module('module_name',
  ...
  include_directories: inc_np,
  dependencies: [npymath_lib],

半精度函数#

头文件 <numpy/halffloat.h> 提供了处理 IEEE 754-2008 16位浮点值的函数.虽然这种格式通常不用于数值计算,但对于需要浮点但不需高精度的存储值非常有用.它还可以作为教育工具,帮助理解浮点舍入误差的性质.

与其他类型类似,NumPy 包含一个 16 位浮点数的 typedef npy_half.与大多数其他类型不同,你不能在 C 语言中将它作为普通类型使用,因为它是对 npy_uint16 的 typedef.例如,1.0 在 C 语言中看起来像 0x3c00,如果你对不同的符号零进行相等比较,你会得到 -0.0 != 0.0（0x8000 != 0x0000）,这是不正确的.

基于这些原因,NumPy 提供了一个 API 来处理 npy_half 值,通过包含 <numpy/halffloat.h> 并链接到 npymath 来访问.对于未直接提供的函数,例如算术运算,首选方法是转换为 float 或 double 并再次转换回来,如下例所示.

npy_half sum(int n, npy_half *array) {
    float ret = 0;
    while(n--) {
        ret += npy_half_to_float(*array++);
    }
    return npy_float_to_half(ret);
}

外部链接:

754-2008 IEEE 浮点运算标准
半精度浮点数维基百科文章
OpenGL Half Float Pixel 支持
The OpenEXR 图像格式.

NPY_HALF_ZERO#: 此宏定义为正零.

NPY_HALF_PZERO#: 此宏定义为正零.

NPY_HALF_NZERO#: 这个宏定义为负零.

NPY_HALF_ONE#: 这个宏定义为 1.0.

NPY_HALF_NEGONE#: 这个宏定义为 -1.0.

NPY_HALF_PINF#: 此宏定义为 +inf.

NPY_HALF_NINF#: 这个宏定义为 -inf.

NPY_HALF_NAN#: 这个宏定义为一个NaN值,保证其符号位未设置.

float npy_half_to_float(npy_half h)#: 将半精度浮点数转换为单精度浮点数.

double npy_half_to_double(npy_half h)#: 将半精度浮点数转换为双精度浮点数.

npy_half npy_float_to_half(float f)#: 将单精度浮点数转换为半精度浮点数.该值四舍五入到最接近的可表示半精度值,平局时取最近的偶数.如果值太小或太大,系统将设置浮点数下溢或上溢位.

npy_half npy_double_to_half(double d)#: 将双精度浮点数转换为半精度浮点数.该值四舍五入到最接近的可表示半精度值,平局时取最近的偶数.如果值太小或太大,系统将设置浮点下溢或上溢位.

int npy_half_eq(npy_half h1, npy_half h2)#: 比较两个半精度浮点数（h1 == h2）.

int npy_half_ne(npy_half h1, npy_half h2)#: 比较两个半精度浮点数（h1 != h2）.

int npy_half_le(npy_half h1, npy_half h2)#: 比较两个半精度浮点数（h1 <= h2）.

int npy_half_lt(npy_half h1, npy_half h2)#: 比较两个半精度浮点数（h1 < h2）.

int npy_half_ge(npy_half h1, npy_half h2)#: 比较两个半精度浮点数（h1 >= h2）.

int npy_half_gt(npy_half h1, npy_half h2)#: 比较两个半精度浮点数（h1 > h2）.

int npy_half_eq_nonan(npy_half h1, npy_half h2)#: 比较两个已知不是 NaN 的半精度浮点数（h1 == h2）.如果一个值是 NaN,结果是未定义的.

int npy_half_lt_nonan(npy_half h1, npy_half h2)#: 比较两个已知不是 NaN 的半精度浮点数（h1 < h2）.如果一个值是 NaN,结果是未定义的.

int npy_half_le_nonan(npy_half h1, npy_half h2)#: 比较两个已知不是 NaN 的半精度浮点数（h1 <= h2）.如果一个值是 NaN,结果是未定义的.

int npy_half_iszero(npy_half h)#: 测试半精度浮点数是否等于零.这可能比调用 npy_half_eq(h, NPY_ZERO) 稍微快一些.

int npy_half_isnan(npy_half h)#: 测试半精度浮点数是否为 NaN.

int npy_half_isinf(npy_half h)#: 测试半精度浮点数是否为正无穷或负无穷.

int npy_half_isfinite(npy_half h)#: 测试半精度浮点数是否有限（不是 NaN 或 Inf）.

int npy_half_signbit(npy_half h)#: 如果 h 为负数,则返回 1,否则返回 0.

npy_half npy_half_copysign(npy_half x, npy_half y)#: 返回x的值,其符号位从y复制.适用于任何值,包括Inf和NaN.

npy_half npy_half_spacing(npy_half h)#: 对于半精度浮点数,这与低级浮点数部分中描述的 npy_spacing 和 npy_spacingf 相同.

npy_half npy_half_nextafter(npy_half x, npy_half y)#: 这对于半精度浮点数与低级浮点数部分中描述的 npy_nextafter 和 npy_nextafterf 相同.

npy_uint16 npy_floatbits_to_halfbits(npy_uint32 f)#: 低级函数,将存储为 uint32 的 32 位单精度浮点数转换为 16 位半精度浮点数.

npy_uint16 npy_doublebits_to_halfbits(npy_uint64 d)#: 低级函数,将一个以 uint64 存储的 64 位双精度浮点数转换为 16 位半精度浮点数.

npy_uint32 npy_halfbits_to_floatbits(npy_uint16 h)#: 低级函数,将16位半精度浮点数转换为32位单精度浮点数,存储为uint32.

npy_uint64 npy_halfbits_to_doublebits(npy_uint16 h)#: 将16位半精度浮点数转换为64位双精度浮点数的低级函数,存储为uint64.