Actual source code: ex77f.F90

petsc-3.13.4 2020-08-01
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:       Vec             cx,cb
 11:       Mat             A
 12:       KSP             ksp
 13:       PetscScalar, pointer :: x_v(:),b_v(:)
 14:       PetscInt        m,n,L,K
 15:       PetscViewer     viewer
 16:       character*(128) dir,name
 17:       PetscBool       flg
 18:       PetscErrorCode  ierr

 20:       call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
 21:       if (ierr .ne. 0) then
 22:         print *,'Unable to initialize PETSc'
 23:         stop
 24:       endif
 25:       dir = '.'
 26:       call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-load_dir',dir,flg,ierr);CHKERRA(ierr)
 27:       K = 5
 28:       call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-n',K,flg,ierr);CHKERRA(ierr)
 29:       call MatCreate(PETSC_COMM_WORLD,A,ierr);CHKERRA(ierr)
 30:       call KSPCreate(PETSC_COMM_WORLD,ksp,ierr);CHKERRA(ierr)
 31:       call KSPSetOperators(ksp,A,A,ierr);CHKERRA(ierr)
 32:       write (name,'(a)')trim(dir)//'/A_400.dat'
 33:       call PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
 34:       call MatLoad(A,viewer,ierr);CHKERRA(ierr)
 35:       call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
 36:       call MatGetLocalSize(A,m,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)
 37:       call MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,K,PETSC_NULL_SCALAR,B,ierr);CHKERRA(ierr)
 38:       call MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,K,PETSC_NULL_SCALAR,X,ierr);CHKERRA(ierr)
 39:       call MatSetRandom(B,PETSC_NULL_RANDOM,ierr);CHKERRA(ierr)
 40:       call KSPSetFromOptions(ksp,ierr);CHKERRA(ierr)
 41:       call KSPSetUp(ksp,ierr);CHKERRA(ierr)
 42:       call PetscObjectTypeCompare(ksp,KSPHPDDM,flg,ierr);CHKERRA(ierr)
 43: #if defined(PETSC_HAVE_HPDDM)
 44:       if (flg) then
 45:         call KSPHPDDMMatSolve(ksp,B,X,ierr);CHKERRA(ierr)
 46:       else
 47: #endif
 48:         call MatGetSize(A,L,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)
 49:         do 50 n=0,K-1
 50:           call MatDenseGetColumnF90(B,n,b_v,ierr);CHKERRA(ierr)
 51:           call MatDenseGetColumnF90(X,n,x_v,ierr);CHKERRA(ierr)
 52:           call VecCreateMPIWithArray(PETSC_COMM_WORLD,1,m,L,b_v,cb,ierr);CHKERRA(ierr)
 53:           call VecCreateMPIWithArray(PETSC_COMM_WORLD,1,m,L,x_v,cx,ierr);CHKERRA(ierr)
 54:           call KSPSolve(ksp,cb,cx,ierr);CHKERRA(ierr)
 55:           call VecDestroy(cx,ierr);CHKERRA(ierr)
 56:           call VecDestroy(cb,ierr);CHKERRA(ierr)
 57:           call MatDenseRestoreColumnF90(X,x_v,ierr);CHKERRA(ierr)
 58:           call MatDenseRestoreColumnF90(B,b_v,ierr);CHKERRA(ierr)
 59:   50    continue
 60: #if defined(PETSC_HAVE_HPDDM)
 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: hpddm 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: !         output_file: output/ex77_2_ksp_hpddm_type-gmres.out
 83: !         args: -ksp_type hpddm -pc_type asm -ksp_hpddm_type gmres
 84: !      test:
 85: !         suffix: 2b
 86: !         output_file: output/ex77_2_ksp_hpddm_type-bgmres.out
 87: !         args: -ksp_type hpddm -pc_type asm -ksp_hpddm_type bgmres
 88: !      test:
 89: !         suffix: 3a
 90: !         output_file: output/ex77_3_ksp_hpddm_type-gcrodr.out
 91: !         args: -ksp_type hpddm -ksp_hpddm_recycle 10 -ksp_hpddm_type gcrodr
 92: !      test:
 93: !         suffix: 3b
 94: !         output_file: output/ex77_3_ksp_hpddm_type-bgcrodr.out
 95: !         args: -ksp_type hpddm -ksp_hpddm_recycle 10 -ksp_hpddm_type bgcrodr
 96: !
 97: !TEST*/