TABLE OF CONTENTS
ABINIT/m_predict_pimd [ Modules ]
NAME
m_predict_pimd
FUNCTION
COPYRIGHT
Copyright (C) 2010-2024 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 .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_predict_pimd 23 24 use defs_basis 25 use m_abicore 26 use m_pimd 27 use m_xmpi 28 use m_results_img 29 30 use defs_abitypes, only : MPI_type 31 use m_geometry, only : mkradim, mkrdim 32 use m_pimd_langevin, only : pimd_langevin_npt, pimd_langevin_nvt 33 use m_pimd_nosehoover, only : pimd_nosehoover_npt, pimd_nosehoover_nvt 34 35 implicit none 36 37 private
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.
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
SOURCE
96 subroutine predict_pimd(imgmov,itimimage,itimimage_eff,mpi_enreg,natom,nimage,nimage_tot,& 97 & ntimimage_stored,pimd_param,prtvolimg,results_img) 98 99 !Arguments ------------------------------------ 100 !scalars 101 integer,intent(in) :: imgmov,itimimage,itimimage_eff,natom 102 integer,intent(in) :: nimage,nimage_tot,ntimimage_stored,prtvolimg 103 type(MPI_type),intent(in) :: mpi_enreg 104 type(pimd_type),intent(inout) :: pimd_param 105 !arrays 106 type(results_img_type) :: results_img(nimage,ntimimage_stored) 107 108 !Local variables------------------------------- 109 !scalars 110 integer :: ierr,ii,itime,itime_next,itime_prev 111 real(dp) :: volume 112 !arrays 113 real(dp) :: rprimd(3,3),rprimd_next(3,3),rprimd_prev(3,3),vel_cell(3,3) 114 real(dp),allocatable :: mpibuf(:),mpibuffer(:,:,:),mpibuffer_all(:,:,:) 115 real(dp),allocatable :: etotal(:),forces(:,:,:),stressin(:,:,:),vel(:,:,:) 116 real(dp),allocatable :: xred(:,:,:),xred_next(:,:,:),xred_prev(:,:,:) 117 118 ! ************************************************************************* 119 120 !############# Parallelism stuff 1 ####################### 121 122 !Parallelism over image: only one process per image of the cell 123 if (mpi_enreg%me_cell==0) then 124 125 itime=itimimage_eff 126 itime_prev=itime-1;if (itime_prev<1) itime_prev=ntimimage_stored 127 128 if (mpi_enreg%paral_img==0.or.mpi_enreg%me_img==0) then 129 ABI_MALLOC(xred,(3,natom,nimage_tot)) 130 ABI_MALLOC(xred_prev,(3,natom,nimage_tot)) 131 ABI_MALLOC(xred_next,(3,natom,nimage_tot)) 132 ABI_MALLOC(etotal,(nimage_tot)) 133 ABI_MALLOC(forces,(3,natom,nimage_tot)) 134 ABI_MALLOC(stressin,(3,3,nimage_tot)) 135 ABI_MALLOC(vel,(3,natom,nimage_tot)) 136 end if 137 138 ! Parallelism: Gather positions/forces/velocities/stresses/energy from all images 139 if (mpi_enreg%paral_img==1) then 140 ABI_MALLOC(mpibuffer,(12,natom+1,nimage)) 141 do ii=1,nimage 142 mpibuffer(1:3 ,1:natom,ii)=results_img(ii,itime)%xred(1:3,1:natom) 143 mpibuffer(4:6 ,1:natom,ii)=results_img(ii,itime_prev)%xred(1:3,1:natom) 144 mpibuffer(7:9 ,1:natom,ii)=results_img(ii,itime)%results_gs%fcart(1:3,1:natom) 145 mpibuffer(10:12,1:natom,ii)=results_img(ii,itime)%vel(1:3,1:natom) 146 mpibuffer(1:6 ,natom+1,ii)=results_img(ii,itime)%results_gs%strten(1:6) 147 mpibuffer(7:12 ,natom+1,ii)=zero 148 end do 149 if (mpi_enreg%me_img==0) then 150 ABI_MALLOC(mpibuffer_all,(12,natom+1,nimage_tot)) 151 end if 152 call gather_array_img(mpibuffer,mpibuffer_all,mpi_enreg,only_one_per_img=.true.,allgather=.false.) 153 ABI_FREE(mpibuffer) 154 if (mpi_enreg%me_img==0) then 155 do ii=1,nimage_tot 156 xred (1:3,1:natom,ii)=mpibuffer_all(1:3 ,1:natom,ii) 157 xred_prev(1:3,1:natom,ii)=mpibuffer_all(4:6 ,1:natom,ii) 158 forces (1:3,1:natom,ii)=mpibuffer_all(7:9 ,1:natom,ii) 159 vel (1:3,1:natom,ii)=mpibuffer_all(10:12,1:natom,ii) 160 stressin (1,1,ii) =mpibuffer_all(1,natom+1,ii) 161 stressin (2,2,ii) =mpibuffer_all(2,natom+1,ii) 162 stressin (3,3,ii) =mpibuffer_all(3,natom+1,ii) 163 stressin (3,2,ii) =mpibuffer_all(4,natom+1,ii) 164 stressin (3,1,ii) =mpibuffer_all(5,natom+1,ii) 165 stressin (2,1,ii) =mpibuffer_all(6,natom+1,ii) 166 stressin (2,3,ii)=stressin (3,2,ii) 167 stressin (1,3,ii)=stressin (3,1,ii) 168 stressin (1,2,ii)=stressin (2,1,ii) 169 end do 170 ABI_FREE(mpibuffer_all) 171 end if 172 ABI_MALLOC(mpibuf,(nimage)) 173 if (mpi_enreg%me_img/=0) then 174 ABI_MALLOC(etotal,(0)) 175 end if 176 do ii=1,nimage 177 mpibuf(ii)=results_img(ii,itime)%results_gs%etotal 178 end do 179 call xmpi_gather(mpibuf,nimage,etotal,nimage,0,mpi_enreg%comm_img,ierr) 180 ABI_FREE(mpibuf) 181 if (mpi_enreg%me_img/=0) then 182 ABI_FREE(etotal) 183 end if 184 185 ! No parallelism: simply copy positions/forces/velocities/stresses/energy 186 else 187 do ii=1,nimage 188 xred (:,:,ii)=results_img(ii,itime)%xred(:,:) 189 xred_prev(:,:,ii)=results_img(ii,itime_prev)%xred(:,:) 190 forces (:,:,ii)=results_img(ii,itime)%results_gs%fcart(:,:) 191 vel (:,:,ii)=results_img(ii,itime)%vel(:,:) 192 etotal ( ii)=results_img(ii,itime)%results_gs%etotal 193 stressin (1,1,ii)=results_img(ii,itime)%results_gs%strten(1) 194 stressin (2,2,ii)=results_img(ii,itime)%results_gs%strten(2) 195 stressin (3,3,ii)=results_img(ii,itime)%results_gs%strten(3) 196 stressin (3,2,ii)=results_img(ii,itime)%results_gs%strten(4) 197 stressin (3,1,ii)=results_img(ii,itime)%results_gs%strten(5) 198 stressin (2,1,ii)=results_img(ii,itime)%results_gs%strten(6) 199 stressin (2,3,ii)=stressin (3,2,ii) 200 stressin (1,3,ii)=stressin (3,1,ii) 201 stressin (1,2,ii)=stressin (2,1,ii) 202 end do 203 end if 204 205 ! Parallelism over image: only one process does the job 206 if (mpi_enreg%paral_img==0.or.mpi_enreg%me_img==0) then 207 208 ! ############# PIMD MD algorithm ######################### 209 210 ! Some useful quantities about the cells (common to all images) 211 ! Take acell and rprim from 1st image 212 call mkrdim(results_img(1,itime)%acell,results_img(1,itime)%rprim,rprimd) 213 call mkrdim(results_img(1,itime_prev)%acell,results_img(1,itime_prev)%rprim,rprimd_prev) 214 vel_cell(:,:)=results_img(1,itime)%vel_cell(:,:) 215 216 ! Compute the volume of the supercell 217 volume=rprimd(1,1)*(rprimd(2,2)*rprimd(3,3)-rprimd(3,2)*rprimd(2,3))+& 218 & rprimd(2,1)*(rprimd(3,2)*rprimd(1,3)-rprimd(1,2)*rprimd(3,3))+& 219 & rprimd(3,1)*(rprimd(1,2)*rprimd(2,3)-rprimd(2,2)*rprimd(1,3)) 220 volume=abs(volume) 221 222 select case(imgmov) 223 224 case(9,10) !Langevin 225 226 select case(pimd_param%optcell) 227 case(0) !NVT 228 call pimd_langevin_nvt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,& 229 & rprimd,stressin,nimage_tot,vel,volume,xred,xred_next,xred_prev) 230 case(2) !NPT 231 call pimd_langevin_npt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,& 232 & rprimd,rprimd_next,rprimd_prev,stressin,nimage_tot,vel,vel_cell,& 233 & volume,xred,xred_next,xred_prev) 234 end select 235 236 case(13) !Nose Hoover chains 237 238 select case(pimd_param%optcell) 239 case(0) !NVT 240 call pimd_nosehoover_nvt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,& 241 & rprimd,stressin,nimage_tot,vel,volume,xred,xred_next,xred_prev) 242 case(2) !NPT 243 call pimd_nosehoover_npt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,& 244 & rprimd,rprimd_next,rprimd_prev,stressin,nimage_tot,vel,vel_cell,& 245 & volume,xred,xred_next,xred_prev) 246 end select 247 248 end select 249 250 ! ############# Parallelism stuff 2 ######################## 251 252 end if ! mpi_enreg%me_img==0 253 254 ! Parallelism: dispatch results 255 ! The trick: use (9,natom) to store xred,xred_next,vel for all atoms 256 ! use (9 ) to store rprimd_next 257 ABI_MALLOC(mpibuffer,(9,natom+3,nimage)) 258 if (mpi_enreg%paral_img==1) then 259 if (mpi_enreg%me_img==0) then 260 ABI_MALLOC(mpibuffer_all,(9,natom+3,nimage_tot)) 261 do ii=1,nimage_tot 262 mpibuffer_all(1:3,1:natom,ii)=xred_next(1:3,1:natom,ii) 263 mpibuffer_all(4:6,1:natom,ii)=xred(1:3,1:natom,ii) 264 mpibuffer_all(7:9,1:natom,ii)=vel(1:3,1:natom,ii) 265 mpibuffer_all(1:3,natom+1,ii)=rprimd_next(1:3,1) 266 mpibuffer_all(4:6,natom+1,ii)=rprimd_next(1:3,2) 267 mpibuffer_all(7:9,natom+1,ii)=rprimd_next(1:3,3) 268 mpibuffer_all(1:3,natom+2,ii)=vel_cell(1:3,1) 269 mpibuffer_all(4:6,natom+2,ii)=vel_cell(1:3,2) 270 mpibuffer_all(7:9,natom+2,ii)=vel_cell(1:3,3) 271 mpibuffer_all(1:3,natom+3,ii)=rprimd(1:3,1) 272 mpibuffer_all(4:6,natom+3,ii)=rprimd(1:3,2) 273 mpibuffer_all(7:9,natom+3,ii)=rprimd(1:3,3) 274 end do 275 end if 276 call scatter_array_img(mpibuffer,mpibuffer_all,mpi_enreg) 277 if (mpi_enreg%me_img==0) then 278 ABI_FREE(mpibuffer_all) 279 end if 280 else 281 do ii=1,nimage 282 mpibuffer(1:3,1:natom,ii)=xred_next(1:3,1:natom,ii) 283 mpibuffer(4:6,1:natom,ii)=xred(1:3,1:natom,ii) 284 mpibuffer(7:9,1:natom,ii)=vel(1:3,1:natom,ii) 285 mpibuffer(1:3,natom+1,ii)=rprimd_next(1:3,1) 286 mpibuffer(4:6,natom+1,ii)=rprimd_next(1:3,2) 287 mpibuffer(7:9,natom+1,ii)=rprimd_next(1:3,3) 288 mpibuffer(1:3,natom+2,ii)=vel_cell(1:3,1) 289 mpibuffer(4:6,natom+2,ii)=vel_cell(1:3,2) 290 mpibuffer(7:9,natom+2,ii)=vel_cell(1:3,3) 291 mpibuffer(1:3,natom+3,ii)=rprimd(1:3,1) 292 mpibuffer(4:6,natom+3,ii)=rprimd(1:3,2) 293 mpibuffer(7:9,natom+3,ii)=rprimd(1:3,3) 294 end do 295 end if 296 297 if (mpi_enreg%paral_img==0.or.mpi_enreg%me_img==0) then 298 ABI_FREE(xred) 299 ABI_FREE(xred_prev) 300 ABI_FREE(xred_next) 301 ABI_FREE(etotal) 302 ABI_FREE(forces) 303 ABI_FREE(stressin) 304 ABI_FREE(vel) 305 end if 306 307 else 308 ABI_MALLOC(mpibuffer,(9,natom+3,nimage)) 309 310 end if ! mpi_enreg%me_cell==0 311 312 !Send results to all procs treating the same image 313 call xmpi_bcast(mpibuffer,0,mpi_enreg%comm_cell,ierr) 314 315 !Store results in final place 316 itime=itimimage_eff 317 itime_prev=itime-1;if (itime_prev<1) itime_prev=ntimimage_stored 318 itime_next=itime+1;if (itime_next>ntimimage_stored) itime_next=1 319 do ii=1,nimage 320 results_img(ii,itime_next)%xred(1:3,1:natom)=mpibuffer(1:3,1:natom,ii) 321 results_img(ii,itime)%xred(1:3,1:natom)=mpibuffer(4:6,1:natom,ii) 322 results_img(ii,itime_next)%vel(1:3,1:natom)=mpibuffer(7:9,1:natom,ii) 323 end do 324 if (pimd_param%optcell/=0) then 325 do ii=1,nimage 326 rprimd(1:3,1)=mpibuffer(1:3,natom+1,ii) 327 rprimd(1:3,2)=mpibuffer(4:6,natom+1,ii) 328 rprimd(1:3,3)=mpibuffer(7:9,natom+1,ii) 329 call mkradim(results_img(ii,itime_next)%acell,results_img(ii,itime_next)%rprim,rprimd) 330 rprimd_prev(1:3,1)=mpibuffer(1:3,natom+3,ii) 331 rprimd_prev(1:3,2)=mpibuffer(4:6,natom+3,ii) 332 rprimd_prev(1:3,3)=mpibuffer(7:9,natom+3,ii) 333 call mkradim(results_img(ii,itime)%acell,results_img(ii,itime)%rprim,rprimd_prev) 334 results_img(ii,itime_next)%vel_cell(1:3,1)=mpibuffer(1:3,natom+2,ii) 335 results_img(ii,itime_next)%vel_cell(1:3,2)=mpibuffer(4:6,natom+2,ii) 336 results_img(ii,itime_next)%vel_cell(1:3,3)=mpibuffer(7:9,natom+2,ii) 337 end do 338 end if 339 ABI_FREE(mpibuffer) 340 341 end subroutine predict_pimd