---
title: P?POTRF
description: "计算对称正定矩阵或者Hermite正定矩阵的Cholesky分解。"
url: https://www.hikunpeng.com/document/detail/zh/kunpengboostkithistory/240RC1/accel/kunpengaccel_kml_16_0401.html
sourcePath: /source/zh/kunpengboostkithistory/240RC1/accel/kunpengaccel_kml_16_0401.html
indexId: ce258ac3b02202b324f29fbe87feca62bab5476ec0919c1594394c110b96096b76
---
# P?POTRF

计算对称正定矩阵或者Hermite正定矩阵的Cholesky分解。

即或者，上三角矩阵，下三角矩阵。

#### 接口定义

C Interface：

void pspotrf_(const char *uplo, const int *n, float *a, const int *ia, const int *ja, const int *desca, int *info)

void pdpotrf_(const char *uplo, const int *n, double *a, const int *ia, const int *ja, const int *desca, int *info)

void pcpotrf_(const char *uplo, const int *n, float _Complex *a, const int *ia, const int *ja, const int *desca, int *info)

void pzpotrf_(const char *uplo, const int *n, double _Complex *a, const int *ia, const int *ja, const int *desca, int *info)

Fortran Interface：

PSPOTRF(uplo, n, a, ia, ja, desca, info)

PDPOTRF(uplo, n, a, ia, ja, desca, info)

PCPOTRF(uplo, n, a, ia, ja, desca, info)

PZPOTRF(uplo, n, a, ia, ja, desca, info)


#### 参数

| 参数 | 类型 | 范围 | 说明 | 输入/输出 |
| --- | --- | --- | --- | --- |
| uplo | 字符 | 全局 | U为上三角矩阵，L为下三角矩阵。 | 输入 |
| n | 整型 | 全局 | 要操作的列数，比如子矩阵的列数。 | 输入 |
| a | 在pspotrf中为单精度浮点型数组。 在pdpotrf中为双精度浮点型数组。 在pcpotrf中为单精度复数型数组。 在pzpotrf中为双精度复数型数组。 | 本地 | 调用前保存待分解的分布式矩阵A的本地部分。 当uplo = 'U'，调用后根据在全局矩阵中的位置存放对角线以及对角线以上的部分保存上三角矩阵U。 当uplo = 'L'，调用后根据在全局矩阵中的位置存放对角线以及对角线以下的部分保存下三角矩阵L。 | 输入，输出 |
| ia | 整型 | 全局 | 子矩阵在全局矩阵中的行索引。 | 输入 |
| ja | 整型 | 全局 | 子矩阵在全局矩阵中的列索引。 | 输入 |
| desca | 整型数组 | 本地，全局 | 分布式矩阵A的矩阵描述符。 | 输入 |
| info | 整型 | 全局 | 执行结果： 等于0：退出成功。 小于0：第\-info个参数值不合法。 大于0：A中info大小的主子式非正定，无法完成分解。 | 输出 |


#### 依赖

#include <kscalapack.h>


#### 示例

C Interface：
```
int izero=0;
int ione=1;
int myrank_mpi, nprocs_mpi;
MPI_Init( &argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &myrank_mpi);
MPI_Comm_size(MPI_COMM_WORLD, &nprocs_mpi);
int n = 8;       // (Global) Matrix size
int nprow = 2;   // Number of row procs
int npcol = 2;   // Number of column procs
int nb = 4;      // (Global) Block size
char uplo='L';   // Matrix is lower triangular
char layout='R'; // Block cyclic, Row major processor mapping
printf("Usage: ./test matrix_size block_size nprocs_row nprocs_col\n");
if(argc > 1) {
n = atoi(argv[1]);
}
if(argc > 2) {
nb = atoi(argv[2]);
}
if(argc > 3) {
nprow = atoi(argv[3]);
}
if(argc > 4) {
npcol = atoi(argv[4]);
}
assert(nprow * npcol == nprocs_mpi);
// Initialize BLACS
int iam, nprocs;
int zero = 0;
int ictxt, myrow, mycol;
blacs_pinfo_(&iam, &nprocs) ; // BLACS rank and world size
blacs_get_(&zero, &zero, &ictxt ); // -> Create context
blacs_gridinit_(&ictxt, &layout, &nprow, &npcol ); // Context -> Initialize the grid
blacs_gridinfo_(&ictxt, &nprow, &npcol, &myrow, &mycol ); // Context -> Context grid info (# procs row/col, current procs row/col)
// Compute the size of the local matrices
int mpA    = numroc_( &n, &nb, &myrow, &izero, &nprow ); // My proc -> row of local A
int nqA    = numroc_( &n, &nb, &mycol, &izero, &npcol ); // My proc -> col of local A
//allocate and fill the matrices A and B
double *A;
A = (double *)calloc(mpA*nqA,sizeof(double)) ;
if (A==NULL){ printf("Error of memory allocation A on proc %dx%d\n",myrow,mycol); exit(0); }
int k = 0;
for (int j = 0; j < nqA; j++) { // local col
int l_j = j / nb; // which block
int x_j = j % nb; // where within that block
int J   = (l_j * npcol + mycol) * nb + x_j; // global col
for (int i = 0; i < mpA; i++) { // local row
int l_i = i / nb; // which block
int x_i = i % nb; // where within that block
int I   = (l_i * nprow + myrow) * nb + x_i; // global row
assert(I < n);
assert(J < n);
if(I == J) {
A[k] = nb * nb;
} else {
A[k] = I +  J;
}
//printf("%d %d -> %d %d -> %f\n", i, j, I, J, A[k]);
k++;
}
}
//creat descriptor
int descA[9];
int info=0;
int lddA = mpA > 1 ? mpA : 1;
descinit_( descA,  &n, &n, &nb, &nb, &izero, &izero, &ictxt, &lddA, &info);
if(info != 0) {
printf("Error in descinit, info = %d\n", info);
}
//run pdpotrf_ and time
double MPIt1 = MPI_Wtime();
printf("[%dx%d] Starting \n", myrow, mycol);
pdpotrf_(&uplo, &n, A, &ione, &ione, descA, &info);
if (info != 0) {
printf("Error in caculate, info = %d\n", info);
}
double MPIt2 = MPI_Wtime();
//exit and finanlize
blacs_gridexit_(&ictxt);
MPI_Finalize();
return 0;
/*
* Origin A:
*
*16 1 2 3 4 5 6 7
*1 16 3 4 5 6 7 8
*2 3 16 5 6 7 8 9
*3 4 5 16 7 8 9 10
*4 5 6 7 16 9 10 11
*5 6 7 8 9 16 11 12
*6 7 8 9 10 11 16 13
*7 8 9 10 11 12 13 16
*/
/*After PDPOTRF:
*4.00000 1.00000 2.00000 3.00000 4.00000 5.00000 6.00000 7.00000
*0.250000 3.99218 3.00000 4.00000 5.00000 6.00000 7.00000 8.00000
*0.500000 0.720158 3.90274 5.00000 6.00000 7.00000 8.00000 9.00000
*0.750000 0.954992 1.00884 3.67529 7.00000 8.00000 9.00000 10.0000
*1.00000 1.18983 1.18971 1.06481 3.32191 9.00000 10.0000 11.0000
*1.25000 1.42466 1.37058 1.17522 0.955149 2.86983 11.0000 12.0000
*1.50000 1.65949 1.55145 1.28562 0.996646 0.756689 2.31741 13.0000
*1.75000 1.89433 1.73232 1.39603 1.03814 0.734271 0.500013 1.59132
*/
```
