Actual source code: ex75f.F90

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

      5: #include <petsc/finclude/petscksp.h>
      6: program main
      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:   PetscCallA(PetscInitialize(ierr))
     24:   dir = '.'
     25:   PetscCallA(PetscOptionsGetString(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-load_dir', dir, flg, ierr))
     26:   PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-nmat', nmat, flg, ierr))
     27:   reset = PETSC_FALSE
     28:   PetscCallA(PetscOptionsGetBool(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-reset', reset, flg, ierr))
     29:   PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr))
     30:   PetscCallA(KSPCreate(PETSC_COMM_WORLD, ksp, ierr))
     31:   PetscCallA(KSPSetOperators(ksp, A, A, ierr))
     32:   do i = 0, nmat - 1
     33:     j = i + 400
     34:     fmt = '(I3)'
     35:     write (cmat, fmt) j
     36:     write (name, '(a)') trim(dir)//'/A_'//cmat//'.dat'
     37:     PetscCallA(PetscViewerBinaryOpen(PETSC_COMM_WORLD, name, FILE_MODE_READ, viewer, ierr))
     38:     PetscCallA(MatLoad(A, viewer, ierr))
     39:     PetscCallA(PetscViewerDestroy(viewer, ierr))
     40:     if (i == 0) then
     41:       PetscCallA(MatCreateVecs(A, x, b, ierr))
     42:     end if
     43:     write (name, '(a)') trim(dir)//'/rhs_'//cmat//'.dat'
     44:     PetscCallA(PetscViewerBinaryOpen(PETSC_COMM_WORLD, name, FILE_MODE_READ, viewer, ierr))
     45:     PetscCallA(VecLoad(b, viewer, ierr))
     46:     PetscCallA(PetscViewerDestroy(viewer, ierr))
     47:     PetscCallA(KSPSetFromOptions(ksp, ierr))
     48:     PetscCallA(KSPSolve(ksp, b, x, ierr))
     49:     PetscCallA(PetscObjectTypeCompare(ksp, KSPHPDDM, flg, ierr))
     50: #if defined(PETSC_HAVE_HPDDM)
     51:     if (flg .and. reset) then
     52:       PetscCallA(KSPHPDDMGetDeflationMat(ksp, U, ierr))
     53:       PetscCallA(KSPReset(ksp, ierr))
     54:       PetscCallA(KSPSetOperators(ksp, A, A, ierr))
     55:       PetscCallA(KSPSetFromOptions(ksp, ierr))
     56:       PetscCallA(KSPSetUp(ksp, ierr))
     57:       if (.not. PetscObjectIsNull(U)) then
     58:         PetscCallA(KSPHPDDMSetDeflationMat(ksp, U, ierr))
     59:         PetscCallA(MatDestroy(U, ierr))
     60:       end if
     61:     end if
     62: #endif
     63:   end do
     64:   PetscCallA(VecDestroy(x, ierr))
     65:   PetscCallA(VecDestroy(b, ierr))
     66:   PetscCallA(MatDestroy(A, ierr))
     67:   PetscCallA(KSPDestroy(ksp, ierr))
     68:   PetscCallA(PetscFinalize(ierr))
     69: end

     71: !/*TEST
     72: !
     73: !   test:
     74: !      suffix: 1
     75: !      output_file: output/ex75_1.out
     76: !      nsize: 1
     77: !      requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
     78: !      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
     79: !
     80: !   test:
     81: !      requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
     82: !      suffix: 1_icc
     83: !      output_file: output/ex75_1_icc.out
     84: !      nsize: 1
     85: !      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
     86: !
     87: !   testset:
     88: !      requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
     89: !      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
     90: !      test:
     91: !        nsize: 1
     92: !        suffix: 2_seq
     93: !        output_file: output/ex75_2.out
     94: !      test:
     95: !        nsize: 2
     96: !        suffix: 2_par
     97: !        output_file: output/ex75_2.out
     98: !
     99: !   testset:
    100: !      requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
    101: !      output_file: output/ex75_2_icc.out
    102: !      nsize: 1
    103: !      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
    104: !      test:
    105: !        suffix: 2_icc
    106: !        args:
    107: !      test:
    108: !        suffix: 2_icc_atol
    109: !        output_file: output/ex75_2_icc_atol.out
    110: !        args: -ksp_atol 1e-12
    111: !
    112: !TEST*/