TABLE OF CONTENTS
ABINIT/m_predict_neb [ Modules ]
NAME
m_predict_neb
FUNCTION
COPYRIGHT
Copyright (C) 2012-2024 ABINIT group (MT) 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_predict_neb 23 24 use defs_basis 25 use m_splines 26 use m_mep 27 use m_abicore 28 use m_errors 29 use m_xmpi 30 31 use defs_abitypes, only : MPI_type 32 use m_results_img, only : results_img_type, gather_array_img, scatter_array_img, get_geometry_img 33 34 implicit none 35 36 private
ABINIT/predict_neb [ Functions ]
NAME
predict_neb
FUNCTION
Given the past history of images, predict the new set of images using Nudged Elastic Band method. No change of acell, rprim and vel at present.
INPUTS
itimimage=time index for image propagation (itimimage+1 is to be predicted here)
itimimage_eff=time index in the history
list_dynimage(nimage)=list of dynamical images. The non-dynamical ones will not change.
Example: in the NEB of string method, one expect the two end images to be fixed.
This is quite useful when ground states of the A and B states is known
mpi_enreg=MPI-parallelisation information
natom=dimension of vel_timimage and xred_timimage
ndynimage=number of dynamical images
nimage=number of images (on current proc)
nimage_tot=total number of images
ntimimage_stored=number of time steps stored in the history
OUTPUT
SIDE EFFECTS
mep_param=several parameters for Minimal Energy Path (MEP) search
results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
results_img(:,:)%acell(3)
at input, history of the values of acell for all images
at output, the predicted values of acell for all images
results_img(:,:)%results_gs
at input, history of the values of energies and forces for all images
results_img(:,:)%rprim(3,3)
at input, history of the values of rprim for all images
at output, the predicted values of rprim for all images
results_img(:,:)%vel(3,natom)
at input, history of the values of vel for all images
at output, the predicted values of vel for all images
results_img(:,:)%vel_cell(3,3)
at input, history of the values of vel_cell for all images
at output, the predicted values of vel_cell for all images
results_img(:,:)%xred(3,natom)
at input, history of the values of xred for all images
at output, the predicted values of xred for all images
SOURCE
91 subroutine predict_neb(itimimage,itimimage_eff,list_dynimage,mep_param,mpi_enreg,natom,& 92 & ndynimage,nimage,nimage_tot,ntimimage_stored,results_img) 93 94 implicit none 95 96 !Arguments ------------------------------------ 97 !scalars 98 integer,intent(in) :: itimimage,itimimage_eff,natom,ndynimage 99 integer,intent(in) :: nimage,nimage_tot,ntimimage_stored 100 type(mep_type),intent(inout) :: mep_param 101 type(MPI_type),intent(in) :: mpi_enreg 102 !arrays 103 integer,intent(in) :: list_dynimage(ndynimage) 104 type(results_img_type),intent(inout) :: results_img(nimage,ntimimage_stored) 105 106 !Local variables------------------------------- 107 !scalars 108 integer :: ierr,ii,iimage,iimage_min,iimage_max,next_itimimage 109 real(dp) :: dvmax,dvmin,ecur,emax,emin,eref,f_para1,f_para2 110 logical :: test_minus_one,test_plus_one 111 !arrays 112 real(dp),allocatable :: buffer(:,:),buffer_all(:,:) 113 real(dp),allocatable :: coordif(:,:,:),dimage(:),rprimd(:,:,:),spring(:) 114 real(dp),allocatable :: tangent(:,:,:),vect(:,:) 115 real(dp),allocatable,target :: etotal(:),fcart(:,:,:),neb_forces(:,:,:) 116 real(dp),allocatable,target :: xcart(:,:,:),xred(:,:,:) 117 real(dp),pointer :: etotal_all(:),fcart_all(:,:,:),neb_forces_all(:,:,:) 118 real(dp),pointer :: xcart_all(:,:,:) 119 120 ! ************************************************************************* 121 122 ABI_MALLOC(xred,(3,natom,nimage)) 123 124 !Parallelism over images: only one process per image of the cell 125 if (mpi_enreg%me_cell==0) then 126 127 ! Retrieve positions and forces 128 ABI_MALLOC(etotal,(nimage)) 129 ABI_MALLOC(xcart,(3,natom,nimage)) 130 ABI_MALLOC(fcart,(3,natom,nimage)) 131 ABI_MALLOC(rprimd,(3,3,nimage)) 132 call get_geometry_img(etotal,natom,nimage,results_img(:,itimimage_eff),& 133 & fcart,rprimd,xcart,xred) 134 135 ! Array containing effective NEB forces (F_ortho+F_spring) 136 ABI_MALLOC(neb_forces,(3,natom,nimage)) 137 138 ! Parallelism: gather data of all images 139 if (mpi_enreg%paral_img==1) then 140 ABI_MALLOC(buffer,(6*natom+1,nimage)) 141 ABI_MALLOC(buffer_all,(6*natom+1,nimage_tot)) 142 ABI_MALLOC(xcart_all,(3,natom,nimage_tot)) 143 ABI_MALLOC(fcart_all,(3,natom,nimage_tot)) 144 ABI_MALLOC(etotal_all,(nimage_tot)) 145 ABI_MALLOC(neb_forces_all,(3,natom,nimage_tot)) 146 buffer=zero;ii=0 147 buffer(ii+1:ii+3*natom,1:nimage)=reshape(xcart,(/3*natom,nimage/));ii=3*natom 148 buffer(ii+1:ii+3*natom,1:nimage)=reshape(fcart,(/3*natom,nimage/));ii=6*natom 149 buffer(ii+1 ,1:nimage)=etotal(1:nimage);ii=0 150 call gather_array_img(buffer,buffer_all,mpi_enreg,allgather=.true.) 151 xcart_all(:,:,:)=reshape(buffer_all(ii+1:ii+3*natom,1:nimage_tot),(/3,natom,nimage_tot/));ii=3*natom 152 fcart_all(:,:,:)=reshape(buffer_all(ii+1:ii+3*natom,1:nimage_tot),(/3,natom,nimage_tot/));ii=6*natom 153 etotal_all(:)=buffer_all(ii+1,1:nimage_tot);ii=0 154 ABI_FREE(buffer) 155 ABI_FREE(buffer_all) 156 else 157 xcart_all => xcart 158 fcart_all => fcart 159 etotal_all => etotal 160 neb_forces_all => neb_forces 161 end if 162 163 ! coordif is the vector between two images 164 ! dimage is the distance between two images 165 ! tangent is the tangent at image i 166 ABI_MALLOC(coordif,(3,natom,nimage_tot)) 167 ABI_MALLOC(tangent,(3,natom,nimage_tot)) 168 ABI_MALLOC(dimage,(nimage_tot)) 169 170 ! Compute distances between images 171 coordif(:,:,1)=zero;dimage(1)=zero 172 do iimage=2,nimage_tot 173 coordif(:,:,iimage)=xcart_all(:,:,iimage)-xcart_all(:,:,iimage-1) 174 dimage(iimage)=mep_img_norm(coordif(:,:,iimage)) 175 end do 176 177 ! Compute tangent (not normalized) 178 tangent(:,:,1)=zero 179 tangent(:,:,nimage_tot)=zero 180 ! === Original definition of tangent 181 if (mep_param%neb_algo==0) then 182 do iimage=2,nimage_tot-1 183 ! tangent(:,:,iimage)=xcart_all(:,:,iimage+1)-xcart_all(:,:,iimage-1) 184 tangent(:,:,iimage)=coordif(:,:,iimage )/dimage(iimage) & 185 & +coordif(:,:,iimage+1)/dimage(iimage+1) 186 end do 187 ! === Improved tangent (J. Chem. Phys. 113, 9978 (2000) [[cite:Henkelman2000]]) 188 else 189 do iimage=2,nimage_tot-1 190 test_minus_one=(etotal_all(iimage-1)<etotal_all(iimage)) 191 test_plus_one =(etotal_all(iimage )<etotal_all(iimage+1)) 192 if (test_minus_one.and.test_plus_one)then ! V_i-1 < V_i < V_i+1 193 tangent(:,:,iimage)=coordif(:,:,iimage+1) 194 else if ((.not.test_minus_one).and.(.not.test_plus_one))then ! V_i-1 > V_i > V_i+1 195 tangent(:,:,iimage)=coordif(:,:,iimage) 196 else ! V_i-1 < V_i > V_i+1 OR V_i-1 > V_i < V_i+1 197 dvmax=max(abs(etotal_all(iimage+1)-etotal_all(iimage)),& 198 & abs(etotal_all(iimage-1)-etotal_all(iimage))) 199 dvmin=min(abs(etotal_all(iimage+1)-etotal_all(iimage)),& 200 & abs(etotal_all(iimage-1)-etotal_all(iimage))) 201 if (etotal_all(iimage+1)>etotal_all(iimage-1)) then ! V_i+1 > V_i-1 202 tangent(:,:,iimage)=coordif(:,:,iimage+1)*dvmax & 203 & +coordif(:,:,iimage )*dvmin 204 else ! V_i+1 < V_i-1 205 tangent(:,:,iimage)=coordif(:,:,iimage+1)*dvmin & 206 & +coordif(:,:,iimage )*dvmax 207 end if 208 end if 209 end do 210 end if 211 212 ! Normalize tangent 213 do iimage=2,nimage_tot-1 214 tangent(:,:,iimage)=tangent(:,:,iimage)/mep_img_norm(tangent(:,:,iimage)) 215 end do 216 217 ! Compute spring constant(s) 218 ABI_MALLOC(spring,(nimage_tot)) 219 if (abs(mep_param%neb_spring(2)-mep_param%neb_spring(1))<tol8) then ! Constant spring 220 spring(:)=mep_param%neb_spring(1) 221 else ! Variable spring 222 emax=maxval(etotal_all(:)) 223 eref=max(etotal_all(1),etotal_all(nimage_tot)) 224 spring(:)=mep_param%neb_spring(1) 225 do iimage=2,nimage_tot-1 226 ecur=max(etotal_all(iimage-1),etotal_all(iimage)) 227 if (ecur<eref) then 228 spring(iimage)=mep_param%neb_spring(1) 229 else 230 spring(iimage)=mep_param%neb_spring(2) & 231 & -(mep_param%neb_spring(2)-mep_param%neb_spring(1)) & 232 & *(emax-ecur)/(emax-eref) 233 end if 234 end do 235 end if 236 237 ! CI-NEB: determine image(s) with maximal energy 238 iimage_min=-1;iimage_max=-1 239 if (mep_param%neb_algo==2.and.itimimage>=mep_param%cineb_start) then 240 emin=min(etotal_all(1),etotal_all(nimage_tot)) 241 emax=max(etotal_all(1),etotal_all(nimage_tot)) 242 do iimage=2,nimage_tot-1 243 if (etotal_all(iimage)<emin) then 244 iimage_min=iimage;emin=etotal_all(iimage) 245 end if 246 if (etotal_all(iimage)>emax) then 247 iimage_max=iimage;emax=etotal_all(iimage) 248 end if 249 end do 250 end if 251 252 ! Compute NEB forces 253 neb_forces_all(:,:,1)=fcart_all(:,:,1) 254 neb_forces_all(:,:,nimage_tot)=fcart_all(:,:,nimage_tot) 255 do iimage=2,nimage_tot-1 256 ! === Standard NEB 257 if (iimage/=iimage_max) then 258 ! if (iimage/=iimage_min.and.iimage/=iimage_max) then 259 f_para1=mep_img_dotp(fcart_all(:,:,iimage),tangent(:,:,iimage)) 260 if (mep_param%neb_algo==0) then ! Original NEB algo 261 ABI_MALLOC(vect,(3,natom)) 262 vect(:,:)=spring(iimage+1)*coordif(:,:,iimage+1)-spring(iimage)*coordif(:,:,iimage) 263 f_para2=mep_img_dotp(vect,tangent(:,:,iimage)) 264 ABI_FREE(vect) 265 else ! Modification from J. Chem. Phys. 113, 9978 (2000) [[cite:Henkelman2000]] 266 f_para2=spring(iimage+1)*dimage(iimage+1)-spring(iimage)*dimage(iimage) 267 end if 268 neb_forces_all(:,:,iimage)=fcart_all(:,:,iimage) & ! F_ortho 269 & -f_para1*tangent(:,:,iimage) & ! F_ortho 270 & +f_para2*tangent(:,:,iimage) ! F_spring 271 ! === CI-NEB for the image(s) with maximal energy 272 else 273 f_para1=mep_img_dotp(fcart_all(:,:,iimage),tangent(:,:,iimage)) 274 neb_forces_all(:,:,iimage)=fcart_all(:,:,iimage)-two*f_para1*tangent(:,:,iimage) 275 end if 276 end do 277 278 ! Parallelism: distribute NEB forces 279 if (mpi_enreg%paral_img==1) then 280 do iimage=1,nimage 281 ii=mpi_enreg%my_imgtab(iimage) 282 neb_forces(:,:,iimage)=neb_forces_all(:,:,ii) 283 end do 284 end if 285 286 ! Compute new atomic positions in each cell 287 if (mep_param%mep_solver==0) then ! Steepest-descent 288 call mep_steepest(neb_forces,list_dynimage,mep_param,natom,ndynimage,nimage,rprimd,xcart,xred) 289 else if (mep_param%mep_solver==1) then ! Quick-min 290 call mep_qmin(neb_forces,itimimage,list_dynimage,mep_param,natom,ndynimage, & 291 & nimage,rprimd,xcart,xred) 292 else if (mep_param%mep_solver==2) then ! Local BFGS 293 call mep_lbfgs(neb_forces,itimimage,list_dynimage,mep_param,natom,ndynimage,& 294 & nimage,rprimd,xcart,xred) 295 else if (mep_param%mep_solver==3) then ! Global BFGS 296 call mep_gbfgs(neb_forces,itimimage,list_dynimage,mep_param,mpi_enreg,natom,ndynimage,& 297 & nimage,nimage_tot,rprimd,xcart,xred) 298 else 299 ABI_BUG("Inconsistent solver !") 300 end if 301 302 ! Free memory 303 ABI_FREE(spring) 304 ABI_FREE(coordif) 305 ABI_FREE(tangent) 306 ABI_FREE(dimage) 307 if (mpi_enreg%paral_img==1) then 308 ABI_FREE(xcart_all) 309 ABI_FREE(fcart_all) 310 ABI_FREE(etotal_all) 311 ABI_FREE(neb_forces_all) 312 end if 313 314 ABI_FREE(neb_forces) 315 ABI_FREE(etotal) 316 ABI_FREE(xcart) 317 ABI_FREE(fcart) 318 ABI_FREE(rprimd) 319 end if ! mpi_enreg%me_cell==0 320 321 !Store acell, rprim, xred and vel for the new iteration 322 call xmpi_bcast(xred,0,mpi_enreg%comm_cell,ierr) 323 next_itimimage=itimimage_eff+1 324 if (next_itimimage>ntimimage_stored) next_itimimage=1 325 do iimage=1,nimage 326 results_img(iimage,next_itimimage)%xred(:,:) =xred(:,:,iimage) 327 results_img(iimage,next_itimimage)%acell(:) =results_img(iimage,itimimage_eff)%acell(:) 328 results_img(iimage,next_itimimage)%rprim(:,:) =results_img(iimage,itimimage_eff)%rprim(:,:) 329 results_img(iimage,next_itimimage)%vel(:,:) =results_img(iimage,itimimage_eff)%vel(:,:) 330 results_img(iimage,next_itimimage)%vel_cell(:,:)=results_img(iimage,itimimage_eff)%vel_cell(:,:) 331 end do 332 ABI_FREE(xred) 333 334 end subroutine predict_neb