exampleMPI.F90

Description

This Fortran example propagates a 1 particle continuous-time quantum walk on an infinite line

Amongst the features used, it illustrates:
  • calling the libctqwMPI library from Fortran
  • recieving command line options using PETSc
  • creating PETSc log stages
  • the use of the chebyshev algorithm
  • adding a diagonal defects to various nodes
  • viewing vectors via PetscViewer

Compiling and executing

To compile this example, open a terminal in the top directory of the pyCTQW-X.Y folder, and run

$ make fortran
$ make examples

The example can then be run by changing into the examples/fortran folder, and running

$ mpirun -np X exampleMPI

where X is the number of MPI processes to run.

Output

$ mpirun -np 4 exampleMPI
Vector Object: 4 MPI processes
  type: mpi
Process [0]
1.30318e-15 + 9.41288e-16 i
4.26061e-15 - 5.83926e-15 i
-2.6771e-14 - 1.97231e-14 i
-8.90705e-14 + 1.19677e-13 i
5.22435e-13 + 3.92973e-13 i
1.69261e-12 - 2.2254e-12 i
-9.24365e-12 - 7.11278e-12 i
-2.91412e-11 + 3.74127e-11 i
1.47434e-10 + 1.16317e-10 i
4.51971e-10 - 5.6522e-10 i
-2.10616e-09 - 1.70822e-09 i
-6.27416e-09 + 7.62086e-09 i
2.67489e-08 + 2.23733e-08 i
7.73784e-08 - 9.09727e-08 i
-2.99429e-07 - 2.59261e-07 i
-8.40536e-07 + 9.52525e-07 i
2.92435e-06 + 2.63331e-06 i
7.96061e-06 - 8.65084e-06 i
-2.46145e-05 - 2.31843e-05 i
-6.49345e-05 + 6.72304e-05 i
0.000175879 + 0.000174556 i
0.00044937 - 0.000439577 i
-0.00104657 - 0.00110507 i
-0.00258837 + 0.0023656 i
0.00505632 + 0.00575514 i
Process [1]
0.012099 - 0.0101714 i
-0.0191447 - 0.0239345 i
-0.0442924 + 0.0334726 i
0.053856 + 0.0761079 i
0.12025 - 0.0787511 i
-0.102846 - 0.172366 i
-0.219755 + 0.116923 i
0.111189 + 0.241323 i
0.214786 - 0.0830223 i
-0.0457964 - 0.132734 i
-0.0208922 + 0.0294486 i
0.0624386 - 0.0589731 i
-0.0404472 - 0.137791 i
-0.191971 - 0.0860889 i
-0.229877 + 0.138719 i
-0.0383627 + 0.243064 i
0.0596756 + 0.217395 i
0.21466 + 0.189964 i
0.261424 + 0.00431375 i
0.21683 - 0.0655021 i
0.178224 - 0.166243 i
0.0818483 - 0.179087 i
0.0335323 - 0.183573 i
-0.01719 - 0.145368 i
-0.0340458 - 0.120739 i
Process [2]
-0.0344425 - 0.0737224 i
-0.0328421 - 0.0407838 i
0.0103582 - 0.0348961 i
-0.0119036 + 0.0257429 i
-0.00586558 - 0.0254129 i
0.010251 - 0.000336933 i
-0.0147693 - 0.0158695 i
0.0286581 - 0.0184521 i
-0.00508709 + 0.00712859 i
0.0298955 - 0.0310717 i
0.026666 + 0.0413678 i
-0.011159 - 0.00816756 i
0.0425462 + 0.0359954 i
-0.0610777 + 0.0473539 i
-0.018832 - 0.0447773 i
-0.00248435 + 0.0247191 i
-0.0610373 - 0.0547772 i
0.0902541 - 0.0779196 i
0.0756292 + 0.10105 i
-0.0915928 + 0.061552 i
-0.0438766 - 0.0715585 i
0.049748 - 0.0280784 i
0.0163906 + 0.0313776 i
-0.0181928 + 0.00882673 i
-0.0044223 - 0.00978976 i
Process [3]
0.00492519 - 0.00207491 i
0.000916545 + 0.00233018 i
-0.0010417 + 0.000382831 i
-0.000151757 - 0.000441799 i
0.000178365 - 5.72714e-05 i
2.06318e-05 + 6.875e-05 i
-2.53648e-05 + 7.11154e-06 i
-2.35018e-06 - 8.97787e-06 i
3.05473e-06 - 7.45989e-07 i
2.27794e-07 + 1.00097e-06 i
-3.16392e-07 + 6.70093e-08 i
-1.90127e-08 - 9.66128e-08 i
2.85389e-08 - 5.20876e-09 i
1.37913e-09 + 8.16537e-09 i
-2.26541e-09 + 3.53182e-10 i
-8.75365e-11 - 6.10116e-10 i
1.59661e-10 - 2.10083e-11 i
4.88366e-12 + 4.06353e-11 i
-1.00669e-11 + 1.0998e-12 i
-2.39917e-13 - 2.42956e-12 i
5.71638e-13 - 5.06798e-14 i
1.03592e-14 + 1.31214e-13 i
-2.94057e-14 + 2.04654e-15 i
-3.89668e-16 - 6.42454e-15 i
1.43612e-15 - 7.35772e-17 i

