1 2 static char help[] = "Tests MatTranspose(), MatNorm(), and MatAXPY().\n\n"; 3 4 #include <petscmat.h> 5 6 int main(int argc,char **argv) 7 { 8 Mat mat,tmat = 0; 9 PetscInt m = 4,n,i,j; 10 PetscMPIInt size,rank; 11 PetscInt rstart,rend,rect = 0; 12 PetscBool flg; 13 PetscScalar v; 14 PetscReal normf,normi,norm1; 15 MatInfo info; 16 17 PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 18 PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL)); 19 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 20 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 21 n = m; 22 PetscCall(PetscOptionsHasName(NULL,NULL,"-rect1",&flg)); 23 if (flg) {n += 2; rect = 1;} 24 PetscCall(PetscOptionsHasName(NULL,NULL,"-rect2",&flg)); 25 if (flg) {n -= 2; rect = 1;} 26 27 /* Create and assemble matrix */ 28 PetscCall(MatCreate(PETSC_COMM_WORLD,&mat)); 29 PetscCall(MatSetSizes(mat,PETSC_DECIDE,PETSC_DECIDE,m,n)); 30 PetscCall(MatSetFromOptions(mat)); 31 PetscCall(MatSetUp(mat)); 32 PetscCall(MatGetOwnershipRange(mat,&rstart,&rend)); 33 for (i=rstart; i<rend; i++) { 34 for (j=0; j<n; j++) { 35 v = 10*i+j; 36 PetscCall(MatSetValues(mat,1,&i,1,&j,&v,INSERT_VALUES)); 37 } 38 } 39 PetscCall(MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY)); 40 PetscCall(MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY)); 41 42 /* Print info about original matrix */ 43 PetscCall(MatGetInfo(mat,MAT_GLOBAL_SUM,&info)); 44 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"original matrix nonzeros = %" PetscInt_FMT ", allocated nonzeros = %" PetscInt_FMT "\n",(PetscInt)info.nz_used,(PetscInt)info.nz_allocated)); 45 PetscCall(MatNorm(mat,NORM_FROBENIUS,&normf)); 46 PetscCall(MatNorm(mat,NORM_1,&norm1)); 47 PetscCall(MatNorm(mat,NORM_INFINITY,&normi)); 48 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"original: Frobenious norm = %g, one norm = %g, infinity norm = %g\n",(double)normf,(double)norm1,(double)normi)); 49 PetscCall(MatView(mat,PETSC_VIEWER_STDOUT_WORLD)); 50 51 /* Form matrix transpose */ 52 PetscCall(PetscOptionsHasName(NULL,NULL,"-in_place",&flg)); 53 if (flg) { 54 PetscCall(MatTranspose(mat,MAT_INPLACE_MATRIX,&mat)); /* in-place transpose */ 55 tmat = mat; mat = 0; 56 } else { /* out-of-place transpose */ 57 PetscCall(MatTranspose(mat,MAT_INITIAL_MATRIX,&tmat)); 58 } 59 60 /* Print info about transpose matrix */ 61 PetscCall(MatGetInfo(tmat,MAT_GLOBAL_SUM,&info)); 62 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"transpose matrix nonzeros = %" PetscInt_FMT ", allocated nonzeros = %" PetscInt_FMT "\n",(PetscInt)info.nz_used,(PetscInt)info.nz_allocated)); 63 PetscCall(MatNorm(tmat,NORM_FROBENIUS,&normf)); 64 PetscCall(MatNorm(tmat,NORM_1,&norm1)); 65 PetscCall(MatNorm(tmat,NORM_INFINITY,&normi)); 66 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"transpose: Frobenious norm = %g, one norm = %g, infinity norm = %g\n",(double)normf,(double)norm1,(double)normi)); 67 PetscCall(MatView(tmat,PETSC_VIEWER_STDOUT_WORLD)); 68 69 /* Test MatAXPY */ 70 if (mat && !rect) { 71 PetscScalar alpha = 1.0; 72 PetscCall(PetscOptionsGetScalar(NULL,NULL,"-alpha",&alpha,NULL)); 73 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"matrix addition: B = B + alpha * A\n")); 74 PetscCall(MatAXPY(tmat,alpha,mat,DIFFERENT_NONZERO_PATTERN)); 75 PetscCall(MatView(tmat,PETSC_VIEWER_STDOUT_WORLD)); 76 } 77 78 /* Free data structures */ 79 PetscCall(MatDestroy(&tmat)); 80 if (mat) PetscCall(MatDestroy(&mat)); 81 82 PetscCall(PetscFinalize()); 83 return 0; 84 } 85 86 /*TEST 87 88 test: 89 90 testset: 91 args: -rect1 92 test: 93 suffix: r1 94 output_file: output/ex49_r1.out 95 test: 96 suffix: r1_inplace 97 args: -in_place 98 output_file: output/ex49_r1.out 99 test: 100 suffix: r1_par 101 nsize: 2 102 output_file: output/ex49_r1_par.out 103 test: 104 suffix: r1_par_inplace 105 args: -in_place 106 nsize: 2 107 output_file: output/ex49_r1_par.out 108 109 TEST*/ 110