Cholesky用ScaLAPACK [英] Cholesky with ScaLAPACK

查看:559
本文介绍了Cholesky用ScaLAPACK的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我想通过MKL的,您可以看到函数的参数意义。






代码基于问题。输出:

  gsamaras @ pythagoras:〜/ konstantis / check_examples $ ../../mpich-install/bin/mpic++  - o test minor.cpp -I ../../ intel / mkl / include ../../intel/mkl/lib/intel64/libmkl_scalapack_lp64.a -Wl, -  start-group ../../intel /mkl/lib/intel64/libmkl_intel_lp64.a ../../intel/mkl/lib/intel64/libmkl_core.a ../../intel/mkl/lib/intel64/libmkl_sequential.a -Wl, -  end -group ../../intel/mkl/lib/intel64/libmkl_blacs_intelmpi_lp64.a -lpthread -lm -ldl 
gsamaras @ pythagoras:〜/ konstantis / check_examples $ mpiexec -n 4 ./test
进程网格模式:
0 1
2 3
nrows = 3,ncols = 3
节点0上的A_loc
4 1 2
1 0.5 0
2 0 16

nrows = 3,ncols = 2
节点1上的A_loc
2 0.5
0 0
0 0

nrows = 2,ncols = 3
节点2上的A_loc
2 0 0
0.5 0 0

nrows = 2,ncols = 2
节点3上的A_loc
3 0
0 0.625

说明init sucesss!
矩阵不为正定
矩阵A结果:
2 1 2 0.5 2
0.5 0.5 0 0 0
1 0 1 0 -0.25
0.25 -1 -0.5 0.625 0
1 -1 -2 -0.5 14


解决方案

问题可能来自:

  MPI_Bcast(& lda,1,MPI_INT,0,MPI_COMM_WORLD ); 

在此行之前, lda 每个过程如果矩阵的维度是奇数。两个进程处理2行,两个进程处理3行。但是 MPI_Bcast()后, lda 在任何地方都是一样的(3)。



问题是子程序 lda /dd/d22/descinit_8f_source.htmlrel =nofollow> DESCINIT 必须是本地数组的前导维,即2或3。



通过评论 MPI_Bcast(),我得到:

 描述init sucesss! 
SUCCESS
矩阵A结果:
2 1 2 0.5 2
0.5 0.5 0 0 0
1 -1 1 0 0
0.25 -0.25 -0.5 0.5 0
1 -1 -2 -3 1

最后,该程序适用于均匀尺寸,并且对于奇数维度失败!


I am trying to do a Cholesky decomposition via pdpotrf() of MKL-Intel's library, which uses ScaLAPACK. I am reading the whole matrix in the master node and then distribute it like in this example. Everything works fine when the dimension of the SPD matrix is even. However, when it's odd, pdpotrf() thinks that the matrix is not positive definite.

Could it be because the submatrices are not SPD? I am working with this matrix:

and the submatrices are (with 4 processes and blocks of size 2x2):

A_loc on node 0
  4   1   2
  1 0.5   0
  2   0  16

nrows = 3, ncols = 2
A_loc on node 1
  2 0.5
  0   0
  0   0

nrows = 2, ncols = 3
A_loc on node 2
  2   0   0
0.5   0   0

nrows = 2, ncols = 2
A_loc on node 3
  3   0
  0 0.625

Here, every submatrix is not SPD, however, the overall matrix is SPD (have checked with running with 1 process). What should I do? Or there is nothing I can do and pdpotrf() does not work with matrices of odd size?


Here is how I call the routine:

int iZERO = 0;
int descA[9];
// N, M dimensions of matrix. lda = N
// Nb, Mb dimensions of block
descinit_(descA, &N, &M, &Nb, &Mb, &iZERO, &iZERO, &ctxt, &lda, &info);
...
pdpotrf((char*)"L", &ord, A_loc, &IA, &JA, descA, &info);

I also tried this:

// nrows/ncols is the number of rows/columns a submatrix has
descinit_(descA, &N, &M, &nrows, &ncols, &iZERO, &iZERO, &ctxt, &lda, &info);

but I get an error:

{ 0, 0}: On entry to { 0, 1}: On entry to PDPOTR{ 1, 0}: On entry to PDPOTRF parameter number 605 had an illegal value { 1, 1}: On entry to PDPOTRF parameter number 605 had an illegal value F parameter number 605 had an illegal value

PDPOTRF parameter number 605 had an illegal value info < 0: If the i-th argument is an array and the j-entry had an illegal value, then INFO = -(i*100+j), if the i-th argument is a scalar and had an illegal value, then INFO = -i. info = -605

From my answer, you can see what the arguments of the function mean.


The code is based on this question. Output:

gsamaras@pythagoras:~/konstantis/check_examples$ ../../mpich-install/bin/mpic++ -o test minor.cpp -I../../intel/mkl/include ../../intel/mkl/lib/intel64/libmkl_scalapack_lp64.a       -Wl,--start-group       ../../intel/mkl/lib/intel64/libmkl_intel_lp64.a ../../intel/mkl/lib/intel64/libmkl_core.a  ../../intel/mkl/lib/intel64/libmkl_sequential.a    -Wl,--end-group ../../intel/mkl/lib/intel64/libmkl_blacs_intelmpi_lp64.a -lpthread       -lm     -ldl
gsamaras@pythagoras:~/konstantis/check_examples$ mpiexec -n 4 ./test
Processes grid pattern:
0 1
2 3
nrows = 3, ncols = 3
A_loc on node 0
  4   1   2
  1 0.5   0
  2   0  16

nrows = 3, ncols = 2
A_loc on node 1
  2 0.5
  0   0
  0   0

nrows = 2, ncols = 3
A_loc on node 2
  2   0   0
0.5   0   0

nrows = 2, ncols = 2
A_loc on node 3
  3   0
  0 0.625

Description init sucesss!
matrix is not positive definte
Matrix A result:
  2   1   2 0.5   2
0.5 0.5   0   0   0
  1   0   1   0 -0.25
0.25  -1 -0.5 0.625   0
  1  -1  -2 -0.5  14

解决方案

The issue may come from :

MPI_Bcast(&lda, 1, MPI_INT, 0, MPI_COMM_WORLD);

Before this line, lda is different on each process if the dimension of the matrix is odd. Two processes handle 2 rows and two processes handle 3 rows. But after the MPI_Bcast(), lda is the same everywhere (3).

The problem is that the argument lda of the subroutine DESCINIT must be the leading dimension of the local array, that is either 2 or 3.

By commenting MPI_Bcast(), i got:

Description init sucesss!
SUCCESS
Matrix A result:
2   1   2 0.5   2 
0.5 0.5   0   0   0 
1  -1   1   0   0 
0.25 -0.25 -0.5 0.5   0 
1  -1  -2  -3   1 

At last, it would explain that the program works well for even dimensions and fails for odd dimensions !

这篇关于Cholesky用ScaLAPACK的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