1#include "f_defs.h" 2 3!----------------------------------------------------------------------- 4subroutine input_kernel(crys,gvec,kg,kgq,kp,syms,xct,flagbz, & 5 intwfnv,intwfnc) 6!----------------------------------------------------------------------- 7! 8! Read parameters from file WFN_co 9! Initialize k-points sampling, kg type 10! Initialize G-space, gvec 11! 12! input: xct type 13! 14! output: crys,gvec,syms,kg types 15! INT_VWFN_* and INT_CWFN_* files 16! 17 use global_m 18 use checkbz_m 19 use fullbz_m 20 use input_utils_m 21 use misc_m 22 use sort_m 23 use wfn_rho_vxc_io_m 24 use read_rho_vxc_m 25 use epsread_hdf5_m 26#ifdef HDF5 27 use hdf5 28#endif 29 implicit none 30 31 type (crystal), intent(out) :: crys 32 type (gspace), intent(out) :: gvec 33 type (grid), intent(out) :: kg,kgq 34 type (kpoints), intent(out) :: kp 35 type (symmetry), intent(out) :: syms 36 type (xctinfo), intent(inout) :: xct 37 integer, intent(in) :: flagbz 38 type (int_wavefunction), intent(out) :: intwfnc,intwfnv 39 40 41 type (wavefunction) :: wfnc,wfnv 42 character :: tmpfn*16 43 integer :: iwritev,iwritec,iwritek 44 integer :: ickmem,irk 45 integer :: ii,jj,is,isp,ig,ikq,ik,umk 46 integer :: irks,ivband,icband 47 48 real(DP) :: diffvol,vcell,kt(3),div,tol,delta,qq_temp(3) 49 50 real(DP), allocatable :: ek_tmp(:) 51 integer, allocatable :: index(:),indxk(:),indxkq(:),k_tmp(:,:) 52 integer, allocatable :: isrti(:) 53 SCALAR, allocatable :: cg(:,:), cgarray(:) 54 55 character(len=3) :: sheader 56 integer :: iflavor 57 type(gspace) :: gvec_kpt 58 logical :: skip_checkbz 59 60 61 PUSH_SUB(input_kernel) 62 63 call logit('input_kernel: reading WFN_co') 64 if (peinf%inode == 0) call open_file(25,file='WFN_co',form='unformatted',status='old') 65 66 sheader = 'WFN' 67 iflavor = 0 68 call read_binary_header_type(25, sheader, iflavor, kp, gvec, syms, crys, dont_warn_kgrid=xct%patched_sampling_co) 69 call check_trunc_kpts(xct%icutv, kp) 70 71 if(any(kp%shift > TOL_Zero) .and. peinf%inode == 0) then 72 write(0,'(a)') "WARNING: WFN_co has a shift. This is not recommended." 73 endif 74 75 call logit('input_kernel: reading gvec info') 76 SAFE_ALLOCATE(gvec%components, (3, gvec%ng)) 77 call read_binary_gvectors(25, gvec%ng, gvec%ng, gvec%components) 78 79 call get_volume(vcell,crys%bdot) 80 diffvol=abs(crys%celvol-vcell) 81 if (diffvol.gt.TOL_Small) then 82 call die('volume mismatch', only_root_writes = .true.) 83 endif 84 85 ! there is no fine grid to set the Fermi level in kernel 86 call find_efermi(xct%rfermi, xct%efermi, xct%efermi_input, kp, kp%mnband, 1, & 87 "coarse grid", should_search = .true., should_update = .true., write7 = .false.) 88 89 call assess_degeneracies(kp, kp%el(kp%mnband, :, :), kp%mnband - 1, xct%efermi, TOL_Degeneracy) 90 91 xct%nspin=kp%nspin 92 ! consistency check on status of spinor... 93 if (xct%nspinor .eq. 2 .and. kp%nspinor .ne. 2) & 94 call die("Keyword SPINOR used in absorption.inp but WFN is not of spinor type", only_root_writes = .true.) 95 96 if(any(kp%ifmax(:,:) == 0)) & 97 call die("BSE codes cannot handle a system where some k-points have no occupied bands.", only_root_writes = .true.) 98 99 kp%nvband=minval(kp%ifmax(:,:)-kp%ifmin(:,:))+1 100 kp%ncband=kp%mnband-maxval(kp%ifmax(:,:)) 101 102!---------------------------------------------------------------- 103! (gsm) check whether the requested number of bands 104! is available in the wavefunction file 105 106 if(xct%nvb_co.gt.kp%nvband) then 107 call die("The requested number of valence bands is not available in WFN_co.", & 108 only_root_writes=.true.) 109 endif 110 if(xct%ncb_co.gt.kp%ncband) then 111 call die("The requested number of conduction bands is not available in WFN_co.", & 112 only_root_writes=.true.) 113 endif 114 115! DAS: degenerate subspace check 116 117 if (peinf%inode.eq.0) then 118 if(xct%ncb_co.eq.kp%ncband) then 119 call die("You must provide one more conduction band in WFN_co in order to assess degeneracy.", & 120 only_root_writes=.true.) 121 endif 122 do jj = 1, kp%nspin 123 do ii = 1, kp%nrk 124 if(kp%ifmax(ii, jj) - xct%nvb_co > 0) then 125 ! no need to compare against band 0 if all valence are included 126 if(abs(kp%el(kp%ifmax(ii, jj) - xct%nvb_co + 1, ii, jj) & 127 - kp%el(kp%ifmax(ii, jj) - xct%nvb_co, ii, jj)) .lt. TOL_Degeneracy) then 128 if(xct%degeneracy_check_override) then 129 write(0,'(a)') & 130 "WARNING: Selected number of valence bands breaks degenerate subspace in WFN_co. " // & 131 "Run degeneracy_check.x for allowable numbers." 132 write(0,*) 133 else 134 write(0,'(a)') & 135 "Run degeneracy_check.x for allowable numbers, or use keyword " // & 136 "degeneracy_check_override to run anyway (at your peril!)." 137 call die("Selected number of valence bands breaks degenerate subspace in WFN_co.", & 138 only_root_writes=.true.) 139 endif 140 endif 141 endif 142 if(abs(kp%el(kp%ifmax(ii, jj) + xct%ncb_co, ii, jj) & 143 - kp%el(kp%ifmax(ii, jj) + xct%ncb_co + 1, ii, jj)) .lt. TOL_Degeneracy) then 144 if(xct%degeneracy_check_override) then 145 write(0,'(a)') & 146 "WARNING: Selected number of conduction bands breaks degenerate subspace in WFN_co. " // & 147 "Run degeneracy_check.x for allowable numbers." 148 write(0,*) 149 else 150 write(0,'(a)') & 151 "Run degeneracy_check.x for allowable numbers, or use keyword " // & 152 "degeneracy_check_override to run anyway (at your peril!)." 153 call die("Selected number of conduction bands breaks degenerate subspace in WFN_co.",& 154 only_root_writes=.true.) 155 endif 156 endif 157 enddo 158 enddo 159 endif 160 161!----------------------------------------------------------------------- 162! Read the k-point sampling from kpoints (if it exists) or from 163! WFN_co 164 165 if (xct%read_kpoints) then 166 if (peinf%inode.eq.0) then 167 call open_file(9,file='kpoints_co',form='formatted',status='old') 168 read(9,*) kg%nr 169 SAFE_ALLOCATE(kg%r, (3,kg%nr)) 170 do ii=1,kg%nr 171 read(9,*) (kg%r(jj,ii),jj=1,3),div 172 kg%r(:,ii) = kg%r(:,ii)/div 173 enddo 174 call close_file(9) 175 endif 176#ifdef MPI 177 call MPI_BCAST(kg%nr, 1, MPI_INTEGER,0,MPI_COMM_WORLD,mpierr) 178 if(peinf%inode.ne.0) then 179 SAFE_ALLOCATE(kg%r, (3,kg%nr)) 180 endif 181 call MPI_BCAST(kg%r, 3*kg%nr,MPI_REAL_DP,0,MPI_COMM_WORLD,mpierr) 182#endif 183 184!---------------------------------------------------------------- 185! indxk : stores the correspondence between k-points kg%r and kp%rk 186! (it is used to select the set of wavefunctions to be stored) 187! tol : tolerance in the coordinates of k-points 188 189 tol = 1.d-4 190 SAFE_ALLOCATE(indxk, (kg%nr)) 191 indxk=0 192 do jj=1,kg%nr 193 do ii=1,kp%nrk 194 kt(:) = kg%r(:,jj) - kp%rk(:,ii) 195 if ((abs(kt(1)).lt.tol).and.(abs(kt(2)).lt.tol) & 196 .and.(abs(kt(3)).lt.tol)) then 197 if (indxk(jj).ne.0) write(0,*) 'WARNING: multiple definition of k-point',jj,indxk(jj),kg%r(:,jj) 198 indxk(jj)=ii 199 endif 200 enddo 201 if (indxk(jj).eq.0) write(0,*) 'WARNING: could not find vector ',kg%r(:,jj),' in WFN_co' 202! 203! no need to stop here; if indxk.eq.0, the job will stop in genwf 204! 205 enddo 206 else 207 kg%nr=kp%nrk 208 SAFE_ALLOCATE(kg%r, (3,kg%nr)) 209 kg%r(1:3,1:kg%nr)=kp%rk(1:3,1:kp%nrk) 210 SAFE_ALLOCATE(indxk, (kg%nr)) 211 do ii=1,kg%nr 212 indxk(ii) = ii 213 enddo 214 endif 215 216!----------------------------------------------------------------------- 217! Order g-vectors with respect to their kinetic energy 218! 219 220 call logit('input_kernel: reordering gvecs') 221 222 ! FHJ: Figure out ecute and ecutg 223 call get_ecut() 224 225 SAFE_ALLOCATE(index, (gvec%ng)) 226 SAFE_ALLOCATE(gvec%ekin, (gvec%ng)) 227 call kinetic_energies(gvec, crys%bdot, gvec%ekin) 228 call sortrx(gvec%ng, gvec%ekin, index, gvec = gvec%components) 229 230 xct%ng = gcutoff(gvec%ng, gvec%ekin, index, xct%ecutg) 231 if (xct%theory.eq.0) & 232 xct%neps = gcutoff(gvec%ng, gvec%ekin, index, xct%ecute) ! cuteness energy?? 233 234 if(xct%ecute > xct%ecutg .and. xct%theory.eq.0) then 235 write(0,*) 'ecute = ', xct%ecute, ' ecutg = ', xct%ecutg 236 call die("The screened_coulomb_cutoff cannot be greater than the bare_coulomb_cutoff.", only_root_writes = .true.) 237 endif 238 239 SAFE_ALLOCATE(ek_tmp, (gvec%ng)) 240 ek_tmp = gvec%ekin 241 SAFE_ALLOCATE(k_tmp, (3,gvec%ng)) 242 k_tmp = gvec%components 243 do ii=1,gvec%ng 244 gvec%ekin(ii) = ek_tmp(index(ii)) 245 gvec%components(:,ii) = k_tmp(:,index(ii)) 246 enddo 247 248 call gvec_index(gvec) 249 250 ! If we are not doing just purely TDHF then 251 ! read the charge density/fxc 252 if ((xct%theory == 1) .and. ((1.0d0 - xct%coulomb_mod%long_range_frac_fock > TOL_SMALL) .or. & 253 (1.0d0 - xct%coulomb_mod%short_range_frac_fock > TOL_SMALL))) then 254 xct%coul_mod_flag=.true. 255 SAFE_ALLOCATE(isrti, (gvec%ng)) 256 do ii=1,gvec%ng 257 isrti(index(ii)) = ii 258 enddo 259 call read_rho(xct%wpg, gvec, kp, syms, crys, isrti, index, 'WFN_co') 260 SAFE_DEALLOCATE(isrti) 261 endif 262 263 SAFE_DEALLOCATE(index) 264 SAFE_DEALLOCATE(ek_tmp) 265 SAFE_DEALLOCATE(k_tmp) 266 267 268!----------------------------------------------------------------------- 269! Generate full brillouin zone from irreducible wedge, rk -> fk 270! 271! If flagbz.eq.1, only Identity will be used as 272! symmetry operation. In this case, kg%r (irreducible BZ) and kg%f 273! (full BZ) will be identical. 274! 275 if (flagbz.eq.0.and.peinf%inode.eq.0) write(6,801) 276 if (flagbz.eq.1.and.peinf%inode.eq.0) write(6,802) 277801 format(1x,'Using symmetries to expand the coarse grid sampling') 278802 format(1x,'No symmetries used in the coarse grid sampling') 279! 280 call timacc(7,1) 281 if (flagbz.eq.1) then 282 call fullbz(crys,syms,kg,1,skip_checkbz,wigner_seitz=.true.,paranoid=.true.) 283 else 284 call fullbz(crys,syms,kg,syms%ntran,skip_checkbz,wigner_seitz=.true.,paranoid=.true.) 285 endif 286 tmpfn='WFN_co' 287 if (.not. skip_checkbz.and..not.xct%patched_sampling_co) then 288 call checkbz(kg%nf,kg%f,kp%kgrid,kp%shift,crys%bdot, & 289 tmpfn,'k',.true.,xct%freplacebz,xct%fwritebz) 290 endif 291 call timacc(7,2) 292 xct%nkpt_co=kg%nf 293 294 if (peinf%verb_high .and. peinf%inode==0) then 295 write(6,'(/1x,a6,14x,a7,12x,2(1x,a6),3x,a3)') 'i', 'k-point', 'indr', 'itran', 'kg0' 296 write(6,'(1x,6("-"),1x,32("-"),2(1x,6("-")),1x,8("-"))') 297 do ii=1,kg%nf 298 write(6,'(1x,i6,3(1x,f10.6),2(1x,i6),3(1x,i2))') & 299 ii, kg%f(:,ii), kg%indr(ii), kg%itran(ii), kg%kg0(:,ii) 300 enddo 301 endif 302 303!------------------------------------------------------------------------ 304! If there is a finite center-of-mass momentum, Q, find mapping between k 305! and k+Q 306 307 SAFE_ALLOCATE(xct%indexq,(kg%nf)) 308 if (xct%qflag.eq.1) then 309 do ik=1,kg%nf 310 xct%indexq(ik) = ik 311 enddo 312 endif 313 314 SAFE_ALLOCATE(kgq%f,(3,kg%nf)) 315 SAFE_ALLOCATE(kgq%kg0,(3,kg%nf)) 316 SAFE_ALLOCATE(kgq%indr,(kg%nf)) 317 318!------------------------------------------------ 319! Distribute vcks-quadruplets among the PEs 320 321 call logit('input_kernel: calling distrib_kernel') 322 if (xct%qflag.eq.1) then 323 call distrib_kernel(xct,kp%ngkmax,kg,kg,gvec) 324 endif 325 if (peinf%verb_debug .and. peinf%inode==0) then 326 write(6,'(/1x,a,3(1x,i0)/)') "Allocating isort", peinf%iownwfv(peinf%inode+1), peinf%iownwfc(peinf%inode+1), gvec%ng 327 endif 328 if (xct%qflag.eq.1) then 329 SAFE_ALLOCATE(intwfnv%cg, (kp%ngkmax,peinf%iownwfv(peinf%inode+1),kp%nspin*kp%nspinor)) 330 endif 331 SAFE_ALLOCATE(intwfnc%cg, (kp%ngkmax,peinf%iownwfc(peinf%inode+1),kp%nspin*kp%nspinor)) 332 if (xct%qflag.eq.1) then 333 SAFE_ALLOCATE(intwfnv%isort, (gvec%ng,peinf%iownwfk(peinf%inode+1))) 334 SAFE_ALLOCATE(intwfnv%ng, (peinf%iownwfk(peinf%inode+1))) 335 endif 336 SAFE_ALLOCATE(intwfnc%isort, (gvec%ng,peinf%iownwfk(peinf%inode+1))) 337 SAFE_ALLOCATE(intwfnc%ng, (peinf%iownwfk(peinf%inode+1))) 338 intwfnv%nspin=kp%nspin 339 intwfnv%nspinor=kp%nspinor 340 intwfnc%nspin=kp%nspin 341 intwfnc%nspinor=kp%nspinor 342 343#ifdef DEBUG 344 if (peinf%inode.eq.0) then 345 write(6,*) 'kp%ngkmax',kp%ngkmax 346 endif 347#endif 348 349!------------------------------------------------ 350! Begin loop that distributes wave functions 351 352 SAFE_ALLOCATE(wfnv%isort, (gvec%ng)) 353 SAFE_ALLOCATE(wfnc%isort, (gvec%ng)) 354 wfnv%nspin=kp%nspin 355 wfnv%nspinor=kp%nspinor 356 wfnc%nspin=kp%nspin 357 wfnc%nspinor=kp%nspinor 358 359 do irk=1,kp%nrk 360 irks = 0 361 do ii=1,kg%nr 362 if (irk == indxk(ii)) then 363 irks=ii 364 exit 365 endif 366 enddo 367 368 SAFE_ALLOCATE(gvec_kpt%components, (3, kp%ngk(irk))) 369 call read_binary_gvectors(25, kp%ngk(irk), kp%ngk(irk), gvec_kpt%components) 370 371 SAFE_ALLOCATE(cg, (kp%ngk(irk),kp%nspin*kp%nspinor)) 372 if(irks > 0) then 373 do ii = 1, kp%ngk(irk) 374 call findvector(wfnv%isort(ii), gvec_kpt%components(:, ii), gvec) 375 if(wfnv%isort(ii) == 0) call die('could not find gvec', only_root_writes=.true.) 376 enddo 377 378 wfnv%ng=kp%ngk(irk) 379 wfnc%ng=kp%ngk(irk) 380 SAFE_ALLOCATE(cgarray, (kp%ngk(irk))) 381 wfnc%isort(1:gvec%ng)=wfnv%isort(1:gvec%ng) 382 ickmem=0 383 endif 384 385! write(6,*) peinf%inode, 'loop wfnup', irk 386 387! Loop Over Bands 388 389 do ii=1,kp%mnband 390 391 call read_binary_data(25, kp%ngk(irk), kp%ngk(irk), kp%nspin*kp%nspinor, cg,bcast=.false.) 392 393! If we do not need this band, skip it... 394 if(irks == 0) cycle 395! 396 do is=1, kp%nspin 397 if (ii .ge. (kp%ifmax(irk,is)-xct%nvb_co+1) .and. ii .le. (kp%ifmax(irk,is)+xct%ncb_co)) then 398 do isp=1, kp%nspinor 399 if (peinf%inode.eq.0) then 400 ! Check normalization of this band 401 call checknorm('WFN_co',ii,irks,kp%ngk(irk),kp%nspin,kp%nspinor,cg(:,:)) 402 do ig=1, kp%ngk(irk) 403 cgarray(ig)=cg(ig, is*isp) 404 end do 405 end if 406! write(6,*) peinf%inode, 'Read', irk, ii 407 408!----------------------- 409! If ii is one of the selected valence band... 410 if((ii.le.kp%ifmax(irk,is)).and. & 411 (ii.ge.(kp%ifmax(irk,is)-xct%nvb_co+1))) then 412 413 ! Skip current valence band if valence bands will be read from WFNq_co 414 if (xct%qflag.eq.0) cycle 415 416! write(6,*) peinf%inode, 'in val',irk,ii 417 418#ifdef MPI 419 call MPI_BCAST(cgarray,kp%ngk(irk),MPI_SCALAR,0, & 420 MPI_COMM_WORLD,mpierr) 421#endif 422 423 ivband=kp%ifmax(irk,is)-ii+1 424 iwritev=peinf%ipev(peinf%inode+1,ivband,irk) 425 if (xct%qflag.eq.1) then 426 iwritek=peinf%ipek(peinf%inode+1,irk) 427 else if (xct%qflag.eq.2) then 428 iwritek=peinf%ipekq(peinf%inode+1,irk) 429 endif 430 431! write(6,*) peinf%inode, 'bcast val',irk,ii,ivband, 432! > iwritev 433 if(iwritev.ne.0) then 434 intwfnv%cg(1:wfnv%ng,iwritev,is*isp)=cgarray(1:kp%ngk(irk)) 435 endif 436 if(iwritek.ne.0) then 437 intwfnv%isort(1:gvec%ng,iwritek)=wfnv%isort(1:gvec%ng) 438 intwfnv%ng(iwritek)=kp%ngk(irk) 439 endif 440 441! write(*,'(a, 4i4)') 'valence', peinf%inode+1, ivband, irk, iwritev 442! write(*, '(a, 6i5)') "valence", ivband, ii, is, peinf%inode, irk, iwritev 443! write(6,*) peinf%inode, 'wrote val',irk,ii 444 445 endif !ii is one of the selected valence band 446 447!----------------------- 448! If ii is one of the selected conduction band... 449 450 if((ii.ge.(kp%ifmax(irk,is)+1)) & 451 .and.(ii.le.(kp%ifmax(irk,is)+xct%ncb_co))) then 452 453#ifdef MPI 454 call MPI_BCAST(cgarray,kp%ngk(irk),MPI_SCALAR,0, & 455 MPI_COMM_WORLD,mpierr) 456#endif 457 458 icband=ii-kp%ifmax(irk,is) 459 iwritec=peinf%ipec(peinf%inode+1,icband,irk) 460 iwritek=peinf%ipek(peinf%inode+1,irk) 461! write(6,*) 'icband',ii,icband,iwritec 462 if(iwritec.ne.0) then 463 intwfnc%cg(1:wfnc%ng,iwritec,is*isp)=cgarray(1:kp%ngk(irk)) 464 endif 465 if(iwritek.ne.0) then 466 intwfnc%isort(1:gvec%ng,iwritek)=wfnc%isort(1:gvec%ng) 467 intwfnc%ng(iwritek)=kp%ngk(irk) 468 endif 469 470! write(*,'(a, 4i4)') 'conduction', peinf%inode+1, icband, irk, iwritec 471! write(*, '(a, 6i5)') "conduction", icband, ii, is, peinf%inode, irk, iwritec 472 473 endif ! ii is one of the selected conduction bands 474 end do ! isp 475 end if 476 end do ! is 477 if (ii>maxval(kp%ifmax)+xct%ncb_co .and. irk==kp%nrk) then 478 exit 479 endif 480 enddo ! ii (loop on bands) 481 482 SAFE_DEALLOCATE(cg) 483 SAFE_DEALLOCATE(cgarray) 484 485 486 enddo !end loop over k-points 487 488! Write out info about xtal 489 490 if (peinf%inode.eq.0) then 491 write(6,'(/1x,a)') 'Crystal wavefunctions read from file WFN_co:' 492 write(6,'(1x,a,i0)') '- Number of k-points in WFN_co: ', kp%nrk 493 write(6,'(1x,a,i0)') '- Number of k-points in the full BZ of WFN_co: ', kg%nf 494 if (peinf%verb_high) then 495 write(6,'(1x,a)') '- K-points:' 496 write(6,'(1(2x,3(1x,f10.6)))') kg%r(1:3,1:kg%nr) 497 endif 498 call close_file(25) 499 endif ! node 0 500 501 502 SAFE_DEALLOCATE_P(wfnv%isort) 503 SAFE_DEALLOCATE_P(wfnc%isort) 504 SAFE_DEALLOCATE_P(gvec_kpt%components) 505 SAFE_DEALLOCATE(indxk) 506 SAFE_DEALLOCATE(indxkq) 507 508 POP_SUB(input_kernel) 509 510 return 511 512contains 513 514 ! FHJ: Figure out xct%ecute and xct%ecutg, if the user didn`t specify them. 515 subroutine get_ecut() 516#ifdef HDF5 517 integer :: ng, nq, nfreq, nfreq_imag, nmtx_max, qgrid(3), freq_dep 518#endif 519 real(DP) :: ecuts 520 !real(DP) :: raux 521 522 PUSH_SUB(input_kernel.get_ecut) 523 524 if (xct%theory/=0) then 525 xct%ecute = 0 526 if (xct%ecutg < TOL_ZERO) xct%ecutg = kp%ecutwfc 527 return 528 endif 529 if (peinf%inode==0) then 530#ifdef HDF5 531 if(xct%use_hdf5) then 532 call read_eps_grid_sizes_hdf5(ng, nq, ecuts, nfreq, & 533 nfreq_imag, nmtx_max, qgrid, freq_dep, 'eps0mat.h5') 534 else 535#endif 536 call open_file(10, file='eps0mat', form='unformatted', status='old') 537 read(10) 538 read(10) 539 read(10) 540 read(10) 541 read(10) 542 read(10) 543 read(10) ecuts 544 call close_file(10) 545 endif 546#ifdef HDF5 547 endif 548#endif 549#ifdef MPI 550 call MPI_Bcast(ecuts, 1, MPI_REAL_DP, 0, MPI_COMM_WORLD, mpierr) 551#endif 552 if (xct%ecute < TOL_ZERO) xct%ecute = ecuts 553 !raux = maxval(crys%bdot) 554 !if (xct%ecutg < TOL_ZERO) xct%ecutg = xct%ecute + raux 555 if (xct%ecutg < TOL_ZERO) xct%ecutg = kp%ecutwfc 556 557 POP_SUB(input_kernel.get_ecut) 558 559 end subroutine get_ecut 560 561end subroutine input_kernel 562