Line data Source code
1 : !!****m* ABINIT/m_outqmc
2 : !! NAME
3 : !! m_outqmc
4 : !!
5 : !! FUNCTION
6 : !! Interface with Cambridge quantum Monte Carlo program 'CASINO'
7 : !! See www.tcm.phy.cam.ac.uk/~mdt26/casino.html for more details.
8 : !! M.D.Towler (mdt26 at cam.ac.uk) November 2003
9 : !! N.D.M.Hine (nicholas.hine at imperial.ac.uk) November 2004
10 : !!
11 : !! COPYRIGHT
12 : !! Copyright (C) 1998-2026 ABINIT group (DCA, XG, GMR, JYR, MKV, MT, FJ)
13 : !! This file is distributed under the terms of the
14 : !! GNU General Public License, see ~abinit/COPYING
15 : !! or http://www.gnu.org/copyleft/gpl.txt .
16 : !! For the initials of contributors, see ~abinit/doc/developers/contributors.txt.
17 : !!
18 : !! SOURCE
19 :
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
48 : !!***
49 :
50 : !!****f* m_outqmc/outqmc
51 : !! NAME
52 : !! outqmc
53 : !!
54 : !! FUNCTION
55 : !! Write the wave function to a file in 'pwfn.data' format. This file can be
56 : !! read by the Cambridge quantum Monte Carlo program 'CASINO' and used as
57 : !! trial wave function input for a variational or diffusion Monte Carlo calculation.
58 : !! See www.tcm.phy.cam.ac.uk/~mdt26/casino.html for more details.
59 : !!
60 : !! INPUTS
61 : !! cg(2,mcg)=wavefunction coefficients
62 : !! dtset <type(dataset_type)>=all input variables for this dataset
63 : !! eigen(mband*nkpt*nsppol)=array for holding eigenvalues (hartree)
64 : !! gprimd(3,3)=dimensional primitive translations for reciprocal space
65 : !! hdr <type(hdr_type)>=the header of wf, den and pot files
66 : !! kg(3,mpw*mkmem)=reduced planewave coordinates.
67 : !! mcg=size of wave-functions array (cg) =mpw*nspinor*mband*mkmem*nsppol
68 : !! mpi_enreg=information about MPI parallelization
69 : !! npwarr(nkpt)=number of planewaves in basis and on boundary for each k
70 : !! occ(mband*nkpt*nsppol)=occupation number for each band and k
71 : !! psps <type(pseudopotential_type)>=all the information about psps
72 : !! results_gs <type(results_gs_type)>=results (energy and its components,
73 : !! forces and its components, the stress tensor) of a ground-state computation.
74 : !!
75 : !! OUTPUT
76 : !! Writes the file pwfn.data
77 : !!
78 : !! SOURCE
79 :
80 0 : 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 0 : integer :: kgfull(3,dtset%mpw*dtset%mkmem),kgmap(dtset%mpw*dtset%mkmem)
110 0 : real(dp) :: gcart_qmc(3),kptscart_qmc(3,dtset%nkpt)
111 0 : real(dp),allocatable :: cgfull(:,:),xcart_qmc(:,:)
112 : ! *********************************************************************
113 :
114 : !Go away if I am not the master node.
115 0 : write(message,'(a,a)')ch10,' outqmc: enter '
116 0 : call wrtout(ab_out,message,'PERS')
117 :
118 0 : am_master=.true.
119 : if(xmpi_paral==1 .or. mpi_enreg%paral_kgb==1)then
120 0 : comm=mpi_enreg%comm_cell
121 0 : me=xmpi_comm_rank(comm)
122 0 : if(me/=0) am_master=.false.
123 : end if
124 : if(.not.am_master)return
125 :
126 0 : if(mpi_enreg%paral_spinor==1)then
127 0 : message = ' Parallelization over spinors is not currently supported'
128 0 : 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 0 : if (dtset%ndtset<2)then
134 0 : pwfnfilename='pwfn.data'
135 : else
136 0 : pwfnfilename='pwfn'//trim(i2s(dtset%jdtset))//'.data'
137 : end if
138 0 : call wrtout(ab_out,' outqmc: will open CASINO file: '//trim(pwfnfilename),'PERS')
139 :
140 0 : io = get_unit()
141 0 : open(io,file=pwfnfilename,form='formatted',recl=300,status='unknown',iostat=ierr)
142 :
143 0 : if(ierr/=0)then
144 0 : 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 0 : if (dtset%kptopt==1) then
149 0 : close(io,status='delete')
150 0 : write(message,'(3a)')' outqmc: ERROR - kptopt=1 so k-points have been ',&
151 0 : & 'generated in the irreducible Brillouin Zone only. ',&
152 0 : & 'Set kptopt=2 to obtain full set of k-points.'
153 0 : ABI_ERROR(message)
154 : end if
155 :
156 : !Check if the full set of G vectors has been used in this calculation
157 0 : do ikpt=1,dtset%nkpt
158 0 : if (dtset%istwfk(ikpt)/=1) then
159 0 : close(io,status='delete')
160 : write(message,'(a,i5,a,i2,a,a,a,a,a)')&
161 0 : & ' istwfk(',ikpt,')=',dtset%istwfk(ikpt),' (ie /= 1) so some ',&
162 0 : & 'G-vectors may not be present.',ch10,' Set istwfk=1 for each ',&
163 0 : & 'k-point to obtain full set.'
164 0 : ABI_ERROR(message)
165 : end if
166 : end do !ikpt
167 :
168 0 : write(message,'(a)')' outqmc: QMC trial wave function file for CASINO generated by ABINIT'
169 0 : call wrtout(ab_out,message,'PERS')
170 :
171 : !Write the required quantities to pwfn.data
172 0 : write(io,"('QMC trial wave function file for CASINO generated by ABINIT (www.abinit.org).')")
173 :
174 0 : write(io,fmt="(/'BASIC INFO'/'----------')")
175 0 : write(io,fmt="('Generated by:')")
176 0 : write(io,fmt="(' ABINIT ',a)")trim(hdr%codvsn)
177 0 : write(io,fmt="('Method:'/' DFT')")
178 :
179 0 : write(io,fmt="('DFT Functional:')")
180 0 : select case (dtset%ixc)
181 : case(0)
182 0 : dft_functional='No exchange-correlation.'
183 : case(1)
184 0 : dft_functional='L(S)DA (Teter/Pade parametrization)'
185 : case(2)
186 0 : dft_functional='LDA (Perdew-Zunger-Ceperley-Alder parametrization)'
187 : case(3)
188 0 : dft_functional='LDA (old Teter rational polynomial parametrization)'
189 : case(4)
190 0 : dft_functional='LDA (Wigner)'
191 : case(5)
192 0 : dft_functional='LDA (Hedin-Lundqvist)'
193 : case(6)
194 0 : dft_functional='LDA (X-alpha)'
195 : case(7)
196 0 : dft_functional='L(S)DA (Perdew-Wang 92)'
197 : case(8)
198 0 : dft_functional='L(S)DA (Perdew-Wang 92, exchange-only)'
199 : case(9)
200 0 : dft_functional='L(S)DA (Perdew-Wang 92, exchange- and RPA-correlation)'
201 : case(10)
202 0 : dft_functional='Diff. between ixc=7 and 9; use with accurate RPA corr. energy'
203 : case(11)
204 0 : dft_functional='GGA (Perdew-Burke-Ernzerhof)'
205 : case(12)
206 0 : dft_functional='GGA (Perdew-Burke-Ernzerhof, exchange-only)'
207 : case(13)
208 0 : dft_functional='GGA (potential: van Leeuwen-Baerends ; energy: Perdew-Wang 92)'
209 : case(14)
210 0 : dft_functional='GGA (RPBE of Zhang and Yang)'
211 : case(15)
212 0 : dft_functional='GGA (RPBE of Hammer, Hansen and Norskov)'
213 : case(16)
214 0 : dft_functional='GGA (HTCH)'
215 : case(17)
216 0 : dft_functional='Not defined (as of 11/2003).'
217 : case(18)
218 0 : dft_functional='Not defined (as of 11/2003).'
219 : case(19)
220 0 : dft_functional='Not defined (as of 11/2003).'
221 : case(20)
222 0 : dft_functional='Fermi-Amaldi xc for TDDFT'
223 : case(21)
224 0 : dft_functional='Fermi-Amaldi xc for TDDFT with LDA xc kernel'
225 : case(22)
226 0 : dft_functional='Fermi-Amaldi xc for TDDFT with Burke-Petersilka-Gross hybrid xc kernel'
227 : case default
228 0 : dft_functional='Unknown type.'
229 : end select
230 0 : write(io,"(' ABINIT ixc = ',i2,' : ',a)")dtset%ixc,trim(dft_functional)
231 :
232 0 : write(io,"('Pseudopotential (of first atom type)')")
233 0 : select case(psps%pspcod(1))
234 : case(1)
235 0 : pseudo_type='ABINIT type 1' ; pseudo_name='Troullier-Martins'
236 : case(2)
237 0 : pseudo_type='ABINIT type 2' ; pseudo_name='Goedecker-Teter-Hutter (GTH)'
238 : case(3)
239 0 : pseudo_type='ABINIT type 3' ; pseudo_name='Hartwigsen-Goedecker-Hutter'
240 : case(4)
241 0 : pseudo_type='ABINIT type 4'
242 0 : pseudo_name='Teter pseudo generated using the ATOM code'
243 : case(5)
244 0 : pseudo_type='ABINIT type 5'
245 0 : pseudo_name='"Phoney" pseudo built on a Hamman grid'
246 : case(6)
247 0 : pseudo_type='ABINIT type 6'
248 0 : pseudo_name='Fritz-Haber Institut (Troullier Martins)'
249 : case default
250 0 : pseudo_type='Unknown pseudopotential type (as of 11/2003).' ; pseudo_name=''
251 : end select
252 0 : if(dtset%ixc<10)then
253 0 : write(io,"(1x,a,2x,': ',a)")trim(pseudo_type),trim(pseudo_name)
254 : else
255 0 : write(io,"(1x,a,3x,': ',a)")trim(pseudo_type),trim(pseudo_name)
256 : end if
257 :
258 0 : write(io,"('Plane wave cutoff (au)')")
259 0 : tmpr=r2s(hdr%ecut_eff,'(f12.3)')
260 0 : write(io,'(1x,a)')trim(tmpr)
261 :
262 0 : write(io,"('Spin polarized:')")
263 0 : select case(dtset%nspden)
264 : case(1)
265 0 : write(io,"(' .false.')")
266 : case(2)
267 0 : write(io,"(' .true.')")
268 : case(4)
269 0 : close(io,status='delete')
270 0 : write(message,'(a)')' outqmc: ERROR - nspden=4 but CASINO cannot yet deal with non-collinear spins.'
271 0 : ABI_ERROR(message)
272 : case default
273 0 : close(io,status='delete')
274 0 : ABI_ERROR('Unrecognized value of nspden.')
275 : end select
276 :
277 0 : write(io,"('Total energy (au per primitive cell)')")
278 0 : tmpr=r2s(results_gs%etotal,'(f24.14)')
279 0 : write(io,'(1x,a)')trim(tmpr)
280 0 : write(io,"('Kinetic energy')")
281 0 : tmpr=r2s(results_gs%energies%e_kinetic,'(f24.14)')
282 0 : write(io,'(1x,a)')trim(tmpr)
283 0 : write(io,"('Local potential energy (Local pseudopotential energy eei + pseudopotential core-core energy eii)')")
284 0 : tmpr=r2s((results_gs%energies%e_localpsp+results_gs%energies%e_corepsp),'(f24.14)')
285 0 : write(io,'(1x,a)')trim(tmpr)
286 0 : write(io,"('Non-local potential energy')")
287 0 : tmpr=r2s(results_gs%energies%e_nlpsp_vfock,'(f24.14)')
288 0 : write(io,'(1x,a)')trim(tmpr)
289 0 : write(io,"('Electron-electron energy (Hartree Energy + Exchange-Correlation Energy)')")
290 0 : tmpr=r2s((results_gs%energies%e_hartree+results_gs%energies%e_xc),'(f24.14)')
291 0 : write(io,'(1x,a)')trim(tmpr)
292 0 : write(io,"('Ion-ion energy')")
293 0 : tmpr=r2s(results_gs%energies%e_ewald,'(f24.14)')
294 0 : write(io,'(1x,a)')trim(tmpr)
295 0 : write(io,"('Number of electrons per primitive cell')")
296 :
297 0 : nelecs=0
298 0 : do ii=1,dtset%natom
299 0 : nelecs=nelecs+psps%ziontypat(dtset%typat(ii))
300 : end do
301 :
302 0 : write(io,'(1x,i3)')nelecs
303 0 : write(io,*)
304 0 : write(io,"('GEOMETRY'/'--------')")
305 0 : write(io,"('Number of atoms per primitive cell')")
306 0 : write(io,'(1x,i3)')dtset%natom
307 0 : write(io,"('Atomic numbers and positions of atoms (au)')")
308 :
309 0 : ABI_MALLOC(xcart_qmc,(3,dtset%natom))
310 0 : call xred2xcart(dtset%natom,hdr%rprimd,xcart_qmc,hdr%xred)
311 0 : do ii=1,dtset%natom
312 0 : tmpr=r2s(xcart_qmc(1,ii),'(f24.14)')
313 0 : tmpr2=r2s(xcart_qmc(2,ii),'(f24.14)')
314 0 : tmpr3=r2s(xcart_qmc(3,ii),'(f24.14)')
315 0 : jj=psps%znucltypat(dtset%typat(ii))
316 0 : write(io,'(1x,i2,3(1x,a))')jj,trim(tmpr),trim(tmpr2),trim(tmpr3)
317 : end do
318 0 : ABI_FREE(xcart_qmc)
319 :
320 0 : write(io,"('Primitive lattice vectors (au)')")
321 0 : do ii=1,3
322 0 : tmpr=r2s(hdr%rprimd(1,ii),'(f24.14)')
323 0 : tmpr2=r2s(hdr%rprimd(2,ii),'(f24.14)')
324 0 : tmpr3=r2s(hdr%rprimd(3,ii),'(f24.14)')
325 0 : 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 0 : ikgfull=0
330 0 : do ikg=1,npwarr(1)
331 0 : ikgfull=ikgfull+1
332 0 : kgfull(1:3,ikgfull) = kg(1:3,ikg)
333 0 : kgmap(ikg)=ikgfull
334 : end do
335 0 : 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 0 : do ikpt=2,dtset%nkpt
339 0 : do ikg=ikg_shift,ikg_shift+npwarr(ikpt)
340 0 : foundkg = .false.
341 0 : do ii=1,ikgfull
342 : if(kg(1,ikg)==kgfull(1,ii).and.kg(2,ikg)==kgfull(2,ii) &
343 0 : & .and.kg(3,ikg)==kgfull(3,ii)) then
344 0 : foundkg=.true.
345 0 : kgmap(ikg)=ii
346 : exit
347 : end if
348 : end do
349 0 : if(.not.foundkg)then
350 0 : ikgfull=ikgfull+1
351 0 : kgfull(1:3,ikgfull)=kg(1:3,ikg)
352 0 : kgmap(ikg)=ikgfull
353 : end if
354 : end do
355 0 : ikg_shift=ikg_shift+npwarr(ikpt)
356 : end do
357 0 : nkgfull=ikgfull
358 :
359 0 : write(io,*)
360 0 : write(io,"('G VECTORS'/'---------')")
361 0 : write(io,"('Number of G-vectors')")
362 0 : write(io,'(1x,i7)')nkgfull
363 0 : write(io,"('Gx Gy Gz (au)')")
364 :
365 0 : do ikgfull=1,nkgfull
366 : gcart_qmc=2*pi*(kgfull(1,ikgfull)*gprimd(1:3,1)&
367 0 : & +kgfull(2,ikgfull)*gprimd(1:3,2)+kgfull(3,ikgfull)*gprimd(1:3,3))
368 0 : 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 0 : ABI_MALLOC(cgfull,(2,nkgfull*dtset%nspinor*dtset%nsppol*dtset%mband*dtset%nkpt))
377 0 : cgfull(1:2,1:nkgfull*dtset%nspinor*dtset%nsppol*dtset%mband*dtset%nkpt)=0
378 0 : icg_shift=1
379 0 : do isppol=1,dtset%nsppol
380 0 : ikg_shift=1
381 0 : if(isppol==2)then
382 : ! Go back to the beginning of cgfull but skip the first set of isppol=1 bands
383 0 : icgfull_shift=nkgfull*dtset%nband(1)
384 0 : ikpt_shift=dtset%nkpt
385 : else
386 : icgfull_shift=0 ! Start at the beginning of cgfull
387 : ikpt_shift=0
388 : end if
389 0 : do ikpt=1,dtset%nkpt
390 0 : do iband=1,dtset%nband(ikpt+ikpt_shift)
391 0 : ikg=ikg_shift
392 : ! Find the index in Abinit's coefficient list
393 0 : 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 0 : icgfull = kgmap(ikg)+icgfull_shift
397 0 : cgfull(1:2,icgfull)=cg(1:2,icg)
398 0 : ikg=ikg+1
399 : end do !icg
400 0 : icg_shift=icg_shift+npwarr(ikpt)
401 0 : icgfull_shift=icgfull_shift+nkgfull
402 : end do !iband
403 0 : if(isppol==2)then
404 : ! Skip the isppol==1 bands
405 0 : icgfull_shift=icgfull_shift+nkgfull*dtset%nband(ikpt)
406 : else
407 0 : if(dtset%nsppol==2)then
408 : ! Skip the isppol==2 bands
409 0 : 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 0 : 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 0 : icgfull_shift=1
421 0 : ii=0
422 0 : iocc=1
423 :
424 0 : do ikpt=1,dtset%nkpt
425 0 : do isppol=1,dtset%nsppol
426 0 : do iband=1,dtset%nband(ikpt+(isppol-1)*dtset%nkpt)
427 0 : if(occ(iocc)/=int(occ(iocc)))then
428 : write(message,'(a,i5,a,i1,a,i5,a,f11.8,a)')&
429 0 : & 'Non-integer occupation number for kpt ',ikpt,', sppol ',isppol,', band ',iband,': occ=',occ(iocc),'.'
430 0 : ABI_WARNING(message)
431 : end if
432 0 : iocc=iocc+1
433 0 : norm=0
434 0 : do icgfull=icgfull_shift,icgfull_shift+nkgfull-1
435 0 : norm=norm+cgfull(1,icgfull)**2+cgfull(2,icgfull)**2
436 : end do !icgfull
437 0 : icgfull_shift=icgfull_shift+nkgfull
438 0 : if((norm<0.999).or.(norm>1.001))then
439 : write(message,'(a,i5,a,i1,a,i5,a,f11.8,a)')&
440 0 : & 'Incorrectly normalised orbital for kpt ',ikpt,', sppol ',isppol,', band ',iband,': norm=',norm,'.'
441 0 : ABI_WARNING(message)
442 : end if
443 : end do !iband
444 : end do !isppol
445 : end do !ikpt
446 :
447 0 : write(io,*)
448 0 : write(io,"('WAVE FUNCTION'/'-------------')")
449 0 : write(io,"('Number of k-points')")
450 0 : write(io,'(1x,i5)') dtset%nkpt
451 :
452 0 : do ikpt=1,dtset%nkpt
453 : kptscart_qmc(1:3,ikpt)=2*pi*(dtset%kptns(1,ikpt)*gprimd(1:3,1)&
454 0 : & +dtset%kptns(2,ikpt)*gprimd(1:3,2)+dtset%kptns(3,ikpt)*gprimd(1:3,3))
455 : end do !ikpt
456 :
457 0 : iband_kpt_shift=0
458 0 : icg_shift=0
459 0 : icgfull_shift=0
460 0 : do ikpt=1,dtset%nkpt
461 0 : write(io,"('k-point # ; # of bands (up spin/down spin) ; k-point coords (au)')")
462 0 : if(dtset%nsppol==2)then
463 0 : nband1=dtset%nband(ikpt)
464 0 : nband2=dtset%nband(ikpt+dtset%nkpt)
465 : else
466 0 : nband1=dtset%nband(ikpt)
467 0 : nband2=0
468 : end if
469 0 : tmpr=r2s(kptscart_qmc(1,ikpt),'(f24.14)')
470 0 : tmpr2=r2s(kptscart_qmc(2,ikpt),'(f24.14)')
471 0 : tmpr3=r2s(kptscart_qmc(3,ikpt),'(f24.14)')
472 0 : write(io,'(3(1x,i5),3(1x,a))')ikpt,nband1,nband2,trim(tmpr),trim(tmpr2), &
473 0 : & trim(tmpr3)
474 0 : do isppol=1,dtset%nsppol
475 0 : if (isppol==2) then
476 0 : iband_sppol_shift=sum(dtset%nband(1:dtset%nkpt))
477 0 : iband_kpt_shift=sum(dtset%nband((dtset%nkpt+1):(dtset%nkpt+ikpt-1)))
478 : else
479 0 : iband_sppol_shift=0
480 0 : iband_kpt_shift=sum(dtset%nband(1:(ikpt-1)))
481 : end if
482 0 : do iband=1,dtset%nband(ikpt)
483 0 : write(io,"('Band, spin, eigenvalue (au)')")
484 0 : tmpr=r2s(eigen(iband_kpt_shift+iband_sppol_shift+iband),'(f24.14)')
485 0 : write(io,'(2(1x,i5),1x,a)')iband,isppol,tmpr
486 0 : write(io,"('Eigenvector coefficients')")
487 0 : do icgfull=1,nkgfull
488 0 : write(io,"(1x,'(',e23.16,',',e23.16,')')")cgfull(1:2,icgfull+icgfull_shift)
489 : end do !icgfull
490 0 : icgfull_shift=icgfull_shift+nkgfull
491 : end do !iband
492 : end do !isppol
493 : end do !ikpt
494 :
495 0 : close(io)
496 :
497 0 : write(message,'(a,a)')' outqmc: done with writing of QMC trial wave function file for CASINO',ch10
498 0 : call wrtout(ab_out,message,'PERS')
499 :
500 0 : end subroutine outqmc
501 : !!***
502 :
503 : !!****f* m_outqmc/i2s
504 : !! NAME
505 : !! i2s
506 : !!
507 : !! FUNCTION
508 : !! Convert integers to left justified strings that can be printed in the
509 : !! middle of a sentence without introducing large amounts of white space.
510 : !! Calling routine is intended to include something like:
511 : !! integer i
512 : !! i=12
513 : !! write(std_out,*)'Integer number ',trim(i2s(i)),' with words at the end.'
514 : !!
515 : !! INPUTS
516 : !!
517 : !! OUTPUT
518 : !!
519 : !! SOURCE
520 :
521 0 : 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 0 : if(n==0)then
534 0 : i2s='0' ; return
535 : end if
536 0 : sign=' ' ; if(n<0)sign='-'
537 :
538 0 : do i=1,len(i2s)
539 0 : i2s(i:i)=' '
540 : end do
541 :
542 0 : i=abs(n)
543 0 : do j=1,len(i2s)
544 0 : if(i==0)exit
545 0 : i2s(j:j)=achar(ichar('0')+mod(i,10))
546 0 : i=i/10
547 : end do
548 :
549 0 : i=1 ; j=len_trim(i2s)
550 0 : do
551 0 : if(i>=j)exit
552 0 : tmp=i2s(j:j)
553 0 : i2s(j:j)=i2s(i:i)
554 0 : i2s(i:i)=tmp
555 0 : i=i+1
556 0 : j=j-1
557 : end do
558 :
559 0 : i2s=trim(sign)//i2s
560 :
561 : end function i2s
562 : !!***
563 :
564 : !!****f* ABINIT/r2s
565 : !! NAME
566 : !! r2s
567 : !!
568 : !! FUNCTION
569 : !! Converts real variable with arbitrary format to string that can be
570 : !! trimmed and printed in the middle of a sentence without introducing
571 : !! large amounts of white space, as you would if you did
572 : !! write(std_out,'(f12.6)')12.0 or similar. Note you need to pass through the
573 : !! format string e.g. f12.6.
574 : !!
575 : !! Calling routine is intended to include something like:
576 : !! USE utilities
577 : !! REAL(dp) r
578 : !! r=12._dp
579 : !! tmpr=r2s(r,'(f12.6)')
580 : !! write(std_out,*)'Real number ',trim(tmpr),' with words at the end.'
581 : !!
582 : !! Note : DON'T USE R2S IN A WRITE STATEMENT SINCE THIS IS ILLEGAL
583 : !! IN FORTRAN90 (ALTHOUGH NOT IN FORTRAN200X). IF ANYONE HAS TIME, FEEL
584 : !! FREE TO WRITE A VERSION OF THIS WHICH ISN'T ILLEGAL - SIMILAR TO
585 : !! I2S ABOVE - SO THAT PEOPLE WHO HAVEN'T READ THIS NOTE DON'T FEEL
586 : !! TEMPTED TO CALL R2S IN A WRITE STATEMENT.
587 : !!
588 : !! INPUTS
589 : !!
590 : !! OUTPUT
591 : !!
592 : !! SIDE EFFECTS
593 : !!
594 : !! SOURCE
595 :
596 0 : 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 0 : write(r2s,real_format)r
607 0 : r2s=adjustl(r2s)
608 : end if
609 :
610 0 : end function r2s
611 : !!***
612 :
613 : end module m_outqmc
614 : !!***
|
