Actual source code: ex51.c

petsc-3.12.2 2019-11-22
Report Typos and Errors

  2: static char help[] = "Test PCFailedReason.\n\n";

  4:  #include <petscksp.h>

  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;

 19:   PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
 20:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 21:   if (size != 1) SETERRQ(PETSC_COMM_WORLD,1,"This is a uniprocessor example only!");
 22:   PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
 23:   PetscOptionsGetBool(NULL,NULL,"-reverse",&reverse,NULL);

 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;

 29:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 30:          Create the matrix 
 31:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 32:   MatCreate(PETSC_COMM_WORLD,&A);
 33:   MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
 34:   MatSetFromOptions(A);
 35:   MatSetUp(A);

 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:     MatSetValues(A,1,&i,3,col,value,INSERT_VALUES);
 41:   }
 42:   i    = n - 1; col[0] = n - 2; col[1] = n - 1;
 43:   MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);
 44:   i    = 0; col[0] = 0; col[1] = 1; value[0] = 2.0; value[1] = -1.0;
 45:   MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);
 46:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 47:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);

 49:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 50:                 Create the linear solver and set various options
 51:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 52:   KSPCreate(PETSC_COMM_WORLD,&ksp);
 53:   KSPSetOperators(ksp,A,A);
 54:   MatShift(A,reverse?-alpha:-beta);
 55:   KSPGetPC(ksp,&pc);
 56:   PCSetType(pc,PCLU);
 57:   KSPSetFromOptions(ksp);

 59:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 60:                       Factorize first matrix
 61:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 62:   PetscPrintf(PETSC_COMM_WORLD,"First matrix\n");
 63:   KSPSetUp(ksp);
 64:   KSPGetConvergedReason(ksp,&reason);
 65:   if (reason) {
 66:     PetscPrintf(PETSC_COMM_WORLD,"KSPSetUp() failed due to %s\n",KSPConvergedReasons[reason]);
 67:     PCGetFailedReason(pc,&pcreason);
 68:     PetscPrintf(PETSC_COMM_WORLD,"PC reason is %s\n",PCFailedReasons[pcreason]);
 69:   } else {
 70:     PetscPrintf(PETSC_COMM_WORLD,"Success!\n");
 71:   }

 73:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 74:                       Factorize second matrix
 75:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 76:   MatShift(A,reverse?alpha-beta:beta-alpha);
 77:   KSPSetOperators(ksp,A,A);

 79:   PetscPrintf(PETSC_COMM_WORLD,"Second matrix\n");
 80:   KSPSetUp(ksp);
 81:   KSPGetConvergedReason(ksp,&reason);
 82:   if (reason) {
 83:     PetscPrintf(PETSC_COMM_WORLD,"KSPSetUp() failed due to %s\n",KSPConvergedReasons[reason]);
 84:     PCGetFailedReason(pc,&pcreason);
 85:     PetscPrintf(PETSC_COMM_WORLD,"PC reason is %s\n",PCFailedReasons[pcreason]);
 86:   } else {
 87:     PetscPrintf(PETSC_COMM_WORLD,"Success!\n");
 88:     PCGetFailedReason(pc,&pcreason);
 89:     PetscPrintf(PETSC_COMM_WORLD,"PC reason is %s\n",PCFailedReasons[pcreason]);
 90:   }

 92:   /*
 93:      Free work space.
 94:   */
 95:   MatDestroy(&A);
 96:   KSPDestroy(&ksp);

 98:   PetscFinalize();
 99:   return ierr;
100: }


103: /*TEST

105:    test:
106:       args: -reverse

108:    test:
109:       suffix: 2
110:       args: -reverse -pc_type cholesky

112: TEST*/