1: #include "petsc/finclude/petscksp.h"
  2: program main
  3:   use petscksp
  4:   implicit none

  6:   PetscInt N
  7:   PetscBool draw, flg
  8:   PetscReal rnorm
  9:   PetscScalar, parameter :: one = 1.0, mone = -1.0
 10:   Mat A
 11:   Vec b, x, r
 12:   KSP ksp
 13:   PC pc
 14:   PetscErrorCode ierr
 15:   character(80) mattype

 17:   PetscCallA(PetscInitialize(ierr))
 18:   PetscCallA(PetscPythonInitialize(PETSC_NULL_CHARACTER, PETSC_NULL_CHARACTER, ierr))

 20:   N = 100
 21:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-N', N, flg, ierr))
 22:   draw = .false.
 23:   PetscCallA(PetscOptionsGetBool(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-draw', draw, flg, ierr))

 25:   PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr))
 26:   PetscCallA(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, N, N, ierr))
 27:   PetscCallA(MatSetType(A, 'python', ierr))
 28:   PetscCallA(MatPythonSetType(A, 'example100.py:Laplace1D', ierr))
 29:   PetscCallA(MatPythonGetType(A, mattype, ierr))
 30:   PetscCheckA(mattype == 'example100.py:Laplace1D', PETSC_COMM_WORLD, PETSC_ERR_PLIB, 'Error')
 31:   PetscCallA(MatSetUp(A, ierr))

 33:   PetscCallA(MatCreateVecs(A, x, b, ierr))
 34:   PetscCallA(VecSet(b, one, ierr))

 36:   PetscCallA(KSPCreate(PETSC_COMM_WORLD, ksp, ierr))
 37:   PetscCallA(KSPSetType(ksp, 'python', ierr))
 38:   PetscCallA(KSPPythonSetType(ksp, 'example100.py:ConjGrad', ierr))

 40:   PetscCallA(KSPGetPC(ksp, pc, ierr))
 41:   PetscCallA(PCSetType(pc, 'python', ierr))
 42:   PetscCallA(PCPythonSetType(pc, 'example100.py:Jacobi', ierr))

 44:   PetscCallA(KSPSetOperators(ksp, A, A, ierr))
 45:   PetscCallA(KSPSetFromOptions(ksp, ierr))
 46:   PetscCallA(KSPSolve(ksp, b, x, ierr))

 48:   PetscCallA(VecDuplicate(b, r, ierr))
 49:   PetscCallA(MatMult(A, x, r, ierr))
 50:   PetscCallA(VecAYPX(r, mone, b, ierr))
 51:   PetscCallA(VecNorm(r, NORM_2, rnorm, ierr))
 52:   print *, 'error norm = ', rnorm

 54:   if (draw) then
 55:     PetscCallA(VecView(x, PETSC_VIEWER_DRAW_WORLD, ierr))
 56:     PetscCallA(PetscSleep(2.0_PETSC_REAL_KIND, ierr))
 57:   end if

 59:   PetscCallA(VecDestroy(x, ierr))
 60:   PetscCallA(VecDestroy(b, ierr))
 61:   PetscCallA(VecDestroy(r, ierr))
 62:   PetscCallA(MatDestroy(A, ierr))
 63:   PetscCallA(KSPDestroy(ksp, ierr))

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

 68: !/*TEST
 69: !
 70: !    test:
 71: !      requires: petsc4py
 72: !      localrunfiles: example100.py
 73: !
 74: !TEST*/