TABLE OF CONTENTS


ABINIT/m_predict_steepest [ Modules ]

[ Top ] [ Modules ]

NAME

FUNCTION

COPYRIGHT

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

19 #if defined HAVE_CONFIG_H
20 #include "config.h"
21 #endif
22 
23 #include "abi_common.h"
24 
25 module m_predict_steepest
26 
27  use defs_basis
28  use m_abicore
29  use m_mep
30 
31  use m_results_img, only : results_img_type,get_geometry_img
32 
33  implicit none
34 
35  private

ABINIT/predict_steepest [ Functions ]

[ Top ] [ Functions ]

NAME

 predict_steepest

FUNCTION

 Given the past history of images, predict the new set of images.
 Here, simple steepest descent algorithm, based on the value of the forces on the current timimage step.
 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.
 mep_param=several parameters for Minimal Energy Path (MEP) search
 natom=dimension of vel_timimage and xred_timimage
 ndynimage=number of dynamical images
 nimage=number of images
 ntimimage_stored=number of time steps stored in the history

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

      get_geometry_img,mep_steepest

SOURCE

 94 subroutine predict_steepest(itimimage,itimimage_eff,list_dynimage,mep_param,natom,&
 95 &                           ndynimage,nimage,ntimimage_stored,results_img)
 96 
 97 
 98 !This section has been created automatically by the script Abilint (TD).
 99 !Do not modify the following lines by hand.
100 #undef ABI_FUNC
101 #define ABI_FUNC 'predict_steepest'
102 !End of the abilint section
103 
104  implicit none
105 
106 !Arguments ------------------------------------
107 !scalars
108  integer,intent(in) :: itimimage,itimimage_eff,natom,ndynimage
109  integer,intent(in) :: nimage,ntimimage_stored
110  type(mep_type),intent(inout) :: mep_param
111 !arrays
112  integer,intent(in) :: list_dynimage(ndynimage)
113  type(results_img_type) :: results_img(nimage,ntimimage_stored)
114 
115 !Local variables-------------------------------
116 !scalars
117  integer :: iimage,next_itimimage
118 !arrays
119  real(dp),allocatable :: etotal(:),fcart(:,:,:),rprimd(:,:,:),xcart(:,:,:),xred(:,:,:)
120 
121 ! *************************************************************************
122 
123 !Retrieve positions and forces
124  ABI_ALLOCATE(etotal,(nimage))
125  ABI_ALLOCATE(xred,(3,natom,nimage))
126  ABI_ALLOCATE(xcart,(3,natom,nimage))
127  ABI_ALLOCATE(fcart,(3,natom,nimage))
128  ABI_ALLOCATE(rprimd,(3,3,nimage))
129  call get_geometry_img(etotal,natom,nimage,results_img(:,itimimage_eff),&
130 & fcart,rprimd,xcart,xred)
131 
132 !Compute new atomic positions in each cell
133  call mep_steepest(fcart,list_dynimage,mep_param,natom,ndynimage,nimage,rprimd,xcart,xred)
134 
135 !Store acell, rprim, xred and vel for the new iteration
136  next_itimimage=itimimage+1
137  if (next_itimimage>ntimimage_stored) next_itimimage=1
138  do iimage=1,nimage
139    results_img(iimage,next_itimimage)%xred(:,:)    =xred(:,:,iimage)
140    results_img(iimage,next_itimimage)%acell(:)     =results_img(iimage,itimimage_eff)%acell(:)
141    results_img(iimage,next_itimimage)%rprim(:,:)   =results_img(iimage,itimimage_eff)%rprim(:,:)
142    results_img(iimage,next_itimimage)%vel(:,:)     =results_img(iimage,itimimage_eff)%vel(:,:)
143    results_img(iimage,next_itimimage)%vel_cell(:,:)=results_img(iimage,itimimage_eff)%vel_cell(:,:)
144  end do
145 
146  ABI_DEALLOCATE(etotal)
147  ABI_DEALLOCATE(xred)
148  ABI_DEALLOCATE(xcart)
149  ABI_DEALLOCATE(fcart)
150  ABI_DEALLOCATE(rprimd)
151 
152 end subroutine predict_steepest