NumPy 核心数学库#

NumPy 核心数学库 (npymath) 是朝着这个方向迈出的第一步。这个库包含大多数与数学相关的 C99 功能,可以在不支持 C99 的平台上使用。核心数学函数具有与 C99 函数相同的 API,除了 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 的别名:返回与 y 符号相同的 x。适用于任何值,包括无穷大和 NaN。单精度和扩展精度可以通过添加后缀 f 和 l 来获得。

有用的数学常数#

以下数学常数在 npy_math.h 中可用。单精度和扩展精度可以通过分别添加 fl 后缀来获得。

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\))

NPY_PI_2#

圆周率除以 2 (\(\frac{\pi}{2}\))

NPY_PI_4#

圆周率除以 4 (\(\frac{\pi}{4}\))

NPY_1_PI#

圆周率的倒数 (\(\frac{1}{\pi}\))

NPY_2_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 和 +/- 无穷大的间距返回 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

请注意,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()#

清除浮点状态。返回先前的状态掩码。

请注意,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 头文件添加了内置支持。

为了简化跨版本兼容性,已添加使用新 set 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 编译和链接选项添加到您的扩展中。采取的具体步骤将取决于您使用的构建系统。需要采取的通用步骤是

  1. 将 numpy 包含目录(= np.get_include() 的值)添加到您的包含目录中,

  2. npymath 静态库位于 numpy 包含目录旁边的 lib 目录中(即 pathlib.Path(np.get_include()) / '..' / 'lib')。将其添加到您的库搜索目录中,

  3. 链接 libnpymathlibm

注意

请记住,在交叉编译时,您必须使用用于构建平台的 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 提供了一个用于处理 npy_half 值的 API,可以通过包含 <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);
}

外部链接

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 位半精度浮点数转换为以 uint32 存储的 32 位单精度浮点数。

npy_uint64 npy_halfbits_to_doublebits(npy_uint16 h)#

低级函数,将 16 位半精度浮点数转换为以 uint64 存储的 64 位双精度浮点数。