TABLE OF CONTENTS
ABINIT/critics [ Functions ]
NAME
critics
FUNCTION
Search for critical points starting between atom inxat and its neighbors.
COPYRIGHT
Copyright (C) 2002-2017 ABINIT group (PCasek,FF,XG,MM) 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
aim_dtset= the structured entity containing all input variables dstmax=maximum distance to search for neighbors stwo, sthree, sfour: logical switches (TRUE/FALSE) indicating to search CP starting in the middle point of two, three or four atoms. One of these atoms is inxat.
OUTPUT
(see side effects)
SIDE EFFECTS
This routines acts primarily on the data contained in the aim_prom module WARNING This file does not follow the ABINIT coding rules (yet)
PARENTS
drvaim
CHILDREN
bschg1,critic,sort_dp,vgh_rho
SOURCE
42 #if defined HAVE_CONFIG_H 43 #include "config.h" 44 #endif 45 46 #include "abi_common.h" 47 48 49 subroutine critics(aim_dtset,inxat,stwo,sthree,sfour,dstmax) 50 51 use defs_basis 52 use defs_aimprom 53 use defs_parameters 54 use defs_abitypes 55 use m_sort 56 use m_profiling_abi 57 58 !This section has been created automatically by the script Abilint (TD). 59 !Do not modify the following lines by hand. 60 #undef ABI_FUNC 61 #define ABI_FUNC 'critics' 62 use interfaces_63_bader, except_this_one => critics 63 !End of the abilint section 64 65 implicit none 66 67 !Arguments ------------------------------------ 68 !scalars 69 integer,intent(in) :: inxat 70 real(dp),intent(in) :: dstmax 71 logical,intent(in) :: sfour,sthree,stwo 72 !no_abirules 73 type(aim_dataset_type), intent(in) :: aim_dtset 74 75 !Local variables ------------------------------ 76 !scalars 77 integer :: i1,i2,i3,iat,ii,ipos,ires,jii,jj,kjj,kk,ll,n1,n2,n3,nb,nc 78 integer :: nshell 79 real(dp) :: chg,dif1,dif2,diff,dist,olddist,rr 80 ! real(dp) :: ss,uu 81 logical :: found,inter 82 !arrays 83 integer :: ibat(nnpos*natom),inat(nnpos*natom),ipibat(nnpos*natom) 84 integer :: nnat(nnpos*natom),nr(nnpos*natom) 85 real(dp) :: dif(3),dists(nnpos*natom),ev(3),grho(3),hrho(3,3) 86 real(dp) :: pom(3),v1(3),v2(3),v3(3),v4(3),vi(3),vt(3),zz(3,3) 87 88 !************************************************************************ 89 vi(:)=xatm(:,inxat) 90 91 nc=0 92 do jii=1,nnpos 93 do kjj=1,natom 94 dist=0._dp 95 dif(:)=xatm(:,inxat)-xatm(:,kjj)-atp(:,jii) 96 97 ! do ii=1,3 98 ! dif(ii)=xatm(ii,inxat)-xatm(ii,kjj)-atp(ii,jii) 99 ! end do 100 dist=vnorm(dif,0) 101 if (.not.((dist>dstmax).or.(dist<0.001))) then 102 nc=nc+1 103 dists(nc)=dist 104 nnat(nc)=kjj 105 inat(nc)=jii 106 end if 107 end do 108 end do 109 do n1=1,nc 110 nr(n1)=n1 111 end do 112 call sort_dp(nc,dists,nr,tol14) 113 nb=0 114 olddist=0._dp 115 nshell=0 116 !write(std_out,*) ':ORIAT ', (xatm(ii,inxat),ii=1,3) 117 do n1=1,nc 118 n2=nr(n1) 119 n3=nnat(n2) 120 if (dists(n1)<(2*dists(1))) then 121 if ((dists(n1)-olddist)>aim_dlimit) then 122 nshell=nshell+1 123 olddist=dists(n1) 124 if (nshell==5) exit 125 end if 126 nb=nb+1 127 ibat(nb)=n3 128 ipibat(nb)=inat(n2) 129 write(std_out,*) ':NEIG ',inxat,n3,inat(n2),dists(n1) 130 ! write(std_out,*) ':POSAT',(xatm(ii,ibat(nb))+atp(ii,ipibat(nb)),ii=1,3) 131 else 132 exit 133 end if 134 end do 135 136 npc=0 137 npcm3=0 138 139 ! 140 !.....SEARCH BETWEEN EACH PAIR OF ATOMS 141 ! 142 143 if (stwo) then 144 do jii=1,nb 145 do ii=1,3 146 v1(ii)=xatm(ii,inxat) 147 v2(ii)=xatm(ii,ibat(jii))+atp(ii,ipibat(jii)) 148 vt(ii)=(v1(ii)+v2(ii))/2._dp 149 end do 150 inter=.true. 151 diff=0._dp 152 pom(:)=vt(:) 153 pom(:)=pom(:)-vi(:) 154 diff=vnorm(pom,0) 155 if (diff > maxcpdst) inter=.false. 156 if (inter) then 157 call critic(aim_dtset,vt,ev,zz,aim_dmaxcs,ires,0) 158 if (ires==0) then 159 found=.false. 160 if (npc > 0) then 161 do jj=1,npc 162 pom(:)=vt(:)-pc(:,jj) 163 dist=vnorm(pom,0) 164 if (dist < aim_dtset%dpclim) found=.true. 165 end do 166 end if 167 if (.not.found) then 168 pom(:)=vt(:) 169 call bschg1(pom,-1) 170 pcrb(:,npc+1)=pom(:) 171 pom(:)=pom(:)-vi(:) 172 diff=vnorm(pom,0) 173 if (abs(diff) > maxcpdst) found=.true. 174 end if 175 if (.not.found) then 176 npc=npc+1 177 do jj=1,3 178 pc(jj,npc)=vt(jj) 179 evpc(jj,npc)=ev(jj) 180 do kk=1,3 181 zpc(kk,jj,npc)=zz(kk,jj) 182 end do 183 end do 184 i1=ev(1)/abs(ev(1)) 185 i2=ev(2)/abs(ev(2)) 186 i3=ev(3)/abs(ev(3)) 187 icpc(npc)=i1+i2+i3 188 if (icpc(npc)==-3) then 189 npcm3=npcm3+1 190 end if 191 write(std_out,*) 'New critical point found' 192 write(std_out,'("POS: ",3F16.8)') (pc(ii,npc),ii=1,3) 193 write(std_out,'("POS in base: ",3F16.8)') (pcrb(ii,npc),ii=1,3) 194 write(std_out,'("AUTOVAL: ",3F16.8)') ev 195 write(std_out,'("AUTOVEC: ",/,3F16.8,/,3F16.8,/,3F16.8)') & 196 & ((zpc(ii,jj,npc),ii=1,3),jj=1,3) 197 call vgh_rho(vt,chg,grho,hrho,rr,iat,ipos,0) 198 write(22,'(":PC2",3F10.6,3E12.4,I4,2E12.4)') & 199 & (pc(jj,npc),jj=1,3),(ev(jj),jj=1,3),icpc(npc),ev(1)+ev(2)+ev(3),chg 200 write(std_out,'(":PC2",3F10.6,3E12.4,I4,2E12.4)') & 201 & (pc(jj,npc),jj=1,3),(ev(jj),jj=1,3),icpc(npc),ev(1)+ev(2)+ev(3),chg 202 pom(:)=vt(:)-v1(:) 203 dif1=vnorm(pom,0) 204 pom(:)=vt(:)-v2(:) 205 dif2=vnorm(pom,0) 206 write(std_out,'(":DISPC ",2F12.8)') dif1,dif2 207 end if 208 end if 209 end if 210 end do 211 end if 212 ! 213 !.....SEARCH BETWEEN EACH THREE ATOMS 214 ! 215 if(sthree) then 216 do jii=1,nb 217 do kjj=jii+1,nb 218 do ii=1,3 219 v1(ii)=xatm(ii,inxat) 220 v2(ii)=xatm(ii,ibat(jii))+atp(ii,ipibat(jii)) 221 v3(ii)=xatm(ii,ibat(kjj))+atp(ii,ipibat(kjj)) 222 vt(ii)=(v1(ii)+v2(ii)+v3(ii))/3._dp 223 end do 224 inter=.true. 225 pom(:)=vt(:) 226 pom(:)=pom(:)-vi(:) 227 dist=vnorm(pom,0) 228 if (abs(diff)>maxcpdst) then 229 inter=.false. 230 exit 231 end if 232 if (inter) then 233 do jj=1,npc 234 pom(:)=pc(:,jj)-vt(:) 235 diff=vnorm(pom,0) 236 if (diff<aim_dpc0) then 237 inter=.false. 238 exit 239 end if 240 end do 241 end if 242 if (inter) then 243 call critic(aim_dtset,vt,ev,zz,aim_dmaxcs,ires,0) 244 if (ires==0) then 245 found=.false. 246 if (npc>0) then 247 do jj=1,npc 248 pom(:)=vt(:)-pc(:,jj) 249 dist=vnorm(pom,0) 250 if (dist<aim_dtset%dpclim) then 251 found=.true. 252 exit 253 end if 254 end do 255 end if 256 if (.not.found) then 257 pom(:)=vt(:) 258 call bschg1(pom,-1) 259 pcrb(:,npc+1)=pom(:) 260 pom(:)=pom(:)-vi(:) 261 diff=vnorm(pom,0) 262 if (abs(diff)>maxcpdst) found=.true. 263 end if 264 if (.not.found) then 265 npc=npc+1 266 do jj=1,3 267 pc(jj,npc)=vt(jj) 268 evpc(jj,npc)=ev(jj) 269 do kk=1,3 270 zpc(kk,jj,npc)=zz(kk,jj) 271 end do 272 end do 273 i1=ev(1)/abs(ev(1)) 274 i2=ev(2)/abs(ev(2)) 275 i3=ev(3)/abs(ev(3)) 276 icpc(npc)=i1+i2+i3 277 if (icpc(npc)==-3) then 278 npcm3=npcm3+1 279 end if 280 write(std_out,*) 'New critical point found' 281 write(std_out,'("POS: ",3F16.8)') (pc(ii,npc),ii=1,3) 282 write(std_out,'("POS in base: ",3F16.8)') (pcrb(ii,npc),ii=1,3) 283 write(std_out,'("AUTOVAL: ",3F16.8)') ev 284 write(std_out,'("AUTOVEC: ",/,3F16.8,/,3F16.8,/,3F16.8)') & 285 & ((zpc(ii,jj,npc),ii=1,3),jj=1,3) 286 call vgh_rho(vt,chg,grho,hrho,rr,iat,ipos,0) 287 write(22,'(":PC3",3F10.6,3E12.4,I4,2E12.4)') & 288 & (pc(jj,npc),jj=1,3),(ev(jj),jj=1,3),icpc(npc),ev(1)+ev(2)+ev(3),chg 289 write(std_out,'(":PC3",3F10.6,3E12.4,I4,2E12.4)') & 290 & (pc(jj,npc),jj=1,3),(ev(jj),jj=1,3),icpc(npc),ev(1)+ev(2)+ev(3),chg 291 end if 292 end if 293 end if 294 end do 295 end do 296 end if 297 298 ! 299 !.....SEARCH BETWEEN EACH FOUR ATOMS 300 ! 301 if (sfour) then 302 do jii=1,nb 303 do kjj=jii+1,nb 304 do ll=jii+1,nb 305 do ii=1,3 306 v1(ii)=xatm(ii,inxat) 307 v2(ii)=xatm(ii,ibat(jii))+atp(ii,ipibat(jii)) 308 v3(ii)=xatm(ii,ibat(kjj))+atp(ii,ipibat(kjj)) 309 v4(ii)=xatm(ii,ibat(ll))+atp(ii,ipibat(ll)) 310 vt(ii)=(v1(ii)+v2(ii)+v3(ii)+v4(ii))/4._dp 311 end do 312 inter=.true. 313 pom(:)=vt(:) 314 pom(:)=pom(:)-vi(:) 315 diff=vnorm(pom,0) 316 if (abs(diff)>maxcpdst) then 317 inter=.false. 318 exit 319 end if 320 if (inter) then 321 do jj=1,npc 322 pom(:)=pc(:,jj)-vt(:) 323 diff=vnorm(pom,0) 324 if (diff < aim_dpc0) then 325 inter=.false. 326 exit 327 end if 328 end do 329 end if 330 if (inter) then 331 call critic(aim_dtset,vt,ev,zz,aim_dmaxcs,ires,0) 332 if (ires==0) then 333 found=.false. 334 if (npc>0) then 335 do jj=1,npc 336 pom(:)=vt(:)-pc(:,jj) 337 dist=vnorm(pom,0) 338 if (dist < aim_dtset%dpclim) found=.true. 339 end do 340 end if 341 if (.not.found) then 342 pom(:)=vt(:) 343 pcrb(:,npc+1)=pom(:) 344 pom(:)=pom(:)-vi(:) 345 diff=vnorm(pom,0) 346 if (abs(diff)>maxcpdst) found=.true. 347 end if 348 if (.not.found) then 349 npc=npc+1 350 do jj=1,3 351 pc(jj,npc)=vt(jj) 352 evpc(jj,npc)=ev(jj) 353 do kk=1,3 354 zpc(kk,jj,npc)=zz(kk,jj) 355 end do 356 end do 357 i1=ev(1)/abs(ev(1)) 358 i2=ev(2)/abs(ev(2)) 359 i3=ev(3)/abs(ev(3)) 360 icpc(npc)=i1+i2+i3 361 if (icpc(npc)==-3) then 362 npcm3=npcm3+1 363 end if 364 write(std_out,*) 'New critical point found' 365 write(std_out,'("POS: ",3F16.8)') (pc(ii,npc),ii=1,3) 366 write(std_out,'("POS in base: ",3F16.8)') (pcrb(ii,npc),ii=1,3) 367 write(std_out,'("AUTOVAL: ",3F16.8)') ev 368 write(std_out,'("AUTOVEC: ",/,3F16.8,/,3F16.8,/,3F16.8)') & 369 & ((zpc(ii,jj,npc),ii=1,3),jj=1,3) 370 call vgh_rho(vt,chg,grho,hrho,rr,iat,ipos,0) 371 write(22,'(":PC4",3F10.6,3E12.4,I4,2E12.4)') & 372 & (pc(jj,npc),jj=1,3), (ev(jj),jj=1,3),icpc(npc),ev(1)+ev(2)+ev(3),chg 373 write(std_out,'(":PC4",3F10.6,3E12.4,I4,2E12.4)') & 374 & (pc(jj,npc),jj=1,3), (ev(jj),jj=1,3),icpc(npc),ev(1)+ev(2)+ev(3),chg 375 end if 376 end if 377 end if 378 end do 379 end do 380 end do 381 end if 382 383 write(std_out,*) npc 384 end subroutine critics