TABLE OF CONTENTS


ABINIT/m_pred_srkhna14 [ Modules ]

[ Top ] [ Modules ]

NAME

  m_pred_srkna14

FUNCTION

COPYRIGHT

  Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR, JCC, 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_pred_srkhna14
28 
29  use defs_basis
30  use m_abicore
31  use m_abimover
32  use m_abihist
33 
34  use m_geometry,    only : xcart2xred, xred2xcart, metric
35 
36  implicit none
37 
38  private

ABINIT/pred_srkna14 [ Functions ]

[ Top ] [ Functions ]

NAME

 pred_srkna14

FUNCTION

 Ionmov predictors (14) Srkna14 molecular dynamics

 IONMOV 14:
 Simple molecular dynamics with a symplectic algorithm proposed
 by S.Blanes and P.C.Moans, called SRKNa14 in Practical symplectic partitioned
 Runge--Kutta and Runge--Kutta--Nystrom methods, Journal of Computational
 and Applied Mathematics archive, volume 142,  issue 2  (May 2002), pages 313 - 330 [[cite:Blanes2002]].
 of the kind first published by H. Yoshida, Construction of higher order symplectic
 integrators, Physics Letters A, volume 150, number 5 to 7, pages 262 - 268 [[cite:Yoshida1990]]
 This algorithm requires at least 14 evaluation of the forces (actually 15 are done
 within Abinit) per time step. At this cost it usually gives much better
 energy conservation than the verlet algorithm (ionmov 6) for a 30 times bigger
 value of <a href="varrlx.html#dtion">dtion</a>. Notice that the potential
 energy of the initial atomic configuration is never evaluated using this
 algorithm.

INPUTS

 ab_mover <type(abimover)> : Datatype with all the information needed by the preditor
 itime  : Index of the present iteration
 ntime  : Maximal number of iterations
 icycle : Index of the present cycle
 zDEBUG : if true print some debugging information

SIDE EFFECTS

 hist <type(abihist)> : History of positions,forces acell, rprimd, stresses

PARENTS

      mover

CHILDREN

      hist2var,metric,var2hist,xcart2xred,xred2xcart

SOURCE

 86 subroutine pred_srkna14(ab_mover,hist,icycle,zDEBUG,iexit,skipcycle)
 87 
 88 
 89 !This section has been created automatically by the script Abilint (TD).
 90 !Do not modify the following lines by hand.
 91 #undef ABI_FUNC
 92 #define ABI_FUNC 'pred_srkna14'
 93 !End of the abilint section
 94 
 95  implicit none
 96 
 97 !Arguments ------------------------------------
 98 !scalars
 99  type(abimover),intent(in)       :: ab_mover
