TABLE OF CONTENTS


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]
 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].
 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.

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 .
 For the initials of contributors,
 see ~abinit/doc/developers/contributors.txt .

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

OUTPUT

SIDE EFFECTS

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

NOTES

PARENTS

      mover

CHILDREN

      hist2var,metric,var2hist,xcart2xred,xred2xcart

SOURCE

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