TABLE OF CONTENTS
ABINIT/prtxvf [ Functions ]
NAME
prtxvf
FUNCTION
Print the values of xcart, vel, and fcart to unit iout. Also compute and print max and rms forces.
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR) 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
fcart(3,natom)=forces (hartree/bohr) iatfix(3,natom)=1 for frozen or fixed atom along specified direction, else 0 iout=unit number for printing natom=number of atoms in unit cell. prtvel=1 to print velocities, else do not print them vel(3,natom)=velocities xcart(3,natom)=cartesian coordinates (bohr)
OUTPUT
(only writing)
PARENTS
constrf,prtimg
CHILDREN
wrtout
SOURCE
40 #if defined HAVE_CONFIG_H 41 #include "config.h" 42 #endif 43 44 #include "abi_common.h" 45 46 47 subroutine prtxvf(fcart,fred,iatfix,iout,natom,prtvel,vel,xcart,xred) 48 49 use defs_basis 50 use m_profiling_abi 51 use m_errors 52 53 !This section has been created automatically by the script Abilint (TD). 54 !Do not modify the following lines by hand. 55 #undef ABI_FUNC 56 #define ABI_FUNC 'prtxvf' 57 use interfaces_14_hidewrite 58 !End of the abilint section 59 60 implicit none 61 62 !Arguments ------------------------------------ 63 !scalars 64 integer,intent(in) :: iout,natom,prtvel 65 !arrays 66 integer,intent(in) :: iatfix(3,natom) 67 real(dp),intent(in) :: fcart(3,natom),fred(3,natom) 68 real(dp),intent(in) :: xcart(3,natom),xred(3,natom) 69 real(dp),intent(in) :: vel(3,natom) 70 !Local variables------------------------------- 71 !scalars 72 integer :: iatom,mu,unfixd 73 real(dp) :: fmax,frms,val_max,val_rms 74 character(len=500) :: message 75 76 ! ***************************************************************** 77 78 write(message, '(a)' ) ' Cartesian coordinates (xcart) [bohr]' 79 call wrtout(iout,message,'COLL') 80 do iatom=1,natom 81 write(message, '(1p,3e22.14)' )xcart(:,iatom) 82 call wrtout(iout,message,'COLL') 83 end do 84 85 write(message, '(a)' ) ' Reduced coordinates (xred)' 86 call wrtout(iout,message,'COLL') 87 do iatom=1,natom 88 write(message, '(1p,3e22.14)' )xred(:,iatom) 89 call wrtout(iout,message,'COLL') 90 end do 91 92 !Compute max |f| and rms f, EXCLUDING the components determined by iatfix 93 94 fmax=0.0_dp 95 frms=0.0_dp 96 unfixd=0 97 do iatom=1,natom 98 do mu=1,3 99 if (iatfix(mu,iatom) /= 1) then 100 unfixd=unfixd+1 101 frms=frms+fcart(mu,iatom)**2 102 fmax=max(fmax,abs(fcart(mu,iatom))) 103 end if 104 end do 105 end do 106 if ( unfixd /= 0 ) frms=sqrt(frms/dble(unfixd)) 107 108 write(message, '(a,1p,2e12.5,a)' ) & 109 & ' Cartesian forces (fcart) [Ha/bohr]; max,rms=',fmax,frms,' (free atoms)' 110 call wrtout(iout,message,'COLL') 111 do iatom=1,natom 112 write(message, '(1p,3e22.14)' )fcart(:,iatom) 113 call wrtout(iout,message,'COLL') 114 end do 115 116 write(message, '(a)' ) ' Reduced forces (fred)' 117 call wrtout(iout,message,'COLL') 118 do iatom=1,natom 119 write(message, '(1p,3e22.14)' )fred(:,iatom) 120 call wrtout(iout,message,'COLL') 121 end do 122 123 if (prtvel == 1) then 124 125 ! Compute max |v| and rms v, 126 ! EXCLUDING the components determined by iatfix 127 val_max=0.0_dp 128 val_rms=0.0_dp 129 unfixd=0 130 do iatom=1,natom 131 do mu=1,3 132 if (iatfix(mu,iatom) /= 1) then 133 unfixd=unfixd+1 134 val_rms=val_rms+vel(mu,iatom)**2 135 val_max=max(val_max,abs(vel(mu,iatom)**2)) 136 end if 137 end do 138 end do 139 if ( unfixd /= 0 ) val_rms=sqrt(val_rms/dble(unfixd)) 140 141 142 write(message, '(a,1p,2e12.5,a)' ) ' Cartesian velocities (vel) [bohr*Ha/hbar]; max,rms=',& 143 & sqrt(val_max),val_rms,' (free atoms)' 144 call wrtout(iout,message,'COLL') 145 do iatom=1,natom 146 write(message, '(1p,3e22.14)' ) vel(:,iatom) 147 call wrtout(iout,message,'COLL') 148 end do 149 end if 150 151 end subroutine prtxvf