Actual source code: ex77f.F90

petsc-master 2020-07-05
Report Typos and Errors
  1: !
  2: !   Description: Solves a linear system with a block of right-hand sides using KSPHPDDM.
  3: !

  5:       program main
  6:  #include <petsc/finclude/petscksp.h>
  7:       use petscksp
  8:       implicit none
  9:       Mat             X,B
 10:       Mat             A
 11:       KSP             ksp
 12:       PC              pc
 13:       Mat             F
 14:       PetscScalar     alpha
 15:       PetscReal       norm
 16:       PetscInt        m,K
 17:       PetscViewer     viewer
 18:       character*(128) dir,name
 19:       PetscBool       flg
 20:       PetscErrorCode  ierr

 22:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 23:       if (ierr .ne. 0) then
 24:         print *,'Unable to initialize PETSc'
 25:         stop
 26:       endif
 27:       dir = '.'
 28:       call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-load_dir',dir,flg,ierr);CHKERRA(ierr)
 29:       K = 5
 30:       call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-n',K,flg,ierr);CHKERRA(ierr)
 31:       call MatCreate(PETSC_COMM_WORLD,A,ierr);CHKERRA(ierr)
 32:       call KSPCreate(PETSC_COMM_WORLD,ksp,ierr);CHKERRA(ierr)
 33:       call KSPSetOperators(ksp,A,A,ierr);CHKERRA(ierr)
 34:       write (name,'(a)')trim(dir)//'/A_400.dat'
 35:       call PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
 36:       call MatLoad(A,viewer,ierr);CHKERRA(ierr)
 37:       call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
 38:       call MatGetLocalSize(A,m,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)
 39:       call MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,K,PETSC_NULL_SCALAR,B,ierr);CHKERRA(ierr)
 40:       call MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,K,PETSC_NULL_SCALAR,X,ierr);CHKERRA(ierr)
 41:       call MatSetRandom(B,PETSC_NULL_RANDOM,ierr);CHKERRA(ierr)
 42:       call KSPSetFromOptions(ksp,ierr);CHKERRA(ierr)
 43:       call KSPSetUp(ksp,ierr);CHKERRA(ierr)
 44:       call KSPMatSolve(ksp,B,X,ierr);CHKERRA(ierr)
 45:       call KSPGetMatSolveBlockSize(ksp,M,ierr);CHKERRA(ierr)
 46:       if (M .ne. PETSC_DECIDE) then
 47:         call KSPSetMatSolveBlockSize(ksp,PETSC_DECIDE,ierr);CHKERRA(ierr)
 48:         call MatZeroEntries(X,ierr);CHKERRA(ierr)
 49:         call KSPMatSolve(ksp,B,X,ierr);CHKERRA(ierr)
 50:       endif
 51:       call KSPGetPC(ksp,pc,ierr);CHKERRA(ierr)
 52:       call PetscObjectTypeCompare(pc,PCLU,flg,ierr);CHKERRA(ierr)
 53:       if (flg) then
 54:         call PCFactorGetMatrix(pc,F,ierr);CHKERRA(ierr)
 55:         call MatMatSolve(F,B,B,ierr);CHKERRA(ierr)
 56:         alpha = -1.0
 57:         call MatAYPX(B,alpha,X,SAME_NONZERO_PATTERN,ierr);CHKERRA(ierr)
 58:         call MatNorm(B,NORM_INFINITY,norm,ierr);CHKERRA(ierr)
 59:         if (norm > 100*PETSC_MACHINE_EPSILON) then
 60:           SETERRA(PETSC_COMM_WORLD,PETSC_ERR_PLIB,'KSPMatSolve() and MatMatSolve() difference has nonzero norm')
 61:         endif
 62:       endif
 63:       call MatDestroy(X,ierr);CHKERRA(ierr)
 64:       call MatDestroy(B,ierr);CHKERRA(ierr)
 65:       call MatDestroy(A,ierr);CHKERRA(ierr)
 66:       call KSPDestroy(ksp,ierr);CHKERRA(ierr)
 67:       call PetscFinalize(ierr)
 68:       end

 70: !/*TEST
 71: !
 72: !   testset:
 73: !      nsize: 2
 74: !      requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
 75: !      args: -ksp_converged_reason -ksp_max_it 1000 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR
 76: !      test:
 77: !         suffix: 1
 78: !         output_file: output/ex77_1.out
 79: !         args:
 80: !      test:
 81: !         suffix: 2a
 82: !         requires: hpddm
 83: !         output_file: output/ex77_2_ksp_hpddm_type-gmres.out
 84: !         args: -ksp_type hpddm -pc_type asm -ksp_hpddm_type gmres
 85: !      test:
 86: !         suffix: 2b
 87: !         requires: hpddm
 88: !         output_file: output/ex77_2_ksp_hpddm_type-bgmres.out
 89: !         args: -ksp_type hpddm -pc_type asm -ksp_hpddm_type bgmres
 90: !      test:
 91: !         suffix: 3a
 92: !         requires: hpddm
 93: !         output_file: output/ex77_3_ksp_hpddm_type-gcrodr.out
 94: !         args: -ksp_type hpddm -ksp_hpddm_recycle 5 -ksp_hpddm_type gcrodr
 95: !      test:
 96: !         suffix: 3b
 97: !         requires: hpddm
 98: !         output_file: output/ex77_3_ksp_hpddm_type-bgcrodr.out
 99: !         args: -ksp_type hpddm -ksp_hpddm_recycle 5 -ksp_hpddm_type bgcrodr
100: !      test:
101: !         nsize: 4
102: !         suffix: 4
103: !         requires: hpddm
104: !         output_file: output/ex77_4.out
105: !         args: -ksp_rtol 1e-4 -ksp_type hpddm -ksp_hpddm_recycle 5 -ksp_hpddm_type bgcrodr -ksp_view_final_residual -N 12 -ksp_matsolve_block_size 5
106: !   test:
107: !      nsize: 1
108: !      suffix: preonly
109: !      requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
110: !      output_file: output/ex77_preonly.out
111: !      args: -N 6 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR -pc_type lu -ksp_type hpddm -ksp_hpddm_type preonly
112: !   test:
113: !      nsize: 4
114: !      suffix: 4_slepc
115: !      requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) slepc
116: !      output_file: output/ex77_4.out
117: !      filter: sed "/^ksp_hpddm_recycle_ Linear eigensolve converged/d"
118: !      args: -ksp_converged_reason -ksp_max_it 500 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR -ksp_rtol 1e-4 -ksp_type hpddm -ksp_hpddm_recycle 5 -ksp_hpddm_type bgcrodr -ksp_view_final_residual -N 12 -ksp_matsolve_block_size 5 -ksp_hpddm_recycle_redistribute 2 -ksp_hpddm_recycle_mat_type dense -ksp_hpddm_recycle_eps_converged_reason -ksp_hpddm_recycle_st_pc_type redundant
119: !
120: !TEST*/