TABLE OF CONTENTS
ABINIT/xfh_recover_new [ Functions ]
NAME
xfh_recover_new
FUNCTION
Update the contents of the history xfhist taking values from xred, acell, rprim, fred_corrected and strten
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR, JCC, SE) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
OUTPUT
SIDE EFFECTS
PARENTS
pred_bfgs,pred_lbfgs
CHILDREN
hessupdt,xfpack_f2vout,xfpack_x2vin
SOURCE
34 #if defined HAVE_CONFIG_H 35 #include "config.h" 36 #endif 37 38 #include "abi_common.h" 39 40 41 subroutine xfh_recover_new(ab_xfh,ab_mover,acell,acell0,cycl_main,fred,& 42 & hessin,ndim,rprim,rprimd0,strten,ucvol,ucvol0,vin,vin_prev,vout,& 43 & vout_prev,xred) 44 45 use defs_basis 46 use m_profiling_abi 47 use m_abimover 48 49 use m_results_gs , only : results_gs_type 50 use m_bfgs, only : hessupdt 51 52 !This section has been created automatically by the script Abilint (TD). 53 !Do not modify the following lines by hand. 54 #undef ABI_FUNC 55 #define ABI_FUNC 'xfh_recover_new' 56 use interfaces_45_geomoptim, except_this_one => xfh_recover_new 57 !End of the abilint section 58 59 implicit none 60 61 !Arguments ------------------------------------ 62 !scalars 63 64 integer,intent(in) :: ndim 65 integer,intent(out) :: cycl_main 66 real(dp),intent(inout) :: ucvol,ucvol0 67 type(ab_xfh_type),intent(inout) :: ab_xfh 68 type(abimover),intent(in) :: ab_mover 69 70 71 !arrays 72 real(dp),intent(inout) :: acell(3) 73 real(dp),intent(in) :: acell0(3) 74 real(dp),intent(inout) :: hessin(:,:) 75 real(dp),intent(inout) :: xred(3,ab_mover%natom) 76 real(dp),intent(inout) :: rprim(3,3) 77 real(dp),intent(inout) :: rprimd0(3,3) 78 real(dp),intent(inout) :: fred(3,ab_mover%natom) 79 real(dp),intent(inout) :: strten(6) 80 real(dp),intent(inout) :: vin(:) 81 real(dp),intent(inout) :: vin_prev(:) 82 real(dp),intent(inout) :: vout(:) 83 real(dp),intent(inout) :: vout_prev(:) 84 85 !Local variables------------------------------- 86 !scalars 87 integer :: ixfh ! kk,jj 88 89 !********************************************************************* 90 91 if(ab_xfh%nxfh/=0)then 92 ! Loop over previous time steps 93 do ixfh=1,ab_xfh%nxfh 94 95 ! For that time step, get new (x,f) from xfhist 96 xred(:,:) =ab_xfh%xfhist(:,1:ab_mover%natom ,1,ixfh) 97 rprim(1:3,1:3)=ab_xfh%xfhist(:,ab_mover%natom+2:ab_mover%natom+4,1,ixfh) 98 acell(:) =ab_xfh%xfhist(:,ab_mover%natom+1,1,ixfh) 99 fred(:,:) =ab_xfh%xfhist(:,1:ab_mover%natom,2,ixfh) 100 ! This use of results_gs is unusual 101 strten(1:3) =ab_xfh%xfhist(:,ab_mover%natom+2,2,ixfh) 102 strten(4:6) =ab_xfh%xfhist(:,ab_mover%natom+3,2,ixfh) 103 104 ! !DEBUG 105 ! write (ab_out,*) '---READED FROM XFHIST---' 106 107 ! write (ab_out,*) 'XRED' 108 ! do kk=1,ab_mover%natom 109 ! write (ab_out,*) xred(:,kk) 110 ! end do 111 ! write (ab_out,*) 'FRED' 112 ! do kk=1,ab_mover%natom 113 ! write (ab_out,*) fred(:,kk) 114 ! end do 115 ! write(ab_out,*) 'RPRIM' 116 ! do kk=1,3 117 ! write(ab_out,*) rprim(:,kk) 118 ! end do 119 ! write(ab_out,*) 'ACELL' 120 ! write(ab_out,*) acell(:) 121 ! !DEBUG 122 123 ! Transfer it in vin, vout 124 call xfpack_x2vin(acell,acell0,ab_mover%natom,& 125 & ndim,ab_mover%nsym,ab_mover%optcell,rprim,rprimd0,& 126 & ab_mover%symrel,ucvol,ucvol0,vin,xred) 127 call xfpack_f2vout(fred,ab_mover%natom,& 128 & ndim,ab_mover%optcell,ab_mover%strtarget,strten,& 129 & ucvol,vout) 130 ! Get old time step, if any, and update inverse hessian 131 if(ixfh/=1)then 132 xred(:,:) =ab_xfh%xfhist(:,1:ab_mover%natom,1,ixfh-1) 133 rprim(1:3,1:3)=& 134 & ab_xfh%xfhist(:,ab_mover%natom+2:ab_mover%natom+4,1,ixfh-1) 135 acell(:)=ab_xfh%xfhist(:,ab_mover%natom+1,1,ixfh-1) 136 fred(:,:)=ab_xfh%xfhist(:,1:ab_mover%natom,2,ixfh-1) 137 ! This use of results_gs is unusual 138 strten(1:3)=ab_xfh%xfhist(:,ab_mover%natom+2,2,ixfh-1) 139 strten(4:6)=ab_xfh%xfhist(:,ab_mover%natom+3,2,ixfh-1) 140 ! Tranfer it in vin_prev, vout_prev 141 call xfpack_x2vin(acell,acell0,ab_mover%natom,& 142 & ndim,ab_mover%nsym,ab_mover%optcell,rprim,rprimd0,& 143 & ab_mover%symrel,ucvol,ucvol0,vin_prev,xred) 144 call xfpack_f2vout(fred,ab_mover%natom,& 145 & ndim,ab_mover%optcell,ab_mover%strtarget,strten,& 146 & ucvol,vout_prev) 147 148 ! write(ab_out,*) 'Hessian matrix before update',ndim,'x',ndim 149 ! write(ab_out,*) 'ixfh=',ixfh 150 ! do kk=1,ndim 151 ! do jj=1,ndim,3 152 ! if (jj+2<=ndim)then 153 ! write(ab_out,*) jj,hessin(jj:jj+2,kk) 154 ! else 155 ! write(ab_out,*) jj,hessin(jj:ndim,kk) 156 ! end if 157 ! end do 158 ! end do 159 160 call hessupdt(hessin,ab_mover%iatfix,ab_mover%natom,ndim,& 161 & vin,vin_prev,vout,vout_prev) 162 163 ! !DEBUG 164 ! write(ab_out,*) 'Hessian matrix after update',ndim,'x',ndim 165 ! do kk=1,ndim 166 ! do jj=1,ndim,3 167 ! if (jj+2<=ndim)then 168 ! write(ab_out,*) jj,hessin(jj:jj+2,kk) 169 ! else 170 ! write(ab_out,*) jj,hessin(jj:ndim,kk) 171 ! end if 172 ! end do 173 ! end do 174 ! !DEBUG 175 176 end if !if(ab_xfh%nxfh/=0) 177 end do ! End loop over previous time steps 178 179 ! The hessian has been generated, 180 ! as well as the latest vin and vout 181 ! so will cycle the main loop 182 cycl_main=1 183 end if 184 185 end subroutine xfh_recover_new