1 2 static char help[] = "Test PCFailedReason.\n\n"; 3 4 #include <petscksp.h> 5 6 int main(int argc,char **args) 7 { 8 Mat A; /* linear system matrix */ 9 KSP ksp; /* linear solver context */ 10 PC pc; /* preconditioner context */ 11 PetscErrorCode ierr; 12 PetscInt i,n = 10,col[3]; 13 PetscMPIInt size; 14 PetscScalar value[3],alpha,beta,sx; 15 PetscBool reverse=PETSC_FALSE; 16 KSPConvergedReason reason; 17 PCFailedReason pcreason; 18 19 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 20 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr); 21 if (size != 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"This is a uniprocessor example only!"); 22 ierr = PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);CHKERRQ(ierr); 23 ierr = PetscOptionsGetBool(NULL,NULL,"-reverse",&reverse,NULL);CHKERRQ(ierr); 24 25 sx = PetscSinReal(n*PETSC_PI/2/(n+1)); 26 alpha = 4.0*sx*sx; /* alpha is the largest eigenvalue of the matrix */ 27 beta = 4.0; 28 29 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 30 Create the matrix 31 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 32 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 33 ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);CHKERRQ(ierr); 34 ierr = MatSetFromOptions(A);CHKERRQ(ierr); 35 ierr = MatSetUp(A);CHKERRQ(ierr); 36 37 value[0] = -1.0; value[1] = 2.0; value[2] = -1.0; 38 for (i=1; i<n-1; i++) { 39 col[0] = i-1; col[1] = i; col[2] = i+1; 40 ierr = MatSetValues(A,1,&i,3,col,value,INSERT_VALUES);CHKERRQ(ierr); 41 } 42 i = n - 1; col[0] = n - 2; col[1] = n - 1; 43 ierr = MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);CHKERRQ(ierr); 44 i = 0; col[0] = 0; col[1] = 1; value[0] = 2.0; value[1] = -1.0; 45 ierr = MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);CHKERRQ(ierr); 46 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 47 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 48 49 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 50 Create the linear solver and set various options 51 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 52 ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr); 53 ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr); 54 ierr = MatShift(A,reverse?-alpha:-beta);CHKERRQ(ierr); 55 ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr); 56 ierr = PCSetType(pc,PCLU);CHKERRQ(ierr); 57 ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr); 58 59 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 60 Factorize first matrix 61 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 62 ierr = PetscPrintf(PETSC_COMM_WORLD,"First matrix\n");CHKERRQ(ierr); 63 ierr = KSPSetUp(ksp);CHKERRQ(ierr); 64 ierr = KSPGetConvergedReason(ksp,&reason);CHKERRQ(ierr); 65 if (reason) { 66 ierr = PetscPrintf(PETSC_COMM_WORLD,"KSPSetUp() failed due to %s\n",KSPConvergedReasons[reason]);CHKERRQ(ierr); 67 ierr = PCGetFailedReason(pc,&pcreason);CHKERRQ(ierr); 68 ierr = PetscPrintf(PETSC_COMM_WORLD,"PC reason is %s\n",PCFailedReasons[pcreason]);CHKERRQ(ierr); 69 } else { 70 ierr = PetscPrintf(PETSC_COMM_WORLD,"Success!\n");CHKERRQ(ierr); 71 } 72 73 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 74 Factorize second matrix 75 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 76 ierr = MatShift(A,reverse?alpha-beta:beta-alpha);CHKERRQ(ierr); 77 ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr); 78 79 ierr = PetscPrintf(PETSC_COMM_WORLD,"Second matrix\n");CHKERRQ(ierr); 80 ierr = KSPSetUp(ksp);CHKERRQ(ierr); 81 ierr = KSPGetConvergedReason(ksp,&reason);CHKERRQ(ierr); 82 if (reason) { 83 ierr = PetscPrintf(PETSC_COMM_WORLD,"KSPSetUp() failed due to %s\n",KSPConvergedReasons[reason]);CHKERRQ(ierr); 84 ierr = PCGetFailedReason(pc,&pcreason);CHKERRQ(ierr); 85 ierr = PetscPrintf(PETSC_COMM_WORLD,"PC reason is %s\n",PCFailedReasons[pcreason]);CHKERRQ(ierr); 86 } else { 87 ierr = PetscPrintf(PETSC_COMM_WORLD,"Success!\n");CHKERRQ(ierr); 88 ierr = PCGetFailedReason(pc,&pcreason);CHKERRQ(ierr); 89 ierr = PetscPrintf(PETSC_COMM_WORLD,"PC reason is %s\n",PCFailedReasons[pcreason]);CHKERRQ(ierr); 90 } 91 92 /* 93 Free work space. 94 */ 95 ierr = MatDestroy(&A);CHKERRQ(ierr); 96 ierr = KSPDestroy(&ksp);CHKERRQ(ierr); 97 98 ierr = PetscFinalize(); 99 return ierr; 100 } 101 102 103 /*TEST 104 105 test: 106 args: -reverse 107 108 test: 109 suffix: 2 110 args: -reverse -pc_type cholesky 111 112 TEST*/ 113