TABLE OF CONTENTS


ABINIT/m_outqmc [ Modules ]

[ Top ] [ 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-2018 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.

PARENTS

CHILDREN

SOURCE

25 #if defined HAVE_CONFIG_H
26 #include "config.h"
27 #endif
28 
29 #include "abi_common.h"
30 
31 module m_outqmc
32 
33  use defs_basis
34  use defs_datatypes
35  use defs_abitypes
36  use m_errors
37  use m_abicore
38  use m_xmpi
39 
40  use m_io_tools,    only : get_unit
41  use m_geometry,    only : xred2xcart
42  use m_results_gs , only : results_gs_type
43 
44  implicit none
45 
46  private
47  public ::  outqmc
48 
49 contains

ABINIT/r2s [ Functions ]

[ Top ] [ 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

622  function r2s(r,real_format)
623 
624 
625 !This section has been created automatically by the script Abilint (TD).
626 !Do not modify the following lines by hand.
627 #undef ABI_FUNC
628 #define ABI_FUNC 'r2s'
629 !End of the abilint section
630 
631  implicit none
632 
633 !Arguments -------------------------
634  real(dp),intent(in) :: r
635  character(len=*),intent(in) :: real_format
636  character(len=80) :: r2s
637 
638 ! *********************************************************************
639 
640  if(len(r2s)>0)then
641    write(r2s,real_format)r
642    r2s=adjustl(r2s)
643  end if
644 
645 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

538  function i2s(n)
539 
540 
541 !This section has been created automatically by the script Abilint (TD).
542 !Do not modify the following lines by hand.
543 #undef ABI_FUNC
544 #define ABI_FUNC 'i2s'
545 !End of the abilint section
546 
547  implicit none
548 
549 !Arguments ----------------------
550  integer, intent(in) :: n
551  character(len=20) :: i2s
552 
553 !Local variables ----------------
554  integer :: i,j
555  character :: tmp,sign
556 
557 ! *********************************************************************
558 
559  if(n==0)then
560    i2s='0' ; return
561  end if
562  sign=' ' ; if(n<0)sign='-'
563 
564  do i=1,len(i2s)
565    i2s(i:i)=' '
566  end do
567 
568  i=abs(n)
569  do j=1,len(i2s)
570    if(i==0)exit
571    i2s(j:j)=achar(ichar('0')+mod(i,10))
572    i=i/10
573  end do
574 
575  i=1 ; j=len_trim(i2s)
576  do
577    if(i>=j)exit
578    tmp=i2s(j:j)
579    i2s(j:j)=i2s(i:i)
580    i2s(i:i)=tmp
581    i=i+1
582    j=j-1
583  end do
584 
585  i2s=trim(sign)//i2s
586 
587 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

PARENTS

      gstate

CHILDREN

      wrtout,xred2xcart

SOURCE

 88 subroutine outqmc(cg,dtset,eigen,gprimd,hdr,kg,mcg,mpi_enreg,npwarr,occ,psps,results_gs)
 89 
 90 
 91 !This section has been created automatically by the script Abilint (TD).
 92 !Do not modify the following lines by hand.
 93 #undef ABI_FUNC
 94 #define ABI_FUNC 'outqmc'
 95 !End of the abilint section
 96 
 97  implicit none
 98 
 99 !Arguments -------------------------------
100 !scalars
101  integer :: mcg
102  type(dataset_type) :: dtset
103  type(hdr_type) :: hdr
104  type(mpi_type) :: mpi_enreg
105  type(pseudopotential_type) :: psps
106  type(results_gs_type) :: results_gs
107 !arrays
108  integer :: kg(3,dtset%mpw*dtset%mkmem),npwarr(dtset%nkpt)
109  real(dp) :: cg(2,mcg)
110  real(dp) :: eigen(dtset%mband*dtset%nkpt*dtset%nsppol)
111  real(dp) :: gprimd(3,3),occ(dtset%mband*dtset%nkpt*dtset%nsppol)
112 
113 !Local variables -------------------------
114 !scalars
115  integer,parameter :: r2s_length=80
116  integer :: comm,iband,iband_kpt_shift,iband_sppol_shift,icg,icg_shift,icgfull
117  integer :: icgfull_shift,ii,ikg,ikg_shift,ikgfull,ikpt,ikpt_shift,io
118  integer :: iocc,isppol,jj,me,nband1,nband2,nelecs,nkgfull,ierr
119  real(dp) :: norm
120  logical :: am_master,foundkg
121  character(50) :: pwfnfilename
122  character(500) :: message
123  character(80) :: dft_functional,pseudo_name,pseudo_type
124  character(r2s_length) :: tmpr,tmpr2,tmpr3
125 !arrays
126  integer :: kgfull(3,dtset%mpw*dtset%mkmem),kgmap(dtset%mpw*dtset%mkmem)
127  real(dp) :: gcart_qmc(3),kptscart_qmc(3,dtset%nkpt)
128  real(dp),allocatable :: cgfull(:,:),xcart_qmc(:,:)
129 ! *********************************************************************
130 
131 !Go away if I am not the master node.
132  write(message,'(a,a)')ch10,' outqmc: enter '
133  call wrtout(ab_out,message,'PERS')
134 
135  am_master=.true.
136  if(xmpi_paral==1 .or. mpi_enreg%paral_kgb==1)then
137    comm=mpi_enreg%comm_cell
138    me=xmpi_comm_rank(comm)
139    if(me/=0) am_master=.false.
140  end if
141  if(.not.am_master)return
142 
143  if(mpi_enreg%paral_spinor==1)then
144    message = ' Parallelization over spinors is not currently supported'
145    MSG_ERROR(message)
146  end if
147 
148 !write(std_out,*)ch10,'outqmc: DEBUG: dtset%ndtset = ',dtset%ndtset,ch10
149 !Open CASINO pwfn.data file
150  if (dtset%ndtset<2)then
151    pwfnfilename='pwfn.data'
152  else
153    pwfnfilename='pwfn'//trim(i2s(dtset%jdtset))//'.data'
154  end if
155  call wrtout(ab_out,' outqmc: will open CASINO file: '//trim(pwfnfilename),'PERS')
156 
157  io = get_unit()
158  open(io,file=pwfnfilename,form='formatted',recl=300,status='unknown',iostat=ierr)
159 
160  if(ierr/=0)then
161    MSG_ERROR("Unable to open file: "//trim(pwfnfilename))
162  end if
163 
164 !Check if the full set of k vectors has been used in this calculation
165  if (dtset%kptopt==1) then
166    close(io,status='delete')
167    write(message,'(3a)')' outqmc: ERROR - kptopt=1 so k-points have been ',&
168 &   'generated in the irreducible Brillouin Zone only. ',&
169 &   'Set kptopt=2 to obtain full set of k-points.'
170    MSG_ERROR(message)
171  end if
172 
173 !Check if the full set of G vectors has been used in this calculation
174  do ikpt=1,dtset%nkpt
175    if (dtset%istwfk(ikpt)/=1) then
176      close(io,status='delete')
177      write(message,'(a,i5,a,i2,a,a,a,a,a)')&
178 &     '  istwfk(',ikpt,')=',dtset%istwfk(ikpt),' (ie /= 1) so some ',&
179 &     'G-vectors may not be present.',ch10,'  Set istwfk=1 for each ',&
180 &     'k-point to obtain full set.'
181      MSG_ERROR(message)
182    end if
183  end do !ikpt
184 
185  write(message,'(a)')' outqmc: QMC trial wave function file for CASINO generated by ABINIT'
186  call wrtout(ab_out,message,'PERS')
187 
188 !Write the required quantities to pwfn.data
189  write(io,"('QMC trial wave function file for CASINO generated by ABINIT (www.abinit.org).')")
190 
191  write(io,fmt="(/'BASIC INFO'/'----------')")
192  write(io,fmt="('Generated by:')")
193  write(io,fmt="(' ABINIT ',a)")trim(hdr%codvsn)
194  write(io,fmt="('Method:'/' DFT')")
195 
196  write(io,fmt="('DFT Functional:')")
197  select case (dtset%ixc)
198  case(0)
199    dft_functional='No exchange-correlation.'
200  case(1)
201    dft_functional='L(S)DA (Teter/Pade parametrization)'
202  case(2)
203    dft_functional='LDA (Perdew-Zunger-Ceperley-Alder parametrization)'
204  case(3)
205    dft_functional='LDA (old Teter rational polynomial parametrization)'
206  case(4)
207    dft_functional='LDA (Wigner)'
208  case(5)
209    dft_functional='LDA (Hedin-Lundqvist)'
210  case(6)
211    dft_functional='LDA (X-alpha)'
212  case(7)
213    dft_functional='L(S)DA (Perdew-Wang 92)'
214  case(8)
215    dft_functional='L(S)DA (Perdew-Wang 92, exchange-only)'
216  case(9)
217    dft_functional='L(S)DA (Perdew-Wang 92, exchange- and RPA-correlation)'
218  case(10)
219    dft_functional='Diff. between ixc=7 and 9; use with accurate RPA corr. energy'
220  case(11)
221    dft_functional='GGA (Perdew-Burke-Ernzerhof)'
222  case(12)
223    dft_functional='GGA (Perdew-Burke-Ernzerhof, exchange-only)'
224  case(13)
225    dft_functional='GGA (potential: van Leeuwen-Baerends ; energy: Perdew-Wang 92)'
226  case(14)
227    dft_functional='GGA (RPBE of Zhang and Yang)'
228  case(15)
229    dft_functional='GGA (RPBE of Hammer, Hansen and Norskov)'
230  case(16)
231    dft_functional='GGA (HTCH)'
232  case(17)
233    dft_functional='Not defined (as of 11/2003).'
234  case(18)
235    dft_functional='Not defined (as of 11/2003).'
236  case(19)
237    dft_functional='Not defined (as of 11/2003).'
238  case(20)
239    dft_functional='Fermi-Amaldi xc for TDDFT'
240  case(21)
241    dft_functional='Fermi-Amaldi xc for TDDFT with LDA xc kernel'
242  case(22)
243    dft_functional='Fermi-Amaldi xc for TDDFT with Burke-Petersilka-Gross hybrid xc kernel'
244  case default
245    dft_functional='Unknown type.'
246  end select
247  write(io,"(' ABINIT ixc = ',i2,' : ',a)")dtset%ixc,trim(dft_functional)
248 
249  write(io,"('Pseudopotential (of first atom type)')")
250  select case(psps%pspcod(1))
251  case(1)
252    pseudo_type='ABINIT type 1' ; pseudo_name='Troullier-Martins'
253  case(2)
254    pseudo_type='ABINIT type 2' ; pseudo_name='Goedecker-Teter-Hutter (GTH)'
255  case(3)
256    pseudo_type='ABINIT type 3' ; pseudo_name='Hartwigsen-Goedecker-Hutter'
257  case(4)
258    pseudo_type='ABINIT type 4'
259    pseudo_name='Teter pseudo generated using the ATOM code'
260  case(5)
261    pseudo_type='ABINIT type 5'
262    pseudo_name='"Phoney" pseudo built on a Hamman grid'
263  case(6)
264    pseudo_type='ABINIT type 6'
265    pseudo_name='Fritz-Haber Institut (Troullier Martins)'
266  case default
267    pseudo_type='Unknown pseudopotential type (as of 11/2003).' ; pseudo_name=''
268  end select
269  if(dtset%ixc<10)then
270    write(io,"(1x,a,2x,': ',a)")trim(pseudo_type),trim(pseudo_name)
271  else
272    write(io,"(1x,a,3x,': ',a)")trim(pseudo_type),trim(pseudo_name)
273  end if
274 
275  write(io,"('Plane wave cutoff (au)')")
276  tmpr=r2s(hdr%ecut_eff,'(f12.3)')
277  write(io,'(1x,a)')trim(tmpr)
278 
279  write(io,"('Spin polarized:')")
280  select case(dtset%nspden)
281  case(1)
282    write(io,"(' .false.')")
283  case(2)
284    write(io,"(' .true.')")
285  case(4)
286    close(io,status='delete')
287    write(message,'(a)')' outqmc: ERROR - nspden=4 but CASINO cannot yet deal with non-collinear spins.'
288    MSG_ERROR(message)
289  case default
290    close(io,status='delete')
291    MSG_ERROR('Unrecognized value of nspden.')
292  end select
293 
294  write(io,"('Total energy (au per primitive cell)')")
295  tmpr=r2s(results_gs%etotal,'(f24.14)')
296  write(io,'(1x,a)')trim(tmpr)
297  write(io,"('Kinetic energy')")
298  tmpr=r2s(results_gs%energies%e_kinetic,'(f24.14)')
299  write(io,'(1x,a)')trim(tmpr)
300  write(io,"('Local potential energy (Local pseudopotential energy eei + pseudopotential core-core energy eii)')")
301  tmpr=r2s((results_gs%energies%e_localpsp+results_gs%energies%e_corepsp),'(f24.14)')
302  write(io,'(1x,a)')trim(tmpr)
303  write(io,"('Non-local potential energy')")
304  tmpr=r2s(results_gs%energies%e_nonlocalpsp,'(f24.14)')
305  write(io,'(1x,a)')trim(tmpr)
306  write(io,"('Electron-electron energy (Hartree Energy + Exchange-Correlation Energy)')")
307  tmpr=r2s((results_gs%energies%e_hartree+results_gs%energies%e_xc),'(f24.14)')
308  write(io,'(1x,a)')trim(tmpr)
309  write(io,"('Ion-ion energy')")
310  tmpr=r2s(results_gs%energies%e_ewald,'(f24.14)')
311  write(io,'(1x,a)')trim(tmpr)
312  write(io,"('Number of electrons per primitive cell')")
313 
314  nelecs=0
315  do ii=1,dtset%natom
316    nelecs=nelecs+psps%ziontypat(dtset%typat(ii))
317  end do
318 
319  write(io,'(1x,i3)')nelecs
320  write(io,*)
321  write(io,"('GEOMETRY'/'--------')")
322  write(io,"('Number of atoms per primitive cell')")
323  write(io,'(1x,i3)')dtset%natom
324  write(io,"('Atomic numbers and positions of atoms (au)')")
325 
326  ABI_ALLOCATE(xcart_qmc,(3,dtset%natom))
327  call xred2xcart(dtset%natom,hdr%rprimd,xcart_qmc,hdr%xred)
328  do ii=1,dtset%natom
329    tmpr=r2s(xcart_qmc(1,ii),'(f24.14)')
330    tmpr2=r2s(xcart_qmc(2,ii),'(f24.14)')
331    tmpr3=r2s(xcart_qmc(3,ii),'(f24.14)')
332    jj=psps%znucltypat(dtset%typat(ii))
333    write(io,'(1x,i2,3(1x,a))')jj,trim(tmpr),trim(tmpr2),trim(tmpr3)
334  end do
335  ABI_DEALLOCATE(xcart_qmc)
336 
337  write(io,"('Primitive lattice vectors (au)')")
338  do ii=1,3
339    tmpr=r2s(hdr%rprimd(1,ii),'(f24.14)')
340    tmpr2=r2s(hdr%rprimd(2,ii),'(f24.14)')
341    tmpr3=r2s(hdr%rprimd(3,ii),'(f24.14)')
342    write(io,'(3(1x,a))')trim(tmpr),trim(tmpr2),trim(tmpr3)
343  end do
344 
345 !Copy the G vectors for the first k point into kgfull
346  ikgfull=0
347  do ikg=1,npwarr(1)
348    ikgfull=ikgfull+1
349    kgfull(1:3,ikgfull) = kg(1:3,ikg)
350    kgmap(ikg)=ikgfull
351  end do
352  ikg_shift = npwarr(1)
353 !Go through the other k points and look for any G vectors that haven't
354 !already been found and add them to the end of kgfull
355  do ikpt=2,dtset%nkpt
356    do ikg=ikg_shift,ikg_shift+npwarr(ikpt)
357      foundkg = .false.
358      do ii=1,ikgfull
359        if(kg(1,ikg)==kgfull(1,ii).and.kg(2,ikg)==kgfull(2,ii) &
360 &       .and.kg(3,ikg)==kgfull(3,ii)) then
361          foundkg=.true.
362          kgmap(ikg)=ii
363          exit
364        end if
365      end do
366      if(.not.foundkg)then
367        ikgfull=ikgfull+1
368        kgfull(1:3,ikgfull)=kg(1:3,ikg)
369        kgmap(ikg)=ikgfull
370      end if
371    end do
372    ikg_shift=ikg_shift+npwarr(ikpt)
373  end do
374  nkgfull=ikgfull
375 
376  write(io,*)
377  write(io,"('G VECTORS'/'---------')")
378  write(io,"('Number of G-vectors')")
379  write(io,'(1x,i7)')nkgfull
380  write(io,"('Gx Gy Gz (au)')")
381 
382  do ikgfull=1,nkgfull
383    gcart_qmc=2*pi*(kgfull(1,ikgfull)*gprimd(1:3,1)&
384 &   +kgfull(2,ikgfull)*gprimd(1:3,2)+kgfull(3,ikgfull)*gprimd(1:3,3))
385    write(io,*)gcart_qmc(1:3) ! '(3e26.16)'
386  end do
387 
388 !Populate the cgfull array, using the mapping in kgmap between the
389 !coefficients for kg in the per-kpoint list and the ones in the full list
390 !The number of xxx_shift's might seem excessive but the re-ordering of the
391 !list from (spin, kpt, band, kg) to (kpt, spin, band, kgfull) is quite
392 !complicated
393  ABI_ALLOCATE(cgfull,(2,nkgfull*dtset%nspinor*dtset%nsppol*dtset%mband*dtset%nkpt))
394  cgfull(1:2,1:nkgfull*dtset%nspinor*dtset%nsppol*dtset%mband*dtset%nkpt)=0
395  icg_shift=1
396  do isppol=1,dtset%nsppol
397    ikg_shift=1
398    if(isppol==2)then
399 !    Go back to the beginning of cgfull but skip the first set of isppol=1 bands
400      icgfull_shift=nkgfull*dtset%nband(1)
401      ikpt_shift=dtset%nkpt
402    else
403      icgfull_shift=0 ! Start at the beginning of cgfull
404      ikpt_shift=0
405    end if
406    do ikpt=1,dtset%nkpt
407      do iband=1,dtset%nband(ikpt+ikpt_shift)
408        ikg=ikg_shift
409 !      Find the index in Abinit's coefficient list
410        do icg=icg_shift,icg_shift+npwarr(ikpt)-1
411 !        Map it to an index in the full CASINO list with the mapping recorded
412 !        when kgfull was read in
413          icgfull = kgmap(ikg)+icgfull_shift
414          cgfull(1:2,icgfull)=cg(1:2,icg)
415          ikg=ikg+1
416        end do !icg
417        icg_shift=icg_shift+npwarr(ikpt)
418        icgfull_shift=icgfull_shift+nkgfull
419      end do !iband
420      if(isppol==2)then
421 !      Skip the isppol==1 bands
422        icgfull_shift=icgfull_shift+nkgfull*dtset%nband(ikpt)
423      else
424        if(dtset%nsppol==2)then
425 !        Skip the isppol==2 bands
426          icgfull_shift=icgfull_shift+nkgfull*dtset%nband(ikpt+dtset%nkpt)
427        else
428          icgfull_shift=icgfull_shift ! Nothing to be skipped
429        end if
430      end if
431      ikg_shift=ikg_shift+npwarr(ikpt)
432    end do !ikpt
433  end do !isppol
434 
435 !See if each orbital is normalised and check for integer occupancy of orbitals.
436 !These are checked by CASINO and it will complain if they are not as expected.
437  icgfull_shift=1
438  ii=0
439  iocc=1
440 
441  do ikpt=1,dtset%nkpt
442    do isppol=1,dtset%nsppol
443      do iband=1,dtset%nband(ikpt+(isppol-1)*dtset%nkpt)
444        if(occ(iocc)/=int(occ(iocc)))then
445          write(message,'(a,i5,a,i1,a,i5,a,f11.8,a)')&
446 &         'Non-integer occupation number for kpt ',ikpt,', sppol ',isppol,', band ',iband,': occ=',occ(iocc),'.'
447          MSG_WARNING(message)
448        end if
449        iocc=iocc+1
450        norm=0
451        do icgfull=icgfull_shift,icgfull_shift+nkgfull-1
452          norm=norm+cgfull(1,icgfull)**2+cgfull(2,icgfull)**2
453        end do !icgfull
454        icgfull_shift=icgfull_shift+nkgfull
455        if((norm<0.999).or.(norm>1.001))then
456          write(message,'(a,i5,a,i1,a,i5,a,f11.8,a)')&
457 &         'Incorrectly normalised orbital for kpt ',ikpt,', sppol ',isppol,', band ',iband,': norm=',norm,'.'
458          MSG_WARNING(message)
459        end if
460      end do !iband
461    end do !isppol
462  end do !ikpt
463 
464  write(io,*)
465  write(io,"('WAVE FUNCTION'/'-------------')")
466  write(io,"('Number of k-points')")
467  write(io,'(1x,i5)') dtset%nkpt
468 
469  do ikpt=1,dtset%nkpt
470    kptscart_qmc(1:3,ikpt)=2*pi*(dtset%kptns(1,ikpt)*gprimd(1:3,1)&
471 &   +dtset%kptns(2,ikpt)*gprimd(1:3,2)+dtset%kptns(3,ikpt)*gprimd(1:3,3))
472  end do !ikpt
473 
474  iband_kpt_shift=0
475  icg_shift=0
476  icgfull_shift=0
477  do ikpt=1,dtset%nkpt
478    write(io,"('k-point # ; # of bands (up spin/down spin) ; k-point coords (au)')")
479    if(dtset%nsppol==2)then
480      nband1=dtset%nband(ikpt)
481      nband2=dtset%nband(ikpt+dtset%nkpt)
482    else
483      nband1=dtset%nband(ikpt)
484      nband2=0
485    end if
486    tmpr=r2s(kptscart_qmc(1,ikpt),'(f24.14)')
487    tmpr2=r2s(kptscart_qmc(2,ikpt),'(f24.14)')
488    tmpr3=r2s(kptscart_qmc(3,ikpt),'(f24.14)')
489    write(io,'(3(1x,i5),3(1x,a))')ikpt,nband1,nband2,trim(tmpr),trim(tmpr2), &
490 &   trim(tmpr3)
491    do isppol=1,dtset%nsppol
492      if (isppol==2) then
493        iband_sppol_shift=sum(dtset%nband(1:dtset%nkpt))
494        iband_kpt_shift=sum(dtset%nband((dtset%nkpt+1):(dtset%nkpt+ikpt-1)))
495      else
496        iband_sppol_shift=0
497        iband_kpt_shift=sum(dtset%nband(1:(ikpt-1)))
498      end if
499      do iband=1,dtset%nband(ikpt)
500        write(io,"('Band, spin, eigenvalue (au)')")
501        tmpr=r2s(eigen(iband_kpt_shift+iband_sppol_shift+iband),'(f24.14)')
502        write(io,'(2(1x,i5),1x,a)')iband,isppol,tmpr
503        write(io,"('Eigenvector coefficients')")
504        do icgfull=1,nkgfull
505          write(io,"(1x,'(',e23.16,',',e23.16,')')")cgfull(1:2,icgfull+icgfull_shift)
506        end do !icgfull
507        icgfull_shift=icgfull_shift+nkgfull
508      end do !iband
509    end do !isppol
510  end do !ikpt
511 
512  close(io)
513 
514  write(message,'(a,a)')' outqmc: done with writing of QMC trial wave function file for CASINO',ch10
515  call wrtout(ab_out,message,'PERS')
516 
517 end subroutine outqmc