1! 2! Copyright (C) 2009 Quantum ESPRESSO group 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8!------------------------------------------------------------------------- 9SUBROUTINE write_wfcfile(filename,wfc,elaux,num) 10!------------------------------------------------------------------------- 11! 12! This subroutine writes a formatted wavefunction file 13! 14 USE kinds, ONLY : DP 15 USE io_global, ONLY : ionode, ionode_id 16 USE ld1inc, ONLY : grid 17 USE radial_grids, ONLY : ndmx 18 USE mp, ONLY : mp_bcast 19 USE mp_world, ONLY : world_comm 20 21 IMPLICIT NONE 22 INTEGER, INTENT(in) :: num ! the number of wavefunctions to write 23 REAL(DP), INTENT(in) :: wfc(ndmx,num) 24 CHARACTER(len=2), INTENT(in) :: elaux(num) 25 CHARACTER(len=256), INTENT(in) :: filename 26 INTEGER :: ios, n, ns 27 28 IF (filename == ' ') RETURN 29 30 IF (ionode) & 31 OPEN(unit=19,file=filename, status='unknown', iostat=ios, err=80) 3280 CALL mp_bcast(ios, ionode_id,world_comm) 33 CALL errore('write_wfcfile','opening file '//trim(filename),abs(ios)) 34 IF (ionode) THEN 35 WRITE(19,'("#",12x,"r",38(18x,a2))') (elaux(n),n=1,num) 36 DO n=1,grid%mesh 37 WRITE(19,'(38f20.12)') grid%r(n), (wfc(n,ns), ns=1,num) 38 ENDDO 39 CLOSE(19) 40 ENDIF 41 RETURN 42END SUBROUTINE write_wfcfile 43 44!------------------------------------------------------------------------- 45SUBROUTINE write_wfcfile_ft(filename_in,wfc,num) 46!------------------------------------------------------------------------- 47! 48! calculate and print the fourier transform of wfc and the 49! convergence of their norm with cutoff 50! 51USE kinds, ONLY : DP 52USE io_global, ONLY : ionode, ionode_id 53USE ld1inc, ONLY : grid, lls 54USE radial_grids, ONLY : ndmx 55USE mp, ONLY : mp_bcast 56USE mp_world, ONLY : world_comm 57 58IMPLICIT NONE 59INTEGER, INTENT(in) :: num ! the number of wavefunctions to write 60real(DP), INTENT(in) :: wfc(ndmx,num) 61CHARACTER(len=256), INTENT(in) :: filename_in 62 63INTEGER :: ios, n, ns, nmax 64CHARACTER(len=256) :: filename 65real(DP) :: q, fac, pi, wrk(ndmx), jlq(ndmx), norm(num), normr(num), work(num) 66real(DP), EXTERNAL :: int_0_inf_dr ! the function calculating the integral 67 68IF (filename_in == ' ') RETURN 69 70IF (ionode) THEN 71 filename = trim(filename_in)//'.q' 72 OPEN(unit=19,file=filename, status='unknown', iostat=ios) 73 filename = trim(filename_in)//'.norm_q' 74 OPEN(unit=29,file=filename, status='unknown', iostat=ios) 75 pi = 4._dp*atan(1._dp) 76 DO ns=1,num 77 wrk(1:grid%mesh)=(wfc(1:grid%mesh,ns)*exp(-0.04*grid%r2(1:grid%mesh)))**2 78 normr(ns) = int_0_inf_dr ( wrk, grid, grid%mesh, 2*lls(ns)+2 ) 79 ENDDO 80 !write (*,*) normr(1:num) 81 fac = pi /grid%r(grid%mesh) * 2._dp/pi 82 norm(:) = 0._dp 83 nmax = int (10 * grid%r(grid%mesh) /pi) 84 DO n=1,nmax 85 q=n * pi /grid%r(grid%mesh) 86 DO ns=1,num 87 CALL sph_bes ( grid%mesh, grid%r, q, lls(ns), jlq ) 88 wrk(1:grid%mesh)=jlq(1:grid%mesh)*wfc(1:grid%mesh,ns)* & 89 exp(-0.04*grid%r2(1:grid%mesh))*grid%r(1:grid%mesh) 90 work(ns) = int_0_inf_dr ( wrk, grid, grid%mesh, 2*lls(ns)+2 ) 91 ENDDO 92 norm(1:num) = norm(1:num) + work(1:num)*work(1:num)*q*q*fac 93 WRITE (19,'(15f12.6)') q, (work(ns), ns=1,num) 94 WRITE (29,'(15f12.6)') q, (norm(ns)/normr(ns), ns=1,num) 95 ENDDO 96 CLOSE (29) 97 CLOSE (19) 98 ! 99 ! end of Fourier analysis 100 ! 101ENDIF 102RETURN 103END SUBROUTINE write_wfcfile_ft 104! 105!------------------------------------------------------------------------- 106SUBROUTINE write_efun(filename,fun,x,npte,num) 107!------------------------------------------------------------------------- 108! 109! This subroutine writes a functions defined on the energy grid 110! 111USE kinds, ONLY : DP 112USE io_global, ONLY : ionode, ionode_id 113USE mp, ONLY : mp_bcast 114USE mp_world, ONLY : world_comm 115 116IMPLICIT NONE 117INTEGER, INTENT(in) :: num ! the number of wavefunctions to write 118INTEGER, INTENT(in) :: npte ! the number of energies 119real(DP), INTENT(in) :: fun(npte,num), x(npte) 120CHARACTER(len=256), INTENT(in) :: filename 121INTEGER :: ios, n, ns 122 123IF (filename == ' ') RETURN 124 125IF (ionode) & 126 OPEN(unit=19,file=filename, status='unknown', iostat=ios, err=800) 127800 CALL mp_bcast(ios, ionode_id, world_comm) 128CALL errore('write_wfcfile','opening file '//trim(filename),abs(ios)) 129IF (ionode) THEN 130 DO n=1,npte 131 WRITE(19,'(38f20.12)') x(n), (max(min(fun(n,ns),9.d4),-9d4), ns=1,num) 132 ENDDO 133 CLOSE(19) 134ENDIF 135 136RETURN 137END SUBROUTINE write_efun 138 139