TABLE OF CONTENTS
ABINIT/fdtion [ Functions ]
NAME
fdtion
FUNCTION
Compute the apropiated "dtion" from the present values of forces, velocity and viscosity INPUTS (in) hist<type abihist>=Historical record of positions, forces | acell, stresses, and energies, ab_mover<type abimover>=Subset of dtset only related with | movement of ions and acell, contains: | dtion: Time step ! natom: Number of atoms | vis: viscosity | iatfix: Index of atoms and directions fixed | amass: Mass of ions itime: Index of time iteration xcart(3,natom)= cartesian coordinates of atoms fcart(3,natom)= forces in cartesian coordinates vel(3,natom)= velocities OUTPUT (out) fdtion = time step computed
SOURCE
68 function fdtion(ab_mover,itime,xcart,fcart,vel) 69 70 implicit none 71 72 !Arguments --------------------------------------------- 73 !scalars 74 type(abimover),intent(in) :: ab_mover 75 integer,intent(in) :: itime 76 real(dp) :: fdtion 77 !arrays 78 real(dp) :: xcart(:,:),fcart(:,:),vel(:,:) 79 80 !Local variables ------------------------------ 81 !scalars 82 integer :: jj,kk 83 real(dp) :: max,min,val 84 real(dp) :: ff,xc,vv,em 85 86 !************************************************************************ 87 88 max=0 89 min=1e6 90 91 do kk=1,ab_mover%natom 92 em=ab_mover%amass(kk) 93 do jj=1,3 94 ff =fcart(jj,kk) 95 xc =xcart(jj,kk) 96 vv=vel(jj,kk) 97 98 if (vv>1e-8) then 99 val=abs(1.0_dp/vv) 100 write(std_out,*) 'vel',kk,jj,val 101 if (val>max) max=val 102 if (val<min) min=val 103 end if 104 105 if (ff>1e-8) then 106 val=sqrt(abs(2*em/ff)) 107 write(std_out,*) 'forces',kk,jj,val,em,ff 108 if (val>max) max=val 109 if (val<min) min=val 110 end if 111 112 end do 113 114 end do 115 116 write(std_out,*) "DTION max=",max 117 write(std_out,*) "DTION min=",min 118 119 if (itime==1)then 120 fdtion=min/10 121 else 122 fdtion=min/10 123 end if 124 125 end function fdtion
ABINIT/m_predtk [ Modules ]
NAME
m_predtk
FUNCTION
Low-level procedures used by 45_geomoptim routines
COPYRIGHT
Copyright (C) 1998-2024 ABINIT group (DCA, XG, GMR, SE) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_predtk 23 24 use defs_basis 25 use m_abicore 26 use m_abimover 27 28 implicit none 29 30 private
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.
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)
SOURCE
151 subroutine prtxvf(fcart,gred,iatfix,iout,natom,prtvel,vel,xcart,xred) 152 153 implicit none 154 155 !Arguments ------------------------------------ 156 !scalars 157 integer,intent(in) :: iout,natom,prtvel 158 !arrays 159 integer,intent(in) :: iatfix(3,natom) 160 real(dp),intent(in) :: fcart(3,natom),gred(3,natom) 161 real(dp),intent(in) :: xcart(3,natom),xred(3,natom) 162 real(dp),intent(in) :: vel(3,natom) 163 !Local variables------------------------------- 164 !scalars 165 integer :: iatom,mu,unfixd 166 real(dp) :: fmax,frms,val_max,val_rms 167 character(len=500) :: msg 168 169 ! ***************************************************************** 170 171 write(msg, '(a)' ) ' Cartesian coordinates (xcart) [bohr]' 172 call wrtout(iout,msg,'COLL') 173 do iatom=1,natom 174 write(msg, '(1p,3e22.14)' )xcart(:,iatom) 175 call wrtout(iout,msg,'COLL') 176 end do 177 178 write(msg, '(a)' ) ' Reduced coordinates (xred)' 179 call wrtout(iout,msg,'COLL') 180 do iatom=1,natom 181 write(msg, '(1p,3e22.14)' )xred(:,iatom) 182 call wrtout(iout,msg,'COLL') 183 end do 184 185 !Compute max |f| and rms f, EXCLUDING the components determined by iatfix 186 187 fmax=0.0_dp 188 frms=0.0_dp 189 unfixd=0 190 do iatom=1,natom 191 do mu=1,3 192 if (iatfix(mu,iatom) /= 1) then 193 unfixd=unfixd+1 194 frms=frms+fcart(mu,iatom)**2 195 fmax=max(fmax,abs(fcart(mu,iatom))) 196 end if 197 end do 198 end do 199 if ( unfixd /= 0 ) frms=sqrt(frms/dble(unfixd)) 200 201 write(msg, '(a,1p,2e12.5,a)' ) & 202 & ' Cartesian forces (fcart) [Ha/bohr]; max,rms=',fmax,frms,' (free atoms)' 203 call wrtout(iout,msg,'COLL') 204 do iatom=1,natom 205 write(msg, '(1p,3e22.14)' )fcart(:,iatom) 206 call wrtout(iout,msg,'COLL') 207 end do 208 209 write(msg, '(a)' ) ' Gradient of E wrt nuclear positions in reduced coordinates (gred)' 210 call wrtout(iout,msg,'COLL') 211 do iatom=1,natom 212 write(msg, '(1p,3e22.14)' )gred(:,iatom) 213 call wrtout(iout,msg,'COLL') 214 end do 215 216 if (prtvel == 1) then 217 218 ! Compute max |v| and rms v, 219 ! EXCLUDING the components determined by iatfix 220 val_max=0.0_dp 221 val_rms=0.0_dp 222 unfixd=0 223 do iatom=1,natom 224 do mu=1,3 225 if (iatfix(mu,iatom) /= 1) then 226 unfixd=unfixd+1 227 val_rms=val_rms+vel(mu,iatom)**2 228 val_max=max(val_max,abs(vel(mu,iatom)**2)) 229 end if 230 end do 231 end do 232 if ( unfixd /= 0 ) val_rms=sqrt(val_rms/dble(unfixd)) 233 234 235 write(msg, '(a,1p,2e12.5,a)' ) ' Cartesian velocities (vel) [bohr*Ha/hbar]; max,rms=',& 236 & sqrt(val_max),val_rms,' (free atoms)' 237 call wrtout(iout,msg,'COLL') 238 do iatom=1,natom 239 write(msg, '(1p,3e22.14)' ) vel(:,iatom) 240 call wrtout(iout,msg,'COLL') 241 end do 242 end if 243 244 end subroutine prtxvf