100  type(abihist),intent(inout) :: hist
101  integer,intent(in)    :: icycle
102  integer,intent(in)    :: iexit
103  logical,intent(in)    :: zDEBUG
104  logical,intent(out)   :: skipcycle
105 
106 !Local variables-------------------------------
107 !scalars
108  integer  :: ihist_prev,ii,jj,kk
109  real(dp) :: ucvol,ucvol_next
110  real(dp),parameter :: v2tol=tol8
111  real(dp) :: etotal
112  logical  :: jump_end_of_cycle=.FALSE.
113 ! character(len=5000) :: message
114 !arrays
115  real(dp),save :: aa(15),bb(15)
116  real(dp) :: acell(3),acell_next(3)
117  real(dp) :: rprimd(3,3),rprimd_next(3,3)
118  real(dp) :: gprimd(3,3),gmet(3,3),rmet(3,3)
119  real(dp) :: fcart(3,ab_mover%natom),fcart_m(3,ab_mover%natom)
120 !real(dp) :: fred_corrected(3,ab_mover%natom)
121  real(dp) :: xcart(3,ab_mover%natom)
122  real(dp) :: xred(3,ab_mover%natom)
123  real(dp) :: vel(3,ab_mover%natom)
124  real(dp) :: strten(6)
125 
126 !***************************************************************************
127 !Beginning of executable session
128 !***************************************************************************
129 
130  if(iexit/=0)then
131    return
132  end if
133 
134  jump_end_of_cycle=.FALSE.
135  fcart_m(:,:)=zero
136 
137 !write(std_out,*) 'srkna14 03',jump_end_of_cycle
138 !##########################################################
139 !### 03. Obtain the present values from the history
140 
141  call hist2var(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG)
142  call metric(gmet,gprimd,-1,rmet,rprimd,ucvol)
143 
144  call xred2xcart(ab_mover%natom,rprimd,xcart,xred)
145 
146  fcart(:,:)=hist%fcart(:,:,hist%ihist)
147  strten(:) =hist%strten(:,hist%ihist)
148  vel(:,:)  =hist%vel(:,:,hist%ihist)
149  etotal    =hist%etot(hist%ihist)
150 
151  if(zDEBUG)then
152    write (std_out,*) 'fcart:'
153    do kk=1,ab_mover%natom
154      write (std_out,*) fcart(:,kk)
155    end do
156    write (std_out,*) 'vel:'
157    do kk=1,ab_mover%natom
158      write (std_out,*) vel(:,kk)
159    end do
160    write (std_out,*) 'strten:'
161    write (std_out,*) strten(1:3),ch10,strten(4:6)
162    write (std_out,*) 'etotal:'
163    write (std_out,*) etotal
164  end if
165 
166  write(std_out,*) 'RMET'
167  do ii=1,3
168    write(std_out,*) rmet(ii,:)
169  end do
170 
171 !Get rid of mean force on whole unit cell, but only if no
172 !generalized constraints are in effect
173 !  call fcart2fred(fcart,fred_corrected,rprimd,ab_mover%natom)
174 !  if(ab_mover%nconeq==0)then
175 !    amass_tot=sum(ab_mover%amass(:))
176 !    do ii=1,3
177 !      if (ii/=3.or.ab_mover%jellslab==0) then
178 !        favg=sum(fred_corrected(ii,:))/dble(ab_mover%natom)
179 !        fred_corrected(ii,:)=fred_corrected(ii,:)-favg*ab_mover%amass(:)/amass_tot
180 !      end if
181 !    end do
182 !  end if
183 
184 !write(std_out,*) 'srkna14 04',jump_end_of_cycle
185 !##########################################################
186 !### 04. Compute the next values (Only for the first cycle)
187 
188  if (icycle==1) then
189 
190    if(zDEBUG) then
191      write(std_out,*) 'Entering only for first cycle'
192    end if
193 
194    aa(1) =  0.0378593198406116_dp;
195    aa(2) =  0.102635633102435_dp;
196    aa(3) = -0.0258678882665587_dp;
197    aa(4) =  0.314241403071447_dp;
198    aa(5) = -0.130144459517415_dp;
199    aa(6) =  0.106417700369543_dp;
200    aa(7) = -0.00879424312851058_dp;
201    aa(8) =  1._dp -&
202 &   2._dp*(aa(1)+aa(2)+aa(3)+aa(4)+aa(5)+aa(6)+aa(7));
203    aa(9) =  aa(7);
204    aa(10)=  aa(6);
205    aa(11)=  aa(5);
206    aa(12)=  aa(4);
207    aa(13)=  aa(3);
208    aa(14)=  aa(2);
209    aa(15)=  aa(1);
210 
211    bb(1) =  0.0_dp
212    bb(2) =  0.09171915262446165_dp;
213    bb(3) =  0.183983170005006_dp;
214    bb(4) = -0.05653436583288827_dp;
215    bb(5) =  0.004914688774712854_dp;
216    bb(6) =  0.143761127168358_dp;
217    bb(7) =  0.328567693746804_dp;
218    bb(8) =  0.5_dp - (bb(1)+bb(2)+bb(3)+bb(4)+bb(5)+bb(6)+bb(7));
219    bb(9) =  0.5_dp - (bb(1)+bb(2)+bb(3)+bb(4)+bb(5)+bb(6)+bb(7));
220    bb(10)=  bb(7);
221    bb(11)=  bb(6);
222    bb(12)=  bb(5);
223    bb(13)=  bb(4);
224    bb(14)=  bb(3);
225    bb(15)=  bb(2);
226 
227    acell_next(:)=acell(:)
228    ucvol_next=ucvol
229    rprimd_next(:,:)=rprimd(:,:)
230 
231 !  step 1 of 15
232 
233 !  Convert input xred (reduced coordinates) to xcart (cartesian)
234    call xred2xcart(ab_mover%natom,rprimd,xcart,xred)
235 
236    vel(:,:) = vel(:,:) + bb(1) * ab_mover%dtion * fcart_m(:,:)
237 
238    do ii=1,3
239      do jj=1,ab_mover%natom
240        write(std_out,*) xcart(ii,jj), ab_mover%dtion, aa(1), vel(ii,jj)
241        xcart(ii,jj) = xcart(ii,jj) + ab_mover%dtion * aa(1) * vel(ii,jj)
242        write(std_out,*) xcart(ii,jj)
243      end do
244    end do
245 
246 !  xcart(:,:) = xcart(:,:) + ab_mover%dtion * aa(1) * vel(:,:);
247 
248 !  Convert back to xred (reduced coordinates)
249    call xcart2xred(ab_mover%natom,rprimd,xcart,xred)
250 
251  end if ! if (icycle==1)
252 
253 !write(std_out,*) 'srkna14 05',jump_end_of_cycle
254 !##########################################################
255 !### 05. Compute the next values (Only for extra cycles)
256 
257  if (icycle>1) then
258 
259    do ii=1,ab_mover%natom
260      do jj=1,3
261        fcart_m(jj,ii) = fcart(jj,ii)/ab_mover%amass(ii)
262      end do
263    end do
264 
265    if (icycle<16)then
266 
267 !    Update of velocities and positions
268      vel(:,:) = vel(:,:) + bb(icycle) * ab_mover%dtion * fcart_m(:,:)
269      xcart(:,:) = xcart(:,:) +&
270 &     aa(icycle) * ab_mover%dtion * vel(:,:)
271 !    Convert xcart_next to xred_next (reduced coordinates)
272 !    for scfcv
273      call xcart2xred(ab_mover%natom, rprimd, xcart,&
274 &     xred)
275 
276    end if ! (ii<16)
277 
278  end if ! if (icycle>1)
279 
280 !write(std_out,*) 'srkna14 06',jump_end_of_cycle
281 !##########################################################
282 !### 06. Compute the next values (Only for the last cycle)
283 
284  if(jump_end_of_cycle)then
285    skipcycle=.TRUE.
286  else
287    skipcycle=.FALSE.
288  end if
289 
290 !write(std_out,*) 'srkna14 07',jump_end_of_cycle
291 !##########################################################
292 !### 07. Update the history with the prediction
293 
294 !Increase indexes
295  hist%ihist = abihist_findIndex(hist,+1)
296 
297 !Fill the history with the variables
298 !xred, acell, rprimd, vel
299  call var2hist(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG)
300  hist%vel(:,:,hist%ihist)=vel(:,:)
301  ihist_prev = abihist_findIndex(hist,-1)
302  hist%time(hist%ihist)=hist%time(ihist_prev)+ab_mover%dtion
303 
304 end subroutine pred_srkna14