TABLE OF CONTENTS
ABINIT/pred_steepdesc [ Functions ]
NAME
pred_steepdesc
FUNCTION
Ionmov predictor (21) Steepest Descent Algorithm
The update of positions is given by the following
equation:
$$\Delta r_{n,i}=\lambda F_{n,i}$$
r is the position of the 'n' ion along the 'i' direction
F is the force of the 'n' ion along the 'i' direction.
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 .
INPUTS
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 icycle: Index of the internal cycle inside a time step (itime) itime: Index of time iteration zDEBUG : if true print some debugging information
OUTPUT
SIDE EFFECTS
hist<type abihist>=Historical record of positions, forces, stresses, cell and energies. ncycle: Number of cycles of a particular time step
NOTES
* This routine is a predictor, it only produces new positions to be computed in the next iteration, this routine should produce not output at all
PARENTS
mover
CHILDREN
hist2var,mkradim,mkrdim,var2hist,xcart2xred,xred2xcart
SOURCE
56 #if defined HAVE_CONFIG_H 57 #include "config.h" 58 #endif 59 60 #include "abi_common.h" 61 62 63 subroutine pred_steepdesc(ab_mover,forstr,hist,itime,zDEBUG,iexit) 64 65 use defs_basis 66 use m_profiling_abi 67 use m_abimover 68 use m_abihist 69 70 use m_geometry, only : mkradim, mkrdim, xcart2xred, xred2xcart 71 72 !This section has been created automatically by the script Abilint (TD). 73 !Do not modify the following lines by hand. 74 #undef ABI_FUNC 75 #define ABI_FUNC 'pred_steepdesc' 76 use interfaces_45_geomoptim, except_this_one => pred_steepdesc 77 !End of the abilint section 78 79 implicit none 80 81 !Arguments ------------------------------------ 82 !scalars 83 type(abimover),intent(in) :: ab_mover 84 type(abihist),intent(inout),target :: hist 85 type(abiforstr),intent(in) :: forstr 86 integer,intent(in) :: itime,iexit 87 logical,intent(in) :: zDEBUG 88 89 !Local variables------------------------------- 90 !scalars 91 integer :: kk,jj,ihist_prev 92 real(dp) :: em 93 real(dp) :: f_cart 94 real(dp) :: xc,str 95 real(dp) :: xnow,lambda 96 real(dp),save :: hh 97 !arrays 98 real(dp) :: acell(3),strten(6) 99 real(dp) :: rprim(3,3),rprimd(3,3) 100 real(dp) :: xred(3,ab_mover%natom),xcart(3,ab_mover%natom) 101 real(dp) :: residual(3,ab_mover%natom) 102 real(dp), ABI_CONTIGUOUS pointer :: fcart(:,:),vel(:,:) 103 104 !*************************************************************************** 105 !Beginning of executable session 106 !*************************************************************************** 107 108 if(iexit/=0)then 109 return 110 end if 111 112 !write(std_out,*) '01' 113 !########################################################## 114 !### 01. Copy from the history to the variables 115 call hist2var(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG) 116 117 do jj=1,3 118 rprim(jj,1:3)=rprimd(jj,1:3)/acell(1:3) 119 end do 120 121 call xred2xcart(ab_mover%natom,rprimd,xcart,xred) 122 strten(:)=hist%strten(:,hist%ihist) 123 fcart => hist%fcart(:,:,hist%ihist) 124 vel => hist%vel(:,:,hist%ihist) 125 126 !Fill the residual with forces (No preconditioning) 127 !Or the preconditioned forces 128 if (ab_mover%goprecon==0)then 129 residual(:,:)=fcart(:,:) 130 else 131 residual(:,:)= forstr%fcart(:,:) 132 end if 133 134 !write(std_out,*) '01' 135 !########################################################## 136 !### 02. Get or compute de time step dtion 137 138 if (ab_mover%dtion>0)then 139 hh = ab_mover%dtion 140 else 141 hh=fdtion(ab_mover,itime,xcart,fcart,vel) 142 end if 143 144 lambda=hh 145 write(std_out,*) 'Lambda',lambda 146 147 !write(std_out,*) '02' 148 !########################################################## 149 !### 02. For all atoms and directions 150 do kk=1,ab_mover%natom 151 ! Normally this is the mass of the atom 152 ! em=ab_mover%amass(kk) 153 ! But the steepest algorithm is using unitary mass 154 em=1 155 do jj=1,3 156 157 ! write(std_out,*) '03' 158 ! ########################################################## 159 ! ### 03. Filling other values from history (forces and vel) 160 f_cart=residual(jj,kk) 161 xc=xcart(jj,kk) 162 ! This lambda is taken from the kinematical equation 163 ! lambda=hh*hh/(2*em) 164 165 ! write(std_out,*) '04' 166 ! ########################################################## 167 ! ### 04. Take first the atoms that are not allowed to move along 168 ! ### this direction 169 ! ### Warning : implemented in cartesian coordinates 170 if (ab_mover%iatfix(jj,kk)==1) then 171 ! Their positions will be the same as xcart 172 xnow=xc 173 else 174 175 ! This is the main expresion (1) 176 xnow=xc+lambda*f_cart 177 178 end if !if(ab_mover%iatfix(jj,kk)==1) 179 180 ! write(std_out,*) '05' 181 ! ########################################################## 182 ! ### 08. Update history 183 184 xcart(jj,kk)=xnow 185 186 ! write(std_out,*) '06' 187 ! ########################################################## 188 ! ### 09. End loops of atoms and directions 189 end do ! jj=1,3 190 end do ! kk=1,ab_mover%natom 191 192 if (ab_mover%optcell/=0)then 193 194 if (ab_mover%optcell==1)then 195 do jj=1,3 196 acell(jj)=acell(jj)+lambda*strten(jj) 197 end do ! jj=1,3 198 call mkrdim(acell,rprim,rprimd) 199 elseif (ab_mover%optcell==2)then 200 do kk=1,3 201 do jj=1,3 202 if (jj==1 .and. kk==1) str=strten(1) 203 if (jj==2 .and. kk==2) str=strten(2) 204 if (jj==3 .and. kk==3) str=strten(3) 205 if (jj==1 .and. kk==2) str=strten(6) 206 if (jj==1 .and. kk==3) str=strten(5) 207 if (jj==2 .and. kk==1) str=strten(6) 208 if (jj==3 .and. kk==1) str=strten(5) 209 if (jj==2 .and. kk==3) str=strten(4) 210 if (jj==3 .and. kk==2) str=strten(4) 211 rprimd(jj,kk)=rprimd(jj,kk)+lambda*str 212 end do ! jj=1,3 213 end do ! kk=1,3 214 call mkradim(acell,rprim,rprimd) 215 end if 216 217 end if 218 219 220 !write(std_out,*) '08' 221 !########################################################## 222 !### 10. Filling history with the new values 223 224 !Increase indices 225 hist%ihist = abihist_findIndex(hist,+1) 226 227 !Compute xred from xcart, and rprimd 228 call xcart2xred(ab_mover%natom,rprimd,xcart,xred) 229 230 call var2hist(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG) 231 ihist_prev = abihist_findIndex(hist,-1) 232 hist%vel(:,:,hist%ihist)=hist%vel(:,:,ihist_prev) 233 234 end subroutine pred_steepdesc