TABLE OF CONTENTS
ABINIT/xfh_recover_deloc [ Functions ]
NAME
xfh_recover_deloc
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_delocint
CHILDREN
fred2fdeloc,hessupdt,xcart2deloc,xfpack_f2vout,xfpack_x2vin,xred2xcart
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_deloc(ab_xfh,ab_mover,acell,acell0,cycl_main,& 42 & fred,hessin,ndim,rprim,rprimd0,strten,ucvol,ucvol0,vin,vin_prev,& 43 & vout,vout_prev,xred,deloc,deloc_int,deloc_force,btinv,gprimd,prim_int,& 44 & u_matrix) 45 46 use m_profiling_abi 47 use defs_basis 48 use m_abimover 49 50 use m_geometry, only : xred2xcart 51 use m_results_gs , only : results_gs_type 52 use m_bfgs, only : hessupdt 53 54 !This section has been created automatically by the script Abilint (TD). 55 !Do not modify the following lines by hand. 56 #undef ABI_FUNC 57 #define ABI_FUNC 'xfh_recover_deloc' 58 use interfaces_45_geomoptim, except_this_one => xfh_recover_deloc 59 !End of the abilint section 60 61 implicit none 62 63 !Arguments ------------------------------------ 64 !scalars 65 66 integer,intent(in) :: ndim 67 integer,intent(out) :: cycl_main 68 real(dp),intent(inout) :: ucvol,ucvol0 69 type(ab_xfh_type),intent(inout) :: ab_xfh 70 type(abimover),intent(in) :: ab_mover 71 ! DELOCINT specials 72 type(delocint),intent(inout) :: deloc 73 74 !arrays 75 real(dp),intent(inout) :: acell(3) 76 real(dp),intent(in) :: acell0(3) 77 real(dp),intent(inout) :: hessin(:,:) 78 real(dp),intent(inout) :: xred(3,ab_mover%natom) 79 real(dp),intent(inout) :: rprim(3,3) 80 real(dp),intent(inout) :: rprimd0(3,3) 81 real(dp),intent(inout) :: fred(3,ab_mover%natom) 82 real(dp),intent(inout) :: strten(6) 83 real(dp),intent(inout) :: vin(:) 84 real(dp),intent(inout) :: vin_prev(:) 85 real(dp),intent(inout) :: vout(:) 86 real(dp),intent(inout) :: vout_prev(:) 87 ! DELOCINT specials 88 real(dp),intent(inout) :: deloc_force(3*(ab_mover%natom-1)) 89 real(dp),intent(inout) :: deloc_int(3*(ab_mover%natom-1)) 90 real(dp),intent(inout) :: btinv(3*(ab_mover%natom-1),3*ab_mover%natom) 91 real(dp),intent(inout) :: prim_int(:),u_matrix(:,:),gprimd(3,3) 92 93 !Local variables------------------------------- 94 !scalars 95 integer :: ixfh 96 real(dp) :: xcart(3,ab_mover%natom) 97 98 !********************************************************************* 99 100 if(ab_xfh%nxfh/=0)then 101 ! Loop over previous time steps 102 do ixfh=1,ab_xfh%nxfh 103 104 ! For that time step, get new (x,f) from xfhist 105 xred(:,:) =ab_xfh%xfhist(:,1:ab_mover%natom ,1,ixfh) 106 rprim(1:3,1:3)=ab_xfh%xfhist(:,ab_mover%natom+2:ab_mover%natom+4,1,ixfh) 107 acell(:) =ab_xfh%xfhist(:,ab_mover%natom+1,1,ixfh) 108 fred(:,:) =ab_xfh%xfhist(:,1:ab_mover%natom,2,ixfh) 109 ! This use of results_gs is unusual 110 strten(1:3) =ab_xfh%xfhist(:,ab_mover%natom+2,2,ixfh) 111 strten(4:6) =ab_xfh%xfhist(:,ab_mover%natom+3,2,ixfh) 112 113 ! !DEBUG 114 ! write (ab_out,*) '---READ FROM XFHIST---' 115 116 ! write (ab_out,*) 'XRED' 117 ! do kk=1,ab_mover%natom 118 ! write (ab_out,*) xred(:,kk) 119 ! end do 120 ! write (ab_out,*) 'FRED' 121 ! do kk=1,ab_mover%natom 122 ! write (ab_out,*) fred(:,kk) 123 ! end do 124 ! write(ab_out,*) 'RPRIM' 125 ! do kk=1,3 126 ! write(ab_out,*) rprim(:,kk) 127 ! end do 128 ! write(ab_out,*) 'ACELL' 129 ! write(ab_out,*) acell(:) 130 ! !DEBUG 131 132 ! Convert input xred (reduced coordinates) to xcart (cartesian) 133 call xred2xcart(ab_mover%natom,rprimd0,xcart,xred) 134 ! Convert input coordinates in Delocalized internals 135 call xcart2deloc(deloc,ab_mover%natom,rprimd0,xcart,& 136 & btinv,u_matrix,deloc_int,prim_int) 137 ! Convert forces to delocalized coordinates for next step 138 call fred2fdeloc(btinv,deloc_force,fred,ab_mover%natom,gprimd) 139 140 ! Transfer it in vin, vout 141 call xfpack_x2vin(acell,acell0,ab_mover%natom-1,& 142 & ndim,ab_mover%nsym,ab_mover%optcell,rprim,rprimd0,& 143 & ab_mover%symrel,ucvol,ucvol0,vin,deloc_int) 144 call xfpack_f2vout(deloc_force,ab_mover%natom-1,& 145 & ndim,ab_mover%optcell,ab_mover%strtarget,strten,& 146 & ucvol,vout) 147 ! Get old time step, if any, and update inverse hessian 148 if(ixfh/=1)then 149 xred(:,:) =ab_xfh%xfhist(:,1:ab_mover%natom,1,ixfh-1) 150 rprim(1:3,1:3)=& 151 & ab_xfh%xfhist(:,ab_mover%natom+2:ab_mover%natom+4,1,ixfh-1) 152 acell(:)=ab_xfh%xfhist(:,ab_mover%natom+1,1,ixfh-1) 153 fred(:,:)=ab_xfh%xfhist(:,1:ab_mover%natom,2,ixfh-1) 154 ! This use of results_gs is unusual 155 strten(1:3)=ab_xfh%xfhist(:,ab_mover%natom+2,2,ixfh-1) 156 strten(4:6)=ab_xfh%xfhist(:,ab_mover%natom+3,2,ixfh-1) 157 158 ! Convert input xred (reduced coordinates) to xcart (cartesian) 159 call xred2xcart(ab_mover%natom,rprimd0,xcart,xred) 160 ! Convert input coordinates in Delocalized internals 161 call xcart2deloc(deloc,ab_mover%natom,rprimd0,xcart,& 162 & btinv,u_matrix,deloc_int,prim_int) 163 ! Convert forces to delocalized coordinates for next step 164 call fred2fdeloc(btinv,deloc_force,fred,ab_mover%natom,gprimd) 165 166 ! Tranfer it in vin_prev, vout_prev 167 call xfpack_x2vin(acell,acell0,ab_mover%natom-1,& 168 & ndim,ab_mover%nsym,ab_mover%optcell,rprim,rprimd0,& 169 & ab_mover%symrel,ucvol,ucvol0,vin_prev,deloc_int) 170 call xfpack_f2vout(deloc_force,ab_mover%natom-1,& 171 & ndim,ab_mover%optcell,ab_mover%strtarget,strten,& 172 & ucvol,vout_prev) 173 174 ! write(ab_out,*) 'Hessian matrix before update',ndim,'x',ndim 175 ! write(ab_out,*) 'ixfh=',ixfh 176 ! do kk=1,ndim 177 ! do jj=1,ndim,3 178 ! if (jj+2<=ndim)then 179 ! write(ab_out,*) jj,hessin(jj:jj+2,kk) 180 ! else 181 ! write(ab_out,*) jj,hessin(jj:ndim,kk) 182 ! end if 183 ! end do 184 ! end do 185 186 call hessupdt(hessin,ab_mover%iatfix,ab_mover%natom-1,ndim,& 187 & vin,vin_prev,vout,vout_prev) 188 189 ! !DEBUG 190 ! write(ab_out,*) 'Hessian matrix after update',ndim,'x',ndim 191 ! do kk=1,ndim 192 ! do jj=1,ndim,3 193 ! if (jj+2<=ndim)then 194 ! write(ab_out,*) jj,hessin(jj:jj+2,kk) 195 ! else 196 ! write(ab_out,*) jj,hessin(jj:ndim,kk) 197 ! end if 198 ! end do 199 ! end do 200 ! !DEBUG 201 202 end if !if(ab_xfh%nxfh/=0) 203 204 ! End loop over previous time steps 205 end do 206 207 ! The hessian has been generated, 208 ! as well as the latest vin and vout 209 ! so will cycle the main loop 210 cycl_main=1 211 212 end if 213 214 end subroutine xfh_recover_deloc