TABLE OF CONTENTS


ABINIT/m_pred_steepdesc [ Modules ]

[ Top ] [ Modules ]

NAME

  m_pred_steepdesc

FUNCTION

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

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 module m_pred_steepdesc
27 
28  use defs_basis
29  use m_abicore
30  use m_abimover
31  use m_abihist
32 
33  use m_geometry,       only : mkradim, mkrdim, xcart2xred, xred2xcart
34  use m_predtk,         only : fdtion
35 
36  implicit none
37 
38  private

ABINIT/pred_steepdesc [ Functions ]

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

PARENTS

      mover

CHILDREN

      hist2var,mkradim,mkrdim,var2hist,xcart2xred,xred2xcart

SOURCE

 92 subroutine pred_steepdesc(ab_mover,forstr,hist,itime,zDEBUG,iexit)
 93 
 94 
 95 !This section has been created automatically by the script Abilint (TD).
 96 !Do not modify the following lines by hand.
 97 #undef ABI_FUNC
 98 #define ABI_FUNC 'pred_steepdesc'
 99 !End of the abilint section
100 
101 implicit none
102 
103 !Arguments ------------------------------------
104 !scalars
105 type(abimover),intent(in)       :: ab_mover
106 type(abihist),intent(inout),target :: hist
107 type(abiforstr),intent(in) :: forstr
108 integer,intent(in)    :: itime,iexit
109 logical,intent(in)    :: zDEBUG
110 
111 !Local variables-------------------------------
112 !scalars
113 integer  :: kk,jj,ihist_prev
114 real(dp) :: em
115 real(dp) :: f_cart
116 real(dp) :: xc,str
117 real(dp) :: xnow,lambda
118 real(dp),save :: hh
119 !arrays
120 real(dp) :: acell(3),strten(6)
121 real(dp) :: rprim(3,3),rprimd(3,3)
122 real(dp) :: xred(3,ab_mover%natom),xcart(3,ab_mover%natom)
123 real(dp) :: residual(3,ab_mover%natom)
124 real(dp), ABI_CONTIGUOUS pointer :: fcart(:,:),vel(:,:)
125 
126 !***************************************************************************
127 !Beginning of executable session
128 !***************************************************************************
129 
130  if(iexit/=0)then
131    return
132  end if
133 
134 !write(std_out,*) '01'
135 !##########################################################
136 !### 01. Copy from the history to the variables
137  call hist2var(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG)
138 
139  do jj=1,3
140    rprim(jj,1:3)=rprimd(jj,1:3)/acell(1:3)
141  end do
142 
143  call xred2xcart(ab_mover%natom,rprimd,xcart,xred)
144  strten(:)=hist%strten(:,hist%ihist)
145  fcart => hist%fcart(:,:,hist%ihist)
146  vel => hist%vel(:,:,hist%ihist)
147 
148 !Fill the residual with forces (No preconditioning)
149 !Or the preconditioned forces
150  if (ab_mover%goprecon==0)then
151    residual(:,:)=fcart(:,:)
152  else
153    residual(:,:)= forstr%fcart(:,:)
154  end if
155 
156 !write(std_out,*) '01'
157 !##########################################################
158 !### 02. Get or compute de time step dtion
159 
160  if (ab_mover%dtion>0)then
161    hh = ab_mover%dtion
162  else
163    hh=fdtion(ab_mover,itime,xcart,fcart,vel)
164  end if
165 
166  lambda=hh
167  write(std_out,*) 'Lambda',lambda
168 
169 !write(std_out,*) '02'
170 !##########################################################
171 !### 02. For all atoms and directions
172  do kk=1,ab_mover%natom
173 !  Normally this is the mass of the atom
174 !  em=ab_mover%amass(kk)
175 !  But the steepest algorithm is using unitary mass
176    em=1
177    do jj=1,3
178 
179 !    write(std_out,*) '03'
180 !    ##########################################################
181 !    ### 03. Filling other values from history (forces and vel)
182      f_cart=residual(jj,kk)
183      xc=xcart(jj,kk)
184 !    This lambda is taken from the kinematical equation
185 !    lambda=hh*hh/(2*em)
186 
187 !    write(std_out,*) '04'
188 !    ##########################################################
189 !    ### 04. Take first the atoms that are not allowed to move along
190 !    ###     this direction
191 !    ###     Warning : implemented in cartesian coordinates
192      if (ab_mover%iatfix(jj,kk)==1) then
193 !      Their positions will be the same as xcart
194        xnow=xc
195      else
196 
197 !      This is the main expresion (1)
198        xnow=xc+lambda*f_cart
199 
200      end if !if(ab_mover%iatfix(jj,kk)==1)
201 
202 !    write(std_out,*) '05'
203 !    ##########################################################
204 !    ### 08. Update history
205 
206      xcart(jj,kk)=xnow
207 
208 !    write(std_out,*) '06'
209 !    ##########################################################
210 !    ### 09. End loops of atoms and directions
211    end do ! jj=1,3
212  end do ! kk=1,ab_mover%natom
213 
214  if (ab_mover%optcell/=0)then
215 
216    if (ab_mover%optcell==1)then
217      do jj=1,3
218        acell(jj)=acell(jj)+lambda*strten(jj)
219      end do ! jj=1,3
220      call mkrdim(acell,rprim,rprimd)
221    elseif (ab_mover%optcell==2)then
222      do kk=1,3
223        do jj=1,3
224          if (jj==1 .and. kk==1) str=strten(1)
225          if (jj==2 .and. kk==2) str=strten(2)
226          if (jj==3 .and. kk==3) str=strten(3)
227          if (jj==1 .and. kk==2) str=strten(6)
228          if (jj==1 .and. kk==3) str=strten(5)
229          if (jj==2 .and. kk==1) str=strten(6)
230          if (jj==3 .and. kk==1) str=strten(5)
231          if (jj==2 .and. kk==3) str=strten(4)
232          if (jj==3 .and. kk==2) str=strten(4)
233          rprimd(jj,kk)=rprimd(jj,kk)+lambda*str
234        end do ! jj=1,3
235      end do ! kk=1,3
236      call mkradim(acell,rprim,rprimd)
237    end if
238 
239  end if
240 
241 
242 !write(std_out,*) '08'
243 !##########################################################
244 !### 10. Filling history with the new values
245 
246 !Increase indices
247  hist%ihist = abihist_findIndex(hist,+1)
248 
249 !Compute xred from xcart, and rprimd
250  call xcart2xred(ab_mover%natom,rprimd,xcart,xred)
251 
252  call var2hist(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG)
253  ihist_prev = abihist_findIndex(hist,-1)
254  hist%vel(:,:,hist%ihist)=hist%vel(:,:,ihist_prev)
255 
256 end subroutine pred_steepdesc