TABLE OF CONTENTS
ABINIT/m_pred_steepdesc [ Modules ]
NAME
m_pred_steepdesc
FUNCTION
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
15 #if defined HAVE_CONFIG_H 16 #include "config.h" 17 #endif 18 19 #include "abi_common.h" 20 21 module m_pred_steepdesc 22 23 use defs_basis 24 use m_abicore 25 use m_abimover 26 use m_abihist 27 28 use m_geometry, only : mkradim, mkrdim, xcart2xred, xred2xcart 29 use m_predtk, only : fdtion 30 31 implicit none 32 33 private
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.
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
SOURCE
81 subroutine pred_steepdesc(ab_mover,forstr,hist,itime,zDEBUG,iexit) 82 83 implicit none 84 85 !Arguments ------------------------------------ 86 !scalars 87 type(abimover),intent(in) :: ab_mover 88 type(abihist),intent(inout),target :: hist 89 type(abiforstr),intent(in) :: forstr 90 integer,intent(in) :: itime,iexit 91 logical,intent(in) :: zDEBUG 92 93 !Local variables------------------------------- 94 !scalars 95 integer :: kk,jj,ihist_prev 96 real(dp) :: em 97 real(dp) :: f_cart 98 real(dp) :: xc,str 99 real(dp) :: xnow,lambda 100 real(dp),save :: hh 101 !arrays 102 real(dp) :: acell(3),strten(6) 103 real(dp) :: rprim(3,3),rprimd(3,3) 104 real(dp) :: xred(3,ab_mover%natom),xcart(3,ab_mover%natom) 105 real(dp) :: residual(3,ab_mover%natom) 106 real(dp), ABI_CONTIGUOUS pointer :: fcart(:,:),vel(:,:) 107 108 !*************************************************************************** 109 !Beginning of executable session 110 !*************************************************************************** 111 112 if(iexit/=0)then 113 return 114 end if 115 116 !write(std_out,*) '01' 117 !########################################################## 118 !### 01. Copy from the history to the variables 119 call hist2var(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG) 120 121 do jj=1,3 122 rprim(jj,1:3)=rprimd(jj,1:3)/acell(1:3) 123 end do 124 125 call xred2xcart(ab_mover%natom,rprimd,xcart,xred) 126 strten(:)=hist%strten(:,hist%ihist) 127 fcart => hist%fcart(:,:,hist%ihist) 128 vel => hist%vel(:,:,hist%ihist) 129 130 !Fill the residual with forces (No preconditioning) 131 !Or the preconditioned forces 132 if (ab_mover%goprecon==0)then 133 residual(:,:)=fcart(:,:) 134 else 135 residual(:,:)= forstr%fcart(:,:) 136 end if 137 138 !write(std_out,*) '01' 139 !########################################################## 140 !### 02. Get or compute de time step dtion 141 142 if (ab_mover%dtion>0)then 143 hh = ab_mover%dtion 144 else 145 hh=fdtion(ab_mover,itime,xcart,fcart,vel) 146 end if 147 148 lambda=hh 149 write(std_out,*) 'Lambda',lambda 150 151 !write(std_out,*) '02' 152 !########################################################## 153 !### 02. For all atoms and directions 154 do kk=1,ab_mover%natom 155 ! Normally this is the mass of the atom 156 ! em=ab_mover%amass(kk) 157 ! But the steepest algorithm is using unitary mass 158 em=1 159 do jj=1,3 160 161 ! write(std_out,*) '03' 162 ! ########################################################## 163 ! ### 03. Filling other values from history (forces and vel) 164 f_cart=residual(jj,kk) 165 xc=xcart(jj,kk) 166 ! This lambda is taken from the kinematical equation 167 ! lambda=hh*hh/(2*em) 168 169 ! write(std_out,*) '04' 170 ! ########################################################## 171 ! ### 04. Take first the atoms that are not allowed to move along 172 ! ### this direction 173 ! ### Warning : implemented in cartesian coordinates 174 if (ab_mover%iatfix(jj,kk)==1) then 175 ! Their positions will be the same as xcart 176 xnow=xc 177 else 178 179 ! This is the main expresion (1) 180 xnow=xc+lambda*f_cart 181 182 end if !if(ab_mover%iatfix(jj,kk)==1) 183 184 ! write(std_out,*) '05' 185 ! ########################################################## 186 ! ### 08. Update history 187 188 xcart(jj,kk)=xnow 189 190 ! write(std_out,*) '06' 191 ! ########################################################## 192 ! ### 09. End loops of atoms and directions 193 end do ! jj=1,3 194 end do ! kk=1,ab_mover%natom 195 196 if (ab_mover%optcell/=0)then 197 198 if (ab_mover%optcell==1)then 199 do jj=1,3 200 acell(jj)=acell(jj)+lambda*strten(jj) 201 end do ! jj=1,3 202 call mkrdim(acell,rprim,rprimd) 203 elseif (ab_mover%optcell==2)then 204 do kk=1,3 205 do jj=1,3 206 if (jj==1 .and. kk==1) str=strten(1) 207 if (jj==2 .and. kk==2) str=strten(2) 208 if (jj==3 .and. kk==3) str=strten(3) 209 if (jj==1 .and. kk==2) str=strten(6) 210 if (jj==1 .and. kk==3) str=strten(5) 211 if (jj==2 .and. kk==1) str=strten(6) 212 if (jj==3 .and. kk==1) str=strten(5) 213 if (jj==2 .and. kk==3) str=strten(4) 214 if (jj==3 .and. kk==2) str=strten(4) 215 rprimd(jj,kk)=rprimd(jj,kk)+lambda*str 216 end do ! jj=1,3 217 end do ! kk=1,3 218 call mkradim(acell,rprim,rprimd) 219 end if 220 221 end if 222 223 224 !write(std_out,*) '08' 225 !########################################################## 226 !### 10. Filling history with the new values 227 228 !Increase indices 229 hist%ihist = abihist_findIndex(hist,+1) 230 231 !Compute xred from xcart, and rprimd 232 call xcart2xred(ab_mover%natom,rprimd,xcart,xred) 233 234 call var2hist(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG) 235 ihist_prev = abihist_findIndex(hist,-1) 236 hist%vel(:,:,hist%ihist)=hist%vel(:,:,ihist_prev) 237 238 end subroutine pred_steepdesc