Actual source code: ex75f.F90

      1: !
      2: !   Description: Solves a series of linear systems using KSPHPDDM.
      3: !

      5:       program main
      6: #include <petsc/finclude/petscksp.h>
      7:       use petscksp
      8:       implicit none
      9:       Vec                            x,b
     10:       Mat                            A
     11: #if defined(PETSC_HAVE_HPDDM)
     12:       Mat                            U
     13: #endif
     14:       KSP                            ksp
     15:       PetscInt                       i,j,nmat
     16:       PetscViewer                    viewer
     17:       character*(PETSC_MAX_PATH_LEN) dir,name
     18:       character*(8)                  fmt
     19:       character(3)                   cmat
     20:       PetscBool                      flg,reset
     21:       PetscErrorCode                 ierr

     23:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
     24:       if (ierr .ne. 0) then
     25:         print *,'Unable to initialize PETSc'
     26:         stop
     27:       endif
     28:       dir = '.'
     29:       call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-load_dir',dir,flg,ierr);CHKERRA(ierr)
     30:       call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-nmat',nmat,flg,ierr);CHKERRA(ierr)
     31:       reset = PETSC_FALSE
     32:       call PetscOptionsGetBool(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-reset',reset,flg,ierr);CHKERRA(ierr)
     33:       call MatCreate(PETSC_COMM_WORLD,A,ierr);CHKERRA(ierr)
     34:       call KSPCreate(PETSC_COMM_WORLD,ksp,ierr);CHKERRA(ierr)
     35:       call KSPSetOperators(ksp,A,A,ierr);CHKERRA(ierr)
     36:       do 50 i=0,nmat-1
     37:         j = i+400
     38:         fmt = '(I3)'
     39:         write (cmat,fmt) j
     40:         write (name,'(a)')trim(dir)//'/A_'//cmat//'.dat'
     41:         call PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
     42:         call MatLoad(A,viewer,ierr);CHKERRA(ierr)
     43:         call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
     44:         if (i .eq. 0) then
     45:           call MatCreateVecs(A,x,b,ierr);CHKERRA(ierr)
     46:         endif
     47:         write (name,'(a)')trim(dir)//'/rhs_'//cmat//'.dat'
     48:         call PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
     49:         call VecLoad(b,viewer,ierr);CHKERRA(ierr)
     50:         call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
     51:         call KSPSetFromOptions(ksp,ierr);CHKERRA(ierr)
     52:         call KSPSolve(ksp,b,x,ierr);CHKERRA(ierr)
     53:         call PetscObjectTypeCompare(ksp,KSPHPDDM,flg,ierr);CHKERRA(ierr)
     54: #if defined(PETSC_HAVE_HPDDM)
     55:         if (flg .and. reset) then
     56:           call KSPHPDDMGetDeflationSpace(ksp,U,ierr);CHKERRA(ierr)
     57:           call KSPReset(ksp,ierr);CHKERRA(ierr)
     58:           call KSPSetOperators(ksp,A,A,ierr);CHKERRA(ierr)
     59:           call KSPSetFromOptions(ksp,ierr);CHKERRA(ierr)
     60:           call KSPSetUp(ksp,ierr);CHKERRA(ierr)
     61:           if (U .ne. PETSC_NULL_MAT) then
     62:             call KSPHPDDMSetDeflationSpace(ksp,U,ierr);CHKERRA(ierr)
     63:             call MatDestroy(U,ierr);CHKERRA(ierr)
     64:           endif
     65:         endif
     66: #endif
     67:   50  continue
     68:       call VecDestroy(x,ierr);CHKERRA(ierr)
     69:       call VecDestroy(b,ierr);CHKERRA(ierr)
     70:       call MatDestroy(A,ierr);CHKERRA(ierr)
     71:       call KSPDestroy(ksp,ierr);CHKERRA(ierr)
     72:       call PetscFinalize(ierr)
     73:       end

     75: !/*TEST
     76: !
     77: !   test:
     78: !      suffix: 1
     79: !      output_file: output/ex75_1.out
     80: !      nsize: 1
     81: !      requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
     82: !      args: -nmat 1 -pc_type none -ksp_converged_reason -ksp_type {{gmres hpddm}shared output} -ksp_max_it 1000 -ksp_gmres_restart 1000 -ksp_rtol 1e-10 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR
     83: !
     84: !   test:
     85: !      requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
     86: !      suffix: 1_icc
     87: !      output_file: output/ex75_1_icc.out
     88: !      nsize: 1
     89: !      args: -nmat 1 -pc_type icc -ksp_converged_reason -ksp_type {{gmres hpddm}shared output} -ksp_max_it 1000 -ksp_gmres_restart 1000 -ksp_rtol 1e-10 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR
     90: !
     91: !   testset:
     92: !      requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
     93: !      args: -nmat 3 -pc_type none -ksp_converged_reason -ksp_type hpddm -ksp_max_it 1000 -ksp_gmres_restart 40 -ksp_rtol 1e-10 -ksp_hpddm_type gcrodr -ksp_hpddm_recycle 20 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR
     94: !      test:
     95: !        nsize: 1
     96: !        suffix: 2_seq
     97: !        output_file: output/ex75_2.out
     98: !      test:
     99: !        nsize: 2
    100: !        suffix: 2_par
    101: !        output_file: output/ex75_2.out
    102: !
    103: !   testset:
    104: !      requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
    105: !      output_file: output/ex75_2_icc.out
    106: !      nsize: 1
    107: !      args: -nmat 3 -pc_type icc -ksp_converged_reason -ksp_type hpddm -ksp_max_it 1000 -ksp_gmres_restart 40 -ksp_rtol 1e-10 -ksp_hpddm_type gcrodr -ksp_hpddm_recycle 20 -reset {{false true}shared output} -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR
    108: !      test:
    109: !        suffix: 2_icc
    110: !        args:
    111: !      test:
    112: !        suffix: 2_icc_atol
    113: !        output_file: output/ex75_2_icc_atol.out
    114: !        args: -ksp_atol 1e-12
    115: !
    116: !TEST*/