1czora... 2czora...Scale zora eigenvalues and energy 3czora... 4 subroutine dft_zora_scale_so( 5 & rtdb, 6 & g_dens_ini,nexc, 7 & geom, 8 & ao_bas_han, 9 & nbf, 10 & nbf_ao, 11 & nbf_mo, 12 & g_dens, 13 & g_s, 14 & g_moso, 15 & g_zora_scale_sf, 16 & g_zora_scale_so, 17 & evals, 18 & focc, 19 & noc, 20 & nocc, ! FA 21 & ipol, 22 & ener_scal) 23 24 implicit none 25#include "errquit.fh" 26#include "mafdecls.fh" 27#include "stdio.fh" 28#include "global.fh" 29#include "msgids.fh" 30#include "zora.fh" 31#include "rtdb.fh" 32 33 integer ga_create_atom_blocked 34 external ga_create_atom_blocked 35 external get_rhoS_so ! FA 36 external hnd_efgmap_Z4_so ! FA 37 external print_EFGZ4_version ! FA 38 integer g_dens(2) 39 integer g_moso(2) 40 integer g_orb(2) 41 integer g_dens_so(2) 42 integer g_scr(2) 43 integer g_s 44 integer g_zora_scale_sf(2) 45 integer g_zora_scale_so(3) 46 double precision numelecs, ener_kin 47 integer l_vecsre, k_vecsre 48 integer l_vecsim, k_vecsim 49 integer iorb 50 double precision eval_scal 51 double precision ener_scal 52 double precision ener_sf 53 double precision ener_so_x 54 double precision ener_so_y 55 double precision ener_so_z 56 double precision ener_tot 57 integer noc 58 integer ispin 59 integer ipol 60 61 integer geom 62 integer ao_bas_han 63 integer nbf 64 integer nbf_ao 65 integer nbf_mo 66 67 double precision focc(*) ! occupation no. 68 double precision evals(*) ! eigenvalues 69 integer do_replnxn,mem_replnxn,val_in 70 integer l_evals_in, k_evals_in, 71 & l_zora_scale_sf, k_zora_scale_sf, 72 & l_zora_scale_so, k_zora_scale_so, 73 & l_dens_so, k_dens_so, 74 & l_vtmp, k_vtmp 75 double precision zora_denom ! Added by FA 76 integer rtdb ! Added by FA 77 integer nexc,g_dens_ini(2) ! Added by FA 78 integer j,j1,iorb1,nocc(2) ! Added by FA 79 logical do_zgc_old ! Added by FA 80 integer ac_AB(2),iter ! Added by FA 81 character*1 soxyz(3) ! Added by FA 82 double precision ener_sf_1,ener_sf_2 ! Added by FA 83 logical status ! Added by FA 84 integer noelfgz4 ! Added by FA 85 data noelfgZ4/1/ ! Added by FA 86 soxyz(1)='z' 87 soxyz(2)='y' 88 soxyz(3)='x' 89c +++++++++++++++++++++++++++++++++++++++++++++++++++++ 90 status=rtdb_get(rtdb,'prop:efieldgradZ4',MT_INT,1,noelfgZ4) ! FA 91 do_zgc_old=do_zora_get_correction 92 do_zora_get_correction= .true. !-1 ! FORCE-ZORA-CALC-INTS 93c +++++++++++++++++++++++++++++++++++++++++++++++++++++ 94c ----- execute this code ONLY when EFGZ4 is requested -- START 95 if(noelfgZ4.eq.0) then 96 if (ga_nodeid().eq.0) call print_EFGZ4_version() 97 if (.not.ga_create(MT_DBL,2,nocc(1)+nocc(2), ! FA 98 & 'dft_zora_scale: g_Ci',0,0,g_Ci)) ! FA 99 & call errquit('dft_zora_scale_so: g_Ci',0,GA_ERR) ! FA 100 call ga_zero(g_Ci) ! FA 101 endif 102c ----- execute this code ONLY when EFGZ4 is requested -- END 103c allocate memory 104 105 if (.not.MA_Push_Get(MT_Dbl, 2*nbf_mo, 'real vec aux', 106 & l_vecsre, k_vecsre)) 107 & call errquit('dft_zora_scale_so: cannot allocate vec',0, MA_ERR) 108c 109#if 0 110 if (.not.MA_Push_Get(MT_Dbl, nbf_mo, 'imag vec aux', 111 & l_vecsim, k_vecsim)) 112 & call errquit('dft_zora_scale_so: cannot allocate vec',0, MA_ERR) 113#else 114 k_vecsim=k_vecsre+nbf_mo 115#endif 116c 117c spin-orbit vector - real 118 if(.not.ga_create(mt_dbl,nbf_mo,nbf_mo,'orbs re',0,0,g_orb(1))) 119 & call errquit('dft_zora_scale_so: orb real error',0, GA_ERR) 120 call ga_zero(g_orb(1)) 121c 122c spin-orbit vector - imag 123 if(.not.ga_create(mt_dbl,nbf_mo,nbf_mo,'orbs im',0,0,g_orb(2))) 124 & call errquit('dft_zora_scale_so: orb imag 2 error',0, GA_ERR) 125 call ga_zero(g_orb(2)) 126c 127c spin-orbit density matrix - real 128 if(.not.ga_create(mt_dbl,nbf_mo,nbf_mo,'denmxre',0,0, 129 & g_dens_so(1))) 130 & call errquit('dft_zora_scale_so: dens real error',0, GA_ERR) 131 call ga_zero(g_dens_so(1)) 132c 133c spin-orbit density matrix - imag 134 if(.not.ga_create(mt_dbl,nbf_mo,nbf_mo,'denmxim',0,0, 135 & g_dens_so(2))) 136 & call errquit('dft_zora_scale_so: dens imag error',0, GA_ERR) 137 call ga_zero(g_dens_so(2)) 138c 139c scratch array 140 if(.not.ga_duplicate(g_dens(1),g_scr(1),'scratch 1')) 141 & call errquit('dft_zora_scale_so: ga_duplicate failed',1, GA_ERR) 142 call ga_zero(g_scr(1)) 143 if(.not.ga_duplicate(g_dens(2),g_scr(2),'scratch 2')) 144 & call errquit('dft_zora_scale_so: ga_duplicate failed',1, GA_ERR) 145 call ga_zero(g_scr(2)) 146 j=0 ! FA 147 ac_AB(1)=1 148 ac_AB(2)=1 149 ener_scal = 0.d0 150cwhy?? ener_tot = 0.d0 151c 152c check if enough memory is available to allocate a few nxn objects 153c 154 do_replnxn=0 155 mem_replnxn= 156 C nbf_mo*nbf_mo*2 + ! dens_so 157 C nbf_ao*nbf_ao + ! vtmp 158 C nbf_ao*nbf_ao*2 + ! zora_scale_sf 159 C nbf_ao*nbf_ao*3 + ! zora_scale_so 160 C nbf_mo ! evals_in 161 if (MA_inquire_avail(mt_dbl)*8/10.gt.mem_replnxn) then 162 do_replnxn=1 163 else 164 if(ga_nodeid().eq.0) then 165 write(6,*) ' mem needed in MB',mem_replnxn* 166 / 8/1024/1024 167 write(6,*) ' mem available in MB ',MA_inquire_avail(mt_dbl)* 168 / 8/1024/1024 169 endif 170 endif 171 if(noelfgZ4.eq.0) then 172 write(6,*) ' need to compute elfgZ4, do_replnxn=0 ' 173 do_replnxn=0 174 endif 175 if(rtdb_get(rtdb,'zora:do_replnxn',MT_INT,1,val_in)) 176 c do_replnxn=val_in 177 call ga_igop(19481984, do_replnxn, 1, 'min') 178 if(do_replnxn.eq.0) then 179 if(ga_nodeid().eq.0) 180 c write(6,*) ' zoraso: NOT enough mem for repl' 181 do iorb = 1, noc 182 call ga_get(g_moso(1),1,nbf_mo,iorb,iorb,dbl_mb(k_vecsre),1) 183 call ga_zero(g_orb(1)) 184 call ga_put(g_orb(1),1,nbf_mo,iorb,iorb,dbl_mb(k_vecsre),1) 185 call ga_get(g_moso(2),1,nbf_mo,iorb,iorb,dbl_mb(k_vecsim),1) 186 call ga_zero(g_orb(2)) 187 call ga_put(g_orb(2),1,nbf_mo,iorb,iorb,dbl_mb(k_vecsim),1) 188 call dft_densm_so(g_dens_so, g_orb, nbf_ao, noc) 189 call ga_zero(g_scr(1)) 190 call ga_zero(g_scr(2)) 191 call ga_dens_sf(g_scr, g_dens_so(1), nbf_ao) 192 ener_sf_1=ga_ddot(g_scr(1),g_zora_scale_sf(1)) 193 ener_sf_2=ga_ddot(g_scr(2),g_zora_scale_sf(2)) 194 ener_sf=ener_sf_1+ener_sf_2 195 ener_tot=ener_sf 196 do iter=1,3 197 call ga_zero(g_scr(1)) 198 call ga_dens_so(g_scr(1),g_dens_so,nbf_ao,soxyz(iter)) 199 ener_tot=ener_tot+ 200 & ga_ddot(g_scr(1),g_zora_scale_so(iter)) 201 end do 202 zora_denom=1.d0+ener_tot 203 eval_scal = evals(iorb) 204 eval_scal = eval_scal/zora_denom 205 ener_scal = ener_scal - eval_scal*ener_tot 206 evals(iorb) = eval_scal 207c ----- execute this code ONLY when EFGZ4 is requested -- START 208 if(noelfgZ4.eq.0) then 209 j1=j+1 ! FA 210 if (ac_AB(j1).le.nocc(j1)) then ! FA 211 call ga_fill_patch(g_Ci,j1,j1,ac_AB(j1),ac_AB(j1), 212 & zora_denom) ! FA 213 ac_AB(j1)=ac_AB(j1)+1 ! FA 214 end if ! FA 215 j=1-j ! FA 216 endif 217c ----- execute this code ONLY when EFGZ4 is requested -- END 218 enddo ! end-loop-iorb 219 else 220 if(ga_nodeid().eq.0) 221 c write(6,*) ' zoraso: enough mem for repl' 222 if (.not.ma_push_get(mt_dbl, nbf_mo, 'evals_in', 223 & l_evals_in, k_evals_in)) 224 & call errquit('dft_zora_scale_so: cannot allocate evals',0,MA_ERR) 225 if (.not.ma_push_get(mt_dbl, nbf_ao*nbf_ao*2, 'zora_scale_sf', 226 & l_zora_scale_sf, k_zora_scale_sf)) 227 & call errquit('dft_zora_scale_so:cant all zora_scale_sf',0,MA_ERR) 228 if (.not.ma_push_get(mt_dbl, nbf_ao*nbf_ao*3, 'zora_scale_so', 229 & l_zora_scale_so, k_zora_scale_so)) 230 & call errquit('dft_zora_scale_so:cant all zora_scale_so',0,MA_ERR) 231 if (.not.ma_push_get(mt_dbl, nbf_mo*nbf_mo*2, 'dens_so', 232 & l_dens_so, k_dens_so)) 233 & call errquit('dft_zora_scale_so:cant all dens_so',0,MA_ERR) 234 if (.not.ma_push_get(mt_dbl, nbf_ao*nbf_ao, 'vtmp', 235 & l_vtmp, k_vtmp)) 236 & call errquit('dft_zora_scale_so:cant all vtmp',0,MA_ERR) 237c 238c replicate zora_scale_sx 239c 240 call util_mygabcast(g_zora_scale_sf(1), 241 M nbf_ao,nbf_ao, dbl_mb(k_zora_scale_sf),nbf_ao) 242 call util_mygabcast(g_zora_scale_sf(2), 243 M nbf_ao,nbf_ao, dbl_mb(k_zora_scale_sf+nbf_ao*nbf_ao),nbf_ao) 244 245 call util_mygabcast(g_zora_scale_so(1), 246 M nbf_ao,nbf_ao, dbl_mb(k_zora_scale_so),nbf_ao) 247 call util_mygabcast(g_zora_scale_so(2), 248 M nbf_ao,nbf_ao, dbl_mb(k_zora_scale_so+nbf_ao*nbf_ao),nbf_ao) 249 call util_mygabcast(g_zora_scale_so(3), 250 M nbf_ao,nbf_ao, dbl_mb(k_zora_scale_so+2*nbf_ao*nbf_ao),nbf_ao) 251 252c copy input evals 253 call dcopy(nbf_mo,evals,1,dbl_mb(k_evals_in),1) 254c zero evals [output] to get dgop working 255 call dcopy(noc,0d0,0,evals,1) 256 do iorb = ga_nodeid()+1, noc, ga_nnodes() 257 call dft_zora_nogaga(iorb,nbf_ao,nbf_mo,noc,g_moso, 258 C soxyz, 259 D dbl_mb(k_zora_scale_sf),dbl_mb(k_zora_scale_so), 260 E ener_tot, 261 D dbl_mb(k_dens_so), 262 O dbl_mb(k_vtmp),dbl_mb(k_vecsre)) 263 zora_denom=1.d0+ener_tot 264 eval_scal = dbl_mb(k_evals_in+iorb-1)/zora_denom 265 ener_scal = ener_scal - eval_scal*ener_tot 266 evals(iorb) = eval_scal 267#if 0 268NOT WORKING YET 269c ----- execute this code ONLY when EFGZ4 is requested -- START 270 if(noelfgZ4.eq.0) then 271 j1=j+1 ! FA 272 if (ac_AB(j1).le.nocc(j1)) then ! FA 273 call ga_fill_patch(g_Ci,j1,j1,ac_AB(j1),ac_AB(j1), 274 & zora_denom) ! FA 275 ac_AB(j1)=ac_AB(j1)+1 ! FA 276 end if ! FA 277 j=1-j ! FA 278 endif 279#endif 280c ----- execute this code ONLY when EFGZ4 is requested -- END 281 enddo ! end-loop-iorb 282c dgop on evals 283 call ga_dgop(20166102, evals, noc, '+') 284 call ga_dgop(20166103, ener_scal, 1, '+') 285 endif 286c ----- execute this code ONLY when EFGZ4 is requested -- START 287 if(noelfgZ4.eq.0) then 288c if (ga_nodeid().eq.0) write(*,*) "---g_Ci_so -- START" 289c call ga_print(g_Ci) 290c if (ga_nodeid().eq.0) write(*,*) "---g_Ci_so -- END" 291c write(luout,*) "ener_scal:",ener_scal 292 call get_rhoS_so(rtdb,g_dens_ini,nexc, 293 & geom,ao_bas_han, 294 & nbf,nbf_ao,nbf_mo, 295 & g_dens,g_moso, 296 & noc,nocc,ipol) 297c -------- computing EFGs --------- START 298 call hnd_efgmap_Z4_so(rtdb,ao_bas_han,geom, 299 & nbf,nbf_ao,nbf_mo, 300 & g_dens,g_moso, 301 & noc,nocc) 302c -------- computing EFGs --------- END 303 endif 304c ----- execute this code ONLY when EFGZ4 is requested -- START 305 306c deallocate memory 307#if 0 308 if (.not. ma_chop_stack(l_vecsim)) 309 & call errquit('dft_zora_scale_so:l_vecsim', 0, MA_ERR) 310#endif 311 if (.not. ma_chop_stack(l_vecsre)) 312 & call errquit('dft_zora_scale_so:l_vecsre', 0, MA_ERR) 313 314 if (.not. ga_destroy(g_orb(1))) 315 & call errquit('dft_zora_scale_so: ga_destroy failed',0, GA_ERR) 316 if (.not. ga_destroy(g_orb(2))) 317 & call errquit('dft_zora_scale_so: ga_destroy failed',0, GA_ERR) 318 319 if (.not. ga_destroy(g_dens_so(1))) 320 & call errquit('dft_zora_scale_so: ga_destroy failed',0, GA_ERR) 321 if (.not. ga_destroy(g_dens_so(2))) 322 & call errquit('dft_zora_scale_so: ga_destroy failed',0, GA_ERR) 323 324 if (.not. ga_destroy(g_scr(1))) 325 & call errquit('dft_zora_scale_so: ga_destroy failed',0, GA_ERR) 326 if (.not. ga_destroy(g_scr(2))) 327 & call errquit('dft_zora_scale_so: ga_destroy failed',0, GA_ERR) 328c +++++++++++++++ RESTORE VALUE +++++++++++++START 329 do_zora_get_correction=do_zgc_old 330c +++++++++++++++ RESTORE VALUE +++++++++++++END 331 return 332 end 333c 334czora...Inquire if the zora correction file is present 335c 336 logical function dft_zora_inquire_file_so(filename) 337c 338 implicit none 339c 340#include "errquit.fh" 341#include "global.fh" 342#include "tcgmsg.fh" 343#include "msgtypesf.h" 344#include "mafdecls.fh" 345#include "msgids.fh" 346#include "cscfps.fh" 347#include "inp.fh" 348#include "util.fh" 349#include "stdio.fh" 350c 351 character*(*) filename 352 logical found 353c 354 call ga_sync() 355c 356c Inquire if file is present 357 inquire(file=filename,exist=found) 358 dft_zora_inquire_file_so = found 359c 360 call ga_sync() 361c 362 return 363 end 364c 365czora...Read in the zora atomic corrections from disk 366c 367 logical function dft_zora_read_so(filename, nbf, nsets, nmo, 368 & mult, g_zora_sf, g_zora_scale_sf, 369 & g_zora_so, g_zora_scale_so) 370c 371 implicit none 372c 373#include "errquit.fh" 374#include "global.fh" 375#include "tcgmsg.fh" 376#include "msgtypesf.h" 377#include "mafdecls.fh" 378#include "msgids.fh" 379#include "cscfps.fh" 380#include "inp.fh" 381#include "util.fh" 382#include "stdio.fh" 383c 384 character*(*) filename 385 integer iset ! restricted or unrestricted 386c 387 integer g_zora_sf(2) 388 integer g_zora_scale_sf(2) 389 integer g_zora_so(3) 390 integer g_zora_scale_so(3) 391c 392 integer nsets ! restricted or unrestricted 393 integer nbf ! No. of functions in basis 394 integer nmo(nsets) 395 integer mult 396 integer ok, jset, i, j 397 integer l_zora, k_zora 398c 399 integer unitno 400 parameter (unitno = 77) 401 integer inntsize,ddblsize 402c 403 integer nsets_read 404 integer nbf_read 405 integer mult_read 406c 407 l_zora = -1 ! An invalid MA handle 408c 409 inntsize=MA_sizeof(MT_INT,1,MT_BYTE) 410 ddblsize=MA_sizeof(MT_DBL,1,MT_BYTE) 411c 412 call ga_sync() 413 ok = 0 414 if (ga_nodeid() .eq. 0) then 415c 416c Print a message indicating the file being read 417 write(6,22) filename(1:inp_strlen(filename)) 418 22 format(/' Read atomic ZORA corrections from ',a/) 419 call util_flush(luout) 420c 421c Open the file 422 open(unitno, status='old', form='unformatted', file=filename, 423 $ err=1000) 424c 425c Read in some basics to check if they are consistent with the calculation 426 read(unitno, err=1001, end=1001) nsets_read 427 read(unitno, err=1001, end=1001) nbf_read 428 read(unitno, err=1001, end=1001) mult_read 429c 430c Error checks 431 if ((nsets_read .ne. nsets) 432 & .or. (nbf_read .ne. nbf) 433 & .or. (mult_read .ne. mult) ) goto 1003 434c 435c Allocate the temporary buffer 436 if (.not.ma_push_get(mt_dbl,nbf,'dft_zora_read_so',l_zora,k_zora)) 437 $ call errquit('dft_zora_read_so: ma failed', nbf, MA_ERR) 438c 439c Read in g_zora_sf 440 do iset = 1,2 ! 2 components 441 do i = 1,nbf 442 call sread(unitno, dbl_mb(k_zora), nbf) 443 call ga_put(g_zora_sf(iset), 1, nbf, i, i, dbl_mb(k_zora), 1) 444 end do 445 end do 446c 447c Read in g_zora_scale_sf 448 do iset = 1,2 ! 2 components 449 do i = 1,nbf 450 call sread(unitno, dbl_mb(k_zora), nbf) 451 call ga_put(g_zora_scale_sf(iset), 1, nbf, i, i, 452 & dbl_mb(k_zora), 1) 453 end do 454 end do 455c 456c Read in g_zora_so 457 do iset = 1,3 ! 3 components: x,y,z 458 do i = 1,nbf 459 call sread(unitno, dbl_mb(k_zora), nbf) 460 call ga_put(g_zora_so(iset), 1, nbf, i, i, dbl_mb(k_zora), 1) 461 end do 462 end do 463c 464c Read in g_zora_scale_so 465 do iset = 1,3 ! 3 components: x,y,z 466 do i = 1,nbf 467 call sread(unitno, dbl_mb(k_zora), nbf) 468 call ga_put(g_zora_scale_so(iset), 1, nbf, i, i, 469 & dbl_mb(k_zora), 1) 470 end do 471 end do 472c 473c Close the file 474 close(unitno,err=1002) 475 ok = 1 476c 477c Deallocate the temporary buffer 478 if (.not. ma_pop_stack(l_zora)) call errquit 479 $ ('dft_zora_read_so: pop failed', l_zora, MA_ERR) 480c 481 end if 482c 483c Broadcast status to other nodes 484 10 call ga_brdcst(Msg_Vec_Stat+MSGINT, ok, inntsize, 0) ! Propagate status 485 call ga_sync() 486c 487 dft_zora_read_so = ok .eq. 1 488 return 489c 490 1000 write(6,*) 'dft_zora_read_so: failed to open ', 491 $ filename(1:inp_strlen(filename)) 492 call util_flush(luout) 493 ok = 0 494 goto 10 495c 496 1001 write(6,*) 'dft_zora_read_so: failed to read ', 497 $ filename(1:inp_strlen(filename)) 498 call util_flush(luout) 499 ok = 0 500 close(unitno,err=1002) 501 goto 10 502c 503 1003 write(6,*) 'dft_zora_read_so: file inconsistent with calculation', 504 $ filename(1:inp_strlen(filename)) 505 call util_flush(luout) 506 ok = 0 507 close(unitno,err=1002) 508 goto 10 509c 510 1002 write(6,*) 'dft_zora_read_so: failed to close', 511 $ filename(1:inp_strlen(filename)) 512 call util_flush(luout) 513 ok = 0 514 goto 10 515c 516 end 517c 518czora...Write out the zora atomic corrections to disk 519c 520 logical function dft_zora_write_so(rtdb, basis, filename, 521 & nbf, nsets, nmo, mult, 522 & g_zora_sf, g_zora_scale_sf, 523 & g_zora_so, g_zora_scale_so) 524c 525 implicit none 526c 527#include "errquit.fh" 528#include "mafdecls.fh" 529#include "global.fh" 530#include "tcgmsg.fh" 531#include "msgtypesf.h" 532#include "inp.fh" 533#include "msgids.fh" 534#include "cscfps.fh" 535#include "util.fh" 536#include "bas.fh" 537#include "geom.fh" 538#include "rtdb.fh" 539#include "stdio.fh" 540c 541c Temporary routine 542c 543 integer rtdb ! [input] RTDB handle (-1 if not accessible) 544 integer basis ! [input] Basis handle(-1 if not accessible) 545 character*(*) filename ! [input] File to write to 546 integer nbf ! [input] No. of functions in basis 547 integer nsets ! [input] No. of sets of matrices 548 integer nmo(nsets) ! [input] No. of mos in each set 549 integer mult 550c 551 integer g_zora_sf(2) 552 integer g_zora_scale_sf(2) 553 integer g_zora_so(3) 554 integer g_zora_scale_so(3) 555c 556 integer unitno 557 parameter (unitno = 77) 558c 559 integer lentit 560 integer lenbas 561 integer l_zora, k_zora 562 integer ok, iset, i, j 563 integer geom, ma_type, nelem 564 character*26 date 565 character*32 geomsum, basissum, key 566 character*20 scftype20 567 character*128 basis_name, trans_name 568 double precision energy, enrep 569 integer inntsize 570c 571 l_zora = -1 ! An invalid MA handle 572c 573 inntsize=MA_sizeof(MT_INT,1,MT_BYTE) 574 call ga_sync() 575c 576 ok = 0 577c 578c Read routines should be consistent with this 579c 580c Write out the atomic zora corrections 581c 582 if (ga_nodeid() .eq. 0) then 583c 584c Open the file 585 open(unitno, status='unknown', form='unformatted', 586 $ file=filename, err=1000) 587c 588c Write out the number of sets and basis functions 589 write(unitno, err=1001) nsets 590 write(unitno, err=1001) nbf 591 write(unitno, err=1001) mult 592c 593c Allocate the temporary buffer 594 if (.not.ma_push_get(mt_dbl,nbf,'dft_zora_write_so', 595 & l_zora,k_zora)) 596 & call errquit('dft_zora_write_so: ma failed', nbf, MA_ERR) 597c 598c Write out g_zora_sf 599 do iset = 1,2 ! 2 components 600 do i = 1, nbf 601 call ga_get(g_zora_sf(iset), 1, nbf, i, i, dbl_mb(k_zora),1) 602 call swrite(unitno, dbl_mb(k_zora), nbf) 603 end do 604 end do 605c 606c Write out g_zora_scale_sf 607 do iset = 1,2 ! 2 components 608 do i = 1, nbf 609 call ga_get(g_zora_scale_sf(iset), 1, nbf, i, i, 610 & dbl_mb(k_zora),1) 611 call swrite(unitno, dbl_mb(k_zora), nbf) 612 end do 613 end do 614c 615c Write out g_zora_so 616 do iset = 1,3 ! 3 components: x,y,z 617 do i = 1, nbf 618 call ga_get(g_zora_so(iset), 1, nbf, i, i, dbl_mb(k_zora),1) 619 call swrite(unitno, dbl_mb(k_zora), nbf) 620 end do 621 end do 622c 623c Write out g_zora_scale_so 624 do iset = 1,3 ! 3 components: x,y,z 625 do i = 1, nbf 626 call ga_get(g_zora_scale_so(iset), 1, nbf, i, i, 627 & dbl_mb(k_zora),1) 628 call swrite(unitno, dbl_mb(k_zora), nbf) 629 end do 630 end do 631c 632c Deallocate the temporary buffer 633 if (.not. ma_pop_stack(l_zora)) 634 $ call errquit('dft_zora_write_so: ma pop failed', l_zora, MA_ERR) 635c 636c Close the file 637 close(unitno,err=1002) 638c 639 ok = 1 640 end if 641c 642c Broadcast status to other nodes 643 10 call ga_brdcst(Msg_Vec_Stat+MSGINT, ok, inntsize, 0) ! Propagate status 644 call ga_sync() 645c 646 dft_zora_write_so = (ok .eq. 1) 647 if (ga_nodeid() .eq. 0) then 648 write(6,22) filename(1:inp_strlen(filename)) 649 22 format(/' Wrote atomic ZORA corrections to ',a/) 650 call util_flush(luout) 651 endif 652c 653 return 654c 655 1000 write(6,*) 'dft_zora_write_so: failed to open ', 656 $ filename(1:inp_strlen(filename)) 657 call util_flush(luout) 658 ok = 0 659 goto 10 660c 661 1001 write(6,*) 'dft_zora_write_so: failed to write ', 662 $ filename(1:inp_strlen(filename)) 663 call util_flush(luout) 664 ok = 0 665 close(unitno,err=1002) 666 goto 10 667c 668 1002 write(6,*) 'dft_zora_write_so: failed to close', 669 $ filename(1:inp_strlen(filename)) 670 call util_flush(luout) 671 ok = 0 672 goto 10 673c 674 end 675c $Id$ 676 subroutine dft_zora_nogaga(iorb,nbf_ao,nbf_mo,noc,g_moso, 677 C soxyz, 678 D zora_scale_sf,zora_scale_so, 679 E ener_tot, 680 D dens_so, 681 O vtmp,vecs) 682 implicit none 683 integer iorb 684 integer nbf_mo,nbf_ao ! [in] nbf_mo=2*nbf_ao 685 integer noc 686 double precision zora_scale_sf(nbf_ao,nbf_ao,2)! replicated[in] 687 double precision zora_scale_so(nbf_ao,nbf_ao,3)! replicated[in] 688 integer g_moso(2) ! [in] 689 double precision vecs(nbf_mo,2) ! [scratch] 690 double precision dens_so(nbf_mo,nbf_mo,2) ! repl [scratch] 691 double precision vtmp(nbf_ao,nbf_ao) ! repl [scratch] 692 double precision ener_tot ! [out] 693 character*1 soxyz(3) 694 695c 696 integer xyz 697 integer j,i 698c 699 external ddot 700 double precision ddot 701c 702 703 call ga_get(g_moso(1),1,nbf_mo,iorb,iorb,vecs(1,1),1) 704 call ga_get(g_moso(2),1,nbf_mo,iorb,iorb,vecs(1,2),1) 705 call dftv_densm_so(dens_so,vecs,nbf_mo,noc) 706c ener_sf_1 707c strided copy 708 do j=1,nbf_ao 709 call dcopy(nbf_ao,dens_so(1,j,1),1,vtmp(1,j),1) 710 enddo 711 ener_tot=ddot(nbf_ao*nbf_ao, 712 D vtmp,1,zora_scale_sf(1,1,1),1) 713cstrided copy 714 do j=1+nbf_ao,2*nbf_ao 715 call dcopy(nbf_ao,dens_so(1+nbf_ao,j,1),1,vtmp(1,j-nbf_ao),1) 716 enddo 717c ener_sf_2 718 ener_tot=ener_tot+ddot(nbf_ao*nbf_ao, 719 D vtmp,1,zora_scale_sf(1,1,2),1) 720 do xyz=1,3 721 call vec_dens_so(vtmp,dens_so,nbf_ao,soxyz(xyz)) 722 ener_tot=ener_tot+ddot(nbf_ao*nbf_ao, 723 & vtmp,1,zora_scale_so(1,1,xyz),1) 724 end do 725 return 726 end 727