Actual source code: ex100f.F90

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

      5:   PetscInt N
      6:   PetscBool draw, flg
      7:   PetscReal rnorm, rtwo
      8:   PetscScalar one, mone
      9:   Mat A
     10:   Vec b, x, r
     11:   KSP ksp
     12:   PC pc
     13:   PetscErrorCode ierr
     14:   character*80 mattype

     16:   N = 100
     17:   draw = .FALSE.
     18:   one = 1.0
     19:   mone = -1.0
     20:   rtwo = 2.0

     22:   PetscCallA(PetscInitialize(ierr))
     23:   PetscCallA(PetscPythonInitialize(PETSC_NULL_CHARACTER, PETSC_NULL_CHARACTER, ierr))

     25:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-N', N, flg, ierr))
     26:   PetscCallA(PetscOptionsGetBool(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-draw', draw, flg, ierr))

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

     36:   PetscCallA(MatCreateVecs(A, x, b, ierr))
     37:   PetscCallA(VecSet(b, one, ierr))

     39:   PetscCallA(KSPCreate(PETSC_COMM_WORLD, ksp, ierr))
     40:   PetscCallA(KSPSetType(ksp, 'python', ierr))
     41:   PetscCallA(KSPPythonSetType(ksp, 'example100.py:ConjGrad', ierr))

     43:   PetscCallA(KSPGetPC(ksp, pc, ierr))
     44:   PetscCallA(PCSetType(pc, 'python', ierr))
     45:   PetscCallA(PCPythonSetType(pc, 'example100.py:Jacobi', ierr))

     47:   PetscCallA(KSPSetOperators(ksp, A, A, ierr))
     48:   PetscCallA(KSPSetFromOptions(ksp, ierr))
     49:   PetscCallA(KSPSolve(ksp, b, x, ierr))

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

     57:   if (draw) then
     58:     PetscCallA(VecView(x, PETSC_VIEWER_DRAW_WORLD, ierr))
     59:     PetscCallA(PetscSleep(rtwo, ierr))
     60:   end if

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

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

     71: !/*TEST
     72: !
     73: !    test:
     74: !      requires: petsc4py
     75: !      localrunfiles: example100.py
     76: !
     77: !TEST*/