TABLE OF CONTENTS


ABINIT/fdtion [ Functions ]

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

PARENTS

      pred_moldyn

CHILDREN

SOURCE

 78 function fdtion(ab_mover,itime,xcart,fcart,vel)
 79 
 80 
 81 !This section has been created automatically by the script Abilint (TD).
 82 !Do not modify the following lines by hand.
 83 #undef ABI_FUNC
 84 #define ABI_FUNC 'fdtion'
 85 !End of the abilint section
 86 
 87   implicit none
 88 
 89 !Arguments ---------------------------------------------
 90 !scalars
 91   type(abimover),intent(in) :: ab_mover
 92   integer,intent(in) :: itime
 93   real(dp) :: fdtion
 94 !arrays
 95   real(dp) :: xcart(:,:),fcart(:,:),vel(:,:)
 96 
 97 !Local variables ------------------------------
 98 !scalars
 99   integer  :: jj,kk
100   real(dp) :: max,min,val
101   real(dp) :: ff,xc,vv,em
102 
103 !************************************************************************
104 
105  max=0
106  min=1e6
107 
108  do kk=1,ab_mover%natom
109    em=ab_mover%amass(kk)
110    do jj=1,3
111      ff =fcart(jj,kk)
112      xc =xcart(jj,kk)
113      vv=vel(jj,kk)
114 
115      if (vv>1e-8) then
116        val=abs(1.0_dp/vv)
117        write(std_out,*) 'vel',kk,jj,val
118        if (val>max) max=val
119        if (val<min) min=val
120      end if
121 
122      if (ff>1e-8) then
123        val=sqrt(abs(2*em/ff))
124        write(std_out,*) 'forces',kk,jj,val,em,ff
125        if (val>max) max=val
126        if (val<min) min=val
127      end if
128 
129    end do
130 
131  end do
132 
133  write(std_out,*) "DTION max=",max
134  write(std_out,*) "DTION min=",min
135 
136  if (itime==1)then
137    fdtion=min/10
138  else
139    fdtion=min/10
140  end if
141 
142  end function fdtion

ABINIT/m_predtk [ Modules ]

[ Top ] [ Modules ]

NAME

  m_predtk

FUNCTION

  Low-level procedures used by 45_geomoptim routines

COPYRIGHT

  Copyright (C) 1998-2018 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 .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_predtk
28 
29  use defs_basis
30  use m_abicore
31  use m_abimover
32 
33  implicit none
34 
35  private

ABINIT/prtxvf [ Functions ]

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

PARENTS

      constrf,prtimg

CHILDREN

      wrtout

SOURCE

174 subroutine prtxvf(fcart,fred,iatfix,iout,natom,prtvel,vel,xcart,xred)
175 
176 
177 !This section has been created automatically by the script Abilint (TD).
178 !Do not modify the following lines by hand.
179 #undef ABI_FUNC
180 #define ABI_FUNC 'prtxvf'
181 !End of the abilint section
182 
183  implicit none
184 
185 !Arguments ------------------------------------
186 !scalars
187  integer,intent(in) :: iout,natom,prtvel
188 !arrays
189  integer,intent(in) :: iatfix(3,natom)
190  real(dp),intent(in) :: fcart(3,natom),fred(3,natom)
191  real(dp),intent(in) :: xcart(3,natom),xred(3,natom)
192  real(dp),intent(in) :: vel(3,natom)
193 !Local variables-------------------------------
194 !scalars
195  integer :: iatom,mu,unfixd
196  real(dp) :: fmax,frms,val_max,val_rms
197  character(len=500) :: msg
198 
199 ! *****************************************************************
200 
201  write(msg, '(a)' ) ' Cartesian coordinates (xcart) [bohr]'
202  call wrtout(iout,msg,'COLL')
203  do iatom=1,natom
204    write(msg, '(1p,3e22.14)' )xcart(:,iatom)
205    call wrtout(iout,msg,'COLL')
206  end do
207 
208  write(msg, '(a)' ) ' Reduced coordinates (xred)'
209  call wrtout(iout,msg,'COLL')
210  do iatom=1,natom
211    write(msg, '(1p,3e22.14)' )xred(:,iatom)
212    call wrtout(iout,msg,'COLL')
213  end do
214 
215 !Compute max |f| and rms f, EXCLUDING the components determined by iatfix
216 
217  fmax=0.0_dp
218  frms=0.0_dp
219  unfixd=0
220  do iatom=1,natom
221    do mu=1,3
222      if (iatfix(mu,iatom) /= 1) then
223        unfixd=unfixd+1
224        frms=frms+fcart(mu,iatom)**2
225        fmax=max(fmax,abs(fcart(mu,iatom)))
226      end if
227    end do
228  end do
229  if ( unfixd /= 0 ) frms=sqrt(frms/dble(unfixd))
230 
231  write(msg, '(a,1p,2e12.5,a)' ) &
232 & ' Cartesian forces (fcart) [Ha/bohr]; max,rms=',fmax,frms,' (free atoms)'
233  call wrtout(iout,msg,'COLL')
234  do iatom=1,natom
235    write(msg, '(1p,3e22.14)' )fcart(:,iatom)
236    call wrtout(iout,msg,'COLL')
237  end do
238 
239  write(msg, '(a)' ) ' Reduced forces (fred)'
240  call wrtout(iout,msg,'COLL')
241  do iatom=1,natom
242    write(msg, '(1p,3e22.14)' )fred(:,iatom)
243    call wrtout(iout,msg,'COLL')
244  end do
245 
246  if (prtvel == 1) then
247 
248 !  Compute max |v| and rms v,
249 !  EXCLUDING the components determined by iatfix
250    val_max=0.0_dp
251    val_rms=0.0_dp
252    unfixd=0
253    do iatom=1,natom
254      do mu=1,3
255        if (iatfix(mu,iatom) /= 1) then
256          unfixd=unfixd+1
257          val_rms=val_rms+vel(mu,iatom)**2
258          val_max=max(val_max,abs(vel(mu,iatom)**2))
259        end if
260      end do
261    end do
262    if ( unfixd /= 0 ) val_rms=sqrt(val_rms/dble(unfixd))
263 
264 
265    write(msg, '(a,1p,2e12.5,a)' ) ' Cartesian velocities (vel) [bohr*Ha/hbar]; max,rms=',&
266 &   sqrt(val_max),val_rms,' (free atoms)'
267    call wrtout(iout,msg,'COLL')
268    do iatom=1,natom
269      write(msg, '(1p,3e22.14)' ) vel(:,iatom)
270      call wrtout(iout,msg,'COLL')
271    end do
272  end if
273 
274 end subroutine prtxvf