TABLE OF CONTENTS
ABINIT/deloc2xcart [ Functions ]
NAME
deloc2xcart
FUNCTION
Determine the cartesian coordinates which correspond to the given values of the delocalized coordinates. The relationship is non-linear, so use an iterative scheme, as in Baker JCP .105. 192 (1996). Older reference: Pulay and co. JACS 101 2550 (1979)
COPYRIGHT
Copyright (C) 2003-2018 ABINIT group (MVer) 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
deloc <type(delocint)>=Important variables for | pred_delocint | | nang = Number of angles | nbond = Number of bonds | ncart = Number of cartesian directions | (used for constraints) | ndihed = Number of dihedrals | nrshift = Dimension of rshift | ninternal= Number of internal coordinates | ninternal=nbond+nang+ndihed+ncart | | angs(2,3,nang) = Indexes to characterize angles | bonds(2,2,nbond)= For a bond between iatom and jatom | bonds(1,1,nbond) = iatom | bonds(2,1,nbond) = icenter | bonds(1,2,nbond) = jatom | bonds(2,2,nbond) = irshift | carts(2,ncart) = Index of total primitive internal, | and atom (carts(2,:)) | dihedrals(2,4,ndihed)= Indexes to characterize dihedrals | | rshift(3,nrshift)= Shift in xred that must be done to find | all neighbors of a given atom within a | given number of neighboring shells natom = Number of atoms (dtset%natom) rprimd(3,3)=dimensional real space primitive translations (bohr)
OUTPUT
bt_inv_matrix(3*(natom-1),3*natom)=inverse of transpose of B matrix
SIDE EFFECTS
u_matrix(ninternal,3*(natom-1))=eigenvectors of G = BB^T matrix xcart(3,natom)=cartesian coordinates of atoms (bohr)
NOTES
PARENTS
pred_delocint
CHILDREN
dgemv,wrtout,xcart2deloc
SOURCE
66 #if defined HAVE_CONFIG_H 67 #include "config.h" 68 #endif 69 70 #include "abi_common.h" 71 72 73 subroutine deloc2xcart(deloc,natom,rprimd,xcart,deloc_int,btinv,u_matrix) 74 75 use defs_basis 76 use m_abimover 77 use m_errors 78 use m_profiling_abi 79 use m_linalg_interfaces 80 81 !This section has been created automatically by the script Abilint (TD). 82 !Do not modify the following lines by hand. 83 #undef ABI_FUNC 84 #define ABI_FUNC 'deloc2xcart' 85 use interfaces_14_hidewrite 86 use interfaces_45_geomoptim, except_this_one => deloc2xcart 87 !End of the abilint section 88 89 implicit none 90 91 !Arguments ------------------------------------ 92 !scalars 93 integer,intent(in) :: natom 94 type(delocint),intent(inout) :: deloc 95 !arrays 96 real(dp),intent(in) :: deloc_int(3*(natom-1)),rprimd(3,3) 97 real(dp),intent(inout) :: u_matrix(deloc%ninternal,3*(natom-1)) 98 real(dp),intent(inout) :: xcart(3,natom) 99 real(dp),intent(out) :: btinv(3*(natom-1),3*natom) 100 101 !Local variables------------------------------- 102 !scalars 103 integer :: iiter,iprim,niter 104 integer :: ii 105 real(dp) :: minmix, maxmix 106 real(dp) :: mix,tot_diff, toldeloc 107 real(dp) :: lntoldeloc 108 logical :: DEBUG=.FALSE. 109 !arrays 110 real(dp) :: btinv_tmp(3*(natom-1),3*natom) 111 real(dp) :: cgrad(3*natom),cgrad_old(3*natom) 112 real(dp) :: deloc_int_now(3*(natom-1)),prim_int(deloc%ninternal) 113 real(dp) :: tmpxcart(3*natom) 114 real(dp) :: xdeloc_diff(3*(natom-1)) 115 116 character(len=500) :: message 117 118 ! ****************************************************************** 119 120 if (DEBUG) then 121 write(ab_out,*) 'ENTERING DELOC2XCART' 122 123 write (message,*) 'BONDS=',deloc%nbond 124 call wrtout(ab_out,message,'COLL') 125 do ii = 1, deloc%nbond 126 write (message,*) ii, deloc%bonds(:,:,ii) 127 call wrtout(ab_out,message,'COLL') 128 end do 129 130 write (message,*) 'ANGS=',deloc%nang 131 call wrtout(ab_out,message,'COLL') 132 do ii = 1, deloc%nang 133 write (message,*) ii, deloc%angs(:,:,ii) 134 call wrtout(ab_out,message,'COLL') 135 end do 136 137 write (message,*) 'DIHEDRALS=',deloc%ndihed 138 call wrtout(ab_out,message,'COLL') 139 do ii = 1, deloc%ndihed 140 write (message,*) ii, deloc%dihedrals(:,:,ii) 141 call wrtout(ab_out,message,'COLL') 142 end do 143 144 write (message,*) 'CARTS=',deloc%ncart 145 call wrtout(ab_out,message,'COLL') 146 do ii = 1, deloc%ncart 147 write (message,*) ii, deloc%carts(:,ii) 148 call wrtout(ab_out,message,'COLL') 149 end do 150 151 write (ab_out,*) 'xcart (input)' 152 do ii=1,natom 153 write (ab_out,*) xcart(:,ii) 154 end do 155 156 end if 157 158 niter = 200 159 tmpxcart = reshape(xcart,(/3*natom/)) 160 161 cgrad_old(:) = zero 162 cgrad(:) = zero 163 maxmix = 0.9_dp 164 minmix = 0.2_dp 165 toldeloc = tol10 166 lntoldeloc = log(toldeloc) 167 168 do iiter=1,niter 169 if (iiter==1) then 170 mix= minmix 171 else 172 mix = minmix + (maxmix-minmix)*(log(tot_diff)-lntoldeloc) / lntoldeloc 173 end if 174 if (mix < minmix) mix = minmix 175 if (mix > maxmix) mix = maxmix 176 177 tmpxcart(:) = tmpxcart(:) + mix*cgrad(:) 178 xcart = reshape(tmpxcart,(/3,natom/)) 179 call xcart2deloc(deloc,natom,rprimd,xcart,& 180 & btinv_tmp,u_matrix,deloc_int_now,prim_int) 181 ! update the BT^{-1} matrix? 182 btinv(:,:) = btinv_tmp(:,:) 183 184 xdeloc_diff(:) = deloc_int(:) - deloc_int_now(:) 185 186 tot_diff = sum(abs(xdeloc_diff)) 187 if (tot_diff < toldeloc) exit 188 189 cgrad_old(:) = cgrad(:) 190 191 ! gradient vector = btinv^{T} * xdeloc_diff 192 call dgemv('T',3*(natom-1),3*natom,one,& 193 & btinv,3*(natom-1),xdeloc_diff,1,zero,cgrad,1) 194 end do 195 !end iiter do 196 197 call xcart2deloc(deloc,natom,rprimd,xcart,& 198 & btinv,u_matrix,deloc_int_now,prim_int) 199 write (message,'(3a)') 'delocalized internals, after convergence of xcart = ', ch10 200 call wrtout(std_out,message,'COLL') 201 do ii = 1, 3*(natom-1) 202 write (message,'(I6,E20.10,2x)') ii, deloc_int_now(ii) 203 call wrtout(std_out,message,'COLL') 204 end do 205 206 xdeloc_diff(:) = deloc_int(:) - deloc_int_now(:) 207 208 write (message,'(a)') 'Primitive internal coordinate values:' 209 call wrtout(std_out,message,'COLL') 210 do iprim = 1, deloc%nbond 211 write (message,'(i6,E20.10)') iprim, prim_int(iprim) 212 call wrtout(std_out,message,'COLL') 213 end do 214 do iprim = deloc%nbond+1, deloc%nbond+deloc%nang+deloc%ndihed 215 write (message,'(i6,2E20.10)') iprim, prim_int(iprim), prim_int(iprim)/pi*180.0_dp 216 call wrtout(std_out,message,'COLL') 217 end do 218 do iprim = deloc%nbond+deloc%nang+deloc%ndihed+1, deloc%ninternal 219 write (message,'(i6,E20.10)') iprim, prim_int(iprim) 220 call wrtout(std_out,message,'COLL') 221 end do 222 223 if (iiter == niter+1) then 224 write (message,'(a,i6,a,E20.10)') 'deloc2xcart : Error, xcart not converged in ', niter, 'iterations ', tot_diff 225 MSG_ERROR(message) 226 end if 227 228 if(DEBUG)then 229 write (ab_out,*) 'xcart (output)' 230 do ii=1,natom 231 write (ab_out,*) xcart(:,ii) 232 end do 233 write(ab_out,*) 'EXITING DELOC2XCART' 234 end if 235 236 end subroutine deloc2xcart