Source Code

[Download source code]

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
program exampleMPI
    
    use ctqwMPI
    implicit none
#include <finclude/petsc.h>
    
    ! declare variables
    PetscLogStage  :: stage
    PetscLogDouble :: te10, te21,te11, te20, ts0, ts1, tc1, tc0
    PetscMPIInt    :: rank
    PetscErrorCode :: ierr
    PetscBool      :: flag
    PetscInt       :: i, j, its, n, d(2)
    PetscScalar    :: Emin, Emax, t, init_state(2,2), amp(2)
    PetscReal      :: Emin_error, Emax_error
    Mat            :: H, H2
    Vec            :: psi0, psi, psix, psiy
    
    ! initialize SLEPc and PETSc
    call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
    call MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr)
    
    ! get command line arguments
    call PetscOptionsGetInt(PETSC_NULL_CHARACTER,"-n",n,flag,ierr)
    CHKERRQ(ierr)
    if (flag .eqv. .false.) n = 100
    call PetscOptionsGetInt(PETSC_NULL_CHARACTER,"-t",t,flag,ierr)
    CHKERRQ(ierr)
    if (flag .eqv. .false.) t = 10.0
    
    ! create the Hamiltonian
    d = [3,4]
    amp = [2.0,1.5]
    call PetscLogStageRegister('Hamiltonian',stage,ierr)
    call PetscLogStagePush(stage,ierr)
    call MatCreate(PETSC_COMM_WORLD,H,ierr)
    call hamiltonian_p1_line(H,d,amp,size(d),n)
    call PetscBarrier(H,ierr)
    call PetscLogStagePop(ierr)
    ! call MatView(H,PETSC_VIEWER_STDOUT_WORLD,ierr)
    
    ! Eigenvalue solver
    call PetscTime(te10,ierr)
    call PetscLogStageRegister('Emax',stage,ierr)
    call PetscLogStagePush(stage,ierr)
    call min_max_eigs(H,rank,Emax,Emax_error,'max','krylovschur','null',35,1.d-2,0,.false.,ierr)
    call PetscLogStagePop(ierr)
    !if (rank==0 .and. ierr==0) write(*,*)'    ',PetscRealPart(Emax),'+-',Emax_error
    call PetscTime(te11,ierr)
    
    call PetscTime(te20,ierr)
    call PetscLogStageRegister('Emin',stage,ierr)
    call PetscLogStagePush(stage,ierr)
    call min_max_eigs(H,rank,Emin,Emin_error,'min','krylovschur','null',35,1.d-2,0,.false.,ierr)
    call PetscLogStagePop(ierr)
    call PetscBarrier(H,ierr)
   ! if (rank==0 .and. ierr==0) write(*,*)'    ',PetscRealPart(Emin),'+-', Emin_error
    call PetscTime(te21,ierr)
    
    ! create vectors
    call VecCreate(PETSC_COMM_WORLD,psi0,ierr)
    call VecSetSizes(psi0,PETSC_DECIDE,n,ierr)
    call VecSetFromOptions(psi0,ierr)
    
    call VecDuplicate(psi0,psi,ierr)
    
    ! create initial state
    call PetscLogStageRegister('InitState',stage,ierr)
    call PetscLogStagePush(stage,ierr)
    init_state(1,:) = [0.,1.0/sqrt(2.0)]
    init_state(2,:) = [1.,1.0/sqrt(2.0)]
    call p1_init(psi0,init_state,2,n)
    !call VecView(psi0,PETSC_VIEWER_STDOUT_WORLD,ierr)
    call PetscBarrier(psi0,ierr)
    call PetscLogStagePop(ierr)
    
    ! matrix exponential
    call PetscTime(ts0,ierr)
    call PetscLogStageRegister('SLEPc expm',stage,ierr)
    call PetscLogStagePush(stage,ierr)
    call qw_krylov(H,t,psi0,psi)
    !call VecView(psi,PETSC_VIEWER_STDOUT_WORLD,ierr)
    call PetscBarrier(psi,ierr)
    call PetscLogStagePop(ierr)
    call PetscTime(ts1,ierr)
    
    ! QW chebyshev
    call PetscTime(tc0,ierr)
    call PetscLogStageRegister('Chebyshev',stage,ierr)
    call PetscLogStagePush(stage,ierr)
    call qw_cheby(psi0,psi,t,H,0.*PETSc_i,Emax,rank,n)
    call VecView(psi,PETSC_VIEWER_STDOUT_WORLD,ierr)
    call PetscBarrier(psi,ierr)
    call PetscLogStagePop(ierr)
    call PetscTime(tc1,ierr)

    ! destroy matrix/SLEPc
    call MatDestroy(H,ierr)
    call VecDestroy(psi,ierr)
    call VecDestroy(psi0,ierr)
    
    call PetscFinalize(ierr)

end program exampleMPI