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和+/-inf的间距返回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

int npy_clear_floatstatus()#

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

注意,npy_clear_floatstatus_barrier 更为可取,因为它可以防止编译器进行激进的优化,从而避免重新排序相对于设置状态的代码的调用,这可能会导致错误的结果。

int npy_clear_floatstatus_barrier(char*)#

清除浮点状态。传入指向局部变量的指针,以防止编译器进行激进的优化重新排序此函数调用。返回之前的状态掩码。

复数支持#

已添加类似 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 编译和链接选项添加到您的扩展中。要采取的确切步骤取决于您正在使用的构建系统。要采取的通用步骤是

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

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

  3. libnpymathlibm 链接。

注意

请记住,当您进行交叉编译时,必须使用您正在为其构建平台的 numpy,而不是构建机器的原生 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);
}

外部链接

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 位双精度浮点数。