TABLE OF CONTENTS
ABINIT/predict_pimd [ Functions ]
NAME
predict_pimd
FUNCTION
Predicts new positions in Path Integral Molecular Dynamics Given the positions at time t and t-dtion, an estimation of the velocities at time t, the forces and an estimation of the stress at time t, and an estimation of the cell at time t, computes in the Path Integral Molecular Dynamics framework the new positions at time t+dtion, computes self-consistently the velocities, the stress and the cell at time t and produces an estimation of the velocities, stress and new cell at time t+dtion No change of acell and rprim at present.
COPYRIGHT
Copyright (C) 2010-2018 ABINIT group (GG) 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
imgmov=gives the algorithm to be used for prediction of new set of images itimimage=time index for image propagation (itimimage+1 is to be predicted here) itimimage_eff=time index in the history mpi_enreg=MPI-parallelisation information natom=dimension of vel_timimage and xred_timimage nimage=number of images (treated by current proc) nimage_tot=total number of images ntimimage_stored=number of time steps stored in the history results_gs_timimage(ntimimage,nimage)=datastructure that hold all the history of previous computations. pimd_param=datastructure that contains all the parameters necessary to Path-Integral MD prtvolimg=printing volume
OUTPUT
SIDE EFFECTS
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
PARENTS
predictimg
CHILDREN
gather_array_img,mkradim,mkrdim,pimd_langevin_npt,pimd_langevin_nvt pimd_nosehoover_npt,pimd_nosehoover_nvt,scatter_array_img,xmpi_bcast xmpi_gather
SOURCE
67 #if defined HAVE_CONFIG_H 68 #include "config.h" 69 #endif 70 71 #include "abi_common.h" 72 73 subroutine predict_pimd(imgmov,itimimage,itimimage_eff,mpi_enreg,natom,nimage,nimage_tot,& 74 & ntimimage_stored,pimd_param,prtvolimg,results_img) 75 76 use defs_basis 77 use defs_abitypes 78 use m_profiling_abi 79 use m_pimd 80 use m_xmpi 81 use m_results_img 82 83 use m_geometry, only : mkradim, mkrdim 84 85 !This section has been created automatically by the script Abilint (TD). 86 !Do not modify the following lines by hand. 87 #undef ABI_FUNC 88 #define ABI_FUNC 'predict_pimd' 89 use interfaces_45_geomoptim, except_this_one => predict_pimd 90 !End of the abilint section 91 92 implicit none 93 94 !Arguments ------------------------------------ 95 !scalars 96 integer,intent(in) :: imgmov,itimimage,itimimage_eff,natom 97 integer,intent(in) :: nimage,nimage_tot,ntimimage_stored,prtvolimg 98 type(MPI_type),intent(in) :: mpi_enreg 99 type(pimd_type),intent(inout) :: pimd_param 100 !arrays 101 type(results_img_type) :: results_img(nimage,ntimimage_stored) 102 103 !Local variables------------------------------- 104 !scalars 105 integer :: ierr,ii,itime,itime_next,itime_prev 106 real(dp) :: volume 107 !arrays 108 real(dp) :: rprimd(3,3),rprimd_next(3,3),rprimd_prev(3,3),vel_cell(3,3) 109 real(dp),allocatable :: mpibuf(:),mpibuffer(:,:,:),mpibuffer_all(:,:,:) 110 real(dp),allocatable :: etotal(:),forces(:,:,:),stressin(:,:,:),vel(:,:,:) 111 real(dp),allocatable :: xred(:,:,:),xred_next(:,:,:),xred_prev(:,:,:) 112 113 ! ************************************************************************* 114 115 !############# Parallelism stuff 1 ####################### 116 117 !Parallelism over image: only one process per image of the cell 118 if (mpi_enreg%me_cell==0) then 119 120 itime=itimimage_eff 121 itime_prev=itime-1;if (itime_prev<1) itime_prev=ntimimage_stored 122 123 if (mpi_enreg%paral_img==0.or.mpi_enreg%me_img==0) then 124 ABI_ALLOCATE(xred,(3,natom,nimage_tot)) 125 ABI_ALLOCATE(xred_prev,(3,natom,nimage_tot)) 126 ABI_ALLOCATE(xred_next,(3,natom,nimage_tot)) 127 ABI_ALLOCATE(etotal,(nimage_tot)) 128 ABI_ALLOCATE(forces,(3,natom,nimage_tot)) 129 ABI_ALLOCATE(stressin,(3,3,nimage_tot)) 130 ABI_ALLOCATE(vel,(3,natom,nimage_tot)) 131 end if 132 133 ! Parallelism: Gather positions/forces/velocities/stresses/energy from all images 134 if (mpi_enreg%paral_img==1) then 135 ABI_ALLOCATE(mpibuffer,(12,natom+1,nimage)) 136 do ii=1,nimage 137 mpibuffer(1:3 ,1:natom,ii)=results_img(ii,itime)%xred(1:3,1:natom) 138 mpibuffer(4:6 ,1:natom,ii)=results_img(ii,itime_prev)%xred(1:3,1:natom) 139 mpibuffer(7:9 ,1:natom,ii)=results_img(ii,itime)%results_gs%fcart(1:3,1:natom) 140 mpibuffer(10:12,1:natom,ii)=results_img(ii,itime)%vel(1:3,1:natom) 141 mpibuffer(1:6 ,natom+1,ii)=results_img(ii,itime)%results_gs%strten(1:6) 142 mpibuffer(7:12 ,natom+1,ii)=zero 143 end do 144 if (mpi_enreg%me_img==0) then 145 ABI_ALLOCATE(mpibuffer_all,(12,natom+1,nimage_tot)) 146 end if 147 call gather_array_img(mpibuffer,mpibuffer_all,mpi_enreg,only_one_per_img=.true.,allgather=.false.) 148 ABI_DEALLOCATE(mpibuffer) 149 if (mpi_enreg%me_img==0) then 150 do ii=1,nimage_tot 151 xred (1:3,1:natom,ii)=mpibuffer_all(1:3 ,1:natom,ii) 152 xred_prev(1:3,1:natom,ii)=mpibuffer_all(4:6 ,1:natom,ii) 153 forces (1:3,1:natom,ii)=mpibuffer_all(7:9 ,1:natom,ii) 154 vel (1:3,1:natom,ii)=mpibuffer_all(10:12,1:natom,ii) 155 stressin (1,1,ii) =mpibuffer_all(1,natom+1,ii) 156 stressin (2,2,ii) =mpibuffer_all(2,natom+1,ii) 157 stressin (3,3,ii) =mpibuffer_all(3,natom+1,ii) 158 stressin (3,2,ii) =mpibuffer_all(4,natom+1,ii) 159 stressin (3,1,ii) =mpibuffer_all(5,natom+1,ii) 160 stressin (2,1,ii) =mpibuffer_all(6,natom+1,ii) 161 stressin (2,3,ii)=stressin (3,2,ii) 162 stressin (1,3,ii)=stressin (3,1,ii) 163 stressin (1,2,ii)=stressin (2,1,ii) 164 end do 165 ABI_DEALLOCATE(mpibuffer_all) 166 end if 167 ABI_ALLOCATE(mpibuf,(nimage)) 168 if (mpi_enreg%me_img/=0) then 169 ABI_ALLOCATE(etotal,(0)) 170 end if 171 do ii=1,nimage 172 mpibuf(ii)=results_img(ii,itime)%results_gs%etotal 173 end do 174 call xmpi_gather(mpibuf,nimage,etotal,nimage,0,mpi_enreg%comm_img,ierr) 175 ABI_DEALLOCATE(mpibuf) 176 if (mpi_enreg%me_img/=0) then 177 ABI_DEALLOCATE(etotal) 178 end if 179 180 ! No parallelism: simply copy positions/forces/velocities/stresses/energy 181 else 182 do ii=1,nimage 183 xred (:,:,ii)=results_img(ii,itime)%xred(:,:) 184 xred_prev(:,:,ii)=results_img(ii,itime_prev)%xred(:,:) 185 forces (:,:,ii)=results_img(ii,itime)%results_gs%fcart(:,:) 186 vel (:,:,ii)=results_img(ii,itime)%vel(:,:) 187 etotal ( ii)=results_img(ii,itime)%results_gs%etotal 188 stressin (1,1,ii)=results_img(ii,itime)%results_gs%strten(1) 189 stressin (2,2,ii)=results_img(ii,itime)%results_gs%strten(2) 190 stressin (3,3,ii)=results_img(ii,itime)%results_gs%strten(3) 191 stressin (3,2,ii)=results_img(ii,itime)%results_gs%strten(4) 192 stressin (3,1,ii)=results_img(ii,itime)%results_gs%strten(5) 193 stressin (2,1,ii)=results_img(ii,itime)%results_gs%strten(6) 194 stressin (2,3,ii)=stressin (3,2,ii) 195 stressin (1,3,ii)=stressin (3,1,ii) 196 stressin (1,2,ii)=stressin (2,1,ii) 197 end do 198 end if 199 200 ! Parallelism over image: only one process does the job 201 if (mpi_enreg%paral_img==0.or.mpi_enreg%me_img==0) then 202 203 ! ############# PIMD MD algorithm ######################### 204 205 ! Some useful quantities about the cells (common to all images) 206 ! Take acell and rprim from 1st image 207 call mkrdim(results_img(1,itime)%acell,results_img(1,itime)%rprim,rprimd) 208 call mkrdim(results_img(1,itime_prev)%acell,results_img(1,itime_prev)%rprim,rprimd_prev) 209 rprimd_next(:,:)=rprimd_prev(:,:) 210 vel_cell(:,:)=results_img(1,1)%vel_cell(:,:) 211 212 ! Compute the volume of the supercell 213 volume=rprimd(1,1)*(rprimd(2,2)*rprimd(3,3)-rprimd(3,2)*rprimd(2,3))+& 214 & rprimd(2,1)*(rprimd(3,2)*rprimd(1,3)-rprimd(1,2)*rprimd(3,3))+& 215 & rprimd(3,1)*(rprimd(1,2)*rprimd(2,3)-rprimd(2,2)*rprimd(1,3)) 216 volume=abs(volume) 217 218 select case(imgmov) 219 220 case(9,10) !Langevin 221 222 select case(pimd_param%optcell) 223 case(0) !NVT 224 call pimd_langevin_nvt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,& 225 & rprimd,stressin,nimage_tot,vel,volume,xred,xred_next,xred_prev) 226 case(2) !NPT 227 call pimd_langevin_npt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,& 228 & rprimd,rprimd_next,rprimd_prev,stressin,nimage_tot,vel,vel_cell,& 229 & volume,xred,xred_next,xred_prev) 230 end select 231 232 case(13) !Nose Hoover chains 233 234 select case(pimd_param%optcell) 235 case(0) !NVT 236 call pimd_nosehoover_nvt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,& 237 & rprimd,stressin,nimage_tot,vel,volume,xred,xred_next,xred_prev) 238 case(2) !NPT 239 call pimd_nosehoover_npt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,& 240 & rprimd,rprimd_next,rprimd_prev,stressin,nimage_tot,vel,vel_cell,& 241 & volume,xred,xred_next,xred_prev) 242 end select 243 244 end select 245 246 ! ############# Parallelism stuff 2 ######################## 247 248 end if ! mpi_enreg%me_img==0 249 250 ! Parallelism: dispatch results 251 ! The trick: use (9,natom) to store xred,xred_next,vel for all atoms 252 ! use (9 ) to store rprimd_next 253 ABI_ALLOCATE(mpibuffer,(9,natom+2,nimage)) 254 if (mpi_enreg%paral_img==1) then 255 if (mpi_enreg%me_img==0) then 256 ABI_ALLOCATE(mpibuffer_all,(9,natom+2,nimage_tot)) 257 do ii=1,nimage_tot 258 mpibuffer_all(1:3,1:natom,ii)=xred_next(1:3,1:natom,ii) 259 mpibuffer_all(4:6,1:natom,ii)=xred(1:3,1:natom,ii) 260 mpibuffer_all(7:9,1:natom,ii)=vel(1:3,1:natom,ii) 261 mpibuffer_all(1:3,natom+1,ii)=rprimd_next(1:3,1) 262 mpibuffer_all(4:6,natom+1,ii)=rprimd_next(1:3,2) 263 mpibuffer_all(7:9,natom+1,ii)=rprimd_next(1:3,3) 264 mpibuffer_all(1:3,natom+2,ii)=vel_cell(1:3,1) 265 mpibuffer_all(4:6,natom+2,ii)=vel_cell(1:3,2) 266 mpibuffer_all(7:9,natom+2,ii)=vel_cell(1:3,3) 267 end do 268 end if 269 call scatter_array_img(mpibuffer,mpibuffer_all,mpi_enreg) 270 if (mpi_enreg%me_img==0) then 271 ABI_DEALLOCATE(mpibuffer_all) 272 end if 273 else 274 do ii=1,nimage 275 mpibuffer(1:3,1:natom,ii)=xred_next(1:3,1:natom,ii) 276 mpibuffer(4:6,1:natom,ii)=xred(1:3,1:natom,ii) 277 mpibuffer(7:9,1:natom,ii)=vel(1:3,1:natom,ii) 278 mpibuffer(1:3,natom+1,ii)=rprimd_next(1:3,1) 279 mpibuffer(4:6,natom+1,ii)=rprimd_next(1:3,2) 280 mpibuffer(7:9,natom+1,ii)=rprimd_next(1:3,3) 281 mpibuffer(1:3,natom+2,ii)=vel_cell(1:3,1) 282 mpibuffer(4:6,natom+2,ii)=vel_cell(1:3,2) 283 mpibuffer(7:9,natom+2,ii)=vel_cell(1:3,3) 284 end do 285 end if 286 287 if (mpi_enreg%paral_img==0.or.mpi_enreg%me_img==0) then 288 ABI_DEALLOCATE(xred) 289 ABI_DEALLOCATE(xred_prev) 290 ABI_DEALLOCATE(xred_next) 291 ABI_DEALLOCATE(etotal) 292 ABI_DEALLOCATE(forces) 293 ABI_DEALLOCATE(stressin) 294 ABI_DEALLOCATE(vel) 295 end if 296 297 else 298 ABI_ALLOCATE(mpibuffer,(9,natom+2,nimage)) 299 300 end if ! mpi_enreg%me_cell==0 301 302 !Send results to all procs treating the same image 303 call xmpi_bcast(mpibuffer,0,mpi_enreg%comm_cell,ierr) 304 305 !Store results in final place 306 itime=itimimage_eff 307 itime_prev=itime-1;if (itime_prev<1) itime_prev=ntimimage_stored 308 itime_next=itime+1;if (itime_next>ntimimage_stored) itime_next=1 309 do ii=1,nimage 310 results_img(ii,itime_next)%xred(1:3,1:natom)=mpibuffer(1:3,1:natom,ii) 311 results_img(ii,itime)%xred(1:3,1:natom)=mpibuffer(4:6,1:natom,ii) 312 results_img(ii,itime_next)%vel(1:3,1:natom)=mpibuffer(7:9,1:natom,ii) 313 end do 314 if (pimd_param%optcell/=0) then 315 do ii=1,nimage 316 results_img(ii,itime)%acell(1:3)=results_img(ii,itime_prev)%acell(1:3) 317 results_img(ii,itime)%rprim(1:3,1:3)=results_img(ii,itime_prev)%rprim(1:3,1:3) 318 rprimd(1:3,1)=mpibuffer(1:3,natom+1,ii) 319 rprimd(1:3,2)=mpibuffer(4:6,natom+1,ii) 320 rprimd(1:3,3)=mpibuffer(7:9,natom+1,ii) 321 results_img(ii,itime_next)%vel_cell(1:3,1)=mpibuffer(1:3,natom+2,ii) 322 results_img(ii,itime_next)%vel_cell(1:3,2)=mpibuffer(4:6,natom+2,ii) 323 results_img(ii,itime_next)%vel_cell(1:3,3)=mpibuffer(7:9,natom+2,ii) 324 call mkradim(results_img(ii,itime_next)%acell,results_img(ii,itime_next)%rprim,rprimd) 325 end do 326 end if 327 ABI_DEALLOCATE(mpibuffer) 328 329 end subroutine predict_pimd