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

added H2Se test for x2c

Browse source
This commit is contained in:
Jochen Autschbach 2022-03-29 15:27:59 -04:00 • committed by edoapra
parent 91827084a0
commit 44eca95c4f
No known key found for this signature in database
GPG key ID: 48E12DA1EDE9E1B0
4 changed files with 2185 additions and 40 deletions
Download patch file Download diff file Expand all files Collapse all files

3
QA/doafewqmtests.mpi
View file

@ -55,3 +55,6 @@ fi
./runtests.mpi.unix procs $np h2o2-response-uhf-damping
./runtests.mpi.unix procs $np h2o2-prop-notrans
#
#--- relativity tests
./runtests.mpi.unix procs $np x2c-h2se
#

159
QA/tests/x2c-h2se/x2c-h2se.nw Normal file
View file

@ -0,0 +1,159 @@
echo
start testjob
title efgjob
memory total 2000 mb stack 1000 mb
geometry noautoz units bohr nocenter
symmetry c1
Se 0.0000 0.0000 0.0564
H 0.0000 1.0575 -0.9594
H 0.0000 -1.0575 -0.9594
end
BASIS "ao basis" spherical PRINT
# Sapporo-DZP uncontracted
#BASIS SET: (4s,3p) -> [4s,3p]
H S
13.0107003 1.00000000E+00
H S
1.9622571 1.00000000E+00
H S
0.4445380 1.00000000E+00
H S
0.1219496 1.00000000E+00
H P
2.8606730 1.00000000E+00
H P
0.8881900 1.00000000E+00
H P
0.3037760 1.00000000E+00
#BASIS SET: (17s,14p,9d) -> [17s,14p,9d]
Se S
241498.6382409 1.00000000E+00
Se S
36223.7408736 1.00000000E+00
Se S
8243.4547006 1.00000000E+00
Se S
2332.0080828 1.00000000E+00
Se S
760.4702669 1.00000000E+00
Se S
578.5796313 1.00000000E+00
Se S
272.0485848 1.00000000E+00
Se S
179.7486759 1.00000000E+00
Se S
100.7295011 1.00000000E+00
Se S
28.6642790 1.00000000E+00
Se S
19.9031660 1.00000000E+00
Se S
11.8693555 1.00000000E+00
Se S
3.3851284 1.00000000E+00
Se S
2.0728541 1.00000000E+00
Se S
1.4483918 1.00000000E+00
Se S
0.3352235 1.00000000E+00
Se S
0.1251479 1.00000000E+00
Se P
4987.1768459 1.00000000E+00
Se P
1182.1765262 1.00000000E+00
Se P
381.9907152 1.00000000E+00
Se P
144.3524457 1.00000000E+00
Se P
59.8682607 1.00000000E+00
Se P
54.3041566 1.00000000E+00
Se P
26.1507531 1.00000000E+00
Se P
11.1948834 1.00000000E+00
Se P
7.9587015 1.00000000E+00
Se P
3.3433343 1.00000000E+00
Se P
1.3700986 1.00000000E+00
Se P
0.5441071 1.00000000E+00
Se P
0.2107083 1.00000000E+00
Se P
0.0794330 1.00000000E+00
Se D
237.0548228 1.00000000E+00
Se D
70.8504951 1.00000000E+00
Se D
26.4912435 1.00000000E+00
Se D
10.9539115 1.00000000E+00
Se D
4.6411645 1.00000000E+00
Se D
1.9254421 1.00000000E+00
Se D
0.7322816 1.00000000E+00
Se D
0.5653420 1.00000000E+00
Se D
0.2082120 1.00000000E+00
End
charge 0
relativistic
x2c on
#zora on
#douglas-kroll on
zora:cutoff_EFG 1d-10
zora:cutoff 1d-30
#clight 1e6
clight 137.0359895
end
set lindep:n_dep 0
# set dft:no_prune T
set dft:tol_rho 1d-30
set int:cando_txs f
set dft:job_grid_acc 1d-20
dft
grid xfine
direct
#xc slater 1.0 becke88 nonlocal 1.0 lyp 1.0
xc hfexch
convergence energy 1d-10
iterations 50
#print "quadrature"
end
set lindep:n_dep 0
#set sodft:scalar true
#set sodft:spinprojection 0. 0. 1.
#set sodft:collinear true
#set sodft:bfield 0. 1e-4 0.
#set sodft:xcallinone T
#set sodft:alternativexc T
set sodft:debugenergy T
#set sodft:x2c T
set x2c:analyticpvp T
#set x2c:efg T
task dft

