TABLE OF CONTENTS


ABINIT/onestep [ Functions ]

[ Top ] [ Functions ]

NAME

 onestep

FUNCTION

 Advance one step following the gradient from vv(3).
 It returns a new point in vv(3) and the value and gradient of the
 electron density at this point in chg and grho(3)

COPYRIGHT

 Copyright (C) 2002-2017 ABINIT group (PCasek,FF,XG)
 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

  npmax= maximum number of divisions
  hh= determines the initial value of the step (to be multiplied by grho)

OUTPUT

  chg= value of electron density
  deltar= the length of the step thaty was needed
  grho(3)= gradient of electron density
  np= returns the number of divisions that were needed

SIDE EFFECTS

  vv(3)=starting and updated point

 WARNING
 This file does not follow the ABINIT coding rules (yet)

PARENTS

      aim_follow

CHILDREN

      vgh_rho

SOURCE

 42 #if defined HAVE_CONFIG_H
 43 #include "config.h"
 44 #endif
 45 
 46 #include "abi_common.h"
 47 
 48 
 49 subroutine onestep(vv,chg,grho,hh,np,npmax,deltar)
 50 
 51  use m_profiling_abi
 52 
 53  use defs_basis
 54  use defs_parameters
 55  use defs_aimprom
 56 
 57 !This section has been created automatically by the script Abilint (TD).
 58 !Do not modify the following lines by hand.
 59 #undef ABI_FUNC
 60 #define ABI_FUNC 'onestep'
 61  use interfaces_63_bader, except_this_one => onestep
 62 !End of the abilint section
 63 
 64  implicit none
 65 
 66 !Arguments ------------------------------------
 67 !scalars
 68  integer,intent(in) :: npmax
 69  integer,intent(out) :: np
 70  real(dp),intent(in) :: hh
 71  real(dp),intent(out) :: chg,deltar
 72 !arrays
 73  real(dp),intent(inout) :: vv(3)
 74  real(dp),intent(out) :: grho(3)
 75 
 76 !Local variables ------------------------------
 77 !scalars
 78  integer :: iat,ii,ipos,jj
 79  real(dp) :: dt,rr
 80 !arrays
 81  real(dp) :: hrho(3,3),pom(3),vinter(3,200),vk(3),vkold(3)
 82 
 83 !************************************************************************
 84  dt=hh
 85  np=1
 86  deltar=1._dp
 87  vk(1:3)=vv(1:3)
 88  
 89 
 90  do while((np<3).or.((np<=npmax).and.(deltar>aim_deltarmin)))
 91    np=np*2
 92    dt=dt*0.5_dp
 93    vkold(1:3)=vk(1:3)
 94    call vgh_rho(vk,chg,grho,hrho,rr,iat,ipos,0)
 95    vinter(1:3,1)=vv(1:3)+dt*grho(1:3)
 96    do jj=2,np
 97      call vgh_rho(vinter(1,jj-1),chg,grho,hrho,rr,iat,ipos,0)
 98      if(jj.eq.2) then
 99        vinter(1:3,2)=vv(1:3)+2.0*dt*grho(1:3)
100      else
101        vinter(1:3,jj)=vinter(1:3,jj-2)+2.0*dt*grho(1:3)
102      end if
103    end do
104 
105    call vgh_rho(vinter(1,np),chg,grho,hrho,rr,iat,ipos,0)
106    vinter(1:3,np+1)=vinter(1:3,np-1)+dt*grho(1:3)
107    
108    deltar=0._dp
109    do ii=1,3
110      vk(ii)=(vinter(ii,np)+vinter(ii,np+1))*0.5_dp
111      deltar=deltar+(vkold(ii)-vk(ii))*(vkold(ii)-vk(ii))
112    end do
113  end do
114 
115  pom(:)=vk(:)-vv(:)
116  deltar=vnorm(pom,0)
117  vv(1:3)=vk(1:3)
118 
119  call vgh_rho(vv,chg,grho,hrho,rr,iat,ipos,0)
120  if(deb) write(std_out,*) ':VKf ',np,vk
121 
122 end subroutine onestep