Actual source code: ex76f.F90

petsc-3.13.4 2020-08-01
Report Typos and Errors
  1: !
  2: !   Description: Solves a linear systems using PCHPDDM.
  3: !

  5:       program main
  6:  #include <petsc/finclude/petscksp.h>
  7:       use petscksp
  8:       use petscisdef
  9:       implicit none
 10:       Vec             x,b
 11:       Mat             A,aux
 12:       KSP             ksp
 13:       PC              pc
 14:       IS              is,sizes
 15:       PetscScalar     one
 16:       PetscInt, pointer :: idx(:)
 17:       PetscMPIInt     rank,size
 18:       PetscViewer     viewer
 19:       character*(128) dir,name
 20:       character*(8)   fmt
 21:       character(1)    crank,csize
 22:       PetscBool       flg
 23:       PetscErrorCode  ierr

 25:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 26:       if (ierr .ne. 0) then
 27:         print *,'Unable to initialize PETSc'
 28:         stop
 29:       endif
 30:       call MPI_Comm_size(PETSC_COMM_WORLD,size,ierr)
 31:       if (size .ne. 4) then
 32:         print *,'This example requires 4 processes'
 33:         stop
 34:       endif
 35:       call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
 36:       call MatCreate(PETSC_COMM_WORLD,A,ierr);CHKERRA(ierr)
 37:       call MatCreate(PETSC_COMM_SELF,aux,ierr);CHKERRA(ierr)
 38:       call ISCreate(PETSC_COMM_SELF,is,ierr);CHKERRA(ierr)
 39:       dir = '.'
 40:       call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-load_dir',dir,flg,ierr);CHKERRA(ierr)
 41:       fmt = '(I1)'
 42:       write (crank,fmt) rank
 43:       write (csize,fmt) size
 44:       write (name,'(a)')trim(dir)//'/sizes_'//crank//'_'//csize//'.dat'
 45:       call PetscViewerBinaryOpen(PETSC_COMM_SELF,name,FILE_MODE_READ, viewer,ierr);CHKERRA(ierr)
 46:       call ISCreate(PETSC_COMM_SELF,sizes,ierr);CHKERRA(ierr)
 47:       call ISLoad(sizes,viewer,ierr);CHKERRA(ierr)
 48:       call ISGetIndicesF90(sizes,idx,ierr);CHKERRA(ierr)
 49:       call MatSetSizes(A,idx(1),idx(2),idx(3),idx(4),ierr);CHKERRA(ierr)
 50:       call ISRestoreIndicesF90(sizes,idx,ierr);CHKERRA(ierr)
 51:       call ISDestroy(sizes,ierr);CHKERRA(ierr)
 52:       call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
 53:       call MatSetUp(A,ierr);CHKERRA(ierr)
 54:       write (name,'(a)')trim(dir)//'/A.dat'
 55:       call PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
 56:       call MatLoad(A,viewer,ierr);CHKERRA(ierr)
 57:       call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
 58:       write (name,'(a)')trim(dir)//'/is_'//crank//'_'//csize//'.dat'
 59:       call PetscViewerBinaryOpen(PETSC_COMM_SELF,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
 60:       call ISLoad(is,viewer,ierr);CHKERRA(ierr)
 61:       call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
 62:       write (name,'(a)')trim(dir)//'/Neumann_'//crank//'_'//csize//'.dat'
 63:       call PetscViewerBinaryOpen(PETSC_COMM_SELF,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
 64:       call MatSetBlockSizesFromMats(aux,A,A,ierr);CHKERRA(ierr)
 65:       call MatLoad(aux,viewer,ierr);CHKERRA(ierr)
 66:       call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
 67:       call MatSetOption(A,MAT_SYMMETRIC,PETSC_TRUE,ierr);CHKERRA(ierr)
 68:       call MatSetOption(aux,MAT_SYMMETRIC,PETSC_TRUE,ierr);CHKERRA(ierr)
 69:       call KSPCreate(PETSC_COMM_WORLD,ksp,ierr);CHKERRA(ierr)
 70:       call KSPSetOperators(ksp,A,A,ierr);CHKERRA(ierr)
 71:       call KSPGetPC(ksp,pc,ierr);CHKERRA(ierr)
 72:       call PCSetType(pc,PCHPDDM,ierr);CHKERRA(ierr)
 73: #if defined(PETSC_HAVE_HPDDM)
 74:       call PCHPDDMSetAuxiliaryMat(pc,is,aux,PETSC_NULL_FUNCTION,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)
 75:       call PCHPDDMHasNeumannMat(pc,PETSC_FALSE,ierr);CHKERRA(ierr)
 76: #endif
 77:       call ISDestroy(is,ierr);CHKERRA(ierr)
 78:       call MatDestroy(aux,ierr);CHKERRA(ierr)
 79:       call KSPSetFromOptions(ksp,ierr);CHKERRA(ierr)
 80:       call MatCreateVecs(A,x,b,ierr);CHKERRA(ierr)
 81:       one = 1.0
 82:       call VecSet(b,one,ierr);CHKERRA(ierr)
 83:       call KSPSolve(ksp,b,x,ierr);CHKERRA(ierr)
 84:       call VecDestroy(x,ierr);CHKERRA(ierr)
 85:       call VecDestroy(b,ierr);CHKERRA(ierr)
 86:       call KSPDestroy(ksp,ierr);CHKERRA(ierr)
 87:       call MatDestroy(A,ierr);CHKERRA(ierr)
 88:       call PetscFinalize(ierr)
 89:       end

 91: !/*TEST
 92: !
 93: !   test:
 94: !      requires: hpddm slepc datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
 95: !      output_file: output/ex76_1.out
 96: !      nsize: 4
 97: !      args: -ksp_rtol 1e-3 -ksp_converged_reason -pc_type {{bjacobi hpddm}shared output} -pc_hpddm_coarse_sub_pc_type lu -sub_pc_type lu -options_left no -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO
 98: !
 99: !   test:
100: !      requires: hpddm slepc datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
101: !      suffix: geneo
102: !      output_file: output/ex76_geneo_pc_hpddm_levels_1_eps_nev-5.out
103: !      nsize: 4
104: !      args: -ksp_converged_reason -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type cholesky -pc_hpddm_levels_1_eps_nev 5 -pc_hpddm_levels_1_st_pc_type cholesky -pc_hpddm_coarse_p {{1 2}shared output} -pc_hpddm_coarse_pc_type redundant -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO
105: !
106: !   test:
107: !      requires: hpddm slepc datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
108: !      suffix: fgmres_geneo_20_p_2
109: !      output_file: output/ex76_fgmres_geneo_20_p_2.out
110: !      nsize: 4
111: !      args: -ksp_converged_reason -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 20 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_type fgmres -pc_hpddm_coarse_mat_type {{baij sbaij}shared output} -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO
112: !
113: !   test:
114: !      requires: hpddm slepc datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
115: !      suffix: fgmres_geneo_20_p_2_geneo
116: !      output_file: output/ex76_fgmres_geneo_20_p_2.out
117: !      nsize: 4
118: !      args: -ksp_converged_reason -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type cholesky -pc_hpddm_levels_1_eps_nev 20 -pc_hpddm_levels_2_p 2 -pc_hpddm_levels_2_mat_type {{baij sbaij}shared output} -pc_hpddm_levels_2_eps_nev {{5 20}shared output} -pc_hpddm_levels_2_sub_pc_type cholesky -pc_hpddm_levels_2_ksp_type gmres -ksp_type fgmres -pc_hpddm_coarse_mat_type {{baij sbaij}shared output} -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO
119: !
120: !TEST*/