用于使用 igraph C 库
BLAS 是一个高度优化的库,用于基本线性代数运算,例如向量-向量、矩阵-向量和矩阵-矩阵乘积。有关详细信息和 Fortran 中的参考实现,请参阅 http://www.netlib.org/blas/。igraph 包含一些包装函数,可用于以更友好的方式调用 BLAS 例程。并非所有 BLAS 例程都包含在 igraph 中,即使包含的例程也可能没有包装器;包装函数的集合的扩展可能会受到 igraph 内部需求的驱动。包装函数通常会将 BLAS 使用的双精度浮点数组替换为 igraph_vector_t 和 igraph_matrix_t 实例,并且还删除那些可以直接从传递的参数推断的参数(例如行/列数)。
igraph_error_t igraph_blas_ddot(const igraph_vector_t *v1, const igraph_vector_t *v2, igraph_real_t *res);
参数:
|
第一个向量。 |
|
第二个向量。 |
|
指向实数的指针,结果将存储在此处。 |
时间复杂度:O(n),其中 n 是向量的长度。
示例 29.1. 文件 examples/simple/blas.c
#include <igraph.h> int main(void) { igraph_matrix_t m; igraph_vector_t x, y, z; igraph_real_t xz, xx; igraph_vector_init_real(&x, 3, 1.0, 2.0, 3.0); igraph_vector_init_real(&y, 4, 4.0, 5.0, 6.0, 7.0); igraph_vector_init_real(&z, 3, -1.0, 0.0, 0.5); igraph_matrix_init(&m, 4, 3); MATRIX(m, 0, 0) = 1; MATRIX(m, 0, 1) = 2; MATRIX(m, 0, 2) = 3; MATRIX(m, 1, 0) = 2; MATRIX(m, 1, 1) = 3; MATRIX(m, 1, 2) = 4; MATRIX(m, 2, 0) = 3; MATRIX(m, 2, 1) = 4; MATRIX(m, 2, 2) = 5; MATRIX(m, 3, 0) = 4; MATRIX(m, 3, 1) = 5; MATRIX(m, 3, 2) = 6; /* Compute 2 m.x + 3 y and store it in y. */ igraph_blas_dgemv(/* transpose= */ 0, /* alpha= */ 2, &m, &x, /* beta= */ 3, &y); igraph_vector_print(&y); /* Compute the squared norm of x, as well as the dor product of x and z. */ igraph_blas_ddot(&x, &x, &xx); igraph_blas_ddot(&x, &z, &xz); printf("x.x = %g, x.z = %g\n", xx, xz); igraph_matrix_destroy(&m); igraph_vector_destroy(&z); igraph_vector_destroy(&y); igraph_vector_destroy(&x); return 0; }
igraph_real_t igraph_blas_dnrm2(const igraph_vector_t *v);
参数:
|
向量。 |
返回值:
实数值, |
时间复杂度:O(n),其中 n 是向量的长度。
igraph_error_t igraph_blas_dgemv(igraph_bool_t transpose, igraph_real_t alpha, const igraph_matrix_t *a, const igraph_vector_t *x, igraph_real_t beta, igraph_vector_t *y);
此函数是 BLAS 中 dgemv
函数的更友好的接口。dgemv
执行操作 y = alpha*A*x + beta*y
,其中 x
和 y
是向量,A
是适当大小的矩阵(对称或非对称)。
参数:
|
是否转置矩阵 |
|
常数 |
|
矩阵 |
|
向量 |
|
常数 |
|
向量 |
时间复杂度:O(nk),如果矩阵的大小为 n x k
返回值:
|
另请参阅:
如果您有数组而不是向量,请使用 |
示例 29.2. 文件 examples/simple/blas.c
#include <igraph.h> int main(void) { igraph_matrix_t m; igraph_vector_t x, y, z; igraph_real_t xz, xx; igraph_vector_init_real(&x, 3, 1.0, 2.0, 3.0); igraph_vector_init_real(&y, 4, 4.0, 5.0, 6.0, 7.0); igraph_vector_init_real(&z, 3, -1.0, 0.0, 0.5); igraph_matrix_init(&m, 4, 3); MATRIX(m, 0, 0) = 1; MATRIX(m, 0, 1) = 2; MATRIX(m, 0, 2) = 3; MATRIX(m, 1, 0) = 2; MATRIX(m, 1, 1) = 3; MATRIX(m, 1, 2) = 4; MATRIX(m, 2, 0) = 3; MATRIX(m, 2, 1) = 4; MATRIX(m, 2, 2) = 5; MATRIX(m, 3, 0) = 4; MATRIX(m, 3, 1) = 5; MATRIX(m, 3, 2) = 6; /* Compute 2 m.x + 3 y and store it in y. */ igraph_blas_dgemv(/* transpose= */ 0, /* alpha= */ 2, &m, &x, /* beta= */ 3, &y); igraph_vector_print(&y); /* Compute the squared norm of x, as well as the dor product of x and z. */ igraph_blas_ddot(&x, &x, &xx); igraph_blas_ddot(&x, &z, &xz); printf("x.x = %g, x.z = %g\n", xx, xz); igraph_matrix_destroy(&m); igraph_vector_destroy(&z); igraph_vector_destroy(&y); igraph_vector_destroy(&x); return 0; }
igraph_error_t igraph_blas_dgemm(igraph_bool_t transpose_a, igraph_bool_t transpose_b, igraph_real_t alpha, const igraph_matrix_t *a, const igraph_matrix_t *b, igraph_real_t beta, igraph_matrix_t *c);
此函数是 BLAS 中 dgemm
函数的更友好的接口。dgemm
计算 alpha*a*b + beta*c,其中 a、b 和 c 是矩阵,其中 a 和 b 可以转置。
参数:
|
是否转置矩阵 |
|
是否转置矩阵 |
|
常数 |
|
矩阵 |
|
矩阵 |
|
常数 |
|
矩阵 |
时间复杂度:O(n m k),其中矩阵 a 的大小为 n × k,矩阵 b 的大小为 k × m。
返回值:
|
示例 29.3. 文件 examples/simple/blas_dgemm.c
#include <igraph.h> int main(void) { igraph_matrix_t a, b, c; igraph_matrix_init(&a, 2, 2); MATRIX(a, 0, 0) = 1; MATRIX(a, 0, 1) = 2; MATRIX(a, 1, 0) = 3; MATRIX(a, 1, 1) = 4; igraph_matrix_init(&b, 2, 2); MATRIX(b, 0, 0) = 5; MATRIX(b, 0, 1) = 6; MATRIX(b, 1, 0) = 7; MATRIX(b, 1, 1) = 8; igraph_matrix_init(&c, 2, 2); igraph_blas_dgemm(1, 1, 0.5, &a, &b, 0, &c); igraph_matrix_printf(&c, "%g"); igraph_matrix_destroy(&a); igraph_matrix_destroy(&b); igraph_matrix_destroy(&c); return 0; }
igraph_error_t igraph_blas_dgemv_array(igraph_bool_t transpose, igraph_real_t alpha, const igraph_matrix_t* a, const igraph_real_t* x, igraph_real_t beta, igraph_real_t* y);
此函数是 BLAS 中 dgemv
函数的更友好的接口。dgemv
执行操作 y = alpha*A*x + beta*y,其中 x 和 y 是向量,A 是适当大小的矩阵(对称或非对称)。
参数:
|
是否转置矩阵 |
|
常数 |
|
矩阵 |
|
向量 |
|
常数 |
|
向量 |
时间复杂度:O(nk),如果矩阵的大小为 n x k
返回值:
|
另请参阅:
如果您有向量而不是数组,请使用 |
LAPACK 用 Fortran90 编写,并提供用于求解线性方程组、线性方程组的最小二乘解、特征值问题和奇异值问题的例程。还提供了相关的矩阵分解(LU、Cholesky、QR、SVD、Schur、广义 Schur),以及相关的计算,例如 Schur 分解的重新排序和估计条件数。处理密集矩阵和带状矩阵,但不处理一般稀疏矩阵。在所有领域中,都为实矩阵和复矩阵提供类似的功能,包括单精度和双精度。
igraph 提供了一个非常有限的 LAPACK 函数集接口,使用常规 igraph 数据结构。
有关 LAPACK 的更多信息,请参见 http://www.netlib.org/lapack/
igraph_error_t igraph_lapack_dgetrf(igraph_matrix_t *a, igraph_vector_int_t *ipiv, int *info);
分解的形式为 A = P * L * U,其中 P 是置换矩阵,L 是具有单位对角元素的下三角矩阵(如果 m > n 则为下梯形矩阵),U 是上三角矩阵(如果 m < n 则为上梯形矩阵)。
参数:
|
输入/输出矩阵。在输入时,是要分解的 M×N 矩阵。在输出时,是分解 A = P * L * U 中的因子 L 和 U;不存储 L 的单位对角元素。 |
|
一个整数向量,枢轴索引存储在此处,除非它是空指针。矩阵的第 |
|
LAPACK 错误代码。成功退出时为零。如果其值为正数 i,则表示 U(i,i) 恰好为零。分解已完成,但因子 U 恰好是奇异的,如果用于求解方程组,则会发生零除错误。如果 LAPACK 返回错误,即负 info 值,则也会生成 igraph 错误。 |
返回值:
错误代码。 |
时间复杂度:待定。
igraph_error_t igraph_lapack_dgetrs(igraph_bool_t transpose, const igraph_matrix_t *a, const igraph_vector_int_t *ipiv, igraph_matrix_t *b);
此函数调用 LAPACK 来求解线性方程组 A * X = B 或 A' * X = B,其中 A 是一个一般 N×N 矩阵,使用 igraph_lapack_dgetrf
计算的 LU 分解。
参数:
|
布尔值,是否转置输入矩阵。 |
|
包含分解 A = P*L*U 中的因子 L 和 U 的矩阵。L 应为单位三角矩阵,对角线条目是 U 的对角线条目。如果 A 是奇异的,则不会给出警告或错误,并且将返回随机输出。 |
|
一个整数向量, |
|
必须在此处给出右手侧矩阵。解也将放置在此处。 |
返回值:
错误代码。 |
时间复杂度:待定。
igraph_error_t igraph_lapack_dgesv(igraph_matrix_t *a, igraph_vector_int_t *ipiv, igraph_matrix_t *b, int *info);
此函数计算实数线性方程组 A * X = B 的解,其中 A 是一个 N×N 矩阵,X 和 B 是 N×NRHS 矩阵。
使用带部分枢轴和行交换的 LU 分解将 A 分解为 A = P * L * U,其中 P 是置换矩阵,L 是单位下三角矩阵,U 是上三角矩阵。然后使用 A 的分解形式求解方程组 A * X = B。
参数:
|
矩阵。在输入时,是 N×N 系数矩阵,在输出时,是分解 A=P*L*U 中的因子 L 和 U;不存储 L 的单位对角元素。 |
|
一个整数向量或一个空指针。如果不是空指针,则定义置换矩阵 P 的枢轴索引将存储在此处。矩阵的第 i 行与 IPIV(i) 行互换。 |
|
矩阵,在输入时,右手侧矩阵应存储在此处。在输出时,如果没有错误,并且 info 参数为零,则它包含解矩阵 X。 |
|
LAPACK 信息代码。如果它是正数,则 U(info,info) 恰好为零。在这种情况下,分解已完成,但因子 U 恰好是奇异的,因此无法计算解。 |
返回值:
错误代码。 |
时间复杂度:待定。
示例 29.4. 文件 examples/simple/igraph_lapack_dgesv.c
#include <igraph.h> #include <stdio.h> #define DIM 10 void igraph_print_warning(const char *reason, const char *file, int line) { IGRAPH_UNUSED(file); IGRAPH_UNUSED(line); printf("Warning: %s\n", reason); } int main(void) { igraph_matrix_t A, B, RHS; int info; int i, j; /* Identity matrix, you have to start somewhere */ igraph_matrix_init(&A, DIM, DIM); igraph_matrix_init(&B, DIM, 1); for (i = 0; i < DIM; i++) { MATRIX(A, i, i) = 1.0; MATRIX(B, i, 0) = i + 1; } igraph_matrix_init_copy(&RHS, &B); igraph_lapack_dgesv(&A, /*ipiv=*/ 0, &RHS, &info); if (info != 0) { return 1; } if (!igraph_matrix_all_e(&B, &RHS)) { return 2; } igraph_matrix_destroy(&A); igraph_matrix_destroy(&B); igraph_matrix_destroy(&RHS); /* Diagonal matrix */ igraph_matrix_init(&A, DIM, DIM); igraph_matrix_init(&RHS, DIM, 1); for (i = 0; i < DIM; i++) { MATRIX(A, i, i) = i + 1; MATRIX(RHS, i, 0) = i + 1; } igraph_lapack_dgesv(&A, /*ipiv=*/ 0, &RHS, &info); if (info != 0) { return 3; } for (i = 0; i < DIM; i++) { if (MATRIX(RHS, i, 0) != 1.0) { return 4; } } igraph_matrix_destroy(&A); igraph_matrix_destroy(&RHS); /* A general matrix */ igraph_rng_seed(igraph_rng_default(), 42); igraph_matrix_init(&A, DIM, DIM); igraph_matrix_init(&B, DIM, 1); igraph_matrix_init(&RHS, DIM, 1); for (i = 0; i < DIM; i++) { int j; MATRIX(B, i, 0) = igraph_rng_get_integer(igraph_rng_default(), 1, 10); for (j = 0; j < DIM; j++) { MATRIX(A, i, j) = igraph_rng_get_integer(igraph_rng_default(), 1, 10); } } igraph_blas_dgemv_array(/*transpose=*/ 0, /*alpha=*/ 1.0, /*a=*/ &A, /*x-*/ &MATRIX(B, 0, 0), /*beta=*/ 0, /*y=*/ &MATRIX(RHS, 0, 0)); igraph_lapack_dgesv(&A, /*ipiv=*/ 0, &RHS, &info); if (info != 0) { return 5; } for (i = 0; i < DIM; i++) { if (fabs(MATRIX(B, i, 0) - MATRIX(RHS, i, 0)) > 1e-11) { return 6; } } igraph_matrix_destroy(&A); igraph_matrix_destroy(&B); igraph_matrix_destroy(&RHS); /* A singular matrix */ igraph_matrix_init(&A, DIM, DIM); igraph_matrix_init(&B, DIM, 1); igraph_matrix_init(&RHS, DIM, 1); for (i = 0; i < DIM; i++) { MATRIX(B, i, 0) = igraph_rng_get_integer(igraph_rng_default(), 1, 10); for (j = 0; j < DIM; j++) { MATRIX(A, i, j) = i == j ? 1 : 0; } } for (i = 0; i < DIM; i++) { MATRIX(A, DIM - 1, i) = MATRIX(A, 0, i); } igraph_blas_dgemv_array(/*transpose=*/ 0, /*alpha=*/ 1.0, /*a=*/ &A, /*x-*/ &MATRIX(B, 0, 0), /*beta=*/ 0, /*y=*/ &MATRIX(RHS, 0, 0)); igraph_set_warning_handler(igraph_print_warning); igraph_lapack_dgesv(&A, /*ipiv=*/ 0, &RHS, &info); if (info != 10) { printf("LAPACK returned info = %d, should have been 10", info); return 7; } igraph_matrix_destroy(&A); igraph_matrix_destroy(&B); igraph_matrix_destroy(&RHS); return 0; }
igraph_error_t igraph_lapack_dsyevr(const igraph_matrix_t *A, igraph_lapack_dsyev_which_t which, igraph_real_t vl, igraph_real_t vu, int vestimate, int il, int iu, igraph_real_t abstol, igraph_vector_t *values, igraph_matrix_t *vectors, igraph_vector_int_t *support);
调用 DSYEVR LAPACK 函数来计算实对称矩阵 A 的选定特征值,以及可选的特征向量。可以通过指定所需特征值的范围值或范围索引来选择特征值和特征向量。
有关更多信息,请参见 LAPACK 文档。
参数:
|
矩阵,在输入时,它包含对称输入矩阵。只有前导 N×N 上三角部分用于计算。 |
|
给出要计算的特征值(以及可能的相应特征向量)的常量。可能的值为 |
|
如果 |
|
如果 |
|
如果 |
|
如果 |
|
如果 |
|
特征值的绝对误差容限。当近似特征值被确定位于宽度小于或等于 abstol + EPS * max(|a|,|b|) 的区间 [a,b] 中时,它被接受为已收敛,其中 EPS 是机器精度。 |
|
一个已初始化的向量,特征值存储在此处,除非它是一个空指针。它将根据需要调整大小。 |
|
一个已初始化的矩阵,特征向量存储在其列中,除非它是一个空指针。它将根据需要调整大小。 |
|
一个整数向量。如果不是空指针,则会将其调整为 (2*max(1,M))(M 是找到的特征值的总数)。然后, |
返回值:
错误代码。 |
时间复杂度:待定。
示例 29.5. 文件 examples/simple/igraph_lapack_dsyevr.c
#include <igraph.h> int main(void) { igraph_matrix_t A; igraph_matrix_t vectors; igraph_vector_t values; igraph_matrix_init(&A, 2, 2); igraph_matrix_init(&vectors, 0, 0); igraph_vector_init(&values, 0); MATRIX(A, 0, 0) = 2.0; MATRIX(A, 0, 1) = -1.0; MATRIX(A, 1, 0) = -1.0; MATRIX(A, 1, 1) = 3.0; printf("Take a subset:\n"); igraph_lapack_dsyevr(&A, IGRAPH_LAPACK_DSYEV_SELECT, /*vl=*/ 0, /*vu=*/ 0, /*vestimate=*/ 0, /*il=*/ 1, /*iu=*/ 1, /*abstol=*/ 1e-10, &values, &vectors, /*support=*/ 0); printf("eigenvalues:\n"); igraph_vector_print(&values); printf("eigenvectors:\n"); igraph_matrix_print(&vectors); printf("\nTake a subset based on an interval:\n"); igraph_lapack_dsyevr(&A, IGRAPH_LAPACK_DSYEV_INTERVAL, /*vl*/ 3, /*vu*/ 4, /*vestimate=*/ 1, /*il=*/ 0, /*iu=*/ 0, /*abstol=*/ 1e-10, &values, &vectors, /*support=*/ 0); printf("eigenvalues:\n"); igraph_vector_print(&values); printf("eigenvectors:\n"); igraph_matrix_print(&vectors); igraph_vector_destroy(&values); igraph_matrix_destroy(&vectors); igraph_matrix_destroy(&A); return 0; }
igraph_error_t igraph_lapack_dgeev(const igraph_matrix_t *A, igraph_vector_t *valuesreal, igraph_vector_t *valuesimag, igraph_matrix_t *vectorsleft, igraph_matrix_t *vectorsright, int *info);
此函数调用 LAPACK 来计算对于 N×N 实非对称矩阵 A 的特征值,以及可选的左特征向量和/或右特征向量。
A 的右特征向量 v(j) 满足 A * v(j) = lambda(j) * v(j),其中 lambda(j) 是其特征值。A 的左特征向量 u(j) 满足 u(j)^H * A = lambda(j) * u(j)^H,其中 u(j)^H 表示 u(j) 的共轭转置。
计算出的特征向量被归一化为具有等于 1 的欧几里德范数和最大的实数分量。
参数:
|
matrix. 在输入时,它包含 N×N 输入矩阵。 |
|
指向已初始化的向量的指针,或者是指向空指针的指针。如果不是空指针,则特征值的实部存储在此处。向量将根据需要调整大小。 |
|
指向已初始化的向量的指针,或者是指向空指针的指针。如果不是空指针,则特征值的虚部存储在此处。向量将根据需要调整大小。 |
|
指向已初始化的矩阵的指针,或者是指向空指针的指针。如果不是空指针,则左特征向量存储在矩阵的列中。矩阵将根据需要调整大小。 |
|
指向已初始化的矩阵的指针,或者是指向空指针的指针。如果不是空指针,则右特征向量存储在矩阵的列中。矩阵将根据需要调整大小。 |
|
此参数用于两个目的。作为输入参数,它给出了如果 QR 算法无法计算所有特征值是否应生成 igraph 错误。如果 |
返回值:
错误代码。 |
时间复杂度:待定。
示例 29.6. 文件 examples/simple/igraph_lapack_dgeev.c
#include <igraph.h> #include <stdio.h> int main(void) { igraph_matrix_t A; igraph_matrix_t vectors_left, vectors_right; igraph_vector_t values_real, values_imag; igraph_vector_complex_t values; int info = 1; igraph_matrix_init(&A, 2, 2); igraph_matrix_init(&vectors_left, 0, 0); igraph_matrix_init(&vectors_right, 0, 0); igraph_vector_init(&values_real, 0); igraph_vector_init(&values_imag, 0); MATRIX(A, 0, 0) = 1.0; MATRIX(A, 0, 1) = 1.0; MATRIX(A, 1, 0) = -1.0; MATRIX(A, 1, 1) = 1.0; igraph_lapack_dgeev(&A, &values_real, &values_imag, &vectors_left, &vectors_right, &info); igraph_vector_complex_create(&values, &values_real, &values_imag); printf("eigenvalues:\n"); igraph_vector_complex_print(&values); printf("left eigenvectors:\n"); igraph_matrix_print(&vectors_left); printf("right eigenvectors:\n"); igraph_matrix_print(&vectors_right); igraph_vector_destroy(&values_imag); igraph_vector_destroy(&values_real); igraph_vector_complex_destroy(&values); igraph_matrix_destroy(&vectors_right); igraph_matrix_destroy(&vectors_left); igraph_matrix_destroy(&A); return 0; }
igraph_error_t igraph_lapack_dgeevx(igraph_lapack_dgeevx_balance_t balance, const igraph_matrix_t *A, igraph_vector_t *valuesreal, igraph_vector_t *valuesimag, igraph_matrix_t *vectorsleft, igraph_matrix_t *vectorsright, int *ilo, int *ihi, igraph_vector_t *scale, igraph_real_t *abnrm, igraph_vector_t *rconde, igraph_vector_t *rcondv, int *info);
此函数计算非对称 N×N 实矩阵的特征值,以及可选的左特征向量和/或右特征向量。
可选地,它还会计算平衡变换以改善特征值和特征向量的条件(ilo
、ihi
、scale
和 abnrm
)、特征值的倒数条件数(rconde
)和右特征向量的倒数条件数(rcondv
)。
A 的右特征向量 v(j) 满足 A * v(j) = lambda(j) * v(j),其中 lambda(j) 是其特征值。A 的左特征向量 u(j) 满足 u(j)^H * A = lambda(j) * u(j)^H,其中 u(j)^H 表示 u(j) 的共轭转置。
计算出的特征向量被归一化为具有等于 1 的欧几里德范数和最大的实数分量。
平衡矩阵意味着置换行和列以使其更接近上三角矩阵,并应用对角相似变换 D * A * D^(-1),其中 D 是对角矩阵,以使 A 的行和列在范数上更接近,并使特征值和特征向量的条件数更小。计算出的倒数条件数对应于平衡矩阵。置换行和列不会改变条件数(在精确算术中),但对角线缩放会改变条件数。有关平衡的更多说明,请参见 LAPACK 用户指南的第 4.10.2 节。请注意,为平衡矩阵获得的特征向量会反变换为 A
的特征向量。
参数:
|
指示是否应平衡输入矩阵。可能的值
|
||||||||
|
输入矩阵,必须是正方形矩阵。 |
||||||||
|
一个已初始化的向量,或者是一个 |
||||||||
|
一个已初始化的向量,或者是一个 |
||||||||
|
一个已初始化的矩阵,或者是一个 |
||||||||
|
一个已初始化的矩阵,或者是一个 |
||||||||
|
|
||||||||
|
如果不是 NULL, |
||||||||
|
指向已初始化的向量的指针,或者是一个 NULL 指针。如果不是一个 NULL 指针,则平衡
进行互换的顺序是从 N 到 |
||||||||
|
指向实变量的指针,平衡矩阵的 1-范数将存储在此处。(1-范数是任何列中元素绝对值之和的最大值。) |
||||||||
|
一个已初始化的向量,或者是一个 NULL 指针。如果不是一个空指针,则特征值的倒数条件数将存储在此处。 |
||||||||
|
一个已初始化的向量,或者是一个 NULL 指针。如果不是一个空指针,则右特征向量的倒数条件数将存储在此处。 |
||||||||
|
此参数用于两个目的。作为输入参数,它给出了如果 QR 算法无法计算所有特征值是否应生成 igraph 错误。如果 |
返回值:
错误代码。 |
时间复杂度:TODO
示例 29.7. 文件 examples/simple/igraph_lapack_dgeevx.c
#include <igraph.h> #include <stdio.h> void matrix_to_complex_vectors(igraph_vector_complex_t *c1, igraph_vector_complex_t *c2, igraph_matrix_t *m) { IGRAPH_REAL(VECTOR(*c1)[0]) = MATRIX(*m, 0, 0); IGRAPH_REAL(VECTOR(*c1)[1]) = MATRIX(*m, 1, 0); IGRAPH_IMAG(VECTOR(*c1)[0]) = MATRIX(*m, 0, 1); IGRAPH_IMAG(VECTOR(*c1)[1]) = MATRIX(*m, 1, 1); IGRAPH_REAL(VECTOR(*c2)[0]) = MATRIX(*m, 0, 0); IGRAPH_REAL(VECTOR(*c2)[1]) = MATRIX(*m, 1, 0); IGRAPH_IMAG(VECTOR(*c2)[0]) = -MATRIX(*m, 0, 1); IGRAPH_IMAG(VECTOR(*c2)[1]) = -MATRIX(*m, 1, 1); } int main(void) { igraph_matrix_t A; igraph_matrix_t vectors_left, vectors_right; igraph_vector_t values_real, values_imag; igraph_vector_complex_t values; igraph_vector_complex_t eigenvector1; igraph_vector_complex_t eigenvector2; int info = 1; igraph_real_t abnrm; igraph_matrix_init(&A, 2, 2); igraph_matrix_init(&vectors_left, 0, 0); igraph_matrix_init(&vectors_right, 0, 0); igraph_vector_init(&values_real, 0); igraph_vector_init(&values_imag, 0); igraph_vector_complex_init(&eigenvector1, 2); igraph_vector_complex_init(&eigenvector2, 2); MATRIX(A, 0, 0) = 1.0; MATRIX(A, 0, 1) = 1.0; MATRIX(A, 1, 0) = -1.0; MATRIX(A, 1, 1) = 1.0; igraph_lapack_dgeevx(IGRAPH_LAPACK_DGEEVX_BALANCE_BOTH, &A, &values_real, &values_imag, &vectors_left, &vectors_right, NULL, NULL, /*scale=*/ NULL, &abnrm, /*rconde=*/ NULL, /*rcondv=*/ NULL, &info); igraph_vector_complex_create(&values, &values_real, &values_imag); printf("eigenvalues:\n"); igraph_vector_complex_print(&values); printf("\nleft eigenvectors:\n"); /*matrix_to_complex_vectors only works because we have two complex conjugate eigenvalues */ matrix_to_complex_vectors(&eigenvector1, &eigenvector2, &vectors_left); igraph_vector_complex_print(&eigenvector1); igraph_vector_complex_print(&eigenvector2); printf("\nright eigenvectors:\n"); matrix_to_complex_vectors(&eigenvector1, &eigenvector2, &vectors_right); igraph_vector_complex_print(&eigenvector1); igraph_vector_complex_print(&eigenvector2); printf("\nOne-norm of the balanced matrix:\n%g\n", abnrm); igraph_vector_destroy(&values_imag); igraph_vector_destroy(&values_real); igraph_vector_complex_destroy(&values); igraph_vector_complex_destroy(&eigenvector1); igraph_vector_complex_destroy(&eigenvector2); igraph_matrix_destroy(&vectors_right); igraph_matrix_destroy(&vectors_left); igraph_matrix_destroy(&A); return 0; }
ARPACK 是一个用于求解大规模特征值问题的库。该软件包旨在计算一般 n
× n
矩阵 A
的几个特征值和相应的特征向量。它最适用于大型稀疏或结构化矩阵 A
,其中结构化意味着矩阵-向量乘积 w <- Av
需要阶数 n
,而不是通常的阶数 n^2
浮点运算。有关详细信息,请参见 https://github.com/opencollab/arpack-ng。
ARPACK 中的特征值计算(在最简单的情况下)涉及计算 Av
乘积,其中 A
是我们使用的矩阵,v
是一个任意向量。类型为 igraph_arpack_function_t
的用户定义函数应执行此乘积。如果可以有效地进行乘积,例如,如果矩阵是稀疏的,则 ARPACK 通常能够非常快速地计算特征值。
在 igraph 中,特征值/特征向量计算通常涉及以下步骤
使用 igraph_arpack_options_init
初始化 igraph_arpack_options_t
数据结构。
在初始化的 igraph_arpack_options_t
对象中设置一些选项。
定义类型为 igraph_arpack_function_t
的函数。此函数的输入是一个向量,输出应该是输出矩阵乘以输入向量。
调用 igraph_arpack_rssolve()
(如果矩阵是对称的)或 igraph_arpack_rnsolve()
。
igraph_arpack_options_t
对象可以多次使用。
如果我们有许多特征值问题要解决,那么创建 igraph_arpack_storage_t
对象并使用 igraph_arpack_storage_init()
初始化它是值得的。此结构包含 ARPACK 所需的所有内存(具有给定的关于特征值问题大小的上限)。然后可以使用相同的 igraph_arpack_storage_t
对象来解决许多问题,而无需始终重新分配所需的内存。当不再需要 igraph_arpack_storage_t
对象时,需要通过调用 igraph_arpack_storage_destroy()
来销毁它。
igraph 不包含所有 ARPACK 例程,只有处理使用双精度实数的对称和非对称特征值问题的例程。
typedef struct igraph_arpack_options_t { /* INPUT */ char bmat[1]; /* I-standard problem, G-generalized */ int n; /* Dimension of the eigenproblem */ char which[2]; /* LA, SA, LM, SM, BE */ int nev; /* Number of eigenvalues to be computed */ igraph_real_t tol; /* Stopping criterion */ int ncv; /* Number of columns in V */ int ldv; /* Leading dimension of V */ int ishift; /* 0-reverse comm., 1-exact with tridiagonal */ int mxiter; /* Maximum number of update iterations to take */ int nb; /* Block size on the recurrence, only 1 works */ int mode; /* The kind of problem to be solved (1-5) 1: A*x=l*x, A symmetric 2: A*x=l*M*x, A symm. M pos. def. 3: K*x = l*M*x, K symm., M pos. semidef. 4: K*x = l*KG*x, K s. pos. semidef. KG s. indef. 5: A*x = l*M*x, A symm., M symm. pos. semidef. */ int start; /* 0: random, 1: use the supplied vector */ int lworkl; /* Size of temporary storage, default is fine */ igraph_real_t sigma; /* The shift for modes 3,4,5 */ igraph_real_t sigmai; /* The imaginary part of shift for rnsolve */ /* OUTPUT */ int info; /* What happened, see docs */ int ierr; /* What happened in the dseupd call */ int noiter; /* The number of iterations taken */ int nconv; int numop; /* Number of OP*x operations */ int numopb; /* Number of B*x operations if BMAT='G' */ int numreo; /* Number of steps of re-orthogonalizations */ /* INTERNAL */ int iparam[11]; int ipntr[14]; } igraph_arpack_options_t;
此数据结构包含 ARPACK 特征值求解器例程的选项。必须通过在其上调用 igraph_arpack_options_init()
来初始化它。然后它可以用于多个 ARPACK 调用,因为 ARPACK 求解器不会修改它。输入选项
值:
|
Character. 决定是解决标准特征值问题 ('I') 还是广义特征值问题 ('B')。 |
||||||||||||||||||||||
|
特征值问题的维度。 |
||||||||||||||||||||||
|
指定要计算哪些特征值/向量。对称矩阵的可能值
非对称矩阵的可能值
|
||||||||||||||||||||||
|
要计算的特征值的数量。 |
||||||||||||||||||||||
|
停止准则:如果里兹值的误差小于 |
||||||||||||||||||||||
|
要生成的 Lanczos 向量的数量。将其设置为零意味着 |
||||||||||||||||||||||
|
数值标量。在当前的 igraph 实现中应将其设置为零。 |
||||||||||||||||||||||
|
零或一。如果为零,则通过反向通信由用户提供位移。如果为一,则相对于约化三对角矩阵 |
||||||||||||||||||||||
|
允许的最大 Arnoldi 更新迭代次数。 |
||||||||||||||||||||||
|
要在递归中使用的块大小。请始终将其保留为默认值,即 1。 |
||||||||||||||||||||||
|
要解决的特征值问题的类型。如果输入矩阵是对称矩阵,则可能的值
请注意,仅测试了
请注意,仅测试了 |
||||||||||||||||||||||
|
是否使用提供的起始向量 (1),或使用随机起始向量 (0)。起始向量必须在 |
输出选项
值:
|
ARPACK 的错误标志。可能的值
ARPACK 也可以返回其他错误标志,但这些错误标志会转换为 igraph 错误,请参阅 |
||||||
|
第二个 ARPACK 调用的错误标志(通常一次特征值计算涉及两次 ARPACK 调用)。这始终为零,因为其他错误代码会转换为 igraph 错误。 |
||||||
|
Arnoldi 迭代的次数。 |
||||||
|
收敛的里兹值的数量。这表示满足收敛标准的里兹值的数量。 |
||||||
|
矩阵向量乘法的总数。 |
||||||
|
目前未使用。 |
||||||
|
重新正交化的总步数。 |
内部选项
值:
|
不要修改此选项。 |
|
移位反演模式的移位。 |
|
移位的虚部,用于非对称或复数移位反演模式。 |
|
不要修改此选项。 |
|
不要修改此选项。 |
typedef struct igraph_arpack_storage_t { int maxn, maxncv, maxldv; igraph_real_t *v; igraph_real_t *workl; igraph_real_t *workd; igraph_real_t *d; igraph_real_t *resid; igraph_real_t *ax; int *select; /* The following two are only used for non-symmetric problems: */ igraph_real_t *di; igraph_real_t *workev; } igraph_arpack_storage_t;
公共成员,不要直接修改它们,这些被认为是只读的。
值:
|
矩阵的最大秩。 |
|
最大 NCV。 |
|
最大 LDV。 |
这些成员被认为是私有的
值:
|
工作内存。 |
|
工作内存。 |
|
特征值的内存。 |
|
残差的内存。 |
|
工作内存。 |
|
工作内存。 |
|
特征值的内存,仅限非对称情况。 |
|
工作内存,仅限非对称情况。 |
typedef igraph_error_t igraph_arpack_function_t(igraph_real_t *to, const igraph_real_t *from, int n, void *extra);
参数:
|
指向 |
|
指向 |
|
向量的长度(与输入矩阵的阶数相同)。 |
|
矩阵向量计算函数的额外参数。这来自 |
返回值:
错误代码。如果不是 |
void igraph_arpack_options_init(igraph_arpack_options_t *o);
初始化 ARPACK 选项,将它们设置为默认值。您可以始终将初始化的 igraph_arpack_options_t
对象传递给内置的 igraph 函数,而无需进行任何修改。内置的 igraph 函数会修改选项以执行其计算,例如,igraph_pagerank()
始终搜索具有最大模的特征值,而不管提供的值如何。
但是,如果您想使用 ARPACK 实现您自己的涉及特征值计算的函数,您可能需要自己设置字段。
参数:
|
要初始化的 |
时间复杂度:O(1)。
igraph_error_t igraph_arpack_storage_init(igraph_arpack_storage_t *s, igraph_integer_t maxn, igraph_integer_t maxncv, igraph_integer_t maxldv, igraph_bool_t symm);
如果您想使用 ARPACK 运行多个特征值计算,并且想节省每次运行之间的内存分配/释放,则只需要此函数。否则,可以安全地将空指针作为 igraph_arpack_rssolve()
和 igraph_arpack_rnsolve()
的 storage
参数提供,以使内存自动分配和释放。
如果您不再需要它,请不要忘记对存储对象调用 igraph_arpack_storage_destroy()
函数。
参数:
|
要初始化的 |
|
矩阵的最大阶数。 |
|
计划使用的最大 NCV 参数。 |
|
计划使用的最大 LDV 参数。 |
|
是否将使用此 |
返回值:
错误代码。 |
时间复杂度:O(maxncv*(maxldv+maxn))。
void igraph_arpack_storage_destroy(igraph_arpack_storage_t *s);
参数:
|
将释放内存的 |
时间复杂度:操作系统相关。
igraph_error_t igraph_arpack_rssolve(igraph_arpack_function_t *fun, void *extra, igraph_arpack_options_t *options, igraph_arpack_storage_t *storage, igraph_vector_t *values, igraph_matrix_t *vectors);
这是对称矩阵的 ARPACK 求解器。对于非对称矩阵,请使用 igraph_arpack_rnsolve()
。
参数:
|
指向 |
|
要传递给 |
|
一个 |
|
一个 |
|
如果不是空指针,则它应该是指向已初始化向量的指针。特征值将存储在此处。向量将根据需要调整大小。 |
|
如果不是空指针,则它必须是指向已初始化矩阵的指针。特征向量将存储在矩阵的列中。矩阵将根据需要调整大小。如果 ARPACK 迭代是从给定的起始向量开始的,则此参数或 |
返回值:
错误代码。 |
时间复杂度:取决于矩阵向量乘法。通常,少量的迭代就足够了,因此如果矩阵是稀疏的,并且可以在 O(n) 时间(顶点数)内完成矩阵向量乘法,那么特征值也可以在 O(n) 时间内找到。
igraph_error_t igraph_arpack_rnsolve(igraph_arpack_function_t *fun, void *extra, igraph_arpack_options_t *options, igraph_arpack_storage_t *storage, igraph_matrix_t *values, igraph_matrix_t *vectors);
如果您的矩阵是对称的,请始终考虑调用 igraph_arpack_rssolve()
,它会快得多。igraph_arpack_rnsolve()
用于非对称矩阵。
请注意,对于 2x2 矩阵,不会调用 ARPACK,因为在这些情况下存在精确的代数解。
参数:
|
指向 |
|
要传递给 |
|
一个 |
|
一个 |
|
如果不是空指针,则它应该是指向已初始化矩阵的指针。特征值(可能是复数)将存储在此处。矩阵将有两列,第一列包含特征值的实部,第二列包含特征值的虚部。矩阵将根据需要调整大小。 |
|
如果不是空指针,则它必须是指向已初始化矩阵的指针。特征向量将存储在矩阵的列中。矩阵将根据需要调整大小。请注意,实特征值在此矩阵中将在单列中具有实特征向量;但是,复特征值以共轭对出现,结果矩阵将仅存储与具有正虚部的特征值相对应的特征向量。由于在这种情况下特征向量也是复数,因此它将在特征向量矩阵中占据两列(实部和虚部,按此顺序)。注意:如果特征值向量仅返回复共轭特征值对中具有负虚部的特征值,则结果向量将仍然存储与具有正虚部的特征值相对应的特征向量(因为这是 ARPACK 的工作方式)。 |
返回值:
错误代码。 |
时间复杂度:取决于矩阵向量乘法。通常,少量的迭代就足够了,因此如果矩阵是稀疏的,并且可以在 O(n) 时间(顶点数)内完成矩阵向量乘法,那么特征值也可以在 O(n) 时间内找到。
igraph_error_t igraph_arpack_unpack_complex(igraph_matrix_t *vectors, igraph_matrix_t *values, igraph_integer_t nev);
此函数适用于 igraph_arpack_rnsolve
的输出,并稍微润色一下:它只保留 nev
个特征值/向量,并且每个特征向量都存储在 vectors
矩阵的两列中。
非对称 ARPACK 求解器的输出有点难以解析,因为实特征向量在矩阵中仅占据一列,并且根本不存储复共轭特征向量(通常)。另一个问题是求解器可能会返回比请求的更多的特征值。此函数的常见用法是在 igraph_arpack_rnsolve
之后直接调用它,使用它的 vectors
和 values
参数以及 options
->nev 作为 nev
。这将为具有负虚部的特征值添加向量,并将所有向量作为 2 列(实部和虚部)返回。
参数:
|
特征向量矩阵,由 |
|
特征值矩阵,由 |
|
要保留的特征值/向量的数量。可以小于或等于最初从 ARPACK 请求的数量。 |
返回值:
错误代码。 |
时间复杂度:与 vectors
矩阵中的元素数量呈线性关系。
← 第 28 章. 图算子 | 第 30 章. 二分,即双模图 → |