1981
QA/tests/x2c-h2se/x2c-h2se.out Normal file
View file

File diff suppressed because it is too large Load diff

82
src/nwdft/x2c/calc_x2c_1e_scalar_ga.F
View file

@ -105,7 +105,7 @@ c other local variables:
integer dbg
double precision jfac(1),kfac(1)
character*32 pname
character*13 pname
c =================================================================
@ -114,7 +114,7 @@ c dbg = 0: no debug output
dbg=1
me=ga_nodeid()
master = me.eq.0
pname = "calc_x2c_1e_scalar"
pname = "x2c_1e_scalar"
if (dbg>0.and.master) write(luout,*) 'entering '//pname
@ -148,18 +148,18 @@ c make sure the input GAs have the correct dimension
call ga_inquire(g_x2c,type,nc,nc)
if (nc.ne.nbf_ao) then
call errquit('x2c1e: g_x2c has wrong dimension',nc,GA_ERR)
call errquit(pname//': g_x2c has wrong dimension',nc,GA_ERR)
end if
if (.not. bas_numbf(cbas,ntemp))
& call errquit('x2c1e: bas_ao error',ntemp, BASIS_ERR)
& call errquit(pname//': bas_ao error',ntemp, BASIS_ERR)
if (ntemp.ne.nc)
& call errquit('x2c1e: nc.ne.ntemp',ntemp, BASIS_ERR)
& call errquit(pname//': nc.ne.ntemp',ntemp, BASIS_ERR)
if (.not. bas_numbf(ubas,nu))
& call errquit('x2c1e: nu error',nu, BASIS_ERR)
& call errquit(pname//': nu error',nu, BASIS_ERR)
if (master.and.dbg>0) write (luout,*)
& 'x2c: nu, nc =',nu,nc
& pname//': nu, nc =',nu,nc
if (nu.lt.nc) call errquit(
& 'x2c uncontracted basis smaller than contracted', 1,
@ -214,7 +214,7 @@ c in GA, or call relham routines with MA arrays
if (.not.rtdb_get(rtdb, 'x2c:usega', mt_log, 1,
& use_ga)) continue
if (master.and.use_ga) write(luout,*)
& 'x2c: will perform decoupling in GA mode'
& pname//': will perform decoupling in GA mode'
c -----------------------------------------
c allocate GAs for T, V, S, and for W = pVp
@ -257,7 +257,7 @@ c symmetrize, just in case ...
call ga_sync()
if (master.and.dbg>0) write(luout,*) 'x2c: T,V,S done'
if (master.and.dbg>0) write(luout,*) pname//': T,V,S done'
c ------------------------------------------------------
c later we need a transformation from the contracted (c)
@ -283,7 +283,7 @@ c but this is how it is done in the DKH routines)
if (.not. ga_create(MT_DBL,nu,nu,'temp',1,1,g_temp))
& call errquit(pname//': error creating temp',53,GA_ERR)
if (master.and.dbg>0) write(luout,*) 'x2c: u2c mem done'
if (master.and.dbg>0) write(luout,*) pname//': u2c mem done'
call int_1e_ga(ubas, cbas, g_smix, 'overlap' ,oskel)
call ga_copy(g_s,g_temp)
@ -292,12 +292,12 @@ c but this is how it is done in the DKH routines)
& g_temp, g_smix, 0.0d0, g_u2c)
if (.not. ga_destroy(g_temp))
& call errquit('x2c1e: error destroy temp',53, GA_ERR)
& call errquit(pname//': error destroy temp',53, GA_ERR)
if (.not. ga_destroy(g_smix))
& call errquit('x2c1e: error destroy smix',52, GA_ERR)
& call errquit(pname//': error destroy smix',52, GA_ERR)
if (master.and.dbg>0) write(luout,*) 'x2c: u2c done'
if (master.and.dbg>0) write(luout,*) pname//': u2c done'
end if ! needu2c
@ -314,21 +314,21 @@ c Use [(A^T A)**(-1) A^T] A = 1 to define the left-inverse.
if (.not. ga_create(MT_DBL,nc,nc,'temp',1,1,g_temp))
& call errquit(pname//': failed creating temp',55,GA_ERR)
if (master.and.dbg>0) write(luout,*) 'x2c: c2u mem done'
if (master.and.dbg>0) write(luout,*) pname//': c2u mem done'
call ga_dgemm('t','n',nc,nc,nu,1.0d0,g_u2c,g_u2c,
& 0.0d0,g_temp)
call ga_matpow(g_temp, -1.0d0, 1d-12)
if (master.and.dbg>0) write(luout,*) 'x2c: c2u matpow done'
if (master.and.dbg>0) write(luout,*) pname//': c2u matpow done'
call ga_dgemm('n','t',nc,nu,nc,1.0d0,g_temp,g_u2c,
& 0.0d0,g_c2u)
if (.not. ga_destroy(g_temp))
& call errquit('x2c1e: error destroy temp',55, GA_ERR)
& call errquit(pname//': error destroy temp',55, GA_ERR)
if (master.and.dbg>0) write(luout,*) 'x2c: c2u done'
if (master.and.dbg>0) write(luout,*) pname//': c2u done'
c check c2u * u2c
@ -357,9 +357,9 @@ c c2u*u2c must be a unit matrix of dimension nc
end if
if (.not. ga_destroy(g_temp))
& call errquit('x2c1e: error destroy temp',56, GA_ERR)
& call errquit(pname//': error destroy temp',56, GA_ERR)
if (master.and.dbg>0) write(luout,*) 'x2c: c2u done'
if (master.and.dbg>0) write(luout,*) pname//': c2u done'
end if ! needc2u
@ -381,12 +381,12 @@ c P(u) = u2c P(c) u2c^T.
if (.not. ga_create(MT_DBL,nu,nu,'densu',1,1,g_densu))
& call errquit(pname//': error creating densu',58,GA_ERR)
if (master.and.dbg>0) write(luout,*) 'x2c: dens c2u start'
if (master.and.dbg>0) write(luout,*) pname//': dens c2u start'
call ga_zero(g_densu)
call ga_dgemm('n','t',nc,nu,nc, 1.0d0,
& g_dens(1),g_u2c,0.0d0,g_temp)
if (master.and.dbg>0) write(luout,*) 'x2c: dens c2u step1'
if (master.and.dbg>0) write(luout,*) pname//': dens c2u step1'
call ga_dgemm('n','n',nu,nu,nc ,1.0d0,
& g_u2c, g_temp, 0.0d0,g_densu)
call ga_symmetrize(g_densu)
@ -394,7 +394,7 @@ c P(u) = u2c P(c) u2c^T.
if (.not. ga_destroy(g_temp))
& call errquit(pname//': error destroy temp',57, GA_ERR)
if (master.and.dbg>0) write(luout,*) 'x2c: dens c2u done'
if (master.and.dbg>0) write(luout,*) pname//': dens c2u done'
else
@ -419,7 +419,7 @@ c ------------------------------------
& errquit(pname//': ga_create g_temp failed',59,GA_ERR)
call ga_zero(g_temp)
if (master.and.dbg>0) write(luout,*) 'x2c: debug S:'
if (master.and.dbg>0) write(luout,*) pname//': debug S:'
call ga_dgemm('n','n',nu,nu,nu, 1.0d0,
& g_densu,g_s,0.0d0,g_temp)
@ -427,9 +427,9 @@ c ------------------------------------
sum = ga_trace_diag(g_temp)
if (master) write (luout,'(1x,a,1x,e15.7)')
& 'calc_x2c: tr(P * S) : ',sum
& pname//': tr(P * S) : ',sum
if (master.and.dbg>0) write(luout,*) 'x2c: debug T:'
if (master.and.dbg>0) write(luout,*) pname//': debug T:'
call ga_dgemm('n','n',nu,nu,nu, 1.0d0,
& g_densu,g_t,0.0d0,g_temp)
@ -437,10 +437,10 @@ c ------------------------------------
sum = ga_trace_diag(g_temp)
if (master) write (luout,'(1x,a,1x,e15.7)')
& 'calc_x2c: tr(P * T) : ',sum
& pname//': tr(P * T) : ',sum
energy = energy + sum
if (master.and.dbg>0) write(luout,*) 'x2c: debug V:'
if (master.and.dbg>0) write(luout,*) pname//': debug V:'
call ga_dgemm('n','n',nu,nu,nu, 1.0d0,
& g_densu,g_v,0.0d0,g_temp)
@ -448,14 +448,14 @@ c ------------------------------------
sum = ga_trace_diag(g_temp)
if (master) write (luout,'(1x,a,1x,e15.7)')
& 'calc_x2c: tr(P * V) : ',sum
& pname//': tr(P * V) : ',sum
energy = energy + sum
if (.not. ga_destroy(g_temp)) call errquit
& (pname//': temp ga corrupt?',59, GA_ERR)
if (master) write (luout,'(1x,a,1x,e15.7)')
& 'calc_x2c: 1e energy : ',energy
& pname//': 1e energy : ',energy
c ----------------------------------------------------------
c as this routine gets embedded in the scalar DFT framework,
@ -536,7 +536,7 @@ c ---------------------------------------------------------------
if (analytic_pvp) then
if (master) write (luout,*)
& 'x2c: calculating analytic p.Vp integrals'
& pname//': calculating analytic p.Vp integrals'
c add scalar part p.Vp to g_w:
call int_1e_ga(ubas, ubas ,g_w, 'pvp' ,oskel)
@ -552,10 +552,11 @@ c integrals are done analytically
c --------------------------------------------------
if (analytic_pvp .and. .not.modelpot_2e) then
if (master) write(luout,*)'x2c: skipping numerical W integrals'
if (master) write(luout,*)
& pname//': skipping numerical W integrals'
else
if (master) write(luout,*) 'x2c: calculating W numerically'
if (master) write(luout,*) pname//': calculating W numerically'
if (modelpot_2e .and. master)
& write(luout,*) ' (2e model pot. only)'
@ -575,14 +576,14 @@ c g_vxc has the right dimension already
c transform to ubas: tmpv = c2u^T g_tmpv c2u, add to g_w
if(.not.ga_create(mt_dbl,nc,nu,'temp',0,0,g_temp))
& call errquit('x2c1e: error creating temp',148, GA_ERR)
& call errquit(pname//': error creating temp',148, GA_ERR)
call ga_zero(g_temp)
call ga_dgemm('n','n',nc,nu,nc, 1.0d0,
& g_tmpv,g_c2u,0.0d0,g_temp)
if (master.and.dbg>0) write(luout,*) 'x2c: g_w c2u step1'
if (master.and.dbg>0) write(luout,*) pname//': g_w c2u step1'
call ga_dgemm('t','n',nu,nu,nc ,1.0d0,
& g_c2u, g_temp, 1.0d0,g_w)
if (master.and.dbg>0) write(luout,*) 'x2c: g_w c2u step2'
if (master.and.dbg>0) write(luout,*) pname//': g_w c2u step2'
if (.not. ga_destroy(g_tmpv)) call errquit
& (pname//': ga corrupt?',147, GA_ERR)
if (.not. ga_destroy(g_temp)) call errquit
@ -604,7 +605,7 @@ c sanity checks:
if (.not.needu2c) then
if (g_X.ne.g_x2c) call errquit(
& 'x2c1e: g_X/g_x2c inconsistency 1',g_X-g_x2c, UNKNOWN_ERR)
& pname//': g_X/g_x2c inconsistency 1',g_X-g_x2c, UNKNOWN_ERR)
end if
c -----------------------------------------------
@ -918,9 +919,10 @@ c store in g_x2c. Otherwise, we are already done
c ----------------------------------------------
if (needu2c) then
if (dbg>0.and.master) write(luout,*) 'x2c: u2c trafo of result'
if (dbg>0.and.master) write(luout,*)
& pname//': u2c trafo of result'
if (.not. ga_create(MT_DBL,nc,nu,'temp',1,1,g_temp))
& call errquit('x2c1e: ga_create failed temp',77,
& call errquit(pname//': ga_create failed temp',77,
& GA_ERR)
if (g_X.eq.g_x2c) call errquit(
& pname//': g_X/g_x2c inconsistency 3',g_X-g_x2c, UNKNOWN_ERR)
@ -933,7 +935,7 @@ c same. Hence, g_x2c is over-written here
& g_temp,g_u2c, 0.0d0,g_x2c)
if (.not. ga_destroy(g_temp))
& call errquit('x2c1e: ga_destroy failed temp',77,GA_ERR)
& call errquit(pname//': ga_destroy failed temp',77,GA_ERR)
end if ! needu2c
call ga_symmetrize(g_x2c)
@ -942,7 +944,7 @@ c de-allocate g_X, if needed, and related arrays
if (needu2c) then
if (.not. ga_destroy(g_X))
& call errquit('x2c1e: ga_destroy failed g_X',61,GA_ERR)
& call errquit(pname//': ga_destroy failed g_X',61,GA_ERR)
if (.not.keepU) then
if (.not. ga_destroy(g_u2c)) call errquit