1c Subject : Implementation of NMR Hyperfine Tensors 2c using ZORA. 3c Author : Fredy Aquino 4c Date : 03-10-11 5 subroutine get_NMRHFine_ZORA(rtdb, 6 & g_dens_at, 7 & nexc, 8 & geom, 9 & ao_bas_han, 10 & nbf, 11 & focc, 12 & noc, 13 & ipol, 14 & g_densZ4) 15 implicit none 16c 17#include "errquit.fh" 18#include "mafdecls.fh" 19#include "stdio.fh" 20#include "global.fh" 21#include "msgids.fh" 22#include "rtdb.fh" 23#include "geom.fh" 24#include "zora.fh" 25 integer rtdb 26 integer g_densZ4(3) 27 integer g_sdens,g_dens_at(2) 28 integer geom 29 integer ao_bas_han 30 integer nbf,noc(2),ipol,nexc 31 double precision focc(nbf*ipol) ! contains occupation values if modified 32 integer nat_slc,typeprop,nat,l_AtNr,k_AtNr 33 integer alo(3),ahi(3),ld(2) 34 logical status 35 integer i,j,type_nmrdata 36 logical dft_zoraNMR_write 37 integer ga_dia_hfine,ga_para1, 38 & ga_h01_hfine,ga_Fji_hfine 39 logical is_atom,ofinite,Knucl 40 double precision atmass,AtNr_dbl 41 character*255 zorafilename 42 43 external get_Htensor_slow,get_HFine_F1ji, 44 & dft_zoraNMR_write, 45 & get_slctd_atoms 46 logical do_prntNMRCS 47 if(.not.rtdb_get(rtdb,'zora:do_prntNMRCS', ! FA 48 & mt_log,1,do_prntNMRCS)) ! FA 49 & do_prntNMRCS= .false. 50c ------ Read Knucl for including ONLY nuclear part in K ZORA ----- START 51c Note.- stored in rel_input.F(rel_input(rtdb)) 52 Knucl=.false. 53 status=rtdb_get(rtdb,'zora:Knucl',mt_log,1,Knucl) ! Check if gaussian nucl model requested 54 if (ga_nodeid().eq.0) 55 & write(*,*) 'In get_NMRHFine_ZORA:: zora:Knucl=',Knucl 56c ------ Read Knucl for including ONLY nuclear part in K ZORA ----- END 57c ------ Read ofinite to be used by HFine finite calc ---FA-03-21-11-- START 58c Note.- stored in geom_input.F (geom_input(rtdb)) 59 ofinite=.false. 60 status=rtdb_get(rtdb,'prop:ofinite',mt_log,1,ofinite) ! Check if gaussian nucl model requested 61c ------ Read ofinite to be used by HFine finite calc ---FA-03-21-11-- END 62 if (ga_nodeid().eq.0) then 63 write(*,*) 'dft_zora_Hypefine: ofinite=',ofinite 64 endif 65c ---- get spin-densty: g_sdens -------- START 66 if (.not. ga_create(mt_dbl,nbf,nbf, 67 & 'NMRHFine: g_sdens', 68 $ 0,0,g_sdens)) 69 $ call errquit('NMRHFine: g_sdens', 0, 70 & GA_ERR) 71 call ga_add( 1.0d0,g_densZ4(1), 72 & -1.0d0,g_densZ4(2),g_sdens) 73c ---- get spin-densty: g_sdens -------- END 74c === allocate arrays to store diamagnetic tensor === START 75c ------- Read (nat,atmnr) --------- START 76 status=geom_ncent(geom,nat) 77 if (.not.ma_alloc_get( 78 & mt_int,nat,'nmt tmp',l_AtNr,k_AtNr)) 79 & call errquit('dft_zora_Hyperfine: ma failed',0,MA_ERR) 80 typeprop=3 ! =1 EFG =2 Shieldings =3 Hyperfine 81 call get_slctd_atoms(nat_slc, ! out: selected atoms 82 & int_mb(k_AtNr),! out: list of selected atom nr. 83 & nat, ! in : total nr atoms in molecule 84 & rtdb, ! in : rtdb handle 85 & typeprop) ! in : =1,2,3=EFG,Shieldings,Hyperfine 86 if (.not. ga_create(mt_dbl,1,nat_slc, 87 & 'dft_zora_Hyperfine: g_AtNr',0,0,g_AtNr)) 88 $ call errquit('dft_zora_Hyperfine: g_AtNr', 0,GA_ERR) 89 90 do i=1,nat_slc 91 AtNr_dbl=int_mb(k_AtNr+i-1) 92 call ga_put(g_AtNr,1,1,i,i,AtNr_dbl,1) 93 enddo 94 if (ga_nodeid().eq.0) then 95 write(*,*) 'nat_slc=',nat_slc 96 do i=1,nat_slc 97 write(*,7) i,int_mb(k_AtNr+i-1) 98 7 format('In dft_zora_Hyperfine:: atomnr(',i3,')=',i5) 99 enddo 100 endif 101c ------- Read (nat,atmnr) --------- END 102 call get_Htensor_fast( 103 & ga_dia_hfine, ! OUT: dia hyperfine tensor 104 & ga_h01_hfine, ! OUT: h01 hyperfine munu AO matrix 105 & g_sdens, ! IN : spin density 106 & ofinite, ! IN : = .true. if requesting Gaussian Nucl. Model for charge 107 & Knucl, ! in : = .true. for including ONLY nuclear part in K ZORA 108 & nat,nat_slc,int_mb(k_AtNr), 109 & rtdb,g_dens_at,nexc, 110 & geom,ao_bas_han,nbf, 111 & noc,ipol) 112c ---Destroying ga arrays ------- START 113 if (.not. ga_destroy(g_sdens)) call errquit( 114 & 'get_NMRHFine_ZORA: ga_destroy failed ',0, GA_ERR) 115 call get_HFine_F1ji( 116 & ga_Fji_hfine, ! OUT: munu-mat-Fji 117 & rtdb,g_dens_at, 118 & nat, ! in: nr. atoms 119 & ofinite, ! in: = .true. if Gaussian Nucl. Model of charges requested 120 & Knucl, ! in: = .true. if K_ZORA(V=Nuclear pot. only) 121 & nexc, 122 & geom,ao_bas_han,nbf) 123c == get filename for the zora data == 124 type_nmrdata=2 ! =1,2,3=shieldings,hyperfine,gshift 125c Note.- lbl_nmrhyp defined in zora.fh 126 call util_file_name(lbl_nmrhyp,.false.,.false.,zorafilename) 127 if (.not.dft_zoraNMR_write( 128 & zorafilename, 129 & type_nmrdata, ! =1,2,3=shieldings,hyperfine,gshift 130 & nbf, 131 & nat_slc, 132 & g_AtNr, 133 & ga_dia_hfine, 134 & ga_para1, 135 & ga_h01_hfine, 136 & ga_Fji_hfine)) 137 & call errquit('get_NMRHFine_ZORA: dft_zoraNMR_write failed', 138 & 0,DISK_ERR) 139c ---- Destroy stored ga arrays ------ START 140 if (.not. ga_destroy(ga_dia_hfine)) call errquit( 141 & 'get_NMRHFine_ZORA: ga_destroy failed ',0, GA_ERR) 142 if (.not. ga_destroy(ga_h01_hfine)) call errquit( 143 & 'get_NMRHFine_ZORA: ga_destroy failed ',0, GA_ERR) 144 if (.not. ga_destroy(ga_Fji_hfine)) call errquit( 145 & 'get_NMRHFine_ZORA: ga_destroy failed ',0, GA_ERR) 146c if (.not. ga_destroy(g_AtNr)) call errquit( 147c & 'get_NMRHFine_ZORA: ga_destroy failed ',0, GA_ERR) 148c ---- Destroy stored ga arrays ------ END 149 if (.not.ma_free_heap(l_AtNr)) call 150 & errquit('get_NMRHFine_ZORA: ma_free_heap l_AtNr',0,MA_ERR) 151 return 152 end 153 154 subroutine get_Htensor_fast( 155 & gFCSD, ! OUT : dia hyperfine tensor 156 & gPSOSO, ! OUT : h01 para hyperfine munu AO matrix 157 & g_sdens, ! IN : spin density 158 & ofinite, ! IN : = .true. if requesting Gaussian Nucl. Model for charge 159 & Knucl, ! IN : = .true. for including ONLY nuclear part in K ZORA 160 & nat, ! IN : total nr of atoms in molecule 161 & nat_slc, ! IN : nr of selected atoms 162 & atmnr_slc, ! IN : selected atom numbers 163 & rtdb, 164 & g_dens_at, 165 & nexc, 166 & geom, ! IN : geom handle 167 & ao_bas_han, ! IN : basis handle 168 & nbf, ! IN : nr. basis functions 169 & noc, ! IN : nr. occupied MOs 170 & ipol) ! IN : nr. of polarizations 171c Purpose : Compute 172c 1. A_{u,v}^{FC+SD} = 2 g_N B_N /(nA-nB) 173c \sum_{r,s} P_{r,s}^{A-B} < chi_r| h^{(u,v)}| chi_s> 174c GA array gFCSD (dimension 3 x 3 x nlist (nlist, nr of atoms selected) 175c 2. h_{u}^{PSOSO}= < chi_r | h_{u}^{PSOSO}|chi_s 176c h_{u}^{PSOSO}= 2c^2 i . h_{Au,rs}^{ZPSO}= 177c \int dr K/r_A^3 (\vec{r}_A x[chi_r^* \nabla chi_s-\nabla chi_r^* chi_s])_u 178c (Eq. 56 in JA's write-up on Sept 17,2007) 179c 180c Author : Fredy Aquino 181c Date : 03-10-11 182 implicit none 183#include "errquit.fh" 184#include "mafdecls.fh" 185#include "stdio.fh" 186#include "global.fh" 187#include "msgids.fh" 188#include "rtdb.fh" 189#include "geom.fh" 190#include "zora.fh" 191 integer g_sdens,vectors(2) 192 integer gFCSD, ! out 193 & gPSOSO ! out 194 integer rtdb 195 integer g_dens_at(2),g_densZ4(3) 196 integer noc(2),noc1,ipol 197 integer geom,ao_bas_han 198 integer ispin,nexc,iat,iat1, 199 & nat,nat_slc,typeprop 200 integer l_xyzpt,k_xyzpt, 201 & l_zanpt,k_zanpt 202 integer atmnr_slc(nat_slc) 203 logical status,ofinite,Knucl 204 character*16 at_tag 205 integer stat_read 206 integer alo(3),ahi(3),ld(2) 207 integer u,v,i,j,k,t,a,nbf 208 integer dims(3),chunk(3) 209 double precision val,factor 210 double precision xyz_NMRQcoords(3),atmass, 211 & zetanuc_slc 212 character*255 zorafilename 213 integer g_fcsd(3,3),g_scr,g_t1,g_t2,g_t3 214 integer g_zpso(6) 215 integer l_buf,k_buf,cbuf, 216 & l_buf1,k_buf1, 217 & l_zetanuc,k_zetanuc 218c +++++++++ definitions for NLMO analysis ++++++++ START 219 integer ndir,ndir1,n_munu,count,count1,hypfile, 220 & g_munuFCSD,g_munuPSOSO, ! GA contains matrices 221 & g_c1,ndata 222 logical dft_zoraHYP_NLMOAnalysis_write 223c +++++++++ definitions for NLMO analysis ++++++++ END 224 external util_wallsec 225 double precision util_wallsec 226 if(.not.ma_alloc_get(mt_dbl,3*3,'HTensor:buf', 227 & l_buf,k_buf)) 228 & call errquit('HTensor: ma failed',911,MA_ERR) 229 if(.not.ma_alloc_get(MT_DBL,nbf*nbf,'HTensor', 230 & l_buf1,k_buf1)) 231 $ call errquit('HTensor: ma failed',911, MA_ERR) 232c +++++++++++++++creating ga_arrays ++++++++START 233 if (.not. ga_create(mt_dbl,nbf,nbf, 234 & 'HFine: g_t1',0,0,g_t1)) 235 $ call errquit('HFine: g_t1',0,GA_ERR) 236 call ga_zero(g_t1) 237 if (.not. ga_create(mt_dbl,nbf,nbf, 238 & 'HFine: g_t2',0,0,g_t2)) 239 $ call errquit('HFine: g_t2',0,GA_ERR) 240 call ga_zero(g_t2) 241 if (.not. ga_create(mt_dbl,nbf,nbf, 242 & 'HFine: g_t3',0,0,g_t3)) 243 $ call errquit('HFine: g_t3',0,GA_ERR) 244 call ga_zero(g_t3) 245 if (.not. ga_create(mt_dbl,nbf,nbf, 246 & 'HFine: g_scr',0,0,g_scr)) 247 $ call errquit('HFine: g_scr',0,GA_ERR) 248 call ga_zero(g_scr) 249 do i=1,6 250 if (.not. ga_create(mt_dbl,nbf,nbf, 251 & 'HFine: g_zpso',0,0,g_zpso(i))) 252 $ call errquit('HFine: g_zpso',0,GA_ERR) 253 call ga_zero(g_zpso(i)) 254 enddo 255 do i=1,3 256 do j=1,3 257 if (.not. ga_create(mt_dbl,nbf,nbf, 258 & 'HFine: g_fcsd',0,0,g_fcsd(i,j))) 259 $ call errquit('HFine: g_fcsd',0,GA_ERR) 260 call ga_zero(g_fcsd(i,j)) 261 enddo ! end-loop-j 262 enddo ! end-loop-i 263c +++++++++++++++creating ga_arrays ++++++++END 264c +++++ Read Atom Nr for NMR calc ++ 265c----- Allocate memory - FA 266 if (.not. ma_alloc_get(mt_dbl,3*nat_slc, 267 & 'xyz pnt',l_xyzpt,k_xyzpt)) 268 & call errquit('HFine: ma failed',911,MA_ERR) 269 if (.not. ma_alloc_get(mt_dbl,nat_slc, 270 & 'zan pnt',l_zanpt,k_zanpt)) 271 & call errquit('HFine: ma failed',911,MA_ERR) 272c === allocate arrays to store diamagnetic tensor === START 273 alo(1) = nbf 274 alo(2) = -1 275 alo(3) = -1 276 ahi(1) = nbf 277 ahi(2) = nbf 278 ahi(3) = 3*nat_slc 279 if (.not.nga_create(MT_DBL,3,ahi,'H01 matrix',alo,gPSOSO)) call 280 & errquit('get_d2p1: nga_create failed gPSOSO_num',0,GA_ERR) 281 call ga_zero(gPSOSO) 282 alo(1) = 3 283 alo(2) = -1 284 alo(3) = -1 285 ahi(1) = 3 286 ahi(2) = 3 287 ahi(3) = nat_slc ! Total nr. atoms requested 288 if (.not.nga_create(MT_DBL,3,ahi,'gFCSD matrix', 289 & alo,gFCSD)) 290 & call errquit('HFine: nga_create failed gFCSD all', 291 & 0,GA_ERR) 292c === allocate arrays to store diamagnetic tensor === END 293 if (ofinite) then 294 if (.not.ma_alloc_get(mt_dbl,nat, 295 & 'zetanuc1',l_zetanuc,k_zetanuc)) 296 & call errquit('HFine: ma failed',0,MA_ERR) 297 call get_zetanuc_arr(geom,nat,dbl_mb(k_zetanuc)) ! zetanuc_arr(i) i=1,natoms 298 do iat = 1,nat ! == loop over the atoms == 299 if (ga_nodeid().eq.0) 300 & write(*,1) iat,dbl_mb(k_zetanuc+iat-1) 301 1 format('In get_Htensor_fast:: zetanuc_arr(',i3,')=',f35.8) 302 dbl_mb(k_zetanuc+iat-1)=dsqrt(dbl_mb(k_zetanuc+iat-1)) ! Calc sqrt(zetanuc) 303 enddo ! end-loop-iat 304 endif 305c +++++++++ NLMO analysis : create diag matrix ++++++++ START 306 hypfile=0 ! not doing NLMO analysis by default 307 status=rtdb_get(rtdb,'prop:hypfile',mt_int,1,hypfile) ! for NLMO analysis 308 if (hypfile.eq.1) then ! --- if-hypfile--START 309 ndir=6 ! xx,yy,zz,xy,xz,yz 310 n_munu=nbf*(nbf+1)/2*ndir*nat_slc 311 if (.not. ga_create(mt_dbl,1,n_munu, 312 & 'get_Htensor_fast: g_munu',0,0,g_munuFCSD)) 313 $ call errquit('get_Htensor_fast:',0,GA_ERR) 314 call ga_zero(g_munuFCSD) 315 ndir1=3 ! x,y,z 316 n_munu=nbf*(nbf+1)/2*ndir1*nat_slc 317 if (.not. ga_create(mt_dbl,1,n_munu, 318 & 'get_Htensor_fast: g_munu',0,0,g_munuPSOSO)) 319 $ call errquit('get_Htensor_fast:',0,GA_ERR) 320 call ga_zero(g_munuPSOSO) 321 endif ! ------------if-hypfile---------END 322 count=1 323 count1=1 ! for storing g_munu_PSOSO 324c +++++++++ NLMO analysis : create diag matrix ++++++++ END 325 do iat1=1,nat_slc ! nat_slc <= nat 326 iat=atmnr_slc(iat1) 327 status=geom_cent_get(geom,iat,at_tag, 328 & dbl_mb(k_xyzpt+3*(iat1-1)), 329 & dbl_mb(k_zanpt+iat1-1)) 330 if(.not.geom_mass_get(geom, iat, atmass)) call 331 & errquit(' mass_get failed ',iat, GEOM_ERR) 332 call get_znuc(atmass,zetanuc_slc) 333 xyz_NMRQcoords(1)= dbl_mb(k_xyzpt +3*(iat1-1)) 334 xyz_NMRQcoords(2)= dbl_mb(k_xyzpt+1+3*(iat1-1)) 335 xyz_NMRQcoords(3)= dbl_mb(k_xyzpt+2+3*(iat1-1)) 336 337 if (ga_nodeid().eq.0) then 338 write(luout,'(a,f10.2,a)') ' starting at time ', 339 & util_wallsec(),'s' 340 write(*,153) iat,ofinite,atmass,zetanuc_slc 341 153 format('CHECK:(atom,ofinite,atmass,zetanucl_slc)=(', 342 & i4,',',l1,',',f15.8,',',f35.8,')') 343 endif 344 345 call zora_getv_HFine_fast(rtdb,g_dens_at, 346 & ofinite, 347 & dbl_mb(k_zetanuc), 348 & zetanuc_slc, 349 & Knucl, 350 & xyz_NMRQcoords, 351 & g_zpso, ! ZPSO 352 & g_fcsd, ! FC+SD (v,u) term 353 & nexc) 354c ---- prepare g_fcsd-final --------- START 355 call ga_zero(g_scr) 356 do u=1,3 357 call ga_add(1.0d0,g_scr,1.0d0,g_fcsd(u,u),g_scr) 358 enddo 359 do u=1,3 360 call ga_add(-1.0d0,g_scr,+1.0d0,g_fcsd(u,u), 361 & g_fcsd(u,u)) 362 enddo 363c ---- prepare g_fcsd-final --------- END 364 do u=1,3 365c ---- 2nd method (fastest) to calculate g_zpso ---- START 366c ----- g_h01 = < chi_{mu} | h_t^{01}| chi_{nu} > 367c h_t^{01}=K/(2c) (\vec{r} x \vec{p})_t/r_Q^3 + 368c (\vec{r} x \vec{p})_t/r_Q^3 K/(2c) 369 call ga_add(1.0d0,g_zpso(u) , 370 & -1.0d0,g_zpso(u+3),g_t1) ! i.e. u=1 mn-nm=23-32, etc. 371 call ga_transpose(g_t1,g_t2) 372 call ga_zero(g_t3) 373 call ga_add(1.0d0,g_t2,-1.0d0,g_t1,g_t3) ! g_t2=g_h01 374c ---- 2nd method (fastest) to calculate g_zpso ---- END 375 alo(1)=1 376 ahi(1)=nbf 377 alo(2)=1 378 ahi(2)=nbf 379 alo(3)=3*(iat1-1)+u 380 ahi(3)=3*(iat1-1)+u 381 ld(1)=nbf 382 ld(2)=nbf 383 call ga_scale(g_t3,-0.5d0) ! including (-1/2) from h^{u} operator 384c ------------ NLMO analysis para term (PSOSO) ------------- START 385 hypfile=0 ! not doing NLMO analysis by default 386 status=rtdb_get(rtdb,'prop:hypfile',mt_int,1,hypfile) ! for NLMO analysis 387 if (hypfile.eq.1) then 388 call fill_munuPSOSO(g_munuPSOSO, !out : array with matrices (g_t3--> g_munuPSOSO) 389 & count1, !in/out: counting data stored in g_munuPSOSO 390 & g_t3, ! in: nbf 391 & nbf) 392 endif 393c ------------ NLMO analysis para term (PSOSO) ------------- END 394 call ga_get(g_t3,1,nbf,1,nbf,dbl_mb(k_buf1),nbf) 395 call nga_put(gPSOSO,alo,ahi,dbl_mb(k_buf1),ld) ! store gPSOSO_u 396 do v=1,3 397 val=ga_ddot(g_sdens,g_fcsd(u,v)) 398 cbuf=k_buf+(u-1)*3+v-1 399 dbl_mb(cbuf)=val ! missing :2 g_N beta_N /(n_a-n_b) * (1/2) 400 ! Note.- (1/2) is from h^{(u,v)} operator 401 if (ga_nodeid().eq.0) then 402 write(*,10) iat,u,v,val 40310 format('gFCSD(',i3,',',i3,',',i3,')=',f15.8) 404 endif 405 enddo ! end-loop-v 406 enddo ! end-loop-u 407 alo(1)=1 408 ahi(1)=3 409 alo(2)=1 410 ahi(2)=3 411 alo(3)=iat1 412 ahi(3)=iat1 413 ld(1)=3 414 ld(2)=3 415 call nga_put(gFCSD,alo,ahi,dbl_mb(k_buf),ld) 416c +++++++++ NLMO analysis : store diag matrix ++++++++ START 417 hypfile=0 ! not doing NLMO analysis by default 418 status=rtdb_get(rtdb,'prop:hypfile',mt_int,1,hypfile) ! for NLMO analysis 419 if (hypfile.eq.1) then 420c ----- Symmetrize g_fcsd() ------ START 421 do u=1,3 422 do v=u+1,3 423 call ga_add(0.5d0,g_fcsd(u,v), 424 & 0.5d0,g_fcsd(v,u),g_t1) 425 call ga_copy(g_t1,g_fcsd(u,v)) 426 call ga_copy(g_t1,g_fcsd(v,u)) 427 enddo 428 enddo 429c ----- Symmetrize g_fcsd() ------ END 430 call fill_munuFCSD(g_munuFCSD, !out: array with matrices (g_fcsd--> g_munuFCSD) 431 & count, !in/out: counting data stored in g_munuFCSD 432 & g_fcsd, ! in: nbf 433 & nbf) 434 endif 435c +++++++++ NLMO analysis : store diag matrix ++++++++ END 436 end do ! iat loop 437c ---Destroying ga arrays ------- START 438 if (.not. ga_destroy(g_t1)) call errquit( 439 & 'HFine: ga_destroy failed ',0, GA_ERR) 440 if (.not. ga_destroy(g_t2)) call errquit( 441 & 'HFine: ga_destroy failed ',0, GA_ERR) 442 if (.not. ga_destroy(g_t3)) call errquit( 443 & 'HFine: ga_destroy failed ',0, GA_ERR) 444 if (.not. ga_destroy(g_scr)) call errquit( 445 & 'HFine: ga_destroy failed ',0, GA_ERR) 446 do i=1,6 447 if (.not. ga_destroy(g_zpso(i))) call errquit( 448 & 'HFine: ga_destroy failed ',0, GA_ERR) 449 enddo 450 do i=1,3 451 do j=1,3 452 if (.not. ga_destroy(g_fcsd(i,j))) call errquit( 453 & 'HFine: ga_destroy failed ',0, GA_ERR) 454 enddo 455 enddo 456 hypfile=0 ! not doing NLMO analysis by default 457 status=rtdb_get(rtdb,'prop:hypfile',mt_int,1,hypfile) ! for NLMO analysis 458 if (hypfile.eq.1) then 459 ndata=1 ! =1 write FCSD,PSOSO,sdens =2 write g_c1 460 call util_file_name(lbl_nlmohyp,.false.,.false.,zorafilename) 461 if (.not.dft_zoraHYP_NLMOAnalysis_write( 462 & zorafilename, ! in: filename 463 & nbf, ! in: nr basis functions 464 & ndir, ! in: nr of directions: 6 = xx yy zz xy xz yz for g_munuFCSD 465 & ndir1, ! in: nr of directions: 3 = x y z for g_munuPSOSO 466 & nat_slc, ! in: list of selected atoms 467 & noc, ! in: dummy not used yet here 468 & ndata, ! in: =1 write FCSD,PSOSO,sdens =2 write g_c1 469 & g_munuFCSD, ! in: munu dia or Fermi Contact + Spin Dipolar term 470 & g_munuPSOSO, ! in: munu para or PSOSO term 471 & ipol, ! in: nr. of polarizations 472 & vectors, ! in: dummy not used yet here 473 & g_c1, ! in: dummy not used yet here 474 & g_sdens))! in: spin density 475 & call errquit('get_Htensor_fast: dft_zoraHYPNLMO_write failed', 476 & 0,DISK_ERR) 477 if (.not. ga_destroy(g_munuFCSD)) call errquit( 478 & 'HFine: ga_destroy failed ',0, GA_ERR) 479 if (.not. ga_destroy(g_munuPSOSO)) call errquit( 480 & 'HFine: ga_destroy failed ',0, GA_ERR) 481 endif 482c ---Destroying ga arrays ------- END 483c----deallocate memory 484 if (.not.ma_free_heap(l_buf1)) call errquit 485 & ('HFine, ma_free_heap of l_buf failed', 486 & 911,MA_ERR) 487 if (.not.ma_free_heap(l_buf)) call errquit 488 & ('HFine, ma_free_heap of l_buf failed', 489 & 911,MA_ERR) 490 if (.not.ma_free_heap(l_zanpt)) call errquit 491 & ('HFine:, ma_free_heap of l_zanpt failed',911,MA_ERR) 492 if (.not.ma_free_heap(l_xyzpt)) call errquit 493 & ('HFine:, ma_free_heap of l_xyzpt failed',911,MA_ERR) 494 if (ofinite) then 495 if (.not.ma_free_heap(l_zetanuc)) call 496 & errquit('HFine:: ma_free_heap l_zetanuc',0, MA_ERR) 497 endif 498 return 499 end 500 501 subroutine get_Htensor_slow( 502 & gFCSD, ! OUTPUT : Fermi Contact - Spin Dipolar 503 & gPSOSO, ! OUTPUT : 504 & g_sdens, ! IN : spin density 505 & ofinite, ! IN : = .true. if requesting Gaussian Nucl. Model for charge 506 & rtdb,g_dens_at,nexc, 507 & geom, 508 & ao_bas_han, 509 & nbf, 510 & noc, 511 & ipol) 512c Purpose : Compute 513c 1. A_{u,v}^{FC+SD} = 2 g_N B_N /(nA-nB) 514c \sum_{r,s} P_{r,s}^{A-B} < chi_r| h^{(u,v)}| chi_s> 515c GA array gFCSD (dimension 3 x 3 x nlist (nlist, nr of atoms selected) 516c 2. h_{u}^{PSOSO}= < chi_r | h_{u}^{PSOSO}|chi_s 517c h_{u}^{PSOSO}= 2c^2 i . h_{Au,rs}^{ZPSO}= 518c \int dr K/r_A^3 (\vec{r}_A x[chi_r^* \nabla chi_s-\nabla chi_r^* chi_s])_u 519c (Eq. 56 in JA's write-up on Sept 17,2007) 520c 521c Author : Fredy Aquino 522c Date : 11-06-10 523 implicit none 524#include "errquit.fh" 525#include "mafdecls.fh" 526#include "stdio.fh" 527#include "global.fh" 528#include "msgids.fh" 529#include "rtdb.fh" 530#include "geom.fh" 531#include "zora.fh" 532 integer g_sdens 533 integer gFCSD, ! out 534 & gPSOSO ! out 535 integer rtdb 536 integer g_dens_at(2),g_densZ4(3) 537 integer noc(2),noc1,ipol 538 integer geom,ao_bas_han 539 integer ispin,nexc,iat,iat1,nat 540 integer l_xyzpt,k_xyzpt, 541 & l_zanpt,k_zanpt 542 integer l_AtNr,k_AtNr 543 logical status,ofinite 544 character*16 at_tag 545 integer stat_read,read_SLCTD_HFine_Atoms 546 integer alo(3),ahi(3),ld(2) 547 integer u,v,i,j,k,t,a,nbf 548 integer dims(3),chunk(3) 549 double precision val,factor 550 double precision xyz_NMRQcoords(3),atmass 551 integer g_fcsd(3,3),g_scr,g_t1,g_t2,g_t3 552 integer g_zpso(3) 553 integer l_buf,k_buf,cbuf, 554 & l_buf1,k_buf1 555 external zora_getv_HFine, 556 & read_SLCTD_HFine_Atoms, 557 & zora_getv_HFine1 558 if(.not.ma_alloc_get(mt_dbl,3*3,'HTensor:buf', 559 & l_buf,k_buf)) 560 & call errquit('HTensor: ma failed',911,MA_ERR) 561 if(.not.ma_alloc_get(MT_DBL,nbf*nbf,'HTensor', 562 & l_buf1,k_buf1)) 563 $ call errquit('HTensor: ma failed',911, MA_ERR) 564 status=geom_ncent(geom,nat) ! Get nat, # of atoms 565c----- Allocate memory - FA 566 if (.not. ma_alloc_get(mt_dbl,3*nat,'xyz pnt',l_xyzpt,k_xyzpt)) 567 & call errquit('HFine: ma failed',911,MA_ERR) 568 if (.not. ma_alloc_get(mt_dbl,nat,'zan pnt',l_zanpt,k_zanpt)) 569 & call errquit('HFine: ma failed',911,MA_ERR) 570c +++++++++++++++creating ga_arrays ++++++++START 571 if (.not. ga_create(mt_dbl,nbf,nbf, 572 & 'HFine: g_scr',0,0,g_scr)) 573 $ call errquit('HFine: g_scr',0,GA_ERR) 574 call ga_zero(g_scr) 575 do i=1,3 576 if (.not. ga_create(mt_dbl,nbf,nbf, 577 & 'HFine: g_zpso',0,0,g_zpso(i))) 578 $ call errquit('HFine: g_zpso',0,GA_ERR) 579 call ga_zero(g_zpso(i)) 580 enddo 581 do i=1,3 582 do j=1,3 583 if (.not. ga_create(mt_dbl,nbf,nbf, 584 & 'HFine: g_fcsd',0,0,g_fcsd(i,j))) 585 $ call errquit('HFine: g_fcsd',0,GA_ERR) 586 call ga_zero(g_fcsd(i,j)) 587 enddo ! end-loop-j 588 enddo ! end-loop-i 589c +++++++++++++++creating ga_arrays ++++++++END 590c +++++ Read Atom Nr for NMR calc ++ 591 if (.not. ga_create(mt_dbl,1,nat, 592 & 'HFine: g_AtNr',0,0,g_AtNr)) 593 $ call errquit('HFine: g_AtNr',0,GA_ERR) 594 call ga_zero(g_AtNr) 595 stat_read=read_SLCTD_HFine_Atoms(rtdb,nat,nlist,g_AtNr) 596c === allocate arrays to store diamagnetic tensor === START 597 alo(1) = nbf 598 alo(2) = -1 599 alo(3) = -1 600 ahi(1) = nbf 601 ahi(2) = nbf 602 ahi(3) = 3*nlist 603 if (.not.nga_create(MT_DBL,3,ahi,'H01 matrix',alo,gPSOSO)) call 604 & errquit('get_d2p1: nga_create failed gPSOSO_num',0,GA_ERR) 605 call ga_zero(gPSOSO) 606 alo(1) = 3 607 alo(2) = -1 608 alo(3) = -1 609 ahi(1) = 3 610 ahi(2) = 3 611 ahi(3) = nlist ! Total nr. atoms requested 612 if (.not.nga_create(MT_DBL,3,ahi,'gFCSD matrix', 613 & alo,gFCSD)) 614 & call errquit('HFine: nga_create failed gFCSD all', 615 & 0,GA_ERR) 616c === allocate arrays to store diamagnetic tensor === END 617c Allocate memory for l_AtNr,k_AtNr 618 if (.not.ma_alloc_get(mt_dbl,nat, 619 & 'AtNr',l_AtNr,k_AtNr)) 620 & call errquit('HFine: ma failed',0,MA_ERR) 621 call ga_get(g_AtNr,1,1,1,nat,dbl_mb(k_AtNr),1) 622 do iat1=1,nlist ! nlist <= nat 623 iat=dbl_mb(k_AtNr+iat1-1) 624 status=geom_cent_get(geom,iat,at_tag, 625 & dbl_mb(k_xyzpt+3*(iat-1)), 626 & dbl_mb(k_zanpt+iat-1)) 627 if(.not.geom_mass_get(geom, iat, atmass)) call 628 & errquit(' mass_get failed ',iat, GEOM_ERR) 629 xyz_NMRQcoords(1)= dbl_mb(k_xyzpt +3*(iat-1)) 630 xyz_NMRQcoords(2)= dbl_mb(k_xyzpt+1+3*(iat-1)) 631 xyz_NMRQcoords(3)= dbl_mb(k_xyzpt+2+3*(iat-1)) 632 633 if (ga_nodeid().eq.0) then 634 write(*,153) iat,ofinite,atmass 635 153 format('CHECK:(atom,ofinite,atmass)=(', 636 & i4,',',l1,',',f15.8,')') 637 endif 638 639 call zora_getv_HFine_slow(rtdb,g_dens_at, 640 & ofinite, 641 & atmass, 642 & xyz_NMRQcoords, 643 & g_zpso, ! ZPSO 644 & g_fcsd, ! FC+SD (v,u) term 645 & nexc) 646c ---- prepare g_fcsd-final --------- START 647 call ga_zero(g_scr) 648 do u=1,3 649 call ga_add(1.0d0,g_scr,1.0d0,g_fcsd(u,u),g_scr) 650 enddo 651 do u=1,3 652 call ga_add(-1.0d0,g_scr,+1.0d0,g_fcsd(u,u), 653 & g_fcsd(u,u)) 654 enddo 655c ---- prepare g_fcsd-final --------- END 656 do u=1,3 657 alo(1)=1 658 ahi(1)=nbf 659 alo(2)=1 660 ahi(2)=nbf 661 alo(3)=3*(iat1-1)+u 662 ahi(3)=3*(iat1-1)+u 663 ld(1)=nbf 664 ld(2)=nbf 665 call ga_scale(g_zpso(u),0.5d0) 666 call ga_get(g_zpso(u),1,nbf,1,nbf,dbl_mb(k_buf1),nbf) 667 call nga_put(gPSOSO,alo,ahi,dbl_mb(k_buf1),ld) ! store g_h01 668 do v=1,3 669 val=ga_ddot(g_sdens,g_fcsd(u,v)) 670 cbuf=k_buf+(u-1)*3+v-1 671 dbl_mb(cbuf)=val 672 if (ga_nodeid().eq.0) then 673 write(*,10) iat,u,v,val 67410 format('gFCSD(',i3,',',i3,',',i3,')=',f15.8) 675 endif 676 enddo ! end-loop-v 677 enddo ! end-loop-u 678 679 goto 17 680 681 do j=1,3 682 if (ga_nodeid().eq.0) 683 & write(*,1) iat,j 684 1 format('-----ZPSO munu-mat(',i3,',',i3,')----- START') 685 call ga_print(g_zpso(j)) 686 if (ga_nodeid().eq.0) 687 & write(*,2) iat,j 688 2 format('-----ZPSO munu-mat(',i3,',',i3,')----- END') 689 enddo ! end-loop-j (directions xyz) 690 do j=1,3 691 do i=1,3 692 if (ga_nodeid().eq.0) 693 & write(*,3) iat,j,i 694 3 format('-----FC+SD munu-mat(',i3,',',i3,',',i3,')----- START') 695 call ga_print(g_fcsd(j,i)) 696 if (ga_nodeid().eq.0) 697 & write(*,4) iat,j,i 698 4 format('-----FC+SD munu-mat(',i3,',',i3,',',i3,')-----END') 699 enddo 700 enddo 701 702 17 continue 703 704 alo(1)=1 705 ahi(1)=3 706 alo(2)=1 707 ahi(2)=3 708 alo(3)=iat1 709 ahi(3)=iat1 710 ld(1)=3 711 ld(2)=3 712 call nga_put(gFCSD,alo,ahi,dbl_mb(k_buf),ld) 713 714 end do ! iat loop 715c if (ga_nodeid().eq.0) 716c & write(*,*) '------- gFC+SD tensor ------ START' 717c call ga_print(gFCSD) 718c if (ga_nodeid().eq.0) 719c & write(*,*) '------- gFC+SD tensor ------ END' 720c if (ga_nodeid().eq.0) 721c & write(*,*) '------- gPSOSO tensor ------ START' 722c call ga_print(gPSOSO) 723c if (ga_nodeid().eq.0) 724c & write(*,*) '------- gPSOSO tensor ------ END' 725 726c ---Destroying ga arrays ------- START 727 if (.not. ga_destroy(g_scr)) call errquit( 728 & 'HFine: ga_destroy failed ',0, GA_ERR) 729 do i=1,3 730 if (.not. ga_destroy(g_zpso(i))) call errquit( 731 & 'HFine: ga_destroy failed ',0, GA_ERR) 732 enddo 733 do i=1,3 734 do j=1,3 735 if (.not. ga_destroy(g_fcsd(i,j))) call errquit( 736 & 'HFine: ga_destroy failed ',0, GA_ERR) 737 enddo 738 enddo 739c ---Destroying ga arrays ------- END 740c----deallocate memory 741 if (.not.ma_free_heap(l_buf1)) call errquit 742 & ('HFine, ma_free_heap of l_buf failed', 743 & 911,MA_ERR) 744 if (.not.ma_free_heap(l_buf)) call errquit 745 & ('HFine, ma_free_heap of l_buf failed', 746 & 911,MA_ERR) 747 if (.not.ma_free_heap(l_zanpt)) call errquit 748 & ('HFine:, ma_free_heap of l_zanpt failed',911,MA_ERR) 749 if (.not.ma_free_heap(l_xyzpt)) call errquit 750 & ('HFine:, ma_free_heap of l_xyzpt failed',911,MA_ERR) 751 if (.not.ma_free_heap(l_AtNr)) call 752 & errquit('HFine:: ma_free_heap l_AtNr',0, MA_ERR) 753 return 754 end 755 756 subroutine get_HFine_F1ji( 757 & g_Fji, !out: munu-mat-Fji 758 & rtdb,g_dens_at, 759 & nat, ! in: nr. atoms 760 & ofinite, ! in: = .true. if Gaussian Nucl. Model of charges requested 761 & Knucl, ! in: = .true. if K_ZORA(V=Nuclear pot. only) 762 & nexc, 763 & geom,ao_bas_han,nbf) 764 implicit none 765#include "errquit.fh" 766#include "mafdecls.fh" 767#include "stdio.fh" 768#include "global.fh" 769#include "msgids.fh" 770#include "rtdb.fh" 771#include "geom.fh" 772 integer g_Fji ! OUTPUT 773 logical status 774 integer rtdb,g_dens_at(2),geom,ao_bas_han 775 integer nexc 776 integer i,j,k 777 integer nbf 778 integer g_hfine(3) 779 integer l_buf,k_buf 780 integer alo(3),ahi(3),ld(2) 781 integer iat,nat 782 integer l_zetanuc,k_zetanuc 783 logical ofinite,Knucl 784 double precision zetanuc_arr(nat) 785 external zora_getv_NMRHFine_F1ji, 786 & get_zetanuc_arr 787 788 if (ofinite) then 789 if (.not.ma_alloc_get(mt_dbl,nat, 790 & 'zetanuc1',l_zetanuc,k_zetanuc)) 791 & call errquit('HFine: ma failed',0,MA_ERR) 792 call get_zetanuc_arr(geom,nat,dbl_mb(k_zetanuc)) ! zetanuc_arr(i) i=1,natoms 793 do iat = 1,nat ! == loop over the atoms == 794 if (ga_nodeid().eq.0) 795 & write(*,1) iat,dbl_mb(k_zetanuc+iat-1) 796 1 format('In get_Htensor_fast:: zetanuc_arr(',i3,')=',f35.8) 797 dbl_mb(k_zetanuc+iat-1)=dsqrt(dbl_mb(k_zetanuc+iat-1)) ! Calc sqrt(zetanuc) 798 enddo ! end-loop-iat 799 endif 800c +++++++++++++++creating ga_arrays ++++++++START 801 do i=1,3 802 if (.not. ga_create(mt_dbl,nbf,nbf, 803 & 'gFji: g_hfine',0,0,g_hfine(i))) 804 $ call errquit('gFji: g_hfine',0,GA_ERR) 805 call ga_zero(g_hfine(i)) 806 enddo ! end-loop-i 807c +++++++++++++++creating ga_arrays ++++++++END 808 call zora_getv_NMRHFine_F1ji(rtdb,g_dens_at, 809 & g_hfine, 810 & nat, 811 & ofinite, 812 & dbl_mb(k_zetanuc), 813 & Knucl, 814 & nexc) 815 if(.not.ma_alloc_get(mt_dbl,nbf*nbf,'gFji:buf', 816 & l_buf,k_buf)) 817 $ call errquit('gFji: ma failed',911, MA_ERR) 818 alo(1) = nbf 819 alo(2) = -1 820 alo(3) = -1 821 ahi(1) = nbf 822 ahi(2) = nbf 823 ahi(3) = 3 824 if (.not.nga_create(MT_DBL,3,ahi,'Fji matrix',alo,g_Fji)) 825 & call 826 & errquit('gFji: nga_create failed g_Fji', 827 & 0,GA_ERR) 828 call ga_zero(g_Fji) 829 do k=1,3 830 alo(1)=1 831 ahi(1)=nbf 832 alo(2)=1 833 ahi(2)=nbf 834 alo(3)=k 835 ahi(3)=k 836 ld(1) =nbf 837 ld(2) =nbf 838 call ga_get(g_hfine(k),1,nbf,1,nbf,dbl_mb(k_buf),nbf) 839 call nga_put(g_Fji,alo,ahi,dbl_mb(k_buf),ld) 840 enddo ! end-loop-k 841c ---Destroying ga arrays ----- START 842 do i=1,3 843 if (.not. ga_destroy(g_hfine(i))) call errquit( 844 & 'gFij: ga_destroy failed ',0, GA_ERR) 845 enddo 846c ---Destroying ga arrays ----- END 847c----deallocate memory 848 if (.not.ma_free_heap(l_buf)) call errquit 849 & ('gFij, ma_free_heap of l_buf failed',911,MA_ERR) 850 if (ofinite) then 851 if (.not.ma_free_heap(l_zetanuc)) call 852 & errquit('HFine:: ma_free_heap l_zetanuc',0, MA_ERR) 853 endif 854 return 855 end 856c +++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 857c ++++++++++++ NLMO analysis routines +++++++++++++++++ START 858 subroutine fill_munuFCSD(g_munuFCSD, !out: array with matrices 859 & count , !in/out: counting data stored in g_munuFCSD 860 & g_fcsd, ! in: nbf 861 & nbf) ! in: nr. basis functions 862 ! Purpose: g_fcsd(u,v) --> g_munuFCSD 863 ! Note: g_fcsd(u,v) corresponds to ith atom 864 ! g_munuFCSD contains matrices for all atoms 865 ! g_munuFCSD contains unique munu 866 ! in sequence: xx,yy,zz,xy,xz,yz 867 ! each chunk is nbf*(nbf+1)/2 in size 868 ! Author: Fredy W. Aquino 869 ! Date : 06-29-11 870 implicit none 871#include "nwc_const.fh" 872#include "errquit.fh" 873#include "global.fh" 874#include "bas.fh" 875#include "mafdecls.fh" 876#include "geom.fh" 877#include "stdio.fh" 878#include "msgids.fh" 879 880 integer u,v,ii,i,j,k,indx,indx1 881 integer nbf,count 882 double precision val 883 integer g_munuFCSD,g_fcsd(3,3) 884 integer l_fcsd,k_fcsd 885 integer ndir ! in: =6 = xx,yy,zz,xy,xz,yz 886 integer nlist_uv(2,6) ! ndir=6 887 data nlist_uv /1,1, ! xx 888 & 2,2, ! yy 889 & 3,3, ! zz 890 & 1,2, ! xy 891 & 1,3, ! xz 892 & 2,3 / ! yz 893 ndir=6 ! = xx,yy,zz,xy,xz,yz 894 if (.not.ma_alloc_get(mt_dbl,nbf*nbf, 895 & 'fcsd',l_fcsd,k_fcsd)) 896 & call errquit('get_munuFCSD: ma failed',0,MA_ERR) 897 do ii=1,ndir 898 u=nlist_uv(1,ii) 899 v=nlist_uv(2,ii) 900 call ga_get(g_fcsd(u,v),1,nbf,1,nbf, 901 & dbl_mb(k_fcsd),nbf) 902 do i=1,nbf 903 indx=k_fcsd+nbf*(i-1)+i-1 904 val=dbl_mb(indx) 905 call ga_fill_patch(g_munuFCSD,1,1,count,count,val) 906 count=count+1 907 enddo ! end-loop-i 908 do i=2,nbf 909 do j=1,i-1 910 indx=k_fcsd+nbf*(j-1)+i-1 911 val=dbl_mb(indx) 912 call ga_fill_patch(g_munuFCSD,1,1,count,count,val) 913 count=count+1 914 enddo ! end-loop-j 915 enddo ! end-loop-i 916 enddo ! end-loop-ii 917 if (.not.ma_free_heap(l_fcsd)) call 918 & errquit('fill_munuFCSD:: ma_free_heap l_fcsd',0, MA_ERR) 919 return 920 end 921 922 subroutine fill_munuPSOSO(g_munuPSOSO, !out: array with matrices 923 & count , !in/out: counting data stored in g_munuFCSD 924 & g_psoso , ! in: nbf 925 & nbf) ! in: nr. basis functions 926 ! Purpose: g_psoso --> g_munuPSOSO 927 ! Note: g_psoso corresponds to uth component of ith atom (g_psoso_u u=1,2,3) 928 ! g_munuPSOSO contains matrices for all atoms 929 ! g_munuPSOSO contains unique munu elements 930 ! main diag. elements + off diagonal elements 931 ! Author: Fredy W. Aquino 932 ! Date : 07-04-11 933 implicit none 934#include "nwc_const.fh" 935#include "errquit.fh" 936#include "global.fh" 937#include "bas.fh" 938#include "mafdecls.fh" 939#include "geom.fh" 940#include "stdio.fh" 941#include "msgids.fh" 942 943 integer i,j,indx 944 integer nbf,count 945 double precision val 946 integer g_munuPSOSO,g_psoso 947 integer l_psoso,k_psoso 948 949 if (.not.ma_alloc_get(mt_dbl,nbf*nbf, 950 & 'fcsd',l_psoso,k_psoso)) 951 & call errquit('get_munuPSOSO: ma failed',0,MA_ERR) 952 call ga_get(g_psoso,1,nbf,1,nbf, 953 & dbl_mb(k_psoso),nbf) 954 do i=1,nbf 955 indx=k_psoso+nbf*(i-1)+i-1 956 val=dbl_mb(indx) 957 call ga_fill_patch(g_munuPSOSO,1,1,count,count,val) 958 count=count+1 959 enddo ! end-loop-i 960 do i=2,nbf 961 do j=1,i-1 962 indx=k_psoso+nbf*(j-1)+i-1 963 val=dbl_mb(indx) 964 call ga_fill_patch(g_munuPSOSO,1,1,count,count,val) 965 count=count+1 966 enddo ! end-loop-j 967 enddo ! end-loop-i 968 if (.not.ma_free_heap(l_psoso)) call 969 & errquit('fill_munuPSOSO:: ma_free_heap l_psoso',0, MA_ERR) 970 return 971 end 972 973 subroutine fill_munuPSOSO_1( 974 & g_munuPSOSO , ! in: array with matrices 975 & g_munuPSOSO2d , !out: nbf x nbf x 3 munu matrix for ith atom 976 & iat, ! in: atom index = 1, 2, nlist 977 & type_symm , ! in: =1 symmetric =2 antisymmetric 978 & nbf) ! in: nr. basis functions 979 ! Purpose: g_munuPSOSO --> g_munuPSOSO2d 980 ! Note: g_munuPSOSO2d corresponds to uth component of ith atom (g_munuPSOSO2d_u u=1,2,3) 981 ! g_munuPSOSO contains matrices for all atoms 982 ! g_munuPSOSO contains unique munu elements 983 ! main diag. elements + off diagonal elements 984 ! Author: Fredy W. Aquino 985 ! Date : 07-04-11 986 ! from (upper/lower) triang matrix --> full symmetric matrix 987 ! g_munuPSOSO = (11) (22) (33) ...(nn) (21) (31) (32) (41) (42) (43) ...(n (n-1)) 988 ! n_elements(g_munuPSOSO) = n (n-1)/2 n=nbf nr basis functions 989 ! g_munuPSOSO2d is n x n antisymmetric matrix 990 ! g_munuPSOSO ---> g_munuPSOSO2d 991 implicit none 992#include "nwc_const.fh" 993#include "errquit.fh" 994#include "global.fh" 995#include "bas.fh" 996#include "mafdecls.fh" 997#include "geom.fh" 998#include "stdio.fh" 999#include "msgids.fh" 1000 1001 integer i,j,indx,alo(3),ahi(3) 1002 integer nbf,nlst,iat,shift,ndir,ntot,xyz,type_symm 1003 double precision val,valneg 1004 integer g_munuPSOSO,g_munuPSOSO2d 1005 integer l_mat,k_mat 1006 1007 if (.not.(type_symm.eq.1 .or. 1008 & type_symm.eq.2)) then 1009 if (ga_nodeid().eq.0) 1010 & write(*,*) 'Error in fill_munuPSOSO_1: type_symm ne 1 or 2' 1011 stop 1012 endif 1013 1014 nlst=nbf*(nbf+1)/2 1015 ndir=3 ! x, y, z 1016 ntot=nlst*ndir 1017 shift=ntot*(iat-1) 1018 if (.not.ma_alloc_get(mt_dbl,ntot, 1019 & 'fcsd',l_mat,k_mat)) 1020 & call errquit('get_munuPSOSO: ma failed',0,MA_ERR) 1021 call ga_get(g_munuPSOSO,1,1,shift+1,shift+ntot, 1022 & dbl_mb(k_mat),1) 1023 indx=k_mat 1024 do xyz=1,ndir 1025 do i=1,nbf 1026 val=dbl_mb(indx) 1027 alo(1)=i 1028 ahi(1)=i 1029 alo(2)=i 1030 ahi(2)=i 1031 alo(3)=xyz 1032 ahi(3)=xyz 1033 call nga_fill_patch(g_munuPSOSO2d,alo,ahi,val) 1034 indx=indx+1 1035 enddo ! end-loop-i 1036 do i=2,nbf 1037 do j=1,i-1 1038 val=dbl_mb(indx) 1039 alo(1)=i 1040 ahi(1)=i 1041 alo(2)=j 1042 ahi(2)=j 1043 alo(3)=xyz 1044 ahi(3)=xyz 1045 call nga_fill_patch(g_munuPSOSO2d,alo,ahi,val) 1046 if (type_symm.eq.1) then 1047 valneg= val 1048 else if (type_symm.eq.2) then 1049 valneg=-val 1050 endif 1051 alo(1)=j 1052 ahi(1)=j 1053 alo(2)=i 1054 ahi(2)=i 1055 alo(3)=xyz 1056 ahi(3)=xyz 1057 call nga_fill_patch(g_munuPSOSO2d,alo,ahi,valneg) 1058 indx=indx+1 1059 enddo ! end-loop-j 1060 enddo ! end-loop-i 1061 enddo ! end-loop-xyz 1062 if (.not.ma_free_heap(l_mat)) call 1063 & errquit('fill_munuPSOSO_1:: ma_free_heap l_mat',0, MA_ERR) 1064 1065 return 1066 end 1067 1068 logical function dft_zoraHYP_NLMOAnalysis_write( 1069 & filename, ! in: filename 1070 & nbf, ! in: nr basis functions 1071 & ndir, ! in: nr of directions: 6 = xx yy zz xy xz yz for g_munuFCSD 1072 & ndir1, ! in: nr of directions: 3 = x y z for g_munuPSOSO 1073 & nlist, ! in: list of selected atoms 1074 & nocc, ! in: nocc(i) i=1,2 nr. occupations 1075 & ndata, ! in: writing order =1,2 1076 & g_munuFCSD, ! in: munu dia 1077 & g_munuPSOSO, ! in: munu para or PSOSO term 1078 & npol, ! in: nr. of polarizations 1079 & vectors, ! in: MOs 1080 & g_c1, ! in: perturbed MO coeffs 1081 & g_sdens) ! in: spin density 1082c Description: Collecting three matrices to be used 1083c in wefgfile(rtdb) and wnbofile(rtdb) 1084c Those routines are called in prop.F 1085c after hnd_property(rtdb) 1086c The info collected is: 1087c g_munuFCSD, FCSD munu matrix 1088 implicit none 1089#include "errquit.fh" 1090#include "mafdecls.fh" 1091#include "global.fh" 1092#include "tcgmsg.fh" 1093#include "msgtypesf.h" 1094#include "inp.fh" 1095#include "msgids.fh" 1096#include "cscfps.fh" 1097#include "util.fh" 1098#include "stdio.fh" 1099 character*(*) filename ! [input] File to write to 1100 integer nlist, ! = nr. slc atoms 1101 & nmat, ! = nlst*ndir*nlist (ndir=6) 1102 & nmat2, ! = nlst*ndir1*nlist (ndir1=3) 1103 & nmat1, ! = nbf*nbf 1104 & npol 1105 integer g_munuFCSD,g_munuPSOSO,g_sdens, 1106 & vectors(npol),g_c1 1107 integer ndir,ndir1,nbf,nlst,ndata,ntot,nocc(2) 1108 integer unitno 1109 parameter (unitno = 77) 1110 integer l_mat ,k_mat, 1111 & l_mat1,k_mat1, 1112 & l_mat2,k_mat2, 1113 & l_c1,k_c1, 1114 & l_mo,k_mo 1115 integer ok,iset,i,j,alo(3),ahi(3),ld(2) 1116 integer inntsize 1117 1118 inntsize=MA_sizeof(MT_INT,1,MT_BYTE) 1119 call ga_sync() 1120 ok = 0 1121c Read routines should be consistent with this 1122c Write out the atomic zora corrections 1123 if (ga_nodeid() .eq. 0) then 1124c Allocate the temporary buffer 1125 if (ndata.eq.1) then ! First time writing 1126c ++++++++ using ma_alloc_get +++++++++++++++++++ START 1127 nlst=nbf*(nbf+1)/2 1128 nmat=nlst*ndir*nlist 1129 if (.not. ma_alloc_get( 1130 & mt_dbl,nmat,'dft_zoraNLMO_writehyp', 1131 & l_mat,k_mat)) 1132 $ call errquit('dft_zoraNLMO_writehyp: k_mat failed', 1133 & nmat, MA_ERR) 1134 nmat2=nlst*ndir1*nlist 1135 if (.not. ma_alloc_get( 1136 & mt_dbl,nmat2,'dft_zoraNLMO_writehyp', 1137 & l_mat2,k_mat2)) 1138 $ call errquit('dft_zoraNLMO_writehyp: k_mat2 failed', 1139 & nmat2, MA_ERR) 1140 nmat1=nbf*nbf 1141 if (.not. ma_alloc_get( 1142 & mt_dbl,nmat,'dft_zoraNLMO_writehyp', 1143 & l_mat1,k_mat1)) 1144 $ call errquit('dft_zoraNLMO_writehyp: k_mat1 failed', 1145 & nmat1, MA_ERR) 1146c ++++++++ using ma_alloc_get +++++++++++++++++++ END 1147c Open the file - 1st time 1148 open(unitno, status='unknown', form='unformatted', 1149 $ file=filename, err=1000) 1150c Write out the number of sets and basis functions 1151 write(unitno, err=1001) nbf 1152 write(unitno, err=1001) nlst 1153 write(unitno, err=1001) ndir 1154 write(unitno, err=1001) ndir1 1155 write(unitno, err=1001) nlist 1156 call ga_get(g_munuFCSD,1,1,1,nmat, 1157 & dbl_mb(k_mat),1) 1158 call swrite(unitno,dbl_mb(k_mat),nmat) 1159 call ga_get(g_munuPSOSO,1,1,1,nmat2, 1160 & dbl_mb(k_mat2),1) 1161 call swrite(unitno,dbl_mb(k_mat2),nmat2) 1162 call ga_get(g_sdens,1,nbf,1,nbf, 1163 & dbl_mb(k_mat1),nbf) 1164 call swrite(unitno,dbl_mb(k_mat1),nmat1) 1165c Deallocate the temporary buffer 1166c ----- Using ma_free_heap ------------ START 1167 if (.not. ma_free_heap(l_mat)) 1168 $ call errquit('dft_zoraNLMO_writehyp: l_mat free_heap failed', 1169 & 911, MA_ERR) 1170 if (.not. ma_free_heap(l_mat1)) 1171 $ call errquit('dft_zoraNLMO_writehyp: l_mat1 free_heap failed', 1172 & 911, MA_ERR) 1173 if (.not. ma_free_heap(l_mat2)) 1174 $ call errquit('dft_zoraNLMO_writehyp: l_mat2 free_heap failed', 1175 & 911, MA_ERR) 1176c ----- Using ma_free_heap ------------ END 1177 else if (ndata.eq.2) then 1178 if (.not. ma_alloc_get( 1179 & mt_dbl,nbf,'dft_zoraNLMO_writehyp', 1180 & l_mo,k_mo)) 1181 $ call errquit('dft_zoraNLMO_writehyp: l_mo failed', 1182 & nbf,MA_ERR) 1183 ndir = 3 ! x,y,z 1184 ntot = nocc(1)+nocc(2) 1185 nmat = nbf*ndir*ntot 1186 if (.not. ma_alloc_get( 1187 & mt_dbl,nmat,'dft_zoraNLMO_writehyp',l_c1,k_c1)) 1188 $ call errquit('dft_zoraNLMO_writehyp: ma failed', 1189 & nmat, MA_ERR) 1190c Open the file - 2nd time 1191 open(unitno, status='unknown', form='unformatted', 1192 $ file=filename, err=1000,position='append') 1193 write(unitno, err=1001) nocc(1) 1194 write(unitno, err=1001) nocc(2) 1195 write(unitno, err=1001) ntot 1196 write(unitno, err=1001) nmat 1197 write(unitno, err=1001) npol 1198 write(unitno, err=1001) nbf 1199c ----- Add MOs in file ----- START 1200 do i=1,npol 1201 do j=1,nbf 1202 call dcopy(nbf,0.0d0,0,dbl_mb(k_mo),1) ! reset 1203 call ga_get(vectors(i),1,nbf,j,j,dbl_mb(k_mo),1) 1204 call swrite(unitno,dbl_mb(k_mo),nbf) 1205 enddo ! end-loop-j 1206 enddo ! end-loop-i 1207c ----- Add MOs in file ----- END 1208 alo(1)=1 1209 ahi(1)=nbf 1210 alo(2)=1 1211 ahi(2)=ntot 1212 alo(3)=1 1213 ahi(3)=3 1214 ld(1)=nbf 1215 ld(2)=ntot 1216 call nga_get(g_c1,alo,ahi,dbl_mb(k_c1),ld) 1217 call swrite(unitno,dbl_mb(k_c1),nmat) 1218 if (.not. ma_free_heap(l_mo)) 1219 $ call errquit('dft_zoraNLMO_writehyp: ma free_heap failed', 1220 & 911, MA_ERR) 1221 if (.not. ma_free_heap(l_c1)) 1222 $ call errquit('dft_zoraNLMO_writehyp: ma free_heap failed', 1223 & 911, MA_ERR) 1224 endif 1225c Close the file 1226 close(unitno,err=1002) 1227 ok = 1 1228 end if 1229c Broadcast status to other nodes 1230 10 call ga_brdcst(Msg_Vec_Stat+MSGINT, ok, inntsize, 0) ! Propagate status 1231 call ga_sync() 1232 dft_zoraHYP_NLMOAnalysis_write = (ok .eq. 1) 1233 if (ga_nodeid() .eq. 0) then 1234 write(6,22) filename(1:inp_strlen(filename)) 1235 22 format(/' Wrote ZORA NLMO HYP data to ',a/) 1236 call util_flush(luout) 1237 endif 1238 return 1239 1000 write(6,*) 'dft_zoraNLMO_writehyp: failed to open ', 1240 $ filename(1:inp_strlen(filename)) 1241 call util_flush(luout) 1242 ok = 0 1243 goto 10 1244 1001 write(6,*) 'dft_zoraNLMO_writehyp: failed to write ', 1245 $ filename(1:inp_strlen(filename)) 1246 call util_flush(luout) 1247 ok = 0 1248 close(unitno,err=1002) 1249 goto 10 1250 1002 write(6,*) 'dft_zoraNLMO_writehyp: failed to close', 1251 $ filename(1:inp_strlen(filename)) 1252 call util_flush(luout) 1253 ok = 0 1254 goto 10 1255 end 1256 1257 logical function dft_zoraHYP_NLMOAnalysis_read( 1258 & filename, ! in: filename 1259 & nbf, ! in: nr basis functions 1260 & ndir, ! in: nr of directions: 6 = xx yy zz xy xz yz for g_munuFCSD 1261 & ndir1, ! in: nr od directions: 3 = x y z for g_munuPSOSO 1262 & nlist, ! in: list of selected atoms 1263 & nocc, ! in: nocc(i) i=1,2 1264 & npol, ! in: nr polarizations 1265 & g_munuFCSD, ! out: munu dia 1266 & g_munuPSOSO, ! out: munu para 1267 & vectors, ! out: MOs 1268 & g_c1, ! out: perturbed MO coeffs 1269 & g_sdens) ! out: spin density 1270 implicit none 1271#include "errquit.fh" 1272#include "global.fh" 1273#include "tcgmsg.fh" 1274#include "msgtypesf.h" 1275#include "mafdecls.fh" 1276#include "msgids.fh" 1277#include "cscfps.fh" 1278#include "inp.fh" 1279#include "util.fh" 1280#include "stdio.fh" 1281 1282 character*(*) filename ! [input] File to write to 1283 integer nmat,nmat1,nmat2,npol 1284 integer g_munuFCSD,g_munuPSOSO,g_sdens, 1285 & vectors(npol),g_c1 1286 integer nbf,nbf_read,nlst,nlst_read, 1287 & ndir,ndir_read, 1288 & ndir1,ndir1_read, 1289 & nlist,nlist_read, 1290 & ntot,ntot_read, 1291 & nocc(2),nocc_read(2), 1292 & n_c1,n_c1_read, 1293 & npol_read, 1294 & alo(3),ahi(3),ld(2) 1295 integer unitno 1296 parameter (unitno = 77) 1297 integer l_mat ,k_mat, 1298 & l_mat1,k_mat1, 1299 & l_mat2,k_mat2, 1300 & l_c1,k_c1, 1301 & l_mo,k_mo 1302 integer ok,iset,i,j 1303 integer inntsize 1304c Initialise to invalid MA handle 1305 nlst=nbf*(nbf+1)/2 1306 inntsize=MA_sizeof(MT_INT,1,MT_BYTE) 1307 call ga_sync() 1308 ok = 0 1309 nmat=nlst*ndir*nlist 1310 if (.not. ga_create(mt_dbl,1,nmat, 1311 & 'dft_zoraNLMO_read: g_munuFCSD',0,0,g_munuFCSD)) 1312 $ call errquit('dft_zoraNLMO_read: g_munuFCSD',0,GA_ERR) 1313 call ga_zero(g_munuFCSD) 1314 nmat2=nlst*ndir1*nlist 1315 if (.not. ga_create(mt_dbl,1,nmat2, 1316 & 'dft_zoraNLMO_read: g_munuPSOSO',0,0,g_munuPSOSO)) 1317 $ call errquit('dft_zoraNLMO_read: g_munuPSOSO',0,GA_ERR) 1318 call ga_zero(g_munuPSOSO) 1319 if (.not. ga_create(mt_dbl,nbf,nbf, 1320 & 'dft_zoraNLMO_read: g_sdens',0,0,g_sdens)) 1321 $ call errquit('dft_zoraNLMO_read: g_sdens',0,GA_ERR) 1322 call ga_zero(g_sdens) 1323 do i=1,npol 1324 if (.not. ga_create(mt_dbl,nbf,nbf, 1325 & 'dft_zoraNLMO_read: g_sdens',0,0,vectors(i))) 1326 $ call errquit('dft_zoraNLMO_read: vectors',0,GA_ERR) 1327 call ga_zero(vectors(i)) 1328 enddo 1329 ntot=nocc(1)+nocc(2) 1330 n_c1=nbf*3*ntot 1331 alo(1) = nbf 1332 alo(2) = -1 1333 alo(3) = -1 1334 ahi(1) = nbf 1335 ahi(2) = ntot 1336 ahi(3) = 3 1337 if (.not.nga_create(MT_DBL,3,ahi,'c1 matrix',alo,g_c1)) call 1338 & errquit('dft_zoraNMR_read: nga_create failed g_c1', 1339 & 0,GA_ERR) 1340 call ga_zero(g_c1) 1341 1342 if (ga_nodeid() .eq. 0) then 1343c Print a message indicating the file being read 1344 write(6,22) filename(1:inp_strlen(filename)) 1345 22 format(/' Read ZORA HYP NLMO data from ',a/) 1346 call util_flush(luout) 1347c Open the file 1348 open(unitno, status='old', form='unformatted', file=filename, 1349 $ err=1000) 1350c Read in some basics to check if they are consistent with the calculation 1351 read(unitno, err=1001, end=1001) nbf_read 1352 read(unitno, err=1001, end=1001) nlst_read 1353 read(unitno, err=1001, end=1001) ndir_read 1354 read(unitno, err=1001, end=1001) ndir1_read 1355 read(unitno, err=1001, end=1001) nlist_read 1356c Error checks 1357 if ((nbf_read .ne. nbf) .or. 1358 & (nlst_read .ne. nlst) .or. 1359 & (ndir_read .ne. ndir) .or. 1360 & (ndir1_read .ne. ndir1) .or. 1361 & (nlist_read .ne. nlist)) goto 1003 1362c ---> ma_alloc_get: to allocate memory 1363c ---> ma_free_heap: to release allocated memory 1364 nmat=nlst*ndir*nlist 1365 if (.not. ma_alloc_get(mt_dbl,nmat, ! allocate memory 1366 & 'dft_zoraNLMO_read',l_mat,k_mat)) 1367 $ call errquit('dft_zoraNLMO_read: ma failed', 1368 & nmat, MA_ERR) 1369 nmat2=nlst*ndir1*nlist 1370 if (.not. ma_alloc_get(mt_dbl,nmat2, ! allocate memory 1371 & 'dft_zoraNLMO_read',l_mat2,k_mat2)) 1372 $ call errquit('dft_zoraNLMO_read: ma failed', 1373 & nmat2, MA_ERR) 1374 nmat1=nbf*nbf 1375 if (.not. ma_alloc_get(mt_dbl,nmat1, ! allocate memory 1376 & 'dft_zoraNLMO_read',l_mat1,k_mat1)) 1377 $ call errquit('dft_zoraNLMO_read: ma failed', 1378 & nmat1, MA_ERR) 1379 call sread(unitno,dbl_mb(k_mat),nmat) 1380 call ga_put(g_munuFCSD,1,1,1,nmat,dbl_mb(k_mat),1) 1381 call sread(unitno,dbl_mb(k_mat2),nmat2) 1382 call ga_put(g_munuPSOSO,1,1,1,nmat2,dbl_mb(k_mat2),1) 1383 call sread(unitno,dbl_mb(k_mat1),nmat1) 1384 call ga_put(g_sdens,1,nbf,1,nbf,dbl_mb(k_mat1),nbf) 1385 if (.not. ma_free_heap(l_mat)) ! deallocate memory 1386 $ call errquit('dft_zoraNLMO_read: ma free_heap failed', 1387 & 911, MA_ERR) 1388 if (.not. ma_free_heap(l_mat1)) ! deallocate memory 1389 $ call errquit('dft_zoraNLMO_read: ma free_heap failed', 1390 & 911, MA_ERR) 1391 if (.not. ma_free_heap(l_mat2)) ! deallocate memory 1392 $ call errquit('dft_zoraNLMO_read: ma free_heap failed', 1393 & 911, MA_ERR) 1394c -------- Read Perturbed MO coeffs ------------- START 1395 read(unitno, err=1001, end=1001) nocc_read(1) 1396 read(unitno, err=1001, end=1001) nocc_read(2) 1397 read(unitno, err=1001, end=1001) ntot_read 1398 read(unitno, err=1001, end=1001) n_c1_read 1399 read(unitno, err=1001, end=1001) npol_read 1400 read(unitno, err=1001, end=1001) nbf_read 1401c Error checks 1402 if ((nocc_read(1) .ne. nocc(1)) .or. 1403 & (nocc_read(2) .ne. nocc(2)) .or. 1404 & (ntot_read .ne. ntot) .or. 1405 & (nbf_read .ne. nbf) .or. 1406 & (npol_read .ne. npol) .or. 1407 & (n_c1_read .ne. n_c1)) goto 1003 1408c ----- Read MOs ----- START 1409 if (.not. ma_alloc_get( 1410 & mt_dbl,nbf,'dft_zoraNLMO_readhyp', 1411 & l_mo,k_mo)) 1412 $ call errquit('dft_zoraNLMO_readhyp: ma failed', 1413 & nbf,MA_ERR) 1414 do i=1,npol 1415 do j=1,nbf 1416 call dcopy(nbf,0.0d0,0,dbl_mb(k_mo),1) ! reset 1417 call sread(unitno,dbl_mb(k_mo),nbf) 1418 call ga_put(vectors(i),1,nbf,j,j,dbl_mb(k_mo),1) 1419 enddo ! end-loop-j 1420 enddo ! end-loop-i 1421c ----- Read MOs ----- END 1422 if (.not. ma_alloc_get(mt_dbl,n_c1, ! allocate memory 1423 & 'dft_zoraNLMO_read',l_c1,k_c1)) 1424 $ call errquit('dft_zoraNLMO_read: ma failed', 1425 & n_c1, MA_ERR) 1426! Note: n_c1 = nbf*3*ntot [ ntot=nocc(1)+nocc(2) ] 1427 alo(1)=1 1428 ahi(1)=nbf 1429 alo(2)=1 1430 ahi(2)=ntot 1431 alo(3)=1 1432 ahi(3)=3 1433 ld(1) =nbf 1434 ld(2) =ntot 1435 call sread(unitno,dbl_mb(k_c1),n_c1) 1436 call nga_put(g_c1,alo,ahi,dbl_mb(k_c1),ld) 1437 if (.not. ma_free_heap(l_mo)) ! deallocate memory 1438 $ call errquit('dft_zoraNLMO_read: ma free_heap failed', 1439 & 911, MA_ERR) 1440 if (.not. ma_free_heap(l_c1)) ! deallocate memory 1441 $ call errquit('dft_zoraNLMO_read: ma free_heap failed', 1442 & 911, MA_ERR) 1443c -------- Read Perturbed MO coeffs ------------- END 1444c Close the file 1445 close(unitno,err=1002) 1446 ok = 1 1447 end if 1448c 1449c Broadcast status to other nodes 1450 10 call ga_brdcst(Msg_Vec_Stat+MSGINT, ok, inntsize, 0) ! Propagate status 1451 call ga_sync() 1452 dft_zoraHYP_NLMOAnalysis_read = ok .eq. 1 1453 return 1454 1000 write(6,*) 'dft_zoraNLMO_read: failed to open', 1455 $ filename(1:inp_strlen(filename)) 1456 call util_flush(luout) 1457 ok = 0 1458 goto 10 1459 1001 write(6,*) 'dft_zoraNLMO_read: failed to read', 1460 $ filename(1:inp_strlen(filename)) 1461 call util_flush(luout) 1462 ok = 0 1463 close(unitno,err=1002) 1464 goto 10 1465 1003 write(6,*) 'dft_zoraNLMO_read: file inconsistent', 1466 & ' with calculation', 1467 $ filename(1:inp_strlen(filename)) 1468 call util_flush(luout) 1469 ok = 0 1470 close(unitno,err=1002) 1471 goto 10 1472 1002 write(6,*) 'dft_zoraNLMO_read: failed to close', 1473 $ filename(1:inp_strlen(filename)) 1474 call util_flush(luout) 1475 ok = 0 1476 goto 10 1477 end 1478c ++++++++++++ NLMO analysis routines +++++++++++++++++ END 1479c +++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1480 integer function read_SLCTD_HFine_Atoms 1481 & (rtdb,nat,nlist,g_AtNr) 1482c---- GA output: g_AtNr 1483 1484 implicit none 1485#include "errquit.fh" 1486#include "global.fh" 1487#include "tcgmsg.fh" 1488#include "msgtypesf.h" 1489#include "mafdecls.fh" 1490#include "msgids.fh" 1491#include "cscfps.fh" 1492#include "inp.fh" 1493#include "util.fh" 1494#include "stdio.fh" 1495#include "rtdb.fh" 1496#include "context.fh" 1497 1498 integer rtdb,ii,nlist,nat,g_AtNr 1499 integer atomnr(nat) 1500 integer hfineatoms 1501 double precision AtNr_dbl 1502 1503 if (.not. rtdb_get(rtdb, 'hfine:natoms',mt_int, 1504 $ 1,hfineatoms)) 1505 & hfineatoms=0 ! reset 1506 if (hfineatoms.eq.0) then 1507 hfineatoms=nat 1508 nlist=hfineatoms 1509 do ii=1,hfineatoms 1510 AtNr_dbl=ii 1511 call ga_put(g_AtNr,1,1,ii,ii,AtNr_dbl,1) 1512 enddo 1513 else 1514 if (.not. rtdb_get(rtdb, 'hfine:atom list',mt_int, 1515 $ hfineatoms,atomnr)) 1516 $ call errquit('prop_input-hfine: rtdb_get failed', 1517 $ 555, RTDB_ERR) 1518 nlist=hfineatoms 1519 do ii=1,hfineatoms 1520 AtNr_dbl=atomnr(ii) 1521 call ga_put(g_AtNr,1,1,ii,ii,AtNr_dbl,1) 1522 enddo 1523 endif 1524 read_SLCTD_HFine_Atoms = 1 1525 return 1526 end 1527c 1528 subroutine print_NMRHypFine_version() 1529c 1530 implicit none 1531#include "stdio.fh" 1532c 1533 write(LuOut,*) 1534 call util_print_centered(LuOut, 1535 $ 'ZORA NMR Hyperfine Spin-Spin Coupling Constants', 23, .true.) 1536 write(LuOut,*) 1537c 1538 return 1539 end 1540 1541c $Id$ 1542