From 600eb557b26530945104eea9b8af0502f9ead16c Mon Sep 17 00:00:00 2001 From: Pierre-Francois Loos Date: Thu, 20 Jul 2023 21:50:41 +0200 Subject: [PATCH] remove plot_GW --- src/GW/evGW.f90 | 4 -- src/GW/evUGW.f90 | 4 -- src/GW/plot_GW.f90 | 113 --------------------------------------------- 3 files changed, 121 deletions(-) delete mode 100644 src/GW/plot_GW.f90 diff --git a/src/GW/evGW.f90 b/src/GW/evGW.f90 index 54a73fa..9d07238 100644 --- a/src/GW/evGW.f90 +++ b/src/GW/evGW.f90 @@ -212,10 +212,6 @@ subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop ! End main loop !------------------------------------------------------------------------ -! Plot stuff - -! call plot_GW(nBas,nC,nO,nV,nR,nS,eHF,eGW,Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin)) - ! Did it actually converge? if(nSCF == maxSCF+1) then diff --git a/src/GW/evUGW.f90 b/src/GW/evUGW.f90 index 17a5bdc..29f74f8 100644 --- a/src/GW/evUGW.f90 +++ b/src/GW/evUGW.f90 @@ -219,10 +219,6 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE ! End main loop !------------------------------------------------------------------------ -! Plot stuff - -! call plot_GW(nBas,nC,nO,nV,nR,nS,eHF,eGW,Omega(:,ispin),rho(:,:,:,ispin),rhox(:,:,:,ispin)) - ! Did it actually converge? if(nSCF == maxSCF+1) then diff --git a/src/GW/plot_GW.f90 b/src/GW/plot_GW.f90 deleted file mode 100644 index fedc40c..0000000 --- a/src/GW/plot_GW.f90 +++ /dev/null @@ -1,113 +0,0 @@ -subroutine plot_GW(nBas,nC,nO,nV,nR,nS,eHF,eGW,Omega,rho) - -! Dump several GW quantities for external plotting - - implicit none - include 'parameters.h' - -! Input variables - - integer,intent(in) :: nBas,nC,nO,nV,nR,nS - double precision,intent(in) :: eHF(nBas),eGW(nBas),Omega(nS),rho(nBas,nBas,nS) - -! Local variables - - integer :: i,j,a,b,x,jb,g - integer :: nGrid - double precision :: eps,eta,wmin,wmax,dw - double precision,allocatable :: w(:),SigC(:,:),Z(:,:),S(:,:) - -! Infinitesimal - - eta = 1d-3 - -! Construct grid - - nGrid = 1000 - allocate(w(nGrid),SigC(nBas,nGrid),Z(nBas,nGrid),S(nBas,nGrid)) - -! Initialize - - SigC(:,:) = 0d0 - Z(:,:) = 0d0 - -! Minimum and maximum frequency values - - wmin = -5d0 - wmax = +5d0 - dw = (wmax - wmin)/dble(ngrid) - - do g=1,nGrid - w(g) = wmin + dble(g)*dw - enddo - -! Occupied part of the self-energy and renormalization factor - - do g=1,nGrid - do x=nC+1,nBas-nR - do i=nC+1,nO - jb = 0 - do j=nC+1,nO - do b=nO+1,nBas-nR - jb = jb + 1 - eps = w(g) - eHF(i) + Omega(jb) - SigC(x,g) = SigC(x,g) + 2d0*rho(x,i,jb)**2*eps/(eps**2 + eta**2) - Z(x,g) = Z(x,g) + 2d0*rho(x,i,jb)**2/eps**2 - enddo - enddo - enddo - enddo - enddo - -! Virtual part of the self-energy and renormalization factor - - do g=1,nGrid - do x=nC+1,nBas-nR - do a=nO+1,nBas-nR - jb = 0 - do j=nC+1,nO - do b=nO+1,nBas-nR - jb = jb + 1 - eps = w(g) - eHF(a) - Omega(jb) - SigC(x,g) = SigC(x,g) + 2d0*rho(x,a,jb)**2*eps/(eps**2 + eta**2) - Z(x,g) = Z(x,g) + 2d0*rho(x,a,jb)**2/eps**2 - enddo - enddo - enddo - enddo - enddo - - Z(:,:) = 1d0/(1d0 + Z(:,:)) - -! Compute spectral function - - do g=1,nGrid - do x=nC+1,nBas-nR - S(x,g) = eta/((w(g) - eHF(x) - SigC(x,g))**2 + eta**2) - enddo - enddo - - S(:,:) = S(:,:)/pi - -! Dump quantities in files as a function of w - - open(unit=8 ,file='plot/grid.dat') - open(unit=9 ,file='plot/SigC.dat') - open(unit=10 ,file='plot/Z.dat') - open(unit=11 ,file='plot/A.dat') - - do g=1,nGrid - write(8 ,*) w(g)*HaToeV,(SigC(x,g)*HaToeV,x=1,nBas) - write(9 ,*) w(g)*HaToeV,((w(g)-eHF(x))*HaToeV,x=1,nBas) - write(10,*) w(g)*HaToeV,(Z(x,g),x=1,nBas) - write(11,*) w(g)*HaToeV,(S(x,g),x=1,nBas) - enddo - -! Closing files - - close(unit=8) - close(unit=9) - close(unit=10) - close(unit=11) - -end subroutine plot_GW