comp-chem/NWChem
1
0
Fork 0
mirror of https://github.com/nwchemgit/nwchem.git synced 2026-07-21 06:25:21 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 4

Merge pull request #1178 from dmejiar/indicators

Browse source 
Indicators
This commit is contained in:
NWChem: Open Source High-Performance Computational Chemistry 2025-10-21 11:17:20 -07:00 • committed by GitHub
commit 566bd12603
No known key found for this signature in database
GPG key ID: B5690EEEBB952194
9 changed files with 199 additions and 45 deletions
Download patch file Download diff file Expand all files Collapse all files

23
src/bsemol/bse_davidson.F
View file

@ -203,7 +203,7 @@
idavidson = 0
! Restarts should land here
1 continue
c 1 continue
! K-orthogonalization
call bse_davidson_kvec(pars,g_trials,g_wia,g_x,npoles,
@ -252,7 +252,7 @@
call nga_release_block_grid(g_x,(/pars%me,0/))
! Starting vectors
call ga_copy(g_trials,g_x)
1 call ga_copy(g_trials,g_x)
call bse_davidson_kvec(pars,g_x,g_wia,g_y,npoles,
$ ntrials,pars%singlet,1d0)
call bse_davidson_kvec(pars,g_y,g_wia,g_ym,npoles,
@ -459,6 +459,16 @@
$ dbl_mb(k_tmp3),dbl_mb(vkv),dbl_mb(k_tmp4),myld+myld2,
$ ndim,mdim,.false.)
! All trial vectors were linearly-dependent
! restart calculation with current best guess
if (mdim.lt.1) then
nrestart = min(ntrials, pars%nroots*pars%nspace)
call bse_davidson_restart(ntrials,nrestart,npoles,
$ g_x,g_y,g_ym,g_s,g_ks,g_mks,g_trials,l_vkv,vkv,
$ int_mb(k_map),pars%nblocks)
goto 1
endif
if(pars%ipol.eq.1) then
call dcopy(mdim*myld,dbl_mb(k_tmp1),1,dbl_mb(k_w),1)
call dcopy(mdim*myld,dbl_mb(k_tmp2),1,dbl_mb(k_kw),1)
@ -489,15 +499,6 @@
call nga_release_update_block_grid(g_kw,(/pars%me,0/))
call nga_release_update_block_grid(g_mkw,(/pars%me,0/))
! All trial vectors were linearly-dependent
! restart calculation with current best guess
if (mdim.lt.3) then
nrestart = min(ntrials, pars%nroots*pars%nspace)
call bse_davidson_restart(ntrials,nrestart,npoles,
$ g_x,g_y,g_ym,g_s,g_ks,g_mks,g_trials,l_vkv,vkv,
$ int_mb(k_map),pars%nblocks)
goto 1
endif
ndim = mdim

42
src/bsemol/bse_davidson_kortho.F
View file

@ -17,6 +17,7 @@
double precision vkv(n,n),eig(n),work(npoles,n)
double precision sqrtvkv(n,n)
! V**T K V
call dgemm('t','n',n,n,npoles,1d0,v,npoles,kv,npoles,0d0,vkv,n)
call ga_dgop((/14/),vkv,n*n,'+')
@ -36,21 +37,36 @@
$ call errquit('bse_kortho: failed to deallocate work',199,MA_ERR)
! Remove small eigenpairs
m = n
k = 1
100 continue
do i=k,m
if (eig(i).lt.1d-9) then
do j=i+1,n
eig(j-1) = eig(j)
vkv(:,j-1) = vkv(:,j)
enddo
m = m - 1
k = i
goto 100
endif
m = 0
do i=n,1,-1
if(eig(i).lt.1d-10) exit
m = m + 1
enddo
if(m.ne.n) then
j = 1
do i=n-m+1,n
eig(j) = eig(i)
vkv(:,j) = vkv(:,i)
j = j + 1
enddo
endif
c m = n
c k = 1
c 100 continue
c do i=k,m
c if (eig(i).lt.1d-9) then
c do j=i+1,n
c eig(j-1) = eig(j)
c vkv(:,j-1) = vkv(:,j)
c enddo
c m = m - 1
c k = i
c goto 100
c endif
c enddo
! Obtain VKV^(-1/2)
do i=1,m
vkv(:,i) = vkv(:,i)/dsqrt(eig(i))

2
src/bsemol/bse_init.F
View file

@ -121,7 +121,7 @@ c
! Maximum Number of roots
if(.not.rtdb_get(pars%rtdb,'bse:nmax',mt_int,1,pars%nmax))
& call errquit(pname//'could not read from rtdb',123,RTDB_ERR)
if ((pars%nmax.le.5*pars%nroots).and.pars%davidson) then
if ((pars%nmax.lt.5*pars%nroots).and.pars%davidson) then
call errquit(pname//'maxroots should be at least 5*nroots',
& 126,INPUT_ERR)
endif

17
src/bsemol/bse_oscstr.F
View file

@ -48,7 +48,7 @@
parameter (thresh=5.0d-2)
double precision xe(3), xh(3), xe2_exp, xh2_exp, xe_exp_2
double precision xh_exp_2, xhxe_exp, dhtoe, dexc, sigma_h, sigma_e
double precision omega_norm, cov, pcc, dcd
double precision omega_norm, cov, pcc, dcd1, dcd2, dcd3
spin(1)='alpha'
spin(2)='beta '
@ -191,7 +191,7 @@
call ga_transpose(g_td(1), g_tdtot(1))
if(pars%ipol.gt.1) call ga_transpose(g_td(2), g_tdtot(2))
if(pars%davidson) then
if(.not.pars%tda) then
call ga_copy_patch('n',pars%g_amb,1,pars%npoles(1),
I ipole,ipole,
& g_y(1),1,pars%npoles(1),1,1)
@ -278,7 +278,9 @@
dexc = dsqrt(xe2_exp + xh2_exp - 2d0*xhxe_exp)
cov = xhxe_exp - xe(1)*xh(1) - xe(2)*xh(2) - xe(3)*xh(3)
pcc = cov/sigma_e/sigma_h
dcd = dhtoe - 0.5d0*(sigma_e + sigma_h)
dcd1 = dhtoe + dabs(sigma_h - sigma_e)
dcd2 = dhtoe - 0.5d0*(sigma_e + sigma_h)
dcd3 = dhtoe + dexc
if (ltransden) then
call util_file_name('tdens',.false.,.false.,fn_transden)
@ -369,12 +371,13 @@
endif
write(luout,*)
write(luout,9340) omega_norm
write(luout,9300) sigma_e*cau2ang
write(luout,9310) sigma_h*cau2ang
write(luout,9320) dhtoe*cau2ang
write(luout,9330) dexc*cau2ang
write(luout,9360) dcd*cau2ang
write(luout,9360) dcd1*cau2ang
write(luout,9370) dcd2*cau2ang
write(luout,9380) dcd3*cau2ang
write(luout,9350) pcc
@ -506,7 +509,9 @@
9310 format(5x,' sigma_h = ',f8.4, ' Angstrom')
9320 format(5x,' |d_h->e| = ',f8.4, ' Angstrom')
9330 format(5x,' d_exc = ',f8.4, ' Angstrom')
9360 format(5x,' d_CD = ',f8.4, ' Angstrom')
9360 format(5x,' d_CD1 = ',f8.4, ' Angstrom')
9370 format(5x,' d_CD2 = ',f8.4, ' Angstrom')
9380 format(5x,' d_CD3 = ',f8.4, ' Angstrom')
9350 format(5x,'Pearson C.C. = ',f8.4)
9340 format(5x,'Exciton norm = ',f8.4)

2
src/ddscf/int_dip_ga.F
View file

@ -239,7 +239,7 @@ cdbg I ga_nodeid(),'QDISTRIB',g_xx,i0,i1,j0,j1
!
! get the integrals we want
!
call int_mpole(jbas, ishell, ibas, jshell,
call int_mpole(jbas, jshell, ibas, ishell,
$ 2,
$ center,
$ mscratch, dbl_mb(k_scr), max1e, dbl_mb(k_buf))

2
src/libext/libxc/build_libxc.sh
View file

@ -77,7 +77,7 @@ fi
CMAKE_VER_MAJ=$(${CMAKE} --version|cut -d " " -f 3|head -1|cut -d. -f1)
CMAKE_VER_MIN=$(${CMAKE} --version|cut -d " " -f 3|head -1|cut -d. -f2)
echo CMAKE_VER is ${CMAKE_VER_MAJ} ${CMAKE_VER_MIN}
if ((CMAKE_VER_MAJ < 3)) || (((CMAKE_VER_MAJ > 2) && (CMAKE_VER_MIN < 8))); then
if ((CMAKE_VER_MAJ < 3)) || (((CMAKE_VER_MAJ == 2) && (CMAKE_VER_MIN < 24))); then
get_cmake_release $cmake_instdir
status=$?
if [ $status -ne 0 ]; then

2
src/libext/scalapack/build_scalapa.sh
View file

@ -93,7 +93,7 @@ fi
CMAKE_VER_MAJ=$(${CMAKE} --version|cut -d " " -f 3|head -1|cut -d. -f1)
CMAKE_VER_MIN=$(${CMAKE} --version|cut -d " " -f 3|head -1|cut -d. -f2)
echo CMAKE_VER is ${CMAKE_VER_MAJ} ${CMAKE_VER_MIN}
if ((CMAKE_VER_MAJ < 3)) || (((CMAKE_VER_MAJ > 2) && (CMAKE_VER_MIN < 99))); then
if ((CMAKE_VER_MAJ < 3)) || (((CMAKE_VER_MAJ == 3) && (CMAKE_VER_MIN < 24))); then
cmake_instdir=../libext_utils
get_cmake_release $cmake_instdir
status=$?

152
src/nwdft/lr_tddft/tddft_analysis.F
View file

@ -42,6 +42,7 @@ c
#include "rtdb.fh"
#include "sym.fh"
#include "util.fh"
#include "util_params.fh"
#include "msgids.fh"
#include "stdio.fh"
#include "dra.fh"
@ -108,6 +109,8 @@ c
integer g_smat ! GA handle for AO overlap matrix
integer g_corr ! GA handle for alpha-beta MO correlation matrix
integer g_work
integer g_work2, g_work3
integer g_qdrxx, g_qdrxy, g_qdrxz, g_qdryy, g_qdryz, g_qdrzz
c
integer g_x1(2) ! GA handle for X matrix of 1 root
integer g_y1(2) ! GA handle for Y matrix of 1 root
@ -181,6 +184,8 @@ c are hardwired here:
logical debug, fakegiao, lquad, status
logical grid_clinit, xc_gotxc
external grid_clinit, xc_gotxc
integer ga_create_atom_blocked
external ga_create_atom_blocked
c
double precision zero, one, two, three, half, third
parameter (zero = 0d0, one=1d0, two=2d0, three=3d0,
@ -200,6 +205,13 @@ c
logical lstores2
double precision s2_tmp(nroots)
integer ioserr
double precision omega_norm1,omega_norm2,sigma_h,sigma_e,dhtoe
double precision dexc,cov,pcc
double precision dcd1, dcd2, dcd3
double precision xe(3),xe_exp_2,xe2_exp
double precision xh(3),xh_exp_2,xh2_exp
double precision xhxe_exp
c
c MN solvation models -->
c
@ -349,11 +361,33 @@ c == Create and get the dipole matrix element ga: g_dipole ==
$ call errquit(pname//'ga_create failed g_dipole',0, GA_ERR)
call ga_zero(g_dipole(icomp))
end do ! icomp
status = ga_create(MT_DBL, nbf_ao, nbf_ao, 'g_qdrxx',
$ nbf_ao, 0, g_qdrxx)
status = ga_create(MT_DBL, nbf_ao, nbf_ao, 'g_qdrxy',
$ nbf_ao, 0, g_qdrxy)
status = ga_create(MT_DBL, nbf_ao, nbf_ao, 'g_qdrxz',
$ nbf_ao, 0, g_qdrxz)
status = ga_create(MT_DBL, nbf_ao, nbf_ao, 'g_qdryy',
$ nbf_ao, 0, g_qdryy)
status = ga_create(MT_DBL, nbf_ao, nbf_ao, 'g_qdryz',
$ nbf_ao, 0, g_qdryz)
status = ga_create(MT_DBL, nbf_ao, nbf_ao, 'g_qdrzz',
$ nbf_ao, 0, g_qdrzz)
c
c == Compute the dipole integrals for the 3 components ==
call int_dip_ga(ao_bas_han, ao_bas_han,
$ g_dipole(1), g_dipole(2), g_dipole(3))
call int_qdr_ga(ao_bas_han, ao_bas_han, g_qdrxx, g_qdrxy,
$ g_qdrxz, g_qdryy, g_qdryz, g_qdrzz)
c
call ga_add(1d0,g_qdrxx,1d0,g_qdryy,g_qdrxx)
call ga_add(1d0,g_qdrxx,1d0,g_qdrzz,g_qdrxx)
if(.not.(ga_destroy(g_qdrxy).and.ga_destroy(g_qdrxz).and.
$ ga_destroy(g_qdryz).and.ga_destroy(g_qdryy).and.
$ ga_destroy(g_qdrzz)))
$ call errquit('tddft_analysis: deallocate quadrupole',0,GA_ERR)
end if
c
if (cdspectrum .or. oscstr) then
@ -555,14 +589,18 @@ c ----------------------------------
c Get alpha-beta orbital correlation
c ----------------------------------
c
if (.not.ga_create(mt_dbl,nbf_ao,nbf_ao,
1 'AO overlap',-1,-1,g_ovlp))
2 call errquit('tddft_analysis: failed to create g_ovlp',0,
& GA_ERR)
call ga_zero(g_ovlp)
call int_1e_ga(ao_bas_han,ao_bas_han,g_ovlp,'overlap',oskel)
if (oskel) call sym_symmetrize(geom,ao_bas_han,.false.,g_ovlp)
if (ipol.eq.2) then
if (.not.ga_create(mt_dbl,nmo(1),nmo(2),
1 'MO correlation',-1,-1,g_corr))
2 call errquit('tddft_analysis: failed to create g_corr',0,
& GA_ERR)
if (.not.ga_create(mt_dbl,nbf_ao,nbf_ao,
1 'AO overlap',-1,-1,g_ovlp))
2 call errquit('tddft_analysis: failed to create g_ovlp',0,
& GA_ERR)
if (.not.ga_create(mt_dbl,nbf_ao,nmo(2),'work',-1,-1,g_work))
1 call errquit('tddft_analysis: failed to create g_work',0,
@ -570,22 +608,23 @@ c
if (.not.ma_push_get(mt_dbl,nmo(1)*nmo(2),'corr',
1 l_corr,k_corr)) call errquit
2 ('tddft_analysis: failed to allocate corr',0, MA_ERR)
call ga_zero(g_ovlp)
call int_1e_ga(ao_bas_han,ao_bas_han,g_ovlp,'overlap',oskel)
if (oskel) call sym_symmetrize
1 (geom,ao_bas_han,.false.,g_ovlp)
call ga_dgemm('N','N',nbf_ao,nmo(2),nbf_ao,1.0d0,
1 g_ovlp,g_movecs(2),0.0d0,g_work)
call ga_dgemm('T','N',nmo(1),nmo(2),nbf_ao,1.0d0,
1 g_movecs(1),g_work,0.0d0,g_corr)
if (.not.ga_destroy(g_work)) call errquit
1 ('tddft_analysis: failed to destroy g_work',0, GA_ERR)
if (.not.ga_destroy(g_ovlp)) call errquit
1 ('tddft_analysis: failed to destroy g_ovlp',0, GA_ERR)
! if (.not.ga_destroy(g_ovlp)) call errquit
! 1 ('tddft_analysis: failed to destroy g_ovlp',0, GA_ERR)
if (util_print('excited state',print_debug))
1 call ga_print(g_corr)
call ga_get(g_corr,1,nmo(1),1,nmo(2),dbl_mb(k_corr),nmo(1))
endif
if(.not.ga_duplicate(g_ovlp, g_work2, 'g_work'))
$ call errquit('tddft_analysis: failed to allocate work',0,GA_ERR)
if(.not.ga_duplicate(g_ovlp, g_work3, 'g_work'))
$ call errquit('tddft_analysis: failed to allocate work',0,GA_ERR)
c
c ------------
c Tamm-Dancoff
@ -1202,6 +1241,70 @@ c
endif ! .not.tda
c
c excitonic indicators
c
if(ipol.gt.1)
& call ga_add(1d0,g_tdtot(1),1d0,g_tdtot(2),g_tdtot(1))
call ga_dgemm('n','n',nbf_ao,nbf_ao,nbf_ao,1d0,g_tdtot(1),
& g_ovlp, 0d0, g_work2)
call ga_dgemm('n','t',nbf_ao,nbf_ao,nbf_ao,1d0,g_work2,
& g_tdtot(1),0d0,g_work3)
omega_norm1 = ga_ddot(g_work3, g_ovlp)
xe(1) = ga_ddot(g_work3, g_dipole(1))/omega_norm1
xe(2) = ga_ddot(g_work3, g_dipole(2))/omega_norm1
xe(3) = ga_ddot(g_work3, g_dipole(3))/omega_norm1
xe_exp_2 = xe(1)**2 + xe(2)**2 + xe(3)**2
xe2_exp = ga_ddot(g_work3, g_qdrxx)/omega_norm1
call ga_dgemm('t','n',nbf_ao,nbf_ao,nbf_ao,1d0,g_tdtot(1),
& g_ovlp, 0d0, g_work2)
call ga_dgemm('n','n',nbf_ao,nbf_ao,nbf_ao,1d0,g_work2,
& g_tdtot(1),0d0,g_work3)
omega_norm2 = ga_ddot(g_work3, g_ovlp)
xh(1) = ga_ddot(g_work3, g_dipole(1))/omega_norm2
xh(2) = ga_ddot(g_work3, g_dipole(2))/omega_norm2
xh(3) = ga_ddot(g_work3, g_dipole(3))/omega_norm2
xh_exp_2 = xh(1)**2 + xh(2)**2 + xh(3)**2
xh2_exp = ga_ddot(g_work3, g_qdrxx)/omega_norm2
xhxe_exp = 0d0
call ga_dgemm('n','n',nbf_ao,nbf_ao,nbf_ao,1d0,g_tdtot(1),
$ g_dipole(1),0d0,g_work2)
call ga_dgemm('n','t',nbf_ao,nbf_ao,nbf_ao,1d0,g_work2,
$ g_tdtot(1),0d0,g_work3)
xhxe_exp = xhxe_exp + ga_ddot(g_work3,g_dipole(1))/
$ sqrt(omega_norm1*omega_norm2)
call ga_dgemm('n','n',nbf_ao,nbf_ao,nbf_ao,1d0,g_tdtot(1),
$ g_dipole(2),0d0,g_work2)
call ga_dgemm('n','t',nbf_ao,nbf_ao,nbf_ao,1d0,g_work2,
$ g_tdtot(1),0d0,g_work3)
xhxe_exp = xhxe_exp + ga_ddot(g_work3,g_dipole(2))/
$ sqrt(omega_norm1*omega_norm2)
call ga_dgemm('n','n',nbf_ao,nbf_ao,nbf_ao,1d0,g_tdtot(1),
$ g_dipole(3),0d0,g_work2)
call ga_dgemm('n','t',nbf_ao,nbf_ao,nbf_ao,1d0,g_work2,
$ g_tdtot(1),0d0,g_work3)
xhxe_exp = xhxe_exp + ga_ddot(g_work3,g_dipole(3))/
$ sqrt(omega_norm1*omega_norm2)
sigma_h = dsqrt(xh2_exp - xh_exp_2)
sigma_e = dsqrt(xe2_exp - xe_exp_2)
dhtoe = dsqrt((xe(1)-xh(1))**2 +
$ (xe(2)-xh(2))**2 +
$ (xe(3)-xh(3))**2)
dexc = dsqrt(xe2_exp + xh2_exp - 2d0*xhxe_exp)
cov = xhxe_exp - xe(1)*xh(1) - xe(2)*xh(2) - xe(3)*xh(3)
pcc = cov/sigma_e/sigma_h
dcd1 = dhtoe + dabs(sigma_h - sigma_e)
dcd2 = dhtoe - 0.5d0*(sigma_e + sigma_h)
dcd3 = dhtoe + dexc
if(ipol.gt.1)
& call ga_add(1d0,g_tdtot(1),-1d0,g_tdtot(2),g_tdtot(1))
c
c write out the transition density matrix
c
if (ltransden) then
@ -1263,6 +1366,16 @@ c
2 tmom(8),tmom(9),tmom(10),
2 osc_str(1),osc_str(2),osc_str(3),
2 osc_str_tot
c
write(luout,*)
write(luout,9400) sigma_e*cau2ang
write(luout,9410) sigma_h*cau2ang
write(luout,9420) dhtoe*cau2ang
write(luout,9430) dexc*cau2ang
write(luout,9460) dcd1*cau2ang
write(luout,9470) dcd2*cau2ang
write(luout,9480) dcd3*cau2ang
write(luout,9450) pcc
c
c ------------------------
@ -1665,6 +1778,8 @@ c Clean up arrays used for CD spectra
end do
if (.not.ga_destroy(g_dipmag)) call
& errquit(pname//'ga_destroy failed g_dipmag',0,GA_ERR)
if(.not.ga_destroy(g_qdrxx))
& call errquit(pname//'ga_destroy quadrupole',0,GA_ERR)
end if
if (cdspectrum .and. lgiao) then
@ -1737,6 +1852,16 @@ c
& f8.4,' Debye')
9254 format(5x,'Total Ground State Dipole Moment = ',
& f8.4,' Debye')
9400 format(5x,' sigma_e = ',f8.4, ' Angstrom')
9410 format(5x,' sigma_h = ',f8.4, ' Angstrom')
9420 format(5x,' |d_h->e| = ',f8.4, ' Angstrom')
9430 format(5x,' d_exc = ',f8.4, ' Angstrom')
9460 format(5x,' d_CD1 = ',f8.4, ' Angstrom')
9470 format(5x,' d_CD2 = ',f8.4, ' Angstrom')
9480 format(5x,' d_CD3 = ',f8.4, ' Angstrom')
9450 format(5x,'Pearson C.C. = ',f8.4)
9440 format(5x,'Exciton norm = ',f8.4)
c
c -----------
c Target root
@ -1980,6 +2105,13 @@ c
if (.not.ga_destroy(g_tdtot(i))) call errquit
2 ('tddft_analysis: failed to destroy g_tdtot',0, GA_ERR)
enddo
c
if (.not.ga_destroy(g_ovlp)) call errquit
2 ('tddft_analysis: failed to destroy g_ovlp',0, GA_ERR)
if (.not.ga_destroy(g_work2)) call errquit
2 ('tddft_analysis: failed to destroy g_ovlp',0, GA_ERR)
if (.not.ga_destroy(g_work3)) call errquit
2 ('tddft_analysis: failed to destroy g_ovlp',0, GA_ERR)
c
return
end

2
src/tools/GNUmakefile
View file

@ -769,7 +769,7 @@ $(STAMP_DIR):
clean:
@test -e $(BUILDDIR)/Makefile && { cd $(BUILDDIR) && $(MAKE) clean; } || echo "Not configured"
rm -f build install
rm -rf build install
rm -f ../mpi_include.txt ../ga*txt
realclean: