TABLE OF CONTENTS
ABINIT/m_predict_steepest [ Modules ]
NAME
FUNCTION
COPYRIGHT
Copyright (C) 2009-2022 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 .
SOURCE
14 #if defined HAVE_CONFIG_H 15 #include "config.h" 16 #endif 17 18 #include "abi_common.h" 19 20 module m_predict_steepest 21 22 use defs_basis 23 use m_abicore 24 use m_mep 25 26 use m_results_img, only : results_img_type,get_geometry_img 27 28 implicit none 29 30 private
ABINIT/predict_steepest [ 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
SOURCE
83 subroutine predict_steepest(itimimage,itimimage_eff,list_dynimage,mep_param,natom,& 84 & ndynimage,nimage,ntimimage_stored,results_img) 85 86 implicit none 87 88 !Arguments ------------------------------------ 89 !scalars 90 integer,intent(in) :: itimimage,itimimage_eff,natom,ndynimage 91 integer,intent(in) :: nimage,ntimimage_stored 92 type(mep_type),intent(inout) :: mep_param 93 !arrays 94 integer,intent(in) :: list_dynimage(ndynimage) 95 type(results_img_type) :: results_img(nimage,ntimimage_stored) 96 97 !Local variables------------------------------- 98 !scalars 99 integer :: iimage,next_itimimage 100 !arrays 101 real(dp),allocatable :: etotal(:),fcart(:,:,:),rprimd(:,:,:),xcart(:,:,:),xred(:,:,:) 102 103 ! ************************************************************************* 104 105 !Retrieve positions and forces 106 ABI_MALLOC(etotal,(nimage)) 107 ABI_MALLOC(xred,(3,natom,nimage)) 108 ABI_MALLOC(xcart,(3,natom,nimage)) 109 ABI_MALLOC(fcart,(3,natom,nimage)) 110 ABI_MALLOC(rprimd,(3,3,nimage)) 111 call get_geometry_img(etotal,natom,nimage,results_img(:,itimimage_eff),& 112 & fcart,rprimd,xcart,xred) 113 114 !Compute new atomic positions in each cell 115 call mep_steepest(fcart,list_dynimage,mep_param,natom,ndynimage,nimage,rprimd,xcart,xred) 116 117 !Store acell, rprim, xred and vel for the new iteration 118 next_itimimage=itimimage+1 119 if (next_itimimage>ntimimage_stored) next_itimimage=1 120 do iimage=1,nimage 121 results_img(iimage,next_itimimage)%xred(:,:) =xred(:,:,iimage) 122 results_img(iimage,next_itimimage)%acell(:) =results_img(iimage,itimimage_eff)%acell(:) 123 results_img(iimage,next_itimimage)%rprim(:,:) =results_img(iimage,itimimage_eff)%rprim(:,:) 124 results_img(iimage,next_itimimage)%vel(:,:) =results_img(iimage,itimimage_eff)%vel(:,:) 125 results_img(iimage,next_itimimage)%vel_cell(:,:)=results_img(iimage,itimimage_eff)%vel_cell(:,:) 126 end do 127 128 ABI_FREE(etotal) 129 ABI_FREE(xred) 130 ABI_FREE(xcart) 131 ABI_FREE(fcart) 132 ABI_FREE(rprimd) 133 134 end subroutine predict_steepest