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