TABLE OF CONTENTS


ABINIT/m_predict_pimd [ Modules ]

[ Top ] [ Modules ]

NAME

  m_predict_pimd

FUNCTION

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 .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_predict_pimd
28 
29  use defs_basis
30  use defs_abitypes
31  use m_abicore
32  use m_pimd
33  use m_xmpi
34  use m_results_img
35 
36  use m_geometry,       only : mkradim, mkrdim
37  use m_pimd_langevin,  only : pimd_langevin_npt, pimd_langevin_nvt
38  use m_pimd_nosehoover, only : pimd_nosehoover_npt, pimd_nosehoover_nvt
39 
40  implicit none
41 
42  private

ABINIT/predict_pimd [ Functions ]

[ Top ] [ 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

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

109 subroutine predict_pimd(imgmov,itimimage,itimimage_eff,mpi_enreg,natom,nimage,nimage_tot,&
110 &                       ntimimage_stored,pimd_param,prtvolimg,results_img)
111 
112 
113 !This section has been created automatically by the script Abilint (TD).
114 !Do not modify the following lines by hand.
115 #undef ABI_FUNC
116 #define ABI_FUNC 'predict_pimd'
117 !End of the abilint section
118 
119  implicit none
120 
121 !Arguments ------------------------------------
122 !scalars
123  integer,intent(in) :: imgmov,itimimage,itimimage_eff,natom
124  integer,intent(in) :: nimage,nimage_tot,ntimimage_stored,prtvolimg
125  type(MPI_type),intent(in) :: mpi_enreg
126  type(pimd_type),intent(inout) :: pimd_param
127 !arrays
128  type(results_img_type) :: results_img(nimage,ntimimage_stored)
129 
130 !Local variables-------------------------------
131 !scalars
132  integer :: ierr,ii,itime,itime_next,itime_prev
133  real(dp) :: volume
134 !arrays
135  real(dp) :: rprimd(3,3),rprimd_next(3,3),rprimd_prev(3,3),vel_cell(3,3)
136  real(dp),allocatable :: mpibuf(:),mpibuffer(:,:,:),mpibuffer_all(:,:,:)
137  real(dp),allocatable :: etotal(:),forces(:,:,:),stressin(:,:,:),vel(:,:,:)
138  real(dp),allocatable :: xred(:,:,:),xred_next(:,:,:),xred_prev(:,:,:)
139 
140 ! *************************************************************************
141 
142 !############# Parallelism stuff 1 #######################
143 
144 !Parallelism over image: only one process per image of the cell
145  if (mpi_enreg%me_cell==0) then
146 
147    itime=itimimage_eff
148    itime_prev=itime-1;if (itime_prev<1) itime_prev=ntimimage_stored
149 
150    if (mpi_enreg%paral_img==0.or.mpi_enreg%me_img==0) then
151      ABI_ALLOCATE(xred,(3,natom,nimage_tot))
152      ABI_ALLOCATE(xred_prev,(3,natom,nimage_tot))
153      ABI_ALLOCATE(xred_next,(3,natom,nimage_tot))
154      ABI_ALLOCATE(etotal,(nimage_tot))
155      ABI_ALLOCATE(forces,(3,natom,nimage_tot))
156      ABI_ALLOCATE(stressin,(3,3,nimage_tot))
157      ABI_ALLOCATE(vel,(3,natom,nimage_tot))
158    end if
159 
160 !  Parallelism: Gather positions/forces/velocities/stresses/energy from all images
161    if (mpi_enreg%paral_img==1) then
162      ABI_ALLOCATE(mpibuffer,(12,natom+1,nimage))
163      do ii=1,nimage
164        mpibuffer(1:3  ,1:natom,ii)=results_img(ii,itime)%xred(1:3,1:natom)
165        mpibuffer(4:6  ,1:natom,ii)=results_img(ii,itime_prev)%xred(1:3,1:natom)
166        mpibuffer(7:9  ,1:natom,ii)=results_img(ii,itime)%results_gs%fcart(1:3,1:natom)
167        mpibuffer(10:12,1:natom,ii)=results_img(ii,itime)%vel(1:3,1:natom)
168        mpibuffer(1:6  ,natom+1,ii)=results_img(ii,itime)%results_gs%strten(1:6)
169        mpibuffer(7:12 ,natom+1,ii)=zero
170      end do
171      if (mpi_enreg%me_img==0)  then
172        ABI_ALLOCATE(mpibuffer_all,(12,natom+1,nimage_tot))
173      end if
174      call gather_array_img(mpibuffer,mpibuffer_all,mpi_enreg,only_one_per_img=.true.,allgather=.false.)
175      ABI_DEALLOCATE(mpibuffer)
176      if (mpi_enreg%me_img==0) then
177        do ii=1,nimage_tot
178          xred     (1:3,1:natom,ii)=mpibuffer_all(1:3  ,1:natom,ii)
179          xred_prev(1:3,1:natom,ii)=mpibuffer_all(4:6  ,1:natom,ii)
180          forces   (1:3,1:natom,ii)=mpibuffer_all(7:9  ,1:natom,ii)
181          vel      (1:3,1:natom,ii)=mpibuffer_all(10:12,1:natom,ii)
182          stressin (1,1,ii)        =mpibuffer_all(1,natom+1,ii)
183          stressin (2,2,ii)        =mpibuffer_all(2,natom+1,ii)
184          stressin (3,3,ii)        =mpibuffer_all(3,natom+1,ii)
185          stressin (3,2,ii)        =mpibuffer_all(4,natom+1,ii)
186          stressin (3,1,ii)        =mpibuffer_all(5,natom+1,ii)
187          stressin (2,1,ii)        =mpibuffer_all(6,natom+1,ii)
188          stressin (2,3,ii)=stressin (3,2,ii)
189          stressin (1,3,ii)=stressin (3,1,ii)
190          stressin (1,2,ii)=stressin (2,1,ii)
191        end do
192        ABI_DEALLOCATE(mpibuffer_all)
193      end if
194      ABI_ALLOCATE(mpibuf,(nimage))
195      if (mpi_enreg%me_img/=0) then
196        ABI_ALLOCATE(etotal,(0))
197      end if
198      do ii=1,nimage
199        mpibuf(ii)=results_img(ii,itime)%results_gs%etotal
200      end do
201      call xmpi_gather(mpibuf,nimage,etotal,nimage,0,mpi_enreg%comm_img,ierr)
202      ABI_DEALLOCATE(mpibuf)
203      if (mpi_enreg%me_img/=0) then
204        ABI_DEALLOCATE(etotal)
205      end if
206 
207 !    No parallelism: simply copy positions/forces/velocities/stresses/energy
208    else
209      do ii=1,nimage
210        xred     (:,:,ii)=results_img(ii,itime)%xred(:,:)
211        xred_prev(:,:,ii)=results_img(ii,itime_prev)%xred(:,:)
212        forces   (:,:,ii)=results_img(ii,itime)%results_gs%fcart(:,:)
213        vel      (:,:,ii)=results_img(ii,itime)%vel(:,:)
214        etotal   (    ii)=results_img(ii,itime)%results_gs%etotal
215        stressin (1,1,ii)=results_img(ii,itime)%results_gs%strten(1)
216        stressin (2,2,ii)=results_img(ii,itime)%results_gs%strten(2)
217        stressin (3,3,ii)=results_img(ii,itime)%results_gs%strten(3)
218        stressin (3,2,ii)=results_img(ii,itime)%results_gs%strten(4)
219        stressin (3,1,ii)=results_img(ii,itime)%results_gs%strten(5)
220        stressin (2,1,ii)=results_img(ii,itime)%results_gs%strten(6)
221        stressin (2,3,ii)=stressin (3,2,ii)
222        stressin (1,3,ii)=stressin (3,1,ii)
223        stressin (1,2,ii)=stressin (2,1,ii)
224      end do
225    end if
226 
227 !  Parallelism over image: only one process does the job
228    if (mpi_enreg%paral_img==0.or.mpi_enreg%me_img==0) then
229 
230 !    ############# PIMD MD algorithm #########################
231 
232 !    Some useful quantities about the cells (common to all images)
233 !    Take acell and rprim from 1st image
234      call mkrdim(results_img(1,itime)%acell,results_img(1,itime)%rprim,rprimd)
235      call mkrdim(results_img(1,itime_prev)%acell,results_img(1,itime_prev)%rprim,rprimd_prev)
236      rprimd_next(:,:)=rprimd_prev(:,:)
237      vel_cell(:,:)=results_img(1,1)%vel_cell(:,:)
238 
239 !    Compute the volume of the supercell
240      volume=rprimd(1,1)*(rprimd(2,2)*rprimd(3,3)-rprimd(3,2)*rprimd(2,3))+&
241 &     rprimd(2,1)*(rprimd(3,2)*rprimd(1,3)-rprimd(1,2)*rprimd(3,3))+&
242 &     rprimd(3,1)*(rprimd(1,2)*rprimd(2,3)-rprimd(2,2)*rprimd(1,3))
243      volume=abs(volume)
244 
245      select case(imgmov)
246 
247      case(9,10)  !Langevin
248 
249        select case(pimd_param%optcell)
250        case(0)  !NVT
251          call pimd_langevin_nvt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,&
252 &         rprimd,stressin,nimage_tot,vel,volume,xred,xred_next,xred_prev)
253        case(2)  !NPT
254          call pimd_langevin_npt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,&
255 &         rprimd,rprimd_next,rprimd_prev,stressin,nimage_tot,vel,vel_cell,&
256 &         volume,xred,xred_next,xred_prev)
257        end select
258 
259      case(13)  !Nose Hoover chains
260 
261        select case(pimd_param%optcell)
262        case(0)  !NVT
263          call pimd_nosehoover_nvt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,&
264 &         rprimd,stressin,nimage_tot,vel,volume,xred,xred_next,xred_prev)
265        case(2)  !NPT
266          call pimd_nosehoover_npt(etotal,forces,itimimage,natom,pimd_param,prtvolimg,&
267 &         rprimd,rprimd_next,rprimd_prev,stressin,nimage_tot,vel,vel_cell,&
268 &         volume,xred,xred_next,xred_prev)
269        end select
270 
271      end select
272 
273 !    ############# Parallelism stuff 2 ########################
274 
275    end if ! mpi_enreg%me_img==0
276 
277 !  Parallelism: dispatch results
278 !  The trick: use (9,natom) to store xred,xred_next,vel for all atoms
279 !  use (9      ) to store rprimd_next
280    ABI_ALLOCATE(mpibuffer,(9,natom+2,nimage))
281    if (mpi_enreg%paral_img==1) then
282      if (mpi_enreg%me_img==0) then
283        ABI_ALLOCATE(mpibuffer_all,(9,natom+2,nimage_tot))
284        do ii=1,nimage_tot
285          mpibuffer_all(1:3,1:natom,ii)=xred_next(1:3,1:natom,ii)
286          mpibuffer_all(4:6,1:natom,ii)=xred(1:3,1:natom,ii)
287          mpibuffer_all(7:9,1:natom,ii)=vel(1:3,1:natom,ii)
288          mpibuffer_all(1:3,natom+1,ii)=rprimd_next(1:3,1)
289          mpibuffer_all(4:6,natom+1,ii)=rprimd_next(1:3,2)
290          mpibuffer_all(7:9,natom+1,ii)=rprimd_next(1:3,3)
291          mpibuffer_all(1:3,natom+2,ii)=vel_cell(1:3,1)
292          mpibuffer_all(4:6,natom+2,ii)=vel_cell(1:3,2)
293          mpibuffer_all(7:9,natom+2,ii)=vel_cell(1:3,3)
294        end do
295      end if
296      call scatter_array_img(mpibuffer,mpibuffer_all,mpi_enreg)
297      if (mpi_enreg%me_img==0)  then
298        ABI_DEALLOCATE(mpibuffer_all)
299      end if
300    else
301      do ii=1,nimage
302        mpibuffer(1:3,1:natom,ii)=xred_next(1:3,1:natom,ii)
303        mpibuffer(4:6,1:natom,ii)=xred(1:3,1:natom,ii)
304        mpibuffer(7:9,1:natom,ii)=vel(1:3,1:natom,ii)
305        mpibuffer(1:3,natom+1,ii)=rprimd_next(1:3,1)
306        mpibuffer(4:6,natom+1,ii)=rprimd_next(1:3,2)
307        mpibuffer(7:9,natom+1,ii)=rprimd_next(1:3,3)
308        mpibuffer(1:3,natom+2,ii)=vel_cell(1:3,1)
309        mpibuffer(4:6,natom+2,ii)=vel_cell(1:3,2)
310        mpibuffer(7:9,natom+2,ii)=vel_cell(1:3,3)
311      end do
312    end if
313 
314    if (mpi_enreg%paral_img==0.or.mpi_enreg%me_img==0) then
315      ABI_DEALLOCATE(xred)
316      ABI_DEALLOCATE(xred_prev)
317      ABI_DEALLOCATE(xred_next)
318      ABI_DEALLOCATE(etotal)
319      ABI_DEALLOCATE(forces)
320      ABI_DEALLOCATE(stressin)
321      ABI_DEALLOCATE(vel)
322    end if
323 
324  else
325    ABI_ALLOCATE(mpibuffer,(9,natom+2,nimage))
326 
327  end if ! mpi_enreg%me_cell==0
328 
329 !Send results to all procs treating the same image
330  call xmpi_bcast(mpibuffer,0,mpi_enreg%comm_cell,ierr)
331 
332 !Store results in final place
333  itime=itimimage_eff
334  itime_prev=itime-1;if (itime_prev<1) itime_prev=ntimimage_stored
335  itime_next=itime+1;if (itime_next>ntimimage_stored) itime_next=1
336  do ii=1,nimage
337    results_img(ii,itime_next)%xred(1:3,1:natom)=mpibuffer(1:3,1:natom,ii)
338    results_img(ii,itime)%xred(1:3,1:natom)=mpibuffer(4:6,1:natom,ii)
339    results_img(ii,itime_next)%vel(1:3,1:natom)=mpibuffer(7:9,1:natom,ii)
340  end do
341  if (pimd_param%optcell/=0) then
342    do ii=1,nimage
343      results_img(ii,itime)%acell(1:3)=results_img(ii,itime_prev)%acell(1:3)
344      results_img(ii,itime)%rprim(1:3,1:3)=results_img(ii,itime_prev)%rprim(1:3,1:3)
345      rprimd(1:3,1)=mpibuffer(1:3,natom+1,ii)
346      rprimd(1:3,2)=mpibuffer(4:6,natom+1,ii)
347      rprimd(1:3,3)=mpibuffer(7:9,natom+1,ii)
348      results_img(ii,itime_next)%vel_cell(1:3,1)=mpibuffer(1:3,natom+2,ii)
349      results_img(ii,itime_next)%vel_cell(1:3,2)=mpibuffer(4:6,natom+2,ii)
350      results_img(ii,itime_next)%vel_cell(1:3,3)=mpibuffer(7:9,natom+2,ii)
351      call mkradim(results_img(ii,itime_next)%acell,results_img(ii,itime_next)%rprim,rprimd)
352    end do
353  end if
354  ABI_DEALLOCATE(mpibuffer)
355 
356 end subroutine predict_pimd