1!
2      program main
3#include <petsc/finclude/petscksp.h>
4      use petscksp
5      implicit none
6!
7!  Demonstrates using MatFactorGetError() and MatFactorGetErrorZeroPivot()
8!
9
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
21
22      call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
23      if (ierr .ne. 0) then
24        print*,'Unable to initialize PETSc'
25        stop
26      endif
27
28      m = 2
29      n = 2
30      call MatCreate(PETSC_COMM_WORLD,A,ierr)
31      call MatSetSizes(A,m,n,m,n,ierr)
32      call MatSetType(A, MATSEQAIJ,ierr)
33      call MatSetUp(A,ierr)
34      row = 0
35      col = 0
36      two = 2.0
37      one = 1
38      call MatSetValues(A,one,row,one,col,two,INSERT_VALUES,ierr)
39      row = 1
40      col = 1
41      zero = 0.0
42      call MatSetValues(A,one,row,one,col,zero,INSERT_VALUES,ierr)
43      call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)
44      call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)
45
46      call VecCreate(PETSC_COMM_WORLD,b,ierr)
47      call VecSetSizes(b,m,m,ierr)
48      call VecSetType(b,VECSEQ,ierr)
49
50! Set up solution
51      call VecDuplicate(b,x,ierr)
52
53! Solve system
54      call KSPCreate(PETSC_COMM_WORLD,ksp,ierr)
55      call KSPSetOperators(ksp,A,A,ierr)
56      call KSPSetFromOptions(ksp,ierr)
57      call KSPSolve(ksp,b,x,ierr)
58      call KSPGetConvergedReason(ksp,reason,ierr)
59      call KSPGetPC(ksp,pc,ierr)
60      call PCGetFailedReason(pc,pcreason,ierr)
61      call PCFactorGetMatrix(pc,F,ierr)
62      call MatFactorGetError(F,ferr,ierr)
63      call MatFactorGetErrorZeroPivot(F,pivot,row,ierr)
64      write(6,101) ferr,pivot,row
65 101  format('MatFactorError ',i4,' Pivot value ',1pe9.2,' row ',i4)
66
67! Cleanup
68      call KSPDestroy(ksp,ierr)
69      call VecDestroy(b,ierr)
70      call VecDestroy(x,ierr)
71      call MatDestroy(A,ierr)
72
73      call PetscFinalize(ierr)
74      end
75
76!  Nag compiler automatically turns on catching of floating point exceptions and prints message at
77!  end of run about the exceptions seen
78!
79!/*TEST
80!
81!   test:
82!     args: -fp_trap 0
83!     filter: grep -v "Warning: Floating"
84!
85!TEST*/
86