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