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 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr); 26 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 27 28 /* Creat a matrix */ 29 ierr = PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL);CHKERRQ(ierr); 30 ierr = PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL);CHKERRQ(ierr); 31 ierr = MatCreate(PETSC_COMM_WORLD, &A);CHKERRQ(ierr); 32 ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);CHKERRQ(ierr); 33 ierr = MatSetType(A,MATDENSE);CHKERRQ(ierr); 34 ierr = MatSetFromOptions(A);CHKERRQ(ierr); 35 ierr = MatSetUp(A);CHKERRQ(ierr); 36 37 /* Set local matrix entries */ 38 ierr = MatGetOwnershipIS(A,&isrows,&iscols);CHKERRQ(ierr); 39 ierr = ISGetLocalSize(isrows,&nrows);CHKERRQ(ierr); 40 ierr = ISGetIndices(isrows,&rows);CHKERRQ(ierr); 41 ierr = ISGetLocalSize(iscols,&ncols);CHKERRQ(ierr); 42 ierr = ISGetIndices(iscols,&cols);CHKERRQ(ierr); 43 ierr = PetscMalloc1(nrows*ncols,&v);CHKERRQ(ierr); 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 ierr = MatSetValues(A,nrows,rows,ncols,cols,v,INSERT_VALUES);CHKERRQ(ierr); 55 ierr = MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 56 ierr = MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 57 //ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%D] local nrows %D, ncols %D\n",rank,nrows,ncols);CHKERRQ(ierr); 58 //ierr = PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);CHKERRQ(ierr); 59 60 /* Test MatSetValues() by converting A to A_elemental */ 61 ierr = MatGetType(A,&type);CHKERRQ(ierr); 62 if (size == 1) { 63 ierr = PetscObjectTypeCompare((PetscObject)A,MATSEQDENSE,&isDense);CHKERRQ(ierr); 64 ierr = PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&isAIJ);CHKERRQ(ierr); 65 } else { 66 ierr = PetscObjectTypeCompare((PetscObject)A,MATMPIDENSE,&isDense);CHKERRQ(ierr); 67 ierr = PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&isAIJ);CHKERRQ(ierr); 68 } 69 70 if (isDense || isAIJ) { 71 Mat Aexplicit; 72 ierr = MatConvert(A, MATELEMENTAL, MAT_INITIAL_MATRIX, &A_elemental);CHKERRQ(ierr); 73 ierr = MatComputeOperator(A_elemental,isAIJ ? MATAIJ : MATDENSE,&Aexplicit);CHKERRQ(ierr); 74 ierr = MatMultEqual(Aexplicit,A_elemental,5,&flg);CHKERRQ(ierr); 75 if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Aexplicit != A_elemental."); 76 ierr = MatDestroy(&Aexplicit);CHKERRQ(ierr); 77 78 /* Test MAT_REUSE_MATRIX which is only supported for inplace conversion */ 79 ierr = MatConvert(A, MATELEMENTAL, MAT_INPLACE_MATRIX, &A);CHKERRQ(ierr); 80 ierr = MatMultEqual(A_elemental,A,5,&flg);CHKERRQ(ierr); 81 if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"A_elemental != A."); 82 ierr = MatDestroy(&A_elemental);CHKERRQ(ierr); 83 } 84 85 ierr = ISRestoreIndices(isrows,&rows);CHKERRQ(ierr); 86 ierr = ISRestoreIndices(iscols,&cols);CHKERRQ(ierr); 87 ierr = ISDestroy(&isrows);CHKERRQ(ierr); 88 ierr = ISDestroy(&iscols);CHKERRQ(ierr); 89 ierr = PetscFree(v);CHKERRQ(ierr); 90 ierr = MatDestroy(&A);CHKERRQ(ierr); 91 ierr = PetscFinalize(); 92 return ierr; 93 } 94 95 96 /*TEST 97 98 build: 99 requires: elemental 100 101 test: 102 nsize: 6 103 104 test: 105 suffix: 2 106 nsize: 6 107 args: -mat_type aij 108 output_file: output/ex103_1.out 109 110 test: 111 suffix: 3 112 nsize: 6 113 args: -mat_type elemental 114 output_file: output/ex103_1.out 115 116 TEST*/ 117