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

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

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

 37:   PetscCallA(VecCreate(PETSC_COMM_WORLD, b, ierr))
 38:   PetscCallA(VecSetSizes(b, m, m, ierr))
 39:   PetscCallA(VecSetType(b, VECSEQ, ierr))

 41: ! Set up solution
 42:   PetscCallA(VecDuplicate(b, x, ierr))

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

 58: ! Cleanup
 59:   PetscCallA(KSPDestroy(ksp, ierr))
 60:   PetscCallA(VecDestroy(b, ierr))
 61:   PetscCallA(VecDestroy(x, ierr))
 62:   PetscCallA(MatDestroy(A, ierr))

 64:   PetscCallA(PetscFinalize(ierr))
 65: end

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