Actual source code: ex138.c


  2: static char help[] = "Tests MatGetColumnNorms()/Sums()/Means() for matrix read from file.";

  4: #include <petscmat.h>

  6: int main(int argc,char **args)
  7: {
  8:   Mat            A;
  9:   PetscReal      *reductions_real;
 10:   PetscScalar    *reductions_scalar;
 11:   char           file[PETSC_MAX_PATH_LEN];
 12:   PetscBool      flg;
 13:   PetscViewer    fd;
 14:   PetscInt       n;
 15:   PetscMPIInt    rank;

 17:   PetscInitialize(&argc,&args,(char*)0,help);
 18:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 19:   PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg);
 21:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
 22:   MatCreate(PETSC_COMM_WORLD,&A);
 23:   MatSetFromOptions(A);
 24:   MatLoad(A,fd);
 25:   PetscViewerDestroy(&fd);

 27:   MatGetSize(A,NULL,&n);
 28:   PetscMalloc1(n,&reductions_real);
 29:   PetscMalloc1(n,&reductions_scalar);

 31:   MatGetColumnNorms(A,NORM_2,reductions_real);
 32:   if (rank == 0) {
 33:     PetscPrintf(PETSC_COMM_SELF,"NORM_2:\n");
 34:     PetscRealView(n,reductions_real,PETSC_VIEWER_STDOUT_SELF);
 35:   }

 37:   MatGetColumnNorms(A,NORM_1,reductions_real);
 38:   if (rank == 0) {
 39:     PetscPrintf(PETSC_COMM_SELF,"NORM_1:\n");
 40:     PetscRealView(n,reductions_real,PETSC_VIEWER_STDOUT_SELF);
 41:   }

 43:   MatGetColumnNorms(A,NORM_INFINITY,reductions_real);
 44:   if (rank == 0) {
 45:     PetscPrintf(PETSC_COMM_SELF,"NORM_INFINITY:\n");
 46:     PetscRealView(n,reductions_real,PETSC_VIEWER_STDOUT_SELF);
 47:   }

 49:   MatGetColumnSums(A,reductions_scalar);
 50:   if (!rank) {
 51:     PetscPrintf(PETSC_COMM_SELF,"REDUCTION_SUM:\n");
 52:     PetscScalarView(n,reductions_scalar,PETSC_VIEWER_STDOUT_SELF);
 53:   }

 55:   MatGetColumnMeans(A,reductions_scalar);
 56:   if (!rank) {
 57:     PetscPrintf(PETSC_COMM_SELF,"REDUCTION_MEAN:\n");
 58:     PetscScalarView(n,reductions_scalar,PETSC_VIEWER_STDOUT_SELF);
 59:   }

 61:   PetscFree(reductions_real);
 62:   PetscFree(reductions_scalar);
 63:   MatDestroy(&A);
 64:   PetscFinalize();
 65:   return 0;
 66: }

 68: /*TEST

 70:    test:
 71:       suffix: 1
 72:       nsize: 2
 73:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
 74:       args: -f ${DATAFILESPATH}/matrices/small -mat_type aij
 75:       output_file: output/ex138.out

 77:    test:
 78:       suffix: 2
 79:       nsize: {{1 2}}
 80:       requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
 81:       args: -f ${DATAFILESPATH}/matrices/small -mat_type baij -matload_block_size {{2 3}}
 82:       output_file: output/ex138.out

 84:    test:
 85:       suffix: complex
 86:       nsize: 2
 87:       requires: datafilespath complex double !defined(PETSC_USE_64BIT_INDICES)
 88:       args: -f ${DATAFILESPATH}/matrices/nimrod/small_112905 -mat_type aij
 89:       output_file: output/ex138_complex.out
 90:       filter: grep -E "\ 0:|1340:|1344:"

 92: TEST*/