2: static char help[] = "Tests MatTranspose(), MatNorm(), and MatAXPY().\n\n";

  4: #include <petscmat.h>

  8: int main(int argc,char **argv)
  9: {
 10:   Mat            mat,tmat = 0;
 11:   PetscInt       m = 4,n,i,j;
 13:   PetscMPIInt    size,rank;
 14:   PetscInt       rstart,rend,rect = 0;
 15:   PetscBool      flg;
 16:   PetscScalar    v;
 17:   PetscReal      normf,normi,norm1;
 18:   MatInfo        info;

 20:   PetscInitialize(&argc,&argv,(char*)0,help);
 21:   PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);
 22:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 23:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 24:   n    = m;
 25:   PetscOptionsHasName(NULL,NULL,"-rect1",&flg);
 26:   if (flg) {n += 2; rect = 1;}
 27:   PetscOptionsHasName(NULL,NULL,"-rect2",&flg);
 28:   if (flg) {n -= 2; rect = 1;}

 30:   /* Create and assemble matrix */
 31:   MatCreate(PETSC_COMM_WORLD,&mat);
 32:   MatSetSizes(mat,PETSC_DECIDE,PETSC_DECIDE,m,n);
 33:   MatSetFromOptions(mat);
 34:   MatSetUp(mat);
 35:   MatGetOwnershipRange(mat,&rstart,&rend);
 36:   for (i=rstart; i<rend; i++) {
 37:     for (j=0; j<n; j++) {
 38:       v    = 10*i+j;
 39:       MatSetValues(mat,1,&i,1,&j,&v,INSERT_VALUES);
 40:     }
 41:   }
 42:   MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
 43:   MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);

 45:   /* Print info about original matrix */
 46:   MatGetInfo(mat,MAT_GLOBAL_SUM,&info);
 47:   PetscPrintf(PETSC_COMM_WORLD,"original matrix nonzeros = %D, allocated nonzeros = %D\n",
 48:                      (PetscInt)info.nz_used,(PetscInt)info.nz_allocated);
 49:   MatNorm(mat,NORM_FROBENIUS,&normf);
 50:   MatNorm(mat,NORM_1,&norm1);
 51:   MatNorm(mat,NORM_INFINITY,&normi);
 52:   PetscPrintf(PETSC_COMM_WORLD,"original: Frobenious norm = %g, one norm = %g, infinity norm = %g\n",(double)normf,(double)norm1,(double)normi);
 53:   MatView(mat,PETSC_VIEWER_STDOUT_WORLD);

 55:   /* Form matrix transpose */
 56:   PetscOptionsHasName(NULL,NULL,"-in_place",&flg);
 57:   if (!rect && flg) {
 58:     MatTranspose(mat,MAT_REUSE_MATRIX,&mat);   /* in-place transpose */
 59:     tmat = mat; mat = 0;
 60:   } else {      /* out-of-place transpose */
 61:     MatTranspose(mat,MAT_INITIAL_MATRIX,&tmat);
 62:   }

 64:   /* Print info about transpose matrix */
 65:   MatGetInfo(tmat,MAT_GLOBAL_SUM,&info);
 66:   PetscPrintf(PETSC_COMM_WORLD,"transpose matrix nonzeros = %D, allocated nonzeros = %D\n",
 67:                      (PetscInt)info.nz_used,(PetscInt)info.nz_allocated);
 68:   MatNorm(tmat,NORM_FROBENIUS,&normf);
 69:   MatNorm(tmat,NORM_1,&norm1);
 70:   MatNorm(tmat,NORM_INFINITY,&normi);
 71:   PetscPrintf(PETSC_COMM_WORLD,"transpose: Frobenious norm = %g, one norm = %g, infinity norm = %g\n",(double)normf,(double)norm1,(double)normi);
 72:   MatView(tmat,PETSC_VIEWER_STDOUT_WORLD);


 75:   /* Test MatAXPY */
 76:   if (mat && !rect) {
 77:     PetscScalar alpha = 1.0;
 78:     PetscOptionsGetScalar(NULL,NULL,"-alpha",&alpha,NULL);
 79:     PetscPrintf(PETSC_COMM_WORLD,"matrix addition:  B = B + alpha * A\n");
 80:     MatAXPY(tmat,alpha,mat,DIFFERENT_NONZERO_PATTERN);
 81:     MatView(tmat,PETSC_VIEWER_STDOUT_WORLD);
 82:   }

 84:   /* Free data structures */
 85:   MatDestroy(&tmat);
 86:   if (mat) {MatDestroy(&mat);}

 88:   PetscFinalize();
 89:   return 0;
 90: }