1 static char help[] = "Test conversion of ScaLAPACK matrices.\n\n"; 2 3 #include <petscmat.h> 4 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; 15 16 CHKERRQ(PetscInitialize(&argc,&argv,(char*)0,help)); 17 CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 18 CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 19 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL)); 20 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL)); 21 22 /* Create a matrix */ 23 CHKERRQ(MatCreate(PETSC_COMM_WORLD, &A)); 24 mloc = PETSC_DECIDE; 25 CHKERRQ(PetscSplitOwnershipEqual(PETSC_COMM_WORLD,&mloc,&M)); 26 nloc = PETSC_DECIDE; 27 CHKERRQ(PetscSplitOwnershipEqual(PETSC_COMM_WORLD,&nloc,&N)); 28 CHKERRQ(MatSetSizes(A,mloc,nloc,M,N)); 29 CHKERRQ(MatSetType(A,MATDENSE)); 30 CHKERRQ(MatSetFromOptions(A)); 31 CHKERRQ(MatSetUp(A)); 32 33 /* Set local matrix entries */ 34 CHKERRQ(MatGetOwnershipIS(A,&isrows,&iscols)); 35 CHKERRQ(ISGetLocalSize(isrows,&nrows)); 36 CHKERRQ(ISGetIndices(isrows,&rows)); 37 CHKERRQ(ISGetLocalSize(iscols,&ncols)); 38 CHKERRQ(ISGetIndices(iscols,&cols)); 39 CHKERRQ(PetscMalloc1(nrows*ncols,&v)); 40 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 CHKERRQ(MatSetValues(A,nrows,rows,ncols,cols,v,INSERT_VALUES)); 51 CHKERRQ(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 52 CHKERRQ(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 53 54 /* Test MatSetValues() by converting A to A_scalapack */ 55 CHKERRQ(MatGetType(A,&type)); 56 if (size == 1) { 57 CHKERRQ(PetscObjectTypeCompare((PetscObject)A,MATSEQDENSE,&isDense)); 58 CHKERRQ(PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&isAIJ)); 59 } else { 60 CHKERRQ(PetscObjectTypeCompare((PetscObject)A,MATMPIDENSE,&isDense)); 61 CHKERRQ(PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&isAIJ)); 62 } 63 64 if (isDense || isAIJ) { 65 Mat Aexplicit; 66 CHKERRQ(MatConvert(A,MATSCALAPACK,MAT_INITIAL_MATRIX,&A_scalapack)); 67 CHKERRQ(MatComputeOperator(A_scalapack,isAIJ?MATAIJ:MATDENSE,&Aexplicit)); 68 CHKERRQ(MatMultEqual(Aexplicit,A_scalapack,5,&flg)); 69 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Aexplicit != A_scalapack."); 70 CHKERRQ(MatDestroy(&Aexplicit)); 71 72 /* Test MAT_REUSE_MATRIX which is only supported for inplace conversion */ 73 CHKERRQ(MatConvert(A,MATSCALAPACK,MAT_INPLACE_MATRIX,&A)); 74 CHKERRQ(MatMultEqual(A_scalapack,A,5,&flg)); 75 PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"A_scalapack != A."); 76 CHKERRQ(MatDestroy(&A_scalapack)); 77 } 78 79 CHKERRQ(ISRestoreIndices(isrows,&rows)); 80 CHKERRQ(ISRestoreIndices(iscols,&cols)); 81 CHKERRQ(ISDestroy(&isrows)); 82 CHKERRQ(ISDestroy(&iscols)); 83 CHKERRQ(PetscFree(v)); 84 CHKERRQ(MatDestroy(&A)); 85 CHKERRQ(PetscFinalize()); 86 return 0; 87 } 88 89 /*TEST 90 91 build: 92 requires: scalapack 93 94 test: 95 nsize: 6 96 97 test: 98 suffix: 2 99 nsize: 6 100 args: -mat_type aij 101 output_file: output/ex243_1.out 102 103 test: 104 suffix: 3 105 nsize: 6 106 args: -mat_type scalapack 107 output_file: output/ex243_1.out 108 109 TEST*/ 110