TABLE OF CONTENTS
ABINIT/m_pred_srkhna14 [ Modules ]
NAME
m_pred_srkna14
FUNCTION
COPYRIGHT
Copyright (C) 1998-2022 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 .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_pred_srkhna14 23 24 use defs_basis 25 use m_abicore 26 use m_abimover 27 use m_abihist 28 29 use m_geometry, only : xcart2xred, xred2xcart, metric 30 31 implicit none 32 33 private
ABINIT/pred_srkna14 [ 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
SOURCE
75 subroutine pred_srkna14(ab_mover,hist,icycle,zDEBUG,iexit,skipcycle) 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) :: xcart(3,ab_mover%natom) 103 real(dp) :: xred(3,ab_mover%natom) 104 real(dp) :: vel(3,ab_mover%natom) 105 real(dp) :: strten(6) 106 107 !*************************************************************************** 108 !Beginning of executable session 109 !*************************************************************************** 110 111 if(iexit/=0)then 112 return 113 end if 114 115 jump_end_of_cycle=.FALSE. 116 fcart_m(:,:)=zero 117 118 !write(std_out,*) 'srkna14 03',jump_end_of_cycle 119 !########################################################## 120 !### 03. Obtain the present values from the history 121 122 call hist2var(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG) 123 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 124 125 call xred2xcart(ab_mover%natom,rprimd,xcart,xred) 126 127 fcart(:,:)=hist%fcart(:,:,hist%ihist) 128 strten(:) =hist%strten(:,hist%ihist) 129 vel(:,:) =hist%vel(:,:,hist%ihist) 130 etotal =hist%etot(hist%ihist) 131 132 if(zDEBUG)then 133 write (std_out,*) 'fcart:' 134 do kk=1,ab_mover%natom 135 write (std_out,*) fcart(:,kk) 136 end do 137 write (std_out,*) 'vel:' 138 do kk=1,ab_mover%natom 139 write (std_out,*) vel(:,kk) 140 end do 141 write (std_out,*) 'strten:' 142 write (std_out,*) strten(1:3),ch10,strten(4:6) 143 write (std_out,*) 'etotal:' 144 write (std_out,*) etotal 145 end if 146 147 write(std_out,*) 'RMET' 148 do ii=1,3 149 write(std_out,*) rmet(ii,:) 150 end do 151 152 !write(std_out,*) 'srkna14 04',jump_end_of_cycle 153 !########################################################## 154 !### 04. Compute the next values (Only for the first cycle) 155 156 if (icycle==1) then 157 158 if(zDEBUG) then 159 write(std_out,*) 'Entering only for first cycle' 160 end if 161 162 aa(1) = 0.0378593198406116_dp; 163 aa(2) = 0.102635633102435_dp; 164 aa(3) = -0.0258678882665587_dp; 165 aa(4) = 0.314241403071447_dp; 166 aa(5) = -0.130144459517415_dp; 167 aa(6) = 0.106417700369543_dp; 168 aa(7) = -0.00879424312851058_dp; 169 aa(8) = 1._dp -& 170 & 2._dp*(aa(1)+aa(2)+aa(3)+aa(4)+aa(5)+aa(6)+aa(7)); 171 aa(9) = aa(7); 172 aa(10)= aa(6); 173 aa(11)= aa(5); 174 aa(12)= aa(4); 175 aa(13)= aa(3); 176 aa(14)= aa(2); 177 aa(15)= aa(1); 178 179 bb(1) = 0.0_dp 180 bb(2) = 0.09171915262446165_dp; 181 bb(3) = 0.183983170005006_dp; 182 bb(4) = -0.05653436583288827_dp; 183 bb(5) = 0.004914688774712854_dp; 184 bb(6) = 0.143761127168358_dp; 185 bb(7) = 0.328567693746804_dp; 186 bb(8) = 0.5_dp - (bb(1)+bb(2)+bb(3)+bb(4)+bb(5)+bb(6)+bb(7)); 187 bb(9) = 0.5_dp - (bb(1)+bb(2)+bb(3)+bb(4)+bb(5)+bb(6)+bb(7)); 188 bb(10)= bb(7); 189 bb(11)= bb(6); 190 bb(12)= bb(5); 191 bb(13)= bb(4); 192 bb(14)= bb(3); 193 bb(15)= bb(2); 194 195 acell_next(:)=acell(:) 196 ucvol_next=ucvol 197 rprimd_next(:,:)=rprimd(:,:) 198 199 ! step 1 of 15 200 201 ! Convert input xred (reduced coordinates) to xcart (cartesian) 202 call xred2xcart(ab_mover%natom,rprimd,xcart,xred) 203 204 vel(:,:) = vel(:,:) + bb(1) * ab_mover%dtion * fcart_m(:,:) 205 206 do ii=1,3 207 do jj=1,ab_mover%natom 208 write(std_out,*) xcart(ii,jj), ab_mover%dtion, aa(1), vel(ii,jj) 209 xcart(ii,jj) = xcart(ii,jj) + ab_mover%dtion * aa(1) * vel(ii,jj) 210 write(std_out,*) xcart(ii,jj) 211 end do 212 end do 213 214 ! xcart(:,:) = xcart(:,:) + ab_mover%dtion * aa(1) * vel(:,:); 215 216 ! Convert back to xred (reduced coordinates) 217 call xcart2xred(ab_mover%natom,rprimd,xcart,xred) 218 219 end if ! if (icycle==1) 220 221 !write(std_out,*) 'srkna14 05',jump_end_of_cycle 222 !########################################################## 223 !### 05. Compute the next values (Only for extra cycles) 224 225 if (icycle>1) then 226 227 do ii=1,ab_mover%natom 228 do jj=1,3 229 fcart_m(jj,ii) = fcart(jj,ii)/ab_mover%amass(ii) 230 end do 231 end do 232 233 if (icycle<16)then 234 235 ! Update of velocities and positions 236 vel(:,:) = vel(:,:) + bb(icycle) * ab_mover%dtion * fcart_m(:,:) 237 xcart(:,:) = xcart(:,:) +& 238 & aa(icycle) * ab_mover%dtion * vel(:,:) 239 ! Convert xcart_next to xred_next (reduced coordinates) 240 ! for scfcv 241 call xcart2xred(ab_mover%natom, rprimd, xcart,& 242 & xred) 243 244 end if ! (ii<16) 245 246 end if ! if (icycle>1) 247 248 !write(std_out,*) 'srkna14 06',jump_end_of_cycle 249 !########################################################## 250 !### 06. Compute the next values (Only for the last cycle) 251 252 if(jump_end_of_cycle)then 253 skipcycle=.TRUE. 254 else 255 skipcycle=.FALSE. 256 end if 257 258 !write(std_out,*) 'srkna14 07',jump_end_of_cycle 259 !########################################################## 260 !### 07. Update the history with the prediction 261 262 !Increase indexes 263 hist%ihist = abihist_findIndex(hist,+1) 264 265 !Fill the history with the variables 266 !xred, acell, rprimd, vel 267 call var2hist(acell,hist,ab_mover%natom,rprimd,xred,zDEBUG) 268 hist%vel(:,:,hist%ihist)=vel(:,:) 269 ihist_prev = abihist_findIndex(hist,-1) 270 hist%time(hist%ihist)=hist%time(ihist_prev)+ab_mover%dtion 271 272 end subroutine pred_srkna14