TABLE OF CONTENTS


ABINIT/m_predict_string [ Modules ]

[ Top ] [ Modules ]

NAME

  m_predict_string

FUNCTION

COPYRIGHT

  Copyright (C) 2009-2022 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 .

SOURCE

16 #if defined HAVE_CONFIG_H
17 #include "config.h"
18 #endif
19 
20 #include "abi_common.h"
21 
22 module m_predict_string
23 
24  use defs_basis
25  use m_abicore
26  use m_splines
27  use m_mep
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, get_geometry_img
33 
34  implicit none
35 
36  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

SOURCE

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