TABLE OF CONTENTS


ABINIT/m_predict_string [ Modules ]

[ Top ] [ Modules ]

NAME

  m_predict_string

FUNCTION

COPYRIGHT

  Copyright (C) 2009-2018 ABINIT group (XG,ARom,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 .

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_string
28 
29  use defs_basis
30  use defs_abitypes
31  use m_abicore
32  use m_splines
33  use m_mep
34  use m_errors
35  use m_xmpi
36 
37  use m_results_img, only : results_img_type, gather_array_img, get_geometry_img
38 
39  implicit none
40 
41  private

ABINIT/predict_string [ Functions ]

[ Top ] [ Functions ]

NAME

 predict_string

FUNCTION

 Given the past history of images, predict the new set of images using String Method.
 The changes on the geometry and others are predicted by rescaling the path
 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

PARENTS

      predictimg

CHILDREN

      gather_array_img,get_geometry_img,mep_rk4,mep_steepest,spline,splint
      xmpi_bcast

SOURCE

104 subroutine predict_string(itimimage,itimimage_eff,list_dynimage,mep_param,mpi_enreg,natom,&
105 &                         ndynimage,nimage,nimage_tot,ntimimage_stored,results_img)
106 
107 
108 !This section has been created automatically by the script Abilint (TD).
109 !Do not modify the following lines by hand.
110 #undef ABI_FUNC
111 #define ABI_FUNC 'predict_string'
112 !End of the abilint section
113 
114  implicit none
115 
116 !Arguments ------------------------------------
117 !scalars
118  integer,intent(in) :: itimimage,itimimage_eff,natom,ndynimage
119  integer,intent(in) :: nimage,nimage_tot,ntimimage_stored
120  type(mep_type),intent(inout) :: mep_param
121  type(MPI_type),intent(in) :: mpi_enreg
122 !arrays
123  integer,intent(in)     :: list_dynimage(ndynimage)
124  type(results_img_type) :: results_img(nimage,ntimimage_stored)
125 
126 !Local variables-------------------------------
127 !scalars
128  integer :: idynimage,ierr,ii,iimage,iatom,next_itimimage
129  real(dp) :: emax,emin,step
130 !arrays
131  real(dp),allocatable :: buffer(:,:),buffer_all(:,:)
132  real(dp),allocatable :: darc(:),dimage(:),fcart(:,:,:),rprimd(:,:,:),vect(:,:),wimage(:)
133  real(dp),allocatable :: x(:),y(:),z(:),x2(:),y2(:),z2(:)
134  real(dp),allocatable :: xout(:),yout(:),zout(:)
135  real(dp),allocatable,target :: etotal(:),xcart(:,:,:),xred(:,:,:)
136  real(dp),pointer :: etotal_all(:),xcart_all(:,:,:),xred_all(:,:,:)
137 
138 ! *************************************************************************
139 
140  ABI_ALLOCATE(xred,(3,natom,nimage))
141 
142 !Parallelism over images: only one process per image of the cell
143  if (mpi_enreg%me_cell==0) then
144 
145 !  Retrieve positions and forces
146    ABI_ALLOCATE(etotal,(nimage))
147    ABI_ALLOCATE(xcart,(3,natom,nimage))
148    ABI_ALLOCATE(fcart,(3,natom,nimage))
149    ABI_ALLOCATE(rprimd,(3,3,nimage))
150    call get_geometry_img(etotal,natom,nimage,results_img(:,itimimage_eff),&
151 &   fcart,rprimd,xcart,xred)
152 
153 !  EVOLUTION STEP
154 !  ===============================================
155 
156 !  Compute new atomic positions in each cell
157    if      (mep_param%mep_solver==0) then ! Steepest-descent
158      call mep_steepest(fcart,list_dynimage,mep_param,natom,ndynimage,nimage,rprimd,xcart,xred)
159    else if (mep_param%mep_solver==4) then ! 4th-order Runge-Kutta
160      call mep_rk4(fcart,itimimage,list_dynimage,mep_param,natom,ndynimage,nimage,rprimd,xcart,xred)
161    else
162      MSG_BUG("Inconsistent solver !")
163    end if
164 
165 !  REPARAMETRIZATION STEP
166 !  ===============================================
167 
168 !  No reparametrization step in case of Runge-Kutta and mod(istep,4)>0
169    if (mep_param%mep_solver/=4.or.mod(itimimage,4)==0) then
170 
171 !    Parallelism: gather data of all images
172      if (mpi_enreg%paral_img==1) then
173        ABI_ALLOCATE(buffer,(6*natom+1,nimage))
174        ABI_ALLOCATE(buffer_all,(6*natom+1,nimage_tot))
175        ABI_ALLOCATE(xred_all,(3,natom,nimage_tot))
176        ABI_ALLOCATE(xcart_all,(3,natom,nimage_tot))
177        ABI_ALLOCATE(etotal_all,(nimage_tot))
178        buffer=zero;ii=0
179        buffer(ii+1:ii+3*natom,1:nimage)=reshape(xred ,(/3*natom,nimage/));ii=3*natom
180        buffer(ii+1:ii+3*natom,1:nimage)=reshape(xcart,(/3*natom,nimage/));ii=6*natom
181        buffer(ii+1           ,1:nimage)=etotal(1:nimage);ii=0
182        call gather_array_img(buffer,buffer_all,mpi_enreg,allgather=.true.)
183        xred_all(:,:,:) =reshape(buffer_all(ii+1:ii+3*natom,1:nimage_tot),(/3,natom,nimage_tot/));ii=3*natom
184        xcart_all(:,:,:)=reshape(buffer_all(ii+1:ii+3*natom,1:nimage_tot),(/3,natom,nimage_tot/));ii=6*natom
185        etotal_all(:)=buffer_all(ii+1,1:nimage_tot);ii=0
186        ABI_DEALLOCATE(buffer)
187        ABI_DEALLOCATE(buffer_all)
188      else
189        xred_all   => xred
190        xcart_all  => xcart
191        etotal_all => etotal
192      end if
193 
194 !    dimage is the distance between two images
195 !    darc is the parametrization on the string
196 !    wimage is the weight
197      ABI_ALLOCATE(darc,(nimage_tot))
198      ABI_ALLOCATE(dimage,(nimage_tot))
199      ABI_ALLOCATE(wimage,(nimage_tot))
200 
201 !    === Weights for equal arc length
202      if (mep_param%string_algo/=2) then
203        wimage(:)=one
204 
205 !      === Weights for energy-weight arc length
206      else
207        emin=min(etotal_all(1),etotal_all(nimage_tot))
208        emax=maxval(etotal_all(:))
209        wimage(1)=one
210        do iimage=2,nimage_tot
211          wimage(iimage)=exp((half*(etotal_all(iimage)+etotal_all(iimage-1))-emin)/(emax-emin))
212        end do
213      end if
214 
215 !    The distance between images is calculated
216 !    and normalized to a string length of 1.0
217      dimage=zero
218      ABI_ALLOCATE(vect,(3,natom))
219      do iimage=2,nimage_tot
220        dimage(iimage)=dimage(iimage-1)
221 !      MT april 2012: distance must be computed with cartesian coordinates
222 !      vect(:,:)=xred_all(:,:,iimage)-xred_all(:,:,iimage-1)
223        vect(:,:)=xcart_all(:,:,iimage)-xcart_all(:,:,iimage-1)
224        dimage(iimage)=dimage(iimage)+wimage(iimage)*mep_img_norm(vect)
225      end do
226      dimage(:)=dimage(:)/dimage(nimage_tot)
227      ABI_DEALLOCATE(vect)
228 
229 !    Arc lengths
230      darc(1)=zero
231      step=one/dble(nimage_tot-1)
232      do iimage=2,nimage_tot
233        darc(iimage)=darc(iimage-1)+step
234      end do
235 
236 !    New image coordinates are calculated and such that now the mesh is uniform
237      ABI_ALLOCATE(x,(nimage_tot))
238      ABI_ALLOCATE(y,(nimage_tot))
239      ABI_ALLOCATE(z,(nimage_tot))
240      ABI_ALLOCATE(x2,(nimage_tot))
241      ABI_ALLOCATE(y2,(nimage_tot))
242      ABI_ALLOCATE(z2,(nimage_tot))
243      ABI_ALLOCATE(xout,(nimage_tot))
244      ABI_ALLOCATE(yout,(nimage_tot))
245      ABI_ALLOCATE(zout,(nimage_tot))
246      do iatom=1,natom
247        do iimage=1,nimage_tot
248          x(iimage)=xred_all(1,iatom,iimage)
249          y(iimage)=xred_all(2,iatom,iimage)
250          z(iimage)=xred_all(3,iatom,iimage)
251        end do
252        call spline(dimage,x,nimage_tot,greatest_real,greatest_real,x2)
253        call spline(dimage,y,nimage_tot,greatest_real,greatest_real,y2)
254        call spline(dimage,z,nimage_tot,greatest_real,greatest_real,z2)
255        call splint(nimage_tot,dimage,x,x2,nimage_tot,darc,xout)
256        call splint(nimage_tot,dimage,y,y2,nimage_tot,darc,yout)
257        call splint(nimage_tot,dimage,z,z2,nimage_tot,darc,zout)
258 !      After a spline, new image coordinate for that particular
259 !      atom are generated only if they are dynamical
260        do idynimage=1,ndynimage
261          iimage=list_dynimage(idynimage)
262          ii=mpi_enreg%my_imgtab(iimage)
263          xred(1,iatom,iimage)=xout(ii)
264          xred(2,iatom,iimage)=yout(ii)
265          xred(3,iatom,iimage)=zout(ii)
266        end do
267      end do  ! iatom
268 
269 !    Free memory
270      ABI_DEALLOCATE(x)
271      ABI_DEALLOCATE(y)
272      ABI_DEALLOCATE(z)
273      ABI_DEALLOCATE(x2)
274      ABI_DEALLOCATE(y2)
275      ABI_DEALLOCATE(z2)
276      ABI_DEALLOCATE(xout)
277      ABI_DEALLOCATE(yout)
278      ABI_DEALLOCATE(zout)
279      ABI_DEALLOCATE(darc)
280      ABI_DEALLOCATE(dimage)
281      ABI_DEALLOCATE(wimage)
282      if (mpi_enreg%paral_img==1)  then
283        ABI_DEALLOCATE(xred_all)
284        ABI_DEALLOCATE(xcart_all)
285        ABI_DEALLOCATE(etotal_all)
286      end if
287 
288    end if ! Reparametrization
289 
290 !  ===============================================
291 
292    ABI_DEALLOCATE(etotal)
293    ABI_DEALLOCATE(xcart)
294    ABI_DEALLOCATE(fcart)
295    ABI_DEALLOCATE(rprimd)
296  end if ! mpi_enreg%me_cell==0
297 
298 !Store acell, rprim, xred and vel for the new iteration
299  call xmpi_bcast(xred,0,mpi_enreg%comm_cell,ierr)
300  next_itimimage=itimimage_eff+1
301  if (next_itimimage>ntimimage_stored) next_itimimage=1
302  do iimage=1,nimage
303    results_img(iimage,next_itimimage)%xred(:,:)    =xred(:,:,iimage)
304    results_img(iimage,next_itimimage)%acell(:)     =results_img(iimage,itimimage_eff)%acell(:)
305    results_img(iimage,next_itimimage)%rprim(:,:)   =results_img(iimage,itimimage_eff)%rprim(:,:)
306    results_img(iimage,next_itimimage)%vel(:,:)     =results_img(iimage,itimimage_eff)%vel(:,:)
307    results_img(iimage,next_itimimage)%vel_cell(:,:)=results_img(iimage,itimimage_eff)%vel_cell(:,:)
308  end do
309  ABI_DEALLOCATE(xred)
310 
311 end subroutine predict_string