TABLE OF CONTENTS
ABINIT/m_outqmc [ Modules ]
NAME
m_outqmc
FUNCTION
Interface with Cambridge quantum Monte Carlo program 'CASINO' See www.tcm.phy.cam.ac.uk/~mdt26/casino.html for more details. M.D.Towler (mdt26 at cam.ac.uk) November 2003 N.D.M.Hine (nicholas.hine at imperial.ac.uk) November 2004
COPYRIGHT
Copyright (C) 1998-2024 ABINIT group (DCA, XG, GMR, JYR, MKV, MT, FJ) 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.
SOURCE
20 #if defined HAVE_CONFIG_H 21 #include "config.h" 22 #endif 23 24 #include "abi_common.h" 25 26 module m_outqmc 27 28 use defs_basis 29 use m_errors 30 use m_abicore 31 use m_xmpi 32 use m_hdr 33 use m_dtset 34 35 36 use defs_datatypes, only : pseudopotential_type 37 use defs_abitypes, only : mpi_type 38 use m_io_tools, only : get_unit 39 use m_geometry, only : xred2xcart 40 use m_results_gs , only : results_gs_type 41 42 implicit none 43 44 private 45 public :: outqmc 46 47 contains
ABINIT/r2s [ Functions ]
NAME
r2s
FUNCTION
Converts real variable with arbitrary format to string that can be trimmed and printed in the middle of a sentence without introducing large amounts of white space, as you would if you did write(std_out,'(f12.6)')12.0 or similar. Note you need to pass through the format string e.g. f12.6. Calling routine is intended to include something like: USE utilities REAL(dp) r r=12._dp tmpr=r2s(r,'(f12.6)') write(std_out,*)'Real number ',trim(tmpr),' with words at the end.' Note : DON'T USE R2S IN A WRITE STATEMENT SINCE THIS IS ILLEGAL IN FORTRAN90 (ALTHOUGH NOT IN FORTRAN200X). IF ANYONE HAS TIME, FEEL FREE TO WRITE A VERSION OF THIS WHICH ISN'T ILLEGAL - SIMILAR TO I2S ABOVE - SO THAT PEOPLE WHO HAVEN'T READ THIS NOTE DON'T FEEL TEMPTED TO CALL R2S IN A WRITE STATEMENT.
INPUTS
OUTPUT
SIDE EFFECTS
SOURCE
596 function r2s(r,real_format) 597 598 !Arguments ------------------------- 599 real(dp),intent(in) :: r 600 character(len=*),intent(in) :: real_format 601 character(len=80) :: r2s 602 603 ! ********************************************************************* 604 605 if(len(r2s)>0)then 606 write(r2s,real_format)r 607 r2s=adjustl(r2s) 608 end if 609 610 end function r2s
m_outqmc/i2s [ Functions ]
[ Top ] [ m_outqmc ] [ Functions ]
NAME
i2s
FUNCTION
Convert integers to left justified strings that can be printed in the middle of a sentence without introducing large amounts of white space. Calling routine is intended to include something like: integer i i=12 write(std_out,*)'Integer number ',trim(i2s(i)),' with words at the end.'
INPUTS
OUTPUT
SOURCE
521 function i2s(n) 522 523 !Arguments ---------------------- 524 integer, intent(in) :: n 525 character(len=20) :: i2s 526 527 !Local variables ---------------- 528 integer :: i,j 529 character :: tmp,sign 530 531 ! ********************************************************************* 532 533 if(n==0)then 534 i2s='0' ; return 535 end if 536 sign=' ' ; if(n<0)sign='-' 537 538 do i=1,len(i2s) 539 i2s(i:i)=' ' 540 end do 541 542 i=abs(n) 543 do j=1,len(i2s) 544 if(i==0)exit 545 i2s(j:j)=achar(ichar('0')+mod(i,10)) 546 i=i/10 547 end do 548 549 i=1 ; j=len_trim(i2s) 550 do 551 if(i>=j)exit 552 tmp=i2s(j:j) 553 i2s(j:j)=i2s(i:i) 554 i2s(i:i)=tmp 555 i=i+1 556 j=j-1 557 end do 558 559 i2s=trim(sign)//i2s 560 561 end function i2s
m_outqmc/outqmc [ Functions ]
[ Top ] [ m_outqmc ] [ Functions ]
NAME
outqmc
FUNCTION
Write the wave function to a file in 'pwfn.data' format. This file can be read by the Cambridge quantum Monte Carlo program 'CASINO' and used as trial wave function input for a variational or diffusion Monte Carlo calculation. See www.tcm.phy.cam.ac.uk/~mdt26/casino.html for more details.
INPUTS
cg(2,mcg)=wavefunction coefficients dtset <type(dataset_type)>=all input variables for this dataset eigen(mband*nkpt*nsppol)=array for holding eigenvalues (hartree) gprimd(3,3)=dimensional primitive translations for reciprocal space hdr <type(hdr_type)>=the header of wf, den and pot files kg(3,mpw*mkmem)=reduced planewave coordinates. mcg=size of wave-functions array (cg) =mpw*nspinor*mband*mkmem*nsppol mpi_enreg=information about MPI parallelization npwarr(nkpt)=number of planewaves in basis and on boundary for each k occ(mband*nkpt*nsppol)=occupation number for each band and k psps <type(pseudopotential_type)>=all the information about psps results_gs <type(results_gs_type)>=results (energy and its components, forces and its components, the stress tensor) of a ground-state computation.
OUTPUT
Writes the file pwfn.data
SOURCE
80 subroutine outqmc(cg,dtset,eigen,gprimd,hdr,kg,mcg,mpi_enreg,npwarr,occ,psps,results_gs) 81 82 !Arguments ------------------------------- 83 !scalars 84 integer :: mcg 85 type(dataset_type) :: dtset 86 type(hdr_type) :: hdr 87 type(mpi_type) :: mpi_enreg 88 type(pseudopotential_type) :: psps 89 type(results_gs_type) :: results_gs 90 !arrays 91 integer :: kg(3,dtset%mpw*dtset%mkmem),npwarr(dtset%nkpt) 92 real(dp) :: cg(2,mcg) 93 real(dp) :: eigen(dtset%mband*dtset%nkpt*dtset%nsppol) 94 real(dp) :: gprimd(3,3),occ(dtset%mband*dtset%nkpt*dtset%nsppol) 95 96 !Local variables ------------------------- 97 !scalars 98 integer,parameter :: r2s_length=80 99 integer :: comm,iband,iband_kpt_shift,iband_sppol_shift,icg,icg_shift,icgfull 100 integer :: icgfull_shift,ii,ikg,ikg_shift,ikgfull,ikpt,ikpt_shift,io 101 integer :: iocc,isppol,jj,me,nband1,nband2,nelecs,nkgfull,ierr 102 real(dp) :: norm 103 logical :: am_master,foundkg 104 character(50) :: pwfnfilename 105 character(500) :: message 106 character(80) :: dft_functional,pseudo_name,pseudo_type 107 character(r2s_length) :: tmpr,tmpr2,tmpr3 108 !arrays 109 integer :: kgfull(3,dtset%mpw*dtset%mkmem),kgmap(dtset%mpw*dtset%mkmem) 110 real(dp) :: gcart_qmc(3),kptscart_qmc(3,dtset%nkpt) 111 real(dp),allocatable :: cgfull(:,:),xcart_qmc(:,:) 112 ! ********************************************************************* 113 114 !Go away if I am not the master node. 115 write(message,'(a,a)')ch10,' outqmc: enter ' 116 call wrtout(ab_out,message,'PERS') 117 118 am_master=.true. 119 if(xmpi_paral==1 .or. mpi_enreg%paral_kgb==1)then 120 comm=mpi_enreg%comm_cell 121 me=xmpi_comm_rank(comm) 122 if(me/=0) am_master=.false. 123 end if 124 if(.not.am_master)return 125 126 if(mpi_enreg%paral_spinor==1)then 127 message = ' Parallelization over spinors is not currently supported' 128 ABI_ERROR(message) 129 end if 130 131 !write(std_out,*)ch10,'outqmc: DEBUG: dtset%ndtset = ',dtset%ndtset,ch10 132 !Open CASINO pwfn.data file 133 if (dtset%ndtset<2)then 134 pwfnfilename='pwfn.data' 135 else 136 pwfnfilename='pwfn'//trim(i2s(dtset%jdtset))//'.data' 137 end if 138 call wrtout(ab_out,' outqmc: will open CASINO file: '//trim(pwfnfilename),'PERS') 139 140 io = get_unit() 141 open(io,file=pwfnfilename,form='formatted',recl=300,status='unknown',iostat=ierr) 142 143 if(ierr/=0)then 144 ABI_ERROR("Unable to open file: "//trim(pwfnfilename)) 145 end if 146 147 !Check if the full set of k vectors has been used in this calculation 148 if (dtset%kptopt==1) then 149 close(io,status='delete') 150 write(message,'(3a)')' outqmc: ERROR - kptopt=1 so k-points have been ',& 151 & 'generated in the irreducible Brillouin Zone only. ',& 152 & 'Set kptopt=2 to obtain full set of k-points.' 153 ABI_ERROR(message) 154 end if 155 156 !Check if the full set of G vectors has been used in this calculation 157 do ikpt=1,dtset%nkpt 158 if (dtset%istwfk(ikpt)/=1) then 159 close(io,status='delete') 160 write(message,'(a,i5,a,i2,a,a,a,a,a)')& 161 & ' istwfk(',ikpt,')=',dtset%istwfk(ikpt),' (ie /= 1) so some ',& 162 & 'G-vectors may not be present.',ch10,' Set istwfk=1 for each ',& 163 & 'k-point to obtain full set.' 164 ABI_ERROR(message) 165 end if 166 end do !ikpt 167 168 write(message,'(a)')' outqmc: QMC trial wave function file for CASINO generated by ABINIT' 169 call wrtout(ab_out,message,'PERS') 170 171 !Write the required quantities to pwfn.data 172 write(io,"('QMC trial wave function file for CASINO generated by ABINIT (www.abinit.org).')") 173 174 write(io,fmt="(/'BASIC INFO'/'----------')") 175 write(io,fmt="('Generated by:')") 176 write(io,fmt="(' ABINIT ',a)")trim(hdr%codvsn) 177 write(io,fmt="('Method:'/' DFT')") 178 179 write(io,fmt="('DFT Functional:')") 180 select case (dtset%ixc) 181 case(0) 182 dft_functional='No exchange-correlation.' 183 case(1) 184 dft_functional='L(S)DA (Teter/Pade parametrization)' 185 case(2) 186 dft_functional='LDA (Perdew-Zunger-Ceperley-Alder parametrization)' 187 case(3) 188 dft_functional='LDA (old Teter rational polynomial parametrization)' 189 case(4) 190 dft_functional='LDA (Wigner)' 191 case(5) 192 dft_functional='LDA (Hedin-Lundqvist)' 193 case(6) 194 dft_functional='LDA (X-alpha)' 195 case(7) 196 dft_functional='L(S)DA (Perdew-Wang 92)' 197 case(8) 198 dft_functional='L(S)DA (Perdew-Wang 92, exchange-only)' 199 case(9) 200 dft_functional='L(S)DA (Perdew-Wang 92, exchange- and RPA-correlation)' 201 case(10) 202 dft_functional='Diff. between ixc=7 and 9; use with accurate RPA corr. energy' 203 case(11) 204 dft_functional='GGA (Perdew-Burke-Ernzerhof)' 205 case(12) 206 dft_functional='GGA (Perdew-Burke-Ernzerhof, exchange-only)' 207 case(13) 208 dft_functional='GGA (potential: van Leeuwen-Baerends ; energy: Perdew-Wang 92)' 209 case(14) 210 dft_functional='GGA (RPBE of Zhang and Yang)' 211 case(15) 212 dft_functional='GGA (RPBE of Hammer, Hansen and Norskov)' 213 case(16) 214 dft_functional='GGA (HTCH)' 215 case(17) 216 dft_functional='Not defined (as of 11/2003).' 217 case(18) 218 dft_functional='Not defined (as of 11/2003).' 219 case(19) 220 dft_functional='Not defined (as of 11/2003).' 221 case(20) 222 dft_functional='Fermi-Amaldi xc for TDDFT' 223 case(21) 224 dft_functional='Fermi-Amaldi xc for TDDFT with LDA xc kernel' 225 case(22) 226 dft_functional='Fermi-Amaldi xc for TDDFT with Burke-Petersilka-Gross hybrid xc kernel' 227 case default 228 dft_functional='Unknown type.' 229 end select 230 write(io,"(' ABINIT ixc = ',i2,' : ',a)")dtset%ixc,trim(dft_functional) 231 232 write(io,"('Pseudopotential (of first atom type)')") 233 select case(psps%pspcod(1)) 234 case(1) 235 pseudo_type='ABINIT type 1' ; pseudo_name='Troullier-Martins' 236 case(2) 237 pseudo_type='ABINIT type 2' ; pseudo_name='Goedecker-Teter-Hutter (GTH)' 238 case(3) 239 pseudo_type='ABINIT type 3' ; pseudo_name='Hartwigsen-Goedecker-Hutter' 240 case(4) 241 pseudo_type='ABINIT type 4' 242 pseudo_name='Teter pseudo generated using the ATOM code' 243 case(5) 244 pseudo_type='ABINIT type 5' 245 pseudo_name='"Phoney" pseudo built on a Hamman grid' 246 case(6) 247 pseudo_type='ABINIT type 6' 248 pseudo_name='Fritz-Haber Institut (Troullier Martins)' 249 case default 250 pseudo_type='Unknown pseudopotential type (as of 11/2003).' ; pseudo_name='' 251 end select 252 if(dtset%ixc<10)then 253 write(io,"(1x,a,2x,': ',a)")trim(pseudo_type),trim(pseudo_name) 254 else 255 write(io,"(1x,a,3x,': ',a)")trim(pseudo_type),trim(pseudo_name) 256 end if 257 258 write(io,"('Plane wave cutoff (au)')") 259 tmpr=r2s(hdr%ecut_eff,'(f12.3)') 260 write(io,'(1x,a)')trim(tmpr) 261 262 write(io,"('Spin polarized:')") 263 select case(dtset%nspden) 264 case(1) 265 write(io,"(' .false.')") 266 case(2) 267 write(io,"(' .true.')") 268 case(4) 269 close(io,status='delete') 270 write(message,'(a)')' outqmc: ERROR - nspden=4 but CASINO cannot yet deal with non-collinear spins.' 271 ABI_ERROR(message) 272 case default 273 close(io,status='delete') 274 ABI_ERROR('Unrecognized value of nspden.') 275 end select 276 277 write(io,"('Total energy (au per primitive cell)')") 278 tmpr=r2s(results_gs%etotal,'(f24.14)') 279 write(io,'(1x,a)')trim(tmpr) 280 write(io,"('Kinetic energy')") 281 tmpr=r2s(results_gs%energies%e_kinetic,'(f24.14)') 282 write(io,'(1x,a)')trim(tmpr) 283 write(io,"('Local potential energy (Local pseudopotential energy eei + pseudopotential core-core energy eii)')") 284 tmpr=r2s((results_gs%energies%e_localpsp+results_gs%energies%e_corepsp),'(f24.14)') 285 write(io,'(1x,a)')trim(tmpr) 286 write(io,"('Non-local potential energy')") 287 tmpr=r2s(results_gs%energies%e_nlpsp_vfock,'(f24.14)') 288 write(io,'(1x,a)')trim(tmpr) 289 write(io,"('Electron-electron energy (Hartree Energy + Exchange-Correlation Energy)')") 290 tmpr=r2s((results_gs%energies%e_hartree+results_gs%energies%e_xc),'(f24.14)') 291 write(io,'(1x,a)')trim(tmpr) 292 write(io,"('Ion-ion energy')") 293 tmpr=r2s(results_gs%energies%e_ewald,'(f24.14)') 294 write(io,'(1x,a)')trim(tmpr) 295 write(io,"('Number of electrons per primitive cell')") 296 297 nelecs=0 298 do ii=1,dtset%natom 299 nelecs=nelecs+psps%ziontypat(dtset%typat(ii)) 300 end do 301 302 write(io,'(1x,i3)')nelecs 303 write(io,*) 304 write(io,"('GEOMETRY'/'--------')") 305 write(io,"('Number of atoms per primitive cell')") 306 write(io,'(1x,i3)')dtset%natom 307 write(io,"('Atomic numbers and positions of atoms (au)')") 308 309 ABI_MALLOC(xcart_qmc,(3,dtset%natom)) 310 call xred2xcart(dtset%natom,hdr%rprimd,xcart_qmc,hdr%xred) 311 do ii=1,dtset%natom 312 tmpr=r2s(xcart_qmc(1,ii),'(f24.14)') 313 tmpr2=r2s(xcart_qmc(2,ii),'(f24.14)') 314 tmpr3=r2s(xcart_qmc(3,ii),'(f24.14)') 315 jj=psps%znucltypat(dtset%typat(ii)) 316 write(io,'(1x,i2,3(1x,a))')jj,trim(tmpr),trim(tmpr2),trim(tmpr3) 317 end do 318 ABI_FREE(xcart_qmc) 319 320 write(io,"('Primitive lattice vectors (au)')") 321 do ii=1,3 322 tmpr=r2s(hdr%rprimd(1,ii),'(f24.14)') 323 tmpr2=r2s(hdr%rprimd(2,ii),'(f24.14)') 324 tmpr3=r2s(hdr%rprimd(3,ii),'(f24.14)') 325 write(io,'(3(1x,a))')trim(tmpr),trim(tmpr2),trim(tmpr3) 326 end do 327 328 !Copy the G vectors for the first k point into kgfull 329 ikgfull=0 330 do ikg=1,npwarr(1) 331 ikgfull=ikgfull+1 332 kgfull(1:3,ikgfull) = kg(1:3,ikg) 333 kgmap(ikg)=ikgfull 334 end do 335 ikg_shift = npwarr(1) 336 !Go through the other k points and look for any G vectors that haven't 337 !already been found and add them to the end of kgfull 338 do ikpt=2,dtset%nkpt 339 do ikg=ikg_shift,ikg_shift+npwarr(ikpt) 340 foundkg = .false. 341 do ii=1,ikgfull 342 if(kg(1,ikg)==kgfull(1,ii).and.kg(2,ikg)==kgfull(2,ii) & 343 & .and.kg(3,ikg)==kgfull(3,ii)) then 344 foundkg=.true. 345 kgmap(ikg)=ii 346 exit 347 end if 348 end do 349 if(.not.foundkg)then 350 ikgfull=ikgfull+1 351 kgfull(1:3,ikgfull)=kg(1:3,ikg) 352 kgmap(ikg)=ikgfull 353 end if 354 end do 355 ikg_shift=ikg_shift+npwarr(ikpt) 356 end do 357 nkgfull=ikgfull 358 359 write(io,*) 360 write(io,"('G VECTORS'/'---------')") 361 write(io,"('Number of G-vectors')") 362 write(io,'(1x,i7)')nkgfull 363 write(io,"('Gx Gy Gz (au)')") 364 365 do ikgfull=1,nkgfull 366 gcart_qmc=2*pi*(kgfull(1,ikgfull)*gprimd(1:3,1)& 367 & +kgfull(2,ikgfull)*gprimd(1:3,2)+kgfull(3,ikgfull)*gprimd(1:3,3)) 368 write(io,*)gcart_qmc(1:3) ! '(3e26.16)' 369 end do 370 371 !Populate the cgfull array, using the mapping in kgmap between the 372 !coefficients for kg in the per-kpoint list and the ones in the full list 373 !The number of xxx_shift's might seem excessive but the re-ordering of the 374 !list from (spin, kpt, band, kg) to (kpt, spin, band, kgfull) is quite 375 !complicated 376 ABI_MALLOC(cgfull,(2,nkgfull*dtset%nspinor*dtset%nsppol*dtset%mband*dtset%nkpt)) 377 cgfull(1:2,1:nkgfull*dtset%nspinor*dtset%nsppol*dtset%mband*dtset%nkpt)=0 378 icg_shift=1 379 do isppol=1,dtset%nsppol 380 ikg_shift=1 381 if(isppol==2)then 382 ! Go back to the beginning of cgfull but skip the first set of isppol=1 bands 383 icgfull_shift=nkgfull*dtset%nband(1) 384 ikpt_shift=dtset%nkpt 385 else 386 icgfull_shift=0 ! Start at the beginning of cgfull 387 ikpt_shift=0 388 end if 389 do ikpt=1,dtset%nkpt 390 do iband=1,dtset%nband(ikpt+ikpt_shift) 391 ikg=ikg_shift 392 ! Find the index in Abinit's coefficient list 393 do icg=icg_shift,icg_shift+npwarr(ikpt)-1 394 ! Map it to an index in the full CASINO list with the mapping recorded 395 ! when kgfull was read in 396 icgfull = kgmap(ikg)+icgfull_shift 397 cgfull(1:2,icgfull)=cg(1:2,icg) 398 ikg=ikg+1 399 end do !icg 400 icg_shift=icg_shift+npwarr(ikpt) 401 icgfull_shift=icgfull_shift+nkgfull 402 end do !iband 403 if(isppol==2)then 404 ! Skip the isppol==1 bands 405 icgfull_shift=icgfull_shift+nkgfull*dtset%nband(ikpt) 406 else 407 if(dtset%nsppol==2)then 408 ! Skip the isppol==2 bands 409 icgfull_shift=icgfull_shift+nkgfull*dtset%nband(ikpt+dtset%nkpt) 410 else 411 icgfull_shift=icgfull_shift ! Nothing to be skipped 412 end if 413 end if 414 ikg_shift=ikg_shift+npwarr(ikpt) 415 end do !ikpt 416 end do !isppol 417 418 !See if each orbital is normalised and check for integer occupancy of orbitals. 419 !These are checked by CASINO and it will complain if they are not as expected. 420 icgfull_shift=1 421 ii=0 422 iocc=1 423 424 do ikpt=1,dtset%nkpt 425 do isppol=1,dtset%nsppol 426 do iband=1,dtset%nband(ikpt+(isppol-1)*dtset%nkpt) 427 if(occ(iocc)/=int(occ(iocc)))then 428 write(message,'(a,i5,a,i1,a,i5,a,f11.8,a)')& 429 & 'Non-integer occupation number for kpt ',ikpt,', sppol ',isppol,', band ',iband,': occ=',occ(iocc),'.' 430 ABI_WARNING(message) 431 end if 432 iocc=iocc+1 433 norm=0 434 do icgfull=icgfull_shift,icgfull_shift+nkgfull-1 435 norm=norm+cgfull(1,icgfull)**2+cgfull(2,icgfull)**2 436 end do !icgfull 437 icgfull_shift=icgfull_shift+nkgfull 438 if((norm<0.999).or.(norm>1.001))then 439 write(message,'(a,i5,a,i1,a,i5,a,f11.8,a)')& 440 & 'Incorrectly normalised orbital for kpt ',ikpt,', sppol ',isppol,', band ',iband,': norm=',norm,'.' 441 ABI_WARNING(message) 442 end if 443 end do !iband 444 end do !isppol 445 end do !ikpt 446 447 write(io,*) 448 write(io,"('WAVE FUNCTION'/'-------------')") 449 write(io,"('Number of k-points')") 450 write(io,'(1x,i5)') dtset%nkpt 451 452 do ikpt=1,dtset%nkpt 453 kptscart_qmc(1:3,ikpt)=2*pi*(dtset%kptns(1,ikpt)*gprimd(1:3,1)& 454 & +dtset%kptns(2,ikpt)*gprimd(1:3,2)+dtset%kptns(3,ikpt)*gprimd(1:3,3)) 455 end do !ikpt 456 457 iband_kpt_shift=0 458 icg_shift=0 459 icgfull_shift=0 460 do ikpt=1,dtset%nkpt 461 write(io,"('k-point # ; # of bands (up spin/down spin) ; k-point coords (au)')") 462 if(dtset%nsppol==2)then 463 nband1=dtset%nband(ikpt) 464 nband2=dtset%nband(ikpt+dtset%nkpt) 465 else 466 nband1=dtset%nband(ikpt) 467 nband2=0 468 end if 469 tmpr=r2s(kptscart_qmc(1,ikpt),'(f24.14)') 470 tmpr2=r2s(kptscart_qmc(2,ikpt),'(f24.14)') 471 tmpr3=r2s(kptscart_qmc(3,ikpt),'(f24.14)') 472 write(io,'(3(1x,i5),3(1x,a))')ikpt,nband1,nband2,trim(tmpr),trim(tmpr2), & 473 & trim(tmpr3) 474 do isppol=1,dtset%nsppol 475 if (isppol==2) then 476 iband_sppol_shift=sum(dtset%nband(1:dtset%nkpt)) 477 iband_kpt_shift=sum(dtset%nband((dtset%nkpt+1):(dtset%nkpt+ikpt-1))) 478 else 479 iband_sppol_shift=0 480 iband_kpt_shift=sum(dtset%nband(1:(ikpt-1))) 481 end if 482 do iband=1,dtset%nband(ikpt) 483 write(io,"('Band, spin, eigenvalue (au)')") 484 tmpr=r2s(eigen(iband_kpt_shift+iband_sppol_shift+iband),'(f24.14)') 485 write(io,'(2(1x,i5),1x,a)')iband,isppol,tmpr 486 write(io,"('Eigenvector coefficients')") 487 do icgfull=1,nkgfull 488 write(io,"(1x,'(',e23.16,',',e23.16,')')")cgfull(1:2,icgfull+icgfull_shift) 489 end do !icgfull 490 icgfull_shift=icgfull_shift+nkgfull 491 end do !iband 492 end do !isppol 493 end do !ikpt 494 495 close(io) 496 497 write(message,'(a,a)')' outqmc: done with writing of QMC trial wave function file for CASINO',ch10 498 call wrtout(ab_out,message,'PERS') 499 500 end subroutine outqmc