Actual source code: ex26.c


  2: static char help[] = "Tests MatGetRowIJ for SeqAIJ, SeqBAIJ and SeqSBAIJ\n\n";

  4: #include <petscmat.h>

  6: PetscErrorCode DumpCSR(Mat A,PetscInt shift,PetscBool symmetric,PetscBool compressed)
  7: {
  8:   MatType        type;
  9:   PetscInt       i,j,nr,bs = 1;
 10:   const PetscInt *ia,*ja;
 11:   PetscBool      done,isseqbaij,isseqsbaij;

 14:   PetscObjectBaseTypeCompare((PetscObject)A,MATSEQBAIJ,&isseqbaij);
 15:   PetscObjectBaseTypeCompare((PetscObject)A,MATSEQSBAIJ,&isseqsbaij);
 16:   if (isseqbaij || isseqsbaij) {
 17:     MatGetBlockSize(A,&bs);
 18:   }
 19:   MatGetType(A,&type);
 20:   MatGetRowIJ(A,shift,symmetric,compressed,&nr,&ia,&ja,&done);
 21:   PetscPrintf(PETSC_COMM_SELF,"===========================================================\n");
 22:   PetscPrintf(PETSC_COMM_SELF,"CSR for %s: shift %" PetscInt_FMT " symmetric %" PetscInt_FMT " compressed %" PetscInt_FMT "\n",type,shift,(PetscInt)symmetric,(PetscInt)compressed);
 23:   for (i=0;i<nr;i++) {
 24:     PetscPrintf(PETSC_COMM_SELF,"%" PetscInt_FMT ":",i+shift);
 25:     for (j=ia[i];j<ia[i+1];j++) {
 26:       PetscPrintf(PETSC_COMM_SELF," %" PetscInt_FMT,ja[j-shift]);
 27:     }
 28:     PetscPrintf(PETSC_COMM_SELF,"\n");
 29:   }
 30:   MatRestoreRowIJ(A,shift,symmetric,compressed,&nr,&ia,&ja,&done);
 31:   return 0;
 32: }

 34: int main(int argc,char **args)
 35: {
 36:   Mat            A,B,C;
 37:   PetscInt       i,j,k,m = 3,n = 3,bs = 1;
 38:   PetscMPIInt    size;

 40:   PetscInitialize(&argc,&args,(char*)0,help);
 41:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 43:   PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);
 44:   PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
 45:   PetscOptionsGetInt(NULL,NULL,"-bs",&bs,NULL);
 46:   /* adjust sizes by block size */
 47:   if (m%bs) m += bs-m%bs;
 48:   if (n%bs) n += bs-n%bs;

 50:   MatCreate(PETSC_COMM_SELF,&A);
 51:   MatSetSizes(A,m*n,m*n,PETSC_DECIDE,PETSC_DECIDE);
 52:   MatSetBlockSize(A,bs);
 53:   MatSetType(A,MATSEQAIJ);
 54:   MatSetUp(A);
 55:   MatCreate(PETSC_COMM_SELF,&B);
 56:   MatSetSizes(B,m*n,m*n,PETSC_DECIDE,PETSC_DECIDE);
 57:   MatSetBlockSize(B,bs);
 58:   MatSetType(B,MATSEQBAIJ);
 59:   MatSetUp(B);
 60:   MatCreate(PETSC_COMM_SELF,&C);
 61:   MatSetSizes(C,m*n,m*n,PETSC_DECIDE,PETSC_DECIDE);
 62:   MatSetBlockSize(C,bs);
 63:   MatSetType(C,MATSEQSBAIJ);
 64:   MatSetUp(C);
 65:   MatSetOption(C,MAT_IGNORE_LOWER_TRIANGULAR,PETSC_TRUE);

 67:   for (i=0; i<m; i++) {
 68:     for (j=0; j<n; j++) {

 70:       PetscScalar v = -1.0;
 71:       PetscInt    Ii = j + n*i,J;
 72:       J = Ii - n;
 73:       if (J>=0)  {
 74:         MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
 75:         MatSetValues(B,1,&Ii,1,&J,&v,INSERT_VALUES);
 76:         MatSetValues(C,1,&Ii,1,&J,&v,INSERT_VALUES);
 77:       }
 78:       J = Ii + n;
 79:       if (J<m*n) {
 80:         MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
 81:         MatSetValues(B,1,&Ii,1,&J,&v,INSERT_VALUES);
 82:         MatSetValues(C,1,&Ii,1,&J,&v,INSERT_VALUES);
 83:       }
 84:       J = Ii - 1;
 85:       if (J>=0)  {
 86:         MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
 87:         MatSetValues(B,1,&Ii,1,&J,&v,INSERT_VALUES);
 88:         MatSetValues(C,1,&Ii,1,&J,&v,INSERT_VALUES);
 89:       }
 90:       J = Ii + 1;
 91:       if (J<m*n) {
 92:         MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);
 93:         MatSetValues(B,1,&Ii,1,&J,&v,INSERT_VALUES);
 94:         MatSetValues(C,1,&Ii,1,&J,&v,INSERT_VALUES);
 95:       }
 96:       v = 4.0;
 97:       MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
 98:       MatSetValues(B,1,&Ii,1,&Ii,&v,INSERT_VALUES);
 99:       MatSetValues(C,1,&Ii,1,&Ii,&v,INSERT_VALUES);
100:     }
101:   }
102:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
103:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
104:   MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
105:   MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
106:   MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
107:   MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);

109:   /* test MatGetRowIJ for the three Mat types */
110:   MatView(A,NULL);
111:   MatView(B,NULL);
112:   MatView(C,NULL);
113:   for (i=0;i<2;i++) {
114:     PetscInt shift = i;
115:     for (j=0;j<2;j++) {
116:       PetscBool symmetric = ((j>0) ? PETSC_FALSE : PETSC_TRUE);
117:       for (k=0;k<2;k++) {
118:         PetscBool compressed = ((k>0) ? PETSC_FALSE : PETSC_TRUE);
119:         DumpCSR(A,shift,symmetric,compressed);
120:         DumpCSR(B,shift,symmetric,compressed);
121:         DumpCSR(C,shift,symmetric,compressed);
122:       }
123:     }
124:   }
125:   MatDestroy(&A);
126:   MatDestroy(&B);
127:   MatDestroy(&C);
128:   PetscFinalize();
129:   return 0;
130: }

132: /*TEST

134:    test:

136:    test:
137:       suffix: 2
138:       args: -bs 2

140: TEST*/