TABLE OF CONTENTS
ABINIT/exc_build_ham [ Functions ]
NAME
exc_build_ham
FUNCTION
Calculate and write the excitonic Hamiltonian on an external binary file (Fortran file open in random mode) for subsequent treatment in the Bethe-Salpeter code.
COPYRIGHT
Copyright (C) 1992-2009 EXC group (L.Reining, V.Olevano, F.Sottile, S.Albrecht, G.Onida) Copyright (C) 2009-2018 ABINIT group (L.Reining, V.Olevano, F.Sottile, S.Albrecht, G.Onida, M.Giantomassi) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
BSp<excparam>=The parameters for the Bethe-Salpeter calculation. BS_files<excfiles>=File names internally used in the BS code. Cryst<crystal_t>=Info on the crystalline structure. Kmesh<kmesh_t>=The list of k-points in the BZ, IBZ and symmetry tables. Qmesh<kmesh_t>=The list of q-points for epsilon^{-1} and related symmetry tables. ktabr(nfftot_osc,BSp%nkbz)=The FFT index of $(R^{-1}(r-\tau))$ where R is symmetry needed to obtains the k-points from the irreducible image. Used to symmetrize u_Sk where S = \transpose R^{-1} Gsph_x<gsphere_t>=Info on the G-sphere used to describe wavefunctions and W (the largest one is actually stored). Gsph_c<gsphere_t>=Info on the G-sphere used to describe the correlation part. Vcp<vcoul_t>=The Coulomb interaction in reciprocal space. A cutoff can be used W<screen_t>=Data type gathering info and data for W. nfftot_osc=Total Number of FFT points used for the oscillator matrix elements. ngfft_osc(18)=Info on the FFT algorithm used to calculate the oscillator matrix elements. Psps<Pseudopotential_type>=Variables related to pseudopotentials Pawtab(Psps%ntypat)<pawtab_type>=PAW tabulated starting data. Pawang<pawang_type>=PAW angular mesh and related data. Paw_pwff(Cryst%ntypat*Wfd%usepaw)<pawpwff_t>=Form factor used to calculate the onsite mat. elements of a plane wave. Wfd<wfd_t>=Handler for the wavefunctions.
OUTPUT
The excitonic Hamiltonian is saved on an external binary file (see below).
PARENTS
bethe_salpeter
CHILDREN
cwtime,get_bz_item,gsph_fft_tabs,paw_rho_tw_g,pawcprj_alloc pawcprj_free,pawpwij_free,pawpwij_init,rho_tw_g,timab,wfd_change_ngfft wfd_get_cprj,wfd_get_ur,wrtout,xmpi_distab,xmpi_sum
SOURCE
51 #if defined HAVE_CONFIG_H 52 #include "config.h" 53 #endif 54 55 #include "abi_common.h" 56 57 subroutine exc_build_ham(BSp,BS_files,Cryst,Kmesh,Qmesh,ktabr,Gsph_x,Gsph_c,Vcp,& 58 & Wfd,W,Hdr_bse,nfftot_osc,ngfft_osc,Psps,Pawtab,Pawang,Paw_pwff) 59 60 use defs_basis 61 use defs_datatypes 62 use defs_abitypes 63 use m_profiling_abi 64 use m_bs_defs 65 use m_xmpi 66 use m_errors 67 68 use m_gwdefs, only : czero_gw 69 use m_crystal, only : crystal_t 70 use m_gsphere, only : gsphere_t 71 use m_vcoul, only : vcoul_t 72 use m_bz_mesh, only : kmesh_t 73 use m_pawpwij, only : pawpwff_t 74 use m_pawang, only : pawang_type 75 use m_pawtab, only : pawtab_type 76 use m_pawcprj, only : pawcprj_type 77 use m_wfd, only : wfd_t 78 use m_screen, only : screen_t 79 80 !This section has been created automatically by the script Abilint (TD). 81 !Do not modify the following lines by hand. 82 #undef ABI_FUNC 83 #define ABI_FUNC 'exc_build_ham' 84 use interfaces_14_hidewrite 85 use interfaces_18_timing 86 use interfaces_71_bse, except_this_one => exc_build_ham 87 !End of the abilint section 88 89 implicit none 90 91 !Arguments ------------------------------------ 92 !scalars 93 integer,intent(in) :: nfftot_osc 94 type(excparam),intent(in) :: BSp 95 type(excfiles),intent(in) :: BS_files 96 type(screen_t),intent(inout) :: W 97 type(kmesh_t),intent(in) :: Kmesh,Qmesh 98 type(crystal_t),intent(in) :: Cryst 99 type(vcoul_t),intent(in) :: Vcp 100 type(gsphere_t),intent(in) :: Gsph_x,Gsph_c 101 type(Pseudopotential_type),intent(in) :: Psps 102 type(Hdr_type),intent(inout) :: Hdr_bse 103 type(pawang_type),intent(in) :: Pawang 104 type(wfd_t),target,intent(inout) :: Wfd 105 !arrays 106 integer,intent(in) :: ngfft_osc(18) 107 integer,intent(in) :: ktabr(nfftot_osc,Kmesh%nbz) 108 type(Pawtab_type),intent(in) :: Pawtab(Psps%ntypat*Wfd%usepaw) 109 type(pawpwff_t),intent(in) :: Paw_pwff(Psps%ntypat*Wfd%usepaw) 110 111 !Local variables ------------------------------ 112 !scalars 113 logical :: do_resonant,do_coupling 114 !character(len=500) :: msg 115 !arrays 116 real(dp) :: tsec(2) 117 complex(gwpc),allocatable :: all_mgq0(:,:,:,:,:) 118 119 !************************************************************************ 120 121 call timab(670,1,tsec) 122 123 ABI_CHECK(Wfd%nspinor==1,"nspinor==2 not coded") 124 ABI_CHECK(nfftot_osc==PRODUCT(ngfft_osc(1:3)),"mismatch in FFT size") 125 126 if (BSp%have_complex_ene) then 127 MSG_ERROR("Complex energies are not supported yet") 128 end if 129 130 ! Do we have to compute some block? 131 do_resonant = (BS_files%in_hreso == BSE_NOFILE) 132 do_coupling = (BS_files%in_hcoup == BSE_NOFILE) 133 134 if (BSp%use_coupling == 0) then 135 if (.not.do_resonant) then 136 call wrtout(std_out,"Will skip the calculation of resonant block (will use BSR file)","COLL") 137 goto 100 138 end if 139 else 140 if (.not. do_resonant .and. .not. do_coupling) then 141 call wrtout(std_out,"Will skip the calculation of both resonant and coupling block (will use BSR and BSC files)","COLL") 142 goto 100 143 end if 144 end if 145 146 ! Compute M_{k,q=0}^{b,b}(G) for all k-points in the IBZ and each pair b, b' 147 ! used for the exchange part and part of the Coulomb term. 148 call wrtout(std_out," Calculating all matrix elements for q=0 to save CPU time","COLL") 149 150 call wfd_all_mgq0(Wfd,Cryst,Qmesh,Gsph_x,Vcp,Psps,Pawtab,Paw_pwff,& 151 & Bsp%lomo_spin,Bsp%homo_spin,Bsp%humo_spin,nfftot_osc,ngfft_osc,Bsp%npweps,all_mgq0) 152 153 ! ======================== 154 ! ==== Resonant Block ==== 155 ! ======================== 156 if (do_resonant) then 157 call timab(672,1,tsec) 158 call exc_build_block(BSp,Cryst,Kmesh,Qmesh,ktabr,Gsph_x,Gsph_c,Vcp,& 159 & Wfd,W,Hdr_bse,nfftot_osc,ngfft_osc,Psps,Pawtab,Pawang,Paw_pwff,all_mgq0,.TRUE.,BS_files%out_hreso) 160 call timab(672,2,tsec) 161 end if 162 163 ! ======================== 164 ! ==== Coupling Block ==== 165 ! ======================== 166 if (do_coupling.and.BSp%use_coupling>0) then 167 call timab(673,1,tsec) 168 call exc_build_block(BSp,Cryst,Kmesh,Qmesh,ktabr,Gsph_x,Gsph_c,Vcp,& 169 & Wfd,W,Hdr_bse,nfftot_osc,ngfft_osc,Psps,Pawtab,Pawang,Paw_pwff,all_mgq0,.FALSE.,BS_files%out_hcoup) 170 call timab(673,2,tsec) 171 end if 172 ! 173 ! * Free memory. 174 ABI_FREE(all_mgq0) 175 176 100 call timab(670,2,tsec) 177 178 end subroutine exc_build_ham
ABINIT/wfd_all_mgq0 [ Functions ]
NAME
wfd_all_mgq0
FUNCTION
INPUTS
Wfd<wfd_t>=Handler for the wavefunctions. Cryst<crystal_t>=Info on the crystalline structure. Qmesh<kmesh_t>=The list of q-points for epsilon^{-1} and related symmetry tables. Gsph_x<gsphere_t>=G-sphere with the G-vectors in mgq0. Vcp<vcoul_t>=The Coulomb interaction in reciprocal space. A cutoff can be used Psps<Pseudopotential_type>=Variables related to pseudopotentials Pawtab(Psps%ntypat)<pawtab_type>=PAW tabulated starting data. Paw_pwff(Cryst%ntypat*Wfd%usepaw)<pawpwff_t>=Form factor used to calculate the onsite mat. elements of a plane wave. lomo_spin(Wfd%nsppol)=Lowest occupied band for each spin homo_spin(Wfd%nsppol)=Highest occupied band for each spin humo_spin(Wfd%nsppol)=Highest unoccupied band for each spin nfftot_osc=Total Number of FFT points used for the oscillator matrix elements. ngfft_osc(18)=Info on the FFT algorithm used to calculate the oscillator matrix elements. npweps=Number of G-vectors in mgq0.
OUTPUT
mgq0(npweps,lomo_min:humo_max,lomo_min:humo_max,Wfd%nkibz,Wfd%nsppol) Allocated here and filled with the matrix elements on each node.
PARENTS
exc_build_ham
CHILDREN
cwtime,get_bz_item,gsph_fft_tabs,paw_rho_tw_g,pawcprj_alloc pawcprj_free,pawpwij_free,pawpwij_init,rho_tw_g,timab,wfd_change_ngfft wfd_get_cprj,wfd_get_ur,wrtout,xmpi_distab,xmpi_sum
SOURCE
217 subroutine wfd_all_mgq0(Wfd,Cryst,Qmesh,Gsph_x,Vcp,& 218 & Psps,Pawtab,Paw_pwff,lomo_spin,homo_spin,humo_spin,nfftot_osc,ngfft_osc,npweps,mgq0) 219 220 use defs_basis 221 use m_profiling_abi 222 use m_xmpi 223 use m_errors 224 225 use defs_datatypes, only : pseudopotential_type 226 use m_gwdefs, only : czero_gw 227 use m_time, only : cwtime 228 use m_crystal, only : crystal_t 229 use m_gsphere, only : gsphere_t, gsph_fft_tabs 230 use m_vcoul, only : vcoul_t 231 use m_bz_mesh, only : kmesh_t, get_BZ_item 232 use m_pawpwij, only : pawpwff_t, pawpwij_t, pawpwij_init, pawpwij_free, paw_rho_tw_g 233 use m_pawtab, only : pawtab_type 234 use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_free 235 use m_wfd, only : wfd_t, wfd_get_ur, wfd_get_cprj, wfd_change_ngfft, wfd_ihave_ur, wfd_distribute_bbp 236 use m_oscillators, only : rho_tw_g 237 238 !This section has been created automatically by the script Abilint (TD). 239 !Do not modify the following lines by hand. 240 #undef ABI_FUNC 241 #define ABI_FUNC 'wfd_all_mgq0' 242 use interfaces_14_hidewrite 243 use interfaces_18_timing 244 !End of the abilint section 245 246 implicit none 247 248 !Arguments ------------------------------------ 249 !scalars 250 integer,intent(in) :: nfftot_osc,npweps 251 type(kmesh_t),intent(in) :: Qmesh 252 type(crystal_t),intent(in) :: Cryst 253 type(vcoul_t),intent(in) :: Vcp 254 type(gsphere_t),intent(in) :: Gsph_x 255 type(Pseudopotential_type),intent(in) :: Psps 256 type(wfd_t),target,intent(inout) :: Wfd 257 !arrays 258 integer,intent(in) :: lomo_spin(Wfd%nsppol),homo_spin(Wfd%nsppol),humo_spin(Wfd%nsppol) 259 integer,intent(in) :: ngfft_osc(18) 260 complex(gwpc),allocatable,intent(out) :: mgq0(:,:,:,:,:) 261 type(Pawtab_type),intent(in) :: Pawtab(Psps%ntypat) 262 type(pawpwff_t),intent(in) :: Paw_pwff(Psps%ntypat*Wfd%usepaw) 263 264 !Local variables ------------------------------ 265 !scalars 266 integer,parameter :: map2sphere1=1,dim_rtwg1=1,ndat1=1 267 integer :: use_padfft,mgfft_osc,fftalga_osc,ii 268 integer :: ik_ibz,itim_k,isym_k,iq_bz,iq_ibz,isym_q,itim_q,iqbz0 269 integer :: ierr,iv,ic,spin,lomo_min,humo_max !,inv_ipw,ipw 270 real(dp) :: cpu,wall,gflops !q0vol,fcc_const 271 complex(dpc) :: ph_mkt 272 character(len=500) :: msg 273 !arrays 274 integer,allocatable :: igfftg0(:),task_distrib(:,:,:,:) 275 integer,allocatable :: gbound(:,:),id_tab(:) 276 real(dp) :: qbz(3),spinrot_k(4),tsec(2) 277 complex(gwpc),allocatable :: rhotwg1(:) 278 complex(gwpc),target,allocatable :: ur1(:),ur2(:) 279 complex(gwpc),ABI_CONTIGUOUS pointer :: ptr_ur1(:),ptr_ur2(:) 280 type(pawcprj_type),allocatable :: Cp1(:,:),Cp2(:,:) 281 type(pawpwij_t),allocatable :: Pwij_q0(:) 282 283 !************************************************************************ 284 285 call timab(671,1,tsec) 286 287 ABI_CHECK(Wfd%nspinor==1,"nspinor==2 not coded") 288 ABI_CHECK(nfftot_osc==PRODUCT(ngfft_osc(1:3)),"mismatch in FFT size") 289 290 lomo_min = MINVAL(lomo_spin); humo_max = MAXVAL(humo_spin) 291 292 if ( ANY(ngfft_osc(1:3) /= Wfd%ngfft(1:3)) ) then 293 call wfd_change_ngfft(Wfd,Cryst,Psps,ngfft_osc) 294 end if 295 296 mgfft_osc = MAXVAL(ngfft_osc(1:3)) 297 fftalga_osc = ngfft_osc(7)/100 !; fftalgc_osc=MOD(ngfft_osc(7),10) 298 299 ! (temporary) Table used for the wavefunction in the IBZ. 300 ABI_MALLOC(id_tab, (Wfd%nfft)) 301 id_tab = (/(ii, ii=1,Wfd%nfft)/) 302 303 ! Analytic integration of 4pi/q^2 over the volume element: 304 ! $4pi/V \int_V d^3q 1/q^2 =4pi bz_geometric_factor V^(-2/3)$ 305 ! i_sz=4*pi*bz_geometry_factor*q0_vol**(-two_thirds) where q0_vol= V_BZ/N_k 306 ! bz_geometry_factor: sphere=7.79, fcc=7.44, sc=6.188, bcc=6.946, wz=5.255 307 ! (see gwa.pdf, appendix A.4) 308 309 ! If q=0 and C=V then set up rho-twiddle(G=0) to reflect an 310 ! analytic integration of q**-2 over the volume element: 311 ! <q**-2> = 7.44 V**(-2/3) (for fcc cell) 312 313 ! q0vol = (8.0*pi**3) / (Cryst%ucvol*Kmesh%nbz) 314 ! fcc_const = SQRT(7.44*q0vol**(-2.0/3.0)) 315 ! rtw = (6.0*pi**2/(Cryst%ucvol*Kmesh%nkbz))**(1./3.) 316 ! Average of (q+q')**-2 integration for head of Coulomb matrix 317 ! INTRTW(QL) = (2*pi*rtw + pi*(rtw**2/QL-QL)*LOG((QL+rtw)/(QL-rtw))) 318 ! & * (Cryst%ucvol*Kmesh%nbz)/(2*pi)**3. * QL*QL 319 320 if (Wfd%usepaw==1) then 321 ABI_DT_MALLOC(Cp1,(Wfd%natom,Wfd%nspinor)) 322 call pawcprj_alloc(Cp1,0,Wfd%nlmn_atm) 323 ABI_DT_MALLOC(Cp2,(Wfd%natom,Wfd%nspinor)) 324 call pawcprj_alloc(Cp2,0,Wfd%nlmn_atm) 325 end if 326 327 ABI_MALLOC(ur1,(nfftot_osc*Wfd%nspinor)) 328 ABI_MALLOC(ur2,(nfftot_osc*Wfd%nspinor)) 329 330 ! Identify q==0 331 iqbz0=0 332 do iq_bz=1,Qmesh%nbz 333 if (ALL(ABS(Qmesh%bz(:,iq_bz))<tol3)) iqbz0 = iq_bz 334 end do 335 ABI_CHECK(iqbz0/=0,"q=0 not found in q-point list!") 336 337 ! * Get iq_ibz, and symmetries from iqbz0. 338 call get_BZ_item(Qmesh,iqbz0,qbz,iq_ibz,isym_q,itim_q) 339 340 if (Wfd%usepaw==1) then ! Prepare onsite contributions at q==0 341 ABI_DT_MALLOC(Pwij_q0,(Cryst%ntypat)) 342 call pawpwij_init(Pwij_q0,npweps,Qmesh%bz(:,iqbz0),Gsph_x%gvec,Cryst%rprimd,Psps,Pawtab,Paw_pwff) 343 end if 344 ! 345 ! Tables for the FFT of the oscillators. 346 ! a) FFT index of the G sphere (only vertical transitions, unlike cchi0, no need to shift the sphere). 347 ! b) gbound table for the zero-padded FFT performed in rhotwg. 348 ABI_MALLOC(igfftg0,(Gsph_x%ng)) 349 ABI_MALLOC(gbound,(2*mgfft_osc+8,2)) 350 call gsph_fft_tabs(Gsph_x,(/0,0,0/),mgfft_osc,ngfft_osc,use_padfft,gbound,igfftg0) 351 if ( ANY(fftalga_osc == (/2,4/)) ) use_padfft=0 ! Pad-FFT is not coded in rho_tw_g 352 if (use_padfft==0) then 353 ABI_FREE(gbound) 354 ABI_MALLOC(gbound,(2*mgfft_osc+8,2*use_padfft)) 355 end if 356 357 ABI_MALLOC(rhotwg1,(npweps)) 358 359 ABI_STAT_MALLOC(mgq0, (npweps,lomo_min:humo_max,lomo_min:humo_max,Wfd%nkibz,Wfd%nsppol), ierr) 360 ABI_CHECK(ierr==0, "out-of-memory in mgq0") 361 mgq0 = czero 362 363 call cwtime(cpu,wall,gflops,"start") 364 365 do spin=1,Wfd%nsppol 366 ! Distribute the calculation of the matrix elements. 367 ! processors have the entire set of wavefunctions hence we divide the workload 368 ! without checking if the pair of states is available. Last dimension is fake. 369 ABI_MALLOC(task_distrib,(lomo_spin(spin):humo_spin(spin),lomo_spin(spin):humo_spin(spin),Wfd%nkibz,1)) 370 call xmpi_distab(Wfd%nproc,task_distrib) 371 372 ! loop over the k-points in IBZ 373 do ik_ibz=1,Wfd%nkibz 374 if ( ALL(task_distrib(:,:,ik_ibz,1)/= Wfd%my_rank) ) CYCLE 375 376 ! Don't need to symmetrize the wavefunctions. 377 itim_k=1; isym_k=1; ph_mkt=cone; spinrot_k=Cryst%spinrot(:,isym_k) 378 379 do iv=lomo_spin(spin),humo_spin(spin) ! Loop over band V 380 if ( ALL(task_distrib(:,iv,ik_ibz,1)/=Wfd%my_rank) ) CYCLE 381 382 if (wfd_ihave_ur(Wfd,iv,ik_ibz,spin,how="Stored")) then 383 ptr_ur1 => Wfd%Wave(iv,ik_ibz,spin)%ur 384 else 385 call wfd_get_ur(Wfd,iv,ik_ibz,spin,ur1) 386 ptr_ur1 => ur1 387 end if 388 389 if (Wfd%usepaw==1) then 390 call wfd_get_cprj(Wfd,iv,ik_ibz,spin,Cryst,Cp1,sorted=.FALSE.) 391 end if 392 393 ! Loop over band C 394 do ic=lomo_spin(spin),humo_spin(spin) 395 if ( task_distrib(ic,iv,ik_ibz,1)/=Wfd%my_rank ) CYCLE 396 397 if (wfd_ihave_ur(Wfd,ic,ik_ibz,spin,how="Stored")) then 398 ptr_ur2 => Wfd%Wave(ic,ik_ibz,spin)%ur 399 else 400 call wfd_get_ur(Wfd,ic,ik_ibz,spin,ur2) 401 ptr_ur2 => ur2 402 end if 403 404 if (Wfd%usepaw==1) then 405 call wfd_get_cprj(Wfd,ic,ik_ibz,spin,Cryst,Cp2,sorted=.FALSE.) 406 end if 407 408 call rho_tw_g(Wfd%nspinor,npweps,nfftot_osc,ndat1,ngfft_osc,map2sphere1,use_padfft,igfftg0,gbound,& 409 & ptr_ur1,1,id_tab,ph_mkt,spinrot_k,& 410 & ptr_ur2,1,id_tab,ph_mkt,spinrot_k,& 411 & dim_rtwg1,rhotwg1) 412 413 if (Wfd%usepaw==1) then ! Add PAW onsite contribution. 414 call paw_rho_tw_g(npweps,dim_rtwg1,Wfd%nspinor,Cryst%natom,Cryst%ntypat,Cryst%typat,Cryst%xred,Gsph_x%gvec,& 415 & Cp1,Cp2,Pwij_q0,rhotwg1) 416 end if 417 418 ! If q=0 treat Exchange and Coulomb-term independently 419 if (iv <= homo_spin(spin) .and. ic <= homo_spin(spin) .or. & 420 & iv > homo_spin(spin) .and. ic > homo_spin(spin)) then 421 422 if (iv/=ic) then !COULOMB term: C/=V: ignore them 423 rhotwg1(1) = czero_gw 424 else 425 ! If q=0 and C=V then set up rho-twiddle(G=0) to reflect an 426 ! analytic integration of q**-2 over the volume element: 427 ! <q**-2> = 7.44 V**(-2/3) (for fcc cell) 428 !rhotwg1(1) = fcc_const * qpg(1,iqbz0) 429 rhotwg1(1) = SQRT(GWPC_CMPLX(Vcp%i_sz,zero)) / Vcp%vcqlwl_sqrt(1,1) 430 !if (vcut) rhotwg1(1) = 1.0 431 end if 432 433 else 434 ! At present this term is set to zero 435 ! EXCHANGE term: limit value. 436 ! Set up rho-twiddle(G=0) using small vector q instead of zero and k.p perturbation theory (see notes) 437 rhotwg1(1) = czero_gw 438 end if 439 440 mgq0(:,iv,ic,ik_ibz,spin) = rhotwg1(:) 441 end do !ic 442 end do !iv 443 end do !ik_ibz 444 445 ABI_FREE(task_distrib) 446 end do !spin 447 448 ! TODO: One can speedup the calculation by computing the upper triangle of the 449 ! matrix in (b,b') space and then take advantage of the symmetry property: 450 ! 451 ! M_{k,0}{{bb'}(G)^* = M{k,0}{b'b'}(-G) 452 453 #if 0 454 !!!! $OMP PARALLEL DO COLLAPSE(3) PRIVATE(inv_ipw) 455 do spin=1,Wfd%nsppol 456 do ik_ibz=1,Wfd%nkibz 457 do iv=lomo_spin(spin),humo_spin(spin) 458 do ic=1,iv-1 459 do ipw=1,npweps 460 inv_ipw = gsph_x%g2mg(ipw) 461 mgq0(inv_ipw,ic,iv,ik_ibz,spin) = mgq0(ipw,iv,ic,ik_ibz,spin) 462 end do 463 end do 464 end do 465 end do 466 end do 467 #endif 468 ! 469 ! Gather matrix elements on each node. 470 call xmpi_sum(mgq0,Wfd%comm,ierr) 471 472 call cwtime(cpu,wall,gflops,"stop") 473 write(msg,'(2(a,f9.6))')"cpu_time = ",cpu,", wall_time = ",wall 474 call wrtout(std_out,msg,"PERS") 475 476 ABI_FREE(rhotwg1) 477 ABI_FREE(igfftg0) 478 ABI_FREE(gbound) 479 ABI_FREE(ur1) 480 ABI_FREE(ur2) 481 ABI_FREE(id_tab) 482 483 if (Wfd%usepaw==1) then 484 ! Deallocation for PAW. 485 call pawpwij_free(Pwij_q0) 486 ABI_DT_FREE(Pwij_q0) 487 call pawcprj_free(Cp1) 488 ABI_DT_FREE(Cp1) 489 call pawcprj_free(Cp2) 490 ABI_DT_FREE(Cp2) 491 end if 492 493 call timab(671,2,tsec) 494 495 end subroutine wfd_all_mgq0