LAPACK 例程中的 WORK 参数有什么用？

Question

我正在计算对称矩阵的特征值分解scipy.linalg.cython_lapack.syev。从我发现的文档中，我需要传递一个名为 WORK 的数组：

WORK 是 DOUBLE PRECISION 数组，维度 (MAX(1,LWORK)) 在退出时，如果 INFO = 0，WORK(1) 返回最佳 LWORK。

但是，我看不到它的作用（无法理解执行后的值是什么），也看不到它的用途。这个参数的目的是什么？

Answer 1

使用 cython 接口dsyev()from scipy.linalg.cython_lapack 是有道理的：numpy eighwraps dsyevs和scipy eighwraps dsyevr()。但是，遵循 的 Fortran 原型dsyev()，必须提供数组 WORK。

内部使用WORK需要该数组syev（如果LWORK = -1）。

LAPACK 是用Fortran 77编写的，该语言在其标准中不支持堆上的动态分配！动态分配可能依赖于平台或由特定的编译器扩展提供。因此，编写 LAPACK 以便用户可以使用她/他想要的任何东西：静态数组、分配在堆栈上的数组或分配在堆上的数组。

实际上，WORK对库中数组的大小进行硬编码会触发两种逆向情况。要么数组太大，无缘无故增加内存占用，要么数组太小，导致性能不佳或越界错误（分段错误...）。因此，内存管理留给了库的用户。如果提供了数组的最佳大小，则会向用户提供一些帮助LWORK = -1。

如果可以进行动态分配，LAPACK 函数最常见的用途是首先使用 执行工作区查询LWORK = -1，然后使用返回值分配WORK正确大小的数组，最后调用 LAPACK 的例程以获取预期结果。LAPACK 的高端包装器（例如 LAPACKE）的功能就是这样做的：查看LAPACKE 的功能来源LAPACKE_dsyev()！它调用了两次LAPACKE_dsyev_work调用LAPACK_dsyev(wrapping dsyev())的函数。

包装器仍然具有诸如 的功能LAPACKE_dsyev_work()，其中参数work和lwork仍然是必需的。因此，如果通过不在WORK调用之间解除分配，以相似大小多次调用例程，则可以减少分配次数，但用户必须自己执行此操作（请参阅此示例）。此外，LAPACK 的函数的来源ILAENV，用于计算 的优化大小WORK，具有以下文本：

此版本提供了一组参数，这些参数应该在许多当前可用的计算机上提供良好但不是最佳的性能。鼓励用户使用参数中的选项和问题大小信息修改此子例程以设置其特定机器的调整参数。

因此，测试大于工作区查询返回的大小的 WORK 大小可以提高性能。

事实上，LAPACK 中的许多函数都具有WORK和LWORK参数。如果您alloc在文件夹 lapack-3.7.1/SRC 中按 搜索grep -r "alloc" . ，则输出仅包含注释行：

    ./zgejsv.f:*>          Length of CWORK to confirm proper allocation of workspace.
    ./zgejsv.f:*>            In both cases, the allocated CWORK can accommodate blocked runs
    ./zgejsv.f:*>          Length of RWORK to confirm proper allocation of workspace.
    ./zgesdd.f:*       minimal amount of workspace allocated at that point in the code,
    ./zhseqr.f:*     ==== NL allocates some local workspace to help small matrices
    ./dhseqr.f:*     ==== NL allocates some local workspace to help small matrices
    ./dgesdd.f:*       minimal amount of workspace allocated at that point in the code,
    ./shseqr.f:*     ==== NL allocates some local workspace to help small matrices
    ./chseqr.f:*     ==== NL allocates some local workspace to help small matrices
    ./sgesdd.f:*       minimal amount of workspace allocated at that point in the code,
    ./sgejsv.f:*>          Length of WORK to confirm proper allocation of work space.
    ./cgejsv.f:*>          Length of CWORK to confirm proper allocation of workspace.
    ./cgejsv.f:*>            In both cases, the allocated CWORK can accommodate blocked runs
    ./cgejsv.f:*>          Length of RWORK to confirm proper allocation of workspace.
    ./dgejsv.f:*>          Length of WORK to confirm proper allocation of work space.
    ./cgesdd.f:*       minimal amount of workspace allocated at that point in the code,

它表明 LAPACK 的核心不会通过诸如 之类的命令处理堆上的动态内存分配allocate，这对于大型数组很有用：用户必须自己关心。