Actual source code: ex243.c
1: static char help[] = "Test conversion of ScaLAPACK matrices.\n\n";
3: #include <petscmat.h>
5: int main(int argc, char** argv)
6: {
7: Mat A,A_scalapack;
8: PetscInt i,j,M=10,N=5,nloc,mloc,nrows,ncols;
9: PetscMPIInt rank,size;
10: IS isrows,iscols;
11: const PetscInt *rows,*cols;
12: PetscScalar *v;
13: MatType type;
14: PetscBool isDense,isAIJ,flg;
16: PetscInitialize(&argc,&argv,(char*)0,help);
17: MPI_Comm_size(PETSC_COMM_WORLD,&size);
18: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
19: PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL);
20: PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL);
22: /* Create a matrix */
23: MatCreate(PETSC_COMM_WORLD, &A);
24: mloc = PETSC_DECIDE;
25: PetscSplitOwnershipEqual(PETSC_COMM_WORLD,&mloc,&M);
26: nloc = PETSC_DECIDE;
27: PetscSplitOwnershipEqual(PETSC_COMM_WORLD,&nloc,&N);
28: MatSetSizes(A,mloc,nloc,M,N);
29: MatSetType(A,MATDENSE);
30: MatSetFromOptions(A);
31: MatSetUp(A);
33: /* Set local matrix entries */
34: MatGetOwnershipIS(A,&isrows,&iscols);
35: ISGetLocalSize(isrows,&nrows);
36: ISGetIndices(isrows,&rows);
37: ISGetLocalSize(iscols,&ncols);
38: ISGetIndices(iscols,&cols);
39: PetscMalloc1(nrows*ncols,&v);
41: for (i=0; i<nrows; i++) {
42: for (j=0; j<ncols; j++) {
43: if (size == 1) {
44: v[i*ncols+j] = (PetscScalar)(i+j);
45: } else {
46: v[i*ncols+j] = (PetscScalar)rank+j*0.1;
47: }
48: }
49: }
50: MatSetValues(A,nrows,rows,ncols,cols,v,INSERT_VALUES);
51: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
52: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
54: /* Test MatSetValues() by converting A to A_scalapack */
55: MatGetType(A,&type);
56: if (size == 1) {
57: PetscObjectTypeCompare((PetscObject)A,MATSEQDENSE,&isDense);
58: PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&isAIJ);
59: } else {
60: PetscObjectTypeCompare((PetscObject)A,MATMPIDENSE,&isDense);
61: PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&isAIJ);
62: }
64: if (isDense || isAIJ) {
65: Mat Aexplicit;
66: MatConvert(A,MATSCALAPACK,MAT_INITIAL_MATRIX,&A_scalapack);
67: MatComputeOperator(A_scalapack,isAIJ?MATAIJ:MATDENSE,&Aexplicit);
68: MatMultEqual(Aexplicit,A_scalapack,5,&flg);
70: MatDestroy(&Aexplicit);
72: /* Test MAT_REUSE_MATRIX which is only supported for inplace conversion */
73: MatConvert(A,MATSCALAPACK,MAT_INPLACE_MATRIX,&A);
74: MatMultEqual(A_scalapack,A,5,&flg);
76: MatDestroy(&A_scalapack);
77: }
79: ISRestoreIndices(isrows,&rows);
80: ISRestoreIndices(iscols,&cols);
81: ISDestroy(&isrows);
82: ISDestroy(&iscols);
83: PetscFree(v);
84: MatDestroy(&A);
85: PetscFinalize();
86: return 0;
87: }
89: /*TEST
91: build:
92: requires: scalapack
94: test:
95: nsize: 6
97: test:
98: suffix: 2
99: nsize: 6
100: args: -mat_type aij
101: output_file: output/ex243_1.out
103: test:
104: suffix: 3
105: nsize: 6
106: args: -mat_type scalapack
107: output_file: output/ex243_1.out
109: TEST*/