1 static char help[] = "Test MatSetValues() by converting MATDENSE to MATELEMENTAL. \n\ 2 Modified from the code contributed by Yaning Liu @lbl.gov \n\n"; 3 /* 4 Example: 5 mpiexec -n <np> ./ex103 6 mpiexec -n <np> ./ex103 -mat_type elemental -mat_view 7 mpiexec -n <np> ./ex103 -mat_type aij 8 */ 9 10 #include <petscmat.h> 11 12 int main(int argc, char** argv) 13 { 14 Mat A,A_elemental; 15 PetscInt i,j,M=10,N=5,nrows,ncols; 16 PetscErrorCode ierr; 17 PetscMPIInt rank,size; 18 IS isrows,iscols; 19 const PetscInt *rows,*cols; 20 PetscScalar *v; 21 MatType type; 22 PetscBool isDense,isAIJ,flg; 23 24 ierr = PetscInitialize(&argc, &argv, (char*)0, help);if (ierr) return ierr; 25 CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 26 CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 27 28 /* Creat a matrix */ 29 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL)); 30 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL)); 31 CHKERRQ(MatCreate(PETSC_COMM_WORLD, &A)); 32 CHKERRQ(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N)); 33 CHKERRQ(MatSetType(A,MATDENSE)); 34 CHKERRQ(MatSetFromOptions(A)); 35 CHKERRQ(MatSetUp(A)); 36 37 /* Set local matrix entries */ 38 CHKERRQ(MatGetOwnershipIS(A,&isrows,&iscols)); 39 CHKERRQ(ISGetLocalSize(isrows,&nrows)); 40 CHKERRQ(ISGetIndices(isrows,&rows)); 41 CHKERRQ(ISGetLocalSize(iscols,&ncols)); 42 CHKERRQ(ISGetIndices(iscols,&cols)); 43 CHKERRQ(PetscMalloc1(nrows*ncols,&v)); 44 45 for (i=0; i<nrows; i++) { 46 for (j=0; j<ncols; j++) { 47 if (size == 1) { 48 v[i*ncols+j] = (PetscScalar)(i+j); 49 } else { 50 v[i*ncols+j] = (PetscScalar)rank+j*0.1; 51 } 52 } 53 } 54 CHKERRQ(MatSetValues(A,nrows,rows,ncols,cols,v,INSERT_VALUES)); 55 CHKERRQ(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 56 CHKERRQ(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 57 //CHKERRQ(PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%" PetscInt_FMT "] local nrows %" PetscInt_FMT ", ncols %" PetscInt_FMT "\n",rank,nrows,ncols)); 58 //CHKERRQ(PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT)); 59 60 /* Test MatSetValues() by converting A to A_elemental */ 61 CHKERRQ(MatGetType(A,&type)); 62 if (size == 1) { 63 CHKERRQ(PetscObjectTypeCompare((PetscObject)A,MATSEQDENSE,&isDense)); 64 CHKERRQ(PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&isAIJ)); 65 } else { 66 CHKERRQ(PetscObjectTypeCompare((PetscObject)A,MATMPIDENSE,&isDense)); 67 CHKERRQ(PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&isAIJ)); 68 } 69 70 if (isDense || isAIJ) { 71 Mat Aexplicit; 72 CHKERRQ(MatConvert(A, MATELEMENTAL, MAT_INITIAL_MATRIX, &A_elemental)); 73 CHKERRQ(MatComputeOperator(A_elemental,isAIJ ? MATAIJ : MATDENSE,&Aexplicit)); 74 CHKERRQ(MatMultEqual(Aexplicit,A_elemental,5,&flg)); 75 PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Aexplicit != A_elemental."); 76 CHKERRQ(MatDestroy(&Aexplicit)); 77 78 /* Test MAT_REUSE_MATRIX which is only supported for inplace conversion */ 79 CHKERRQ(MatConvert(A, MATELEMENTAL, MAT_INPLACE_MATRIX, &A)); 80 CHKERRQ(MatMultEqual(A_elemental,A,5,&flg)); 81 PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"A_elemental != A."); 82 CHKERRQ(MatDestroy(&A_elemental)); 83 } 84 85 CHKERRQ(ISRestoreIndices(isrows,&rows)); 86 CHKERRQ(ISRestoreIndices(iscols,&cols)); 87 CHKERRQ(ISDestroy(&isrows)); 88 CHKERRQ(ISDestroy(&iscols)); 89 CHKERRQ(PetscFree(v)); 90 CHKERRQ(MatDestroy(&A)); 91 ierr = PetscFinalize(); 92 return ierr; 93 } 94 95 /*TEST 96 97 build: 98 requires: elemental 99 100 test: 101 nsize: 6 102 103 test: 104 suffix: 2 105 nsize: 6 106 args: -mat_type aij 107 output_file: output/ex103_1.out 108 109 test: 110 suffix: 3 111 nsize: 6 112 args: -mat_type elemental 113 output_file: output/ex103_1.out 114 115 TEST*/ 116