Actual source code: ex52f.F90

  1: !
  2:       program main
  3: #include <petsc/finclude/petscksp.h>
  4:       use petscksp
  5:       implicit none
  6: !
  7: !  Demonstrates using MatFactorGetError() and MatFactorGetErrorZeroPivot()
  8: !

 10:       PetscErrorCode  ierr
 11:       PetscInt m,n,one,row,col
 12:       Vec              x,b
 13:       Mat              A,F
 14:       KSP              ksp
 15:       PetscScalar two,zero
 16:       KSPConvergedReason reason
 17:       PCFailedReason pcreason
 18:       PC pc
 19:       MatFactorError ferr
 20:       PetscReal pivot

 22:       PetscCallA(PetscInitialize(ierr))
 23:       m = 2
 24:       n = 2
 25:       PetscCallA(MatCreate(PETSC_COMM_WORLD,A,ierr))
 26:       PetscCallA(MatSetSizes(A,m,n,m,n,ierr))
 27:       PetscCallA(MatSetType(A, MATSEQAIJ,ierr))
 28:       PetscCallA(MatSetUp(A,ierr))
 29:       row = 0
 30:       col = 0
 31:       two = 2.0
 32:       one = 1
 33:       PetscCallA(MatSetValues(A,one,[row],one,[col],[two],INSERT_VALUES,ierr))
 34:       row = 1
 35:       col = 1
 36:       zero = 0.0
 37:       PetscCallA(MatSetValues(A,one,[row],one,[col],[zero],INSERT_VALUES,ierr))
 38:       PetscCallA(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr))
 39:       PetscCallA(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr))

 41:       PetscCallA(VecCreate(PETSC_COMM_WORLD,b,ierr))
 42:       PetscCallA(VecSetSizes(b,m,m,ierr))
 43:       PetscCallA(VecSetType(b,VECSEQ,ierr))

 45: ! Set up solution
 46:       PetscCallA(VecDuplicate(b,x,ierr))

 48: ! Solve system
 49:       PetscCallA(KSPCreate(PETSC_COMM_WORLD,ksp,ierr))
 50:       PetscCallA(KSPSetOperators(ksp,A,A,ierr))
 51:       PetscCallA(KSPSetFromOptions(ksp,ierr))
 52:       PetscCallA(KSPSolve(ksp,b,x,ierr))
 53:       PetscCallA(KSPGetConvergedReason(ksp,reason,ierr))
 54:       PetscCallA(KSPGetPC(ksp,pc,ierr))
 55:       PetscCallA(PCGetFailedReason(pc,pcreason,ierr))
 56:       PetscCallA(PCFactorGetMatrix(pc,F,ierr))
 57:       PetscCallA(MatFactorGetError(F,ferr,ierr))
 58:       PetscCallA(MatFactorGetErrorZeroPivot(F,pivot,row,ierr))
 59:       write(6,101) ferr,pivot,row
 60:  101  format('MatFactorError ',i4,' Pivot value ',1pe9.2,' row ',i4)

 62: ! Cleanup
 63:       PetscCallA(KSPDestroy(ksp,ierr))
 64:       PetscCallA(VecDestroy(b,ierr))
 65:       PetscCallA(VecDestroy(x,ierr))
 66:       PetscCallA(MatDestroy(A,ierr))

 68:       PetscCallA(PetscFinalize(ierr))
 69:       end

 71: !  Nag compiler automatically turns on catching of floating point exceptions and prints message at
 72: !  end of run about the exceptions seen
 73: !
 74: !/*TEST
 75: !
 76: !   test:
 77: !     args: -fp_trap 0
 78: !     filter: grep -v "Warning: Floating"
 79: !
 80: !TEST*/