Actual source code: ex102.c


  2: static char help[] = "Tests MatCreateLRC()\n\n";

  4: #include <petscmat.h>

  6: int main(int argc,char **args)
  7: {
  8:   Vec            x,b,c=NULL;
  9:   Mat            A,U,V,LR,X,LRe;
 10:   PetscInt       M = 5, N = 7;
 11:   PetscBool      flg;

 13:   PetscInitialize(&argc,&args,(char*)0,help);
 14:   PetscOptionsGetInt(NULL,NULL,"-m",&M,NULL);
 15:   PetscOptionsGetInt(NULL,NULL,"-n",&N,NULL);

 17:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 18:          Create the sparse matrix
 19:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 20:   MatCreate(PETSC_COMM_WORLD,&A);
 21:   MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);
 22:   MatSetOptionsPrefix(A,"A_");
 23:   MatSetFromOptions(A);
 24:   MatSeqAIJSetPreallocation(A,5,NULL);
 25:   MatMPIAIJSetPreallocation(A,5,NULL,5,NULL);
 26:   MatSetUp(A);
 27:   MatSetRandom(A,NULL);

 29:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 30:          Create the dense matrices
 31:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 32:   MatCreate(PETSC_COMM_WORLD,&U);
 33:   MatSetSizes(U,PETSC_DECIDE,PETSC_DECIDE,M,3);
 34:   MatSetType(U,MATDENSE);
 35:   MatSetOptionsPrefix(U,"U_");
 36:   MatSetFromOptions(U);
 37:   MatSetUp(U);
 38:   MatSetRandom(U,NULL);

 40:   MatCreate(PETSC_COMM_WORLD,&V);
 41:   MatSetSizes(V,PETSC_DECIDE,PETSC_DECIDE,N,3);
 42:   MatSetType(V,MATDENSE);
 43:   MatSetOptionsPrefix(V,"V_");
 44:   MatSetFromOptions(V);
 45:   MatSetUp(V);
 46:   MatSetRandom(V,NULL);

 48:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 49:          Create a vector to hold the diagonal of C
 50:          A sequential vector can be created as well on each process
 51:          It is user responsibility to ensure the data in the vector
 52:          is consistent across processors
 53:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 54:   PetscOptionsHasName(NULL,NULL,"-use_c",&flg);
 55:   if (flg) {
 56:     VecCreateMPI(PETSC_COMM_WORLD,PETSC_DECIDE,3,&c);
 57:     VecSetRandom(c,NULL);
 58:   }

 60:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 61:          Create low rank correction matrix
 62:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 63:   PetscOptionsHasName(NULL,NULL,"-low_rank",&flg);
 64:   if (flg) {
 65:     /* create a low-rank matrix, with no A-matrix */
 66:     MatCreateLRC(NULL,U,c,V,&LR);
 67:     MatDestroy(&A);
 68:   } else {
 69:     MatCreateLRC(A,U,c,V,&LR);
 70:   }

 72:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 73:          Create the low rank correction matrix explicitly to check for
 74:          correctness
 75:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 76:   MatHermitianTranspose(V,MAT_INITIAL_MATRIX,&X);
 77:   MatDiagonalScale(X,c,NULL);
 78:   MatMatMult(U,X,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&LRe);
 79:   MatDestroy(&X);
 80:   if (A) {
 81:     MatAYPX(LRe,1.0,A,DIFFERENT_NONZERO_PATTERN);
 82:   }

 84:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 85:          Create test vectors
 86:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 87:   MatCreateVecs(LR,&x,&b);
 88:   VecSetRandom(x,NULL);
 89:   MatMult(LR,x,b);
 90:   MatMultTranspose(LR,b,x);
 91:   VecDestroy(&x);
 92:   VecDestroy(&b);

 94:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 95:          Check correctness
 96:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 97:   MatMultEqual(LR,LRe,10,&flg);
 98:   if (!flg) PetscPrintf(PETSC_COMM_WORLD,"Error in MatMult\n");
 99: #if !defined(PETSC_USE_COMPLEX)
100:   MatMultHermitianTransposeEqual(LR,LRe,10,&flg);
101:   if (!flg) PetscPrintf(PETSC_COMM_WORLD,"Error in MatMultTranspose\n");
102: #endif

104:   MatDestroy(&A);
105:   MatDestroy(&LRe);
106:   MatDestroy(&U);
107:   MatDestroy(&V);
108:   VecDestroy(&c);
109:   MatDestroy(&LR);

111:   /*
112:      Always call PetscFinalize() before exiting a program.  This routine
113:        - finalizes the PETSc libraries as well as MPI
114:        - provides summary and diagnostic information if certain runtime
115:          options are chosen (e.g., -log_view).
116:   */
117:   PetscFinalize();
118:   return 0;
119: }

121: /*TEST

123:    testset:
124:       output_file: output/ex102_1.out
125:       nsize: {{1 2}}
126:       args: -low_rank {{0 1}} -use_c {{0 1}}
127:       test:
128:         suffix: standard
129:       test:
130:         suffix: cuda
131:         requires: cuda
132:         args: -A_mat_type aijcusparse -U_mat_type densecuda -V_mat_type densecuda

134: TEST*/