Start BSE branch

This commit is contained in:
Daniel Mejia-Rodriguez 2022-08-11 11:15:06 -07:00 committed by edoapra
parent 554aef3e30
commit e74b2e5920
No known key found for this signature in database
GPG key ID: 9E6A0B70826967BA
8 changed files with 373 additions and 8 deletions

32
src/bsemol/GNUmakefile Normal file
View file

@ -0,0 +1,32 @@
OBJ = bsemol.o bse_input.o
USES_BLAS =
OBJ_OPTIMIZE =
LIBRARY = libbsemol.a
HEADERS = ../nwdft/include/cdft.fh ../nwdft/grid/cgridfile.fh
LIB_DEFINES =
LIB_INCLUDES = -I./
include ../config/makefile.h
ifeq ($(_USE_SCALAPACK),Y)
_USE_ELPA = $(shell ${GA_PATH}/bin/ga-config --use_elpa| awk ' /1/ {print "Y"}')
ifeq ($(_USE_ELPA),Y)
LIB_DEFINES += -DUSE_ELPA
endif
#assuming we have GAPZHEEVD
GOTGAPZHEEVD= Y
ifeq ($(GOTGAPZHEEVD),Y)
LIB_DEFINES += -DUSE_GAPZHEEVD
endif
_PDSYEVR = $(shell [ -e ${NWCHEM_TOP}/src/tools/install/lib/libga.a ] && nm ${NWCHEM_TOP}/src/tools/install/lib/libga.a|egrep T |egrep ga_pdsyevr |awk '/ga_pdsyevr/ {print "Y"}')
ifeq ($(_PDSYEVR),Y)
LIB_DEFINES += -DUSE_PDSYEVR
endif
endif
include ../config/makelib.h

35
src/bsemol/bse.fh Normal file
View file

@ -0,0 +1,35 @@
double precision, parameter :: ha2ev = 27.211386245988d0
type :: bse_params_t
sequence
integer :: rtdb
integer :: me, nprocs, iMaxThreads
integer :: ipol
logical :: main, odft
c
integer :: natoms
integer :: nmo
integer :: nocc(2),nvir(2),npoles(2)
integer :: noqp(2),nvqp(2),nstates(2),mynpoles(2)
integer :: oolo(2),oohi(2)
integer :: ovlo(2),ovhi(2)
integer :: vvlo(2),vvhi(2)
integer :: llmo(2),nqp(2)
integer :: ngl
double precision :: tol2e, eta
c
integer :: l_evals, k_evals
integer :: l_occ, k_occ
c
integer :: fd
c
integer :: g_movecs(2), g_moeris, g_sols
integer :: g_erioo(2), g_eriov(2), g_erivv(2), g_erim
logical :: noio
c
logical :: ncap, davidson, lanczos
double precision :: thresh, fdstep, vxddp(2), vxddm(2)
end type bse_params_t

76
src/bsemol/bse_input.F Normal file
View file

@ -0,0 +1,76 @@
subroutine bse_input(rtdb)
c
c Input parser for molecular BSE module
c
c Daniel Mejia-Rodriguez
c November 2022
c
implicit none
c
#include "errquit.fh"
#include "util.fh"
#include "inp.fh"
#include "rtdb.fh"
#include "mafdecls.fh"
#include "stdio.fh"
c
integer rtdb
character*255 test
character(*), parameter :: pname = "bse_input: "
c
integer, parameter :: num_dirs = 2
character*11 dirs(num_dirs)
data dirs / 'method', 'end' /
c
integer ind, qpiter,eviter,noqp(2),nvqp(2), ngl, first
double precision threshold, eta
double precision, parameter :: ha2ev = 27.211386245988d0
c
c Set default values and put them in RTDB
c
call bse_defaults(rtdb)
c
c Parse the input file
c
10 if (.not. inp_read()) then
call errquit(pname//'failed to read input',0,INPUT_ERR)
endif
c
if (.not. inp_a(test)) then
call errquit(pname//'failed to read keyword',0,INPUT_ERR)
endif
c
if (.not.inp_match(num_dirs, .false., test, dirs, ind)) then
call errquit(pname//'unknown directive', 0, INPUT_ERR)
end if
c
goto (100, 9999) ind
call errquit(pname//'unknown directive',ind, INPUT_ERR)
c
c -----------------------------------------------------------------
c
c method
c
100 if (.not. inp_a(test))
& call errquit(pname//'could not read method', 100, INPUT_ERR)
if (inp_compare(.false.,'analytic',test)) then
if (.not.rtdb_put(rtdb, 'bse:analytic', mt_log, 1, .true.))
& call errquit(pname//'rtdb_put failed', 101, RTDB_ERR)
elseif (inp_compare(.false.,'bse:davidson',test)) then
if (.not.rtdb_put(rtdb, 'bse:davidson', mt_log, 1, .true.))
& call errquit(pname//'rtdb_put failed', 102, RTDB_ERR)
elseif (inp_compare(.false.,'bse:lanczos',test)) then
if (.not.rtdb_put(rtdb, 'bse:davidson', mt_log, 1, .true.))
& call errquit(pname//'rtdb_put failed', 103, RTDB_ERR)
else
call errquit(pname//'unknown BSE method',104,INPUT_ERR)
endif
goto 10
c
c normal termination
c
9999 return
c
end subroutine bse_input

73
src/bsemol/bsemol.F Normal file
View file

@ -0,0 +1,73 @@
logical function bsemol(rtdb)
implicit none
#include "rtdb.fh"
integer rtdb
c
logical bse_doit
external bse doit
c
logical status
c
status = rtdb_cput(rtdb,'dft:theory', 1, 'dft')
if (status) status = bse_doit(rtdb)
bsemol = status
call grid_cleanup(.false.)
c
end function
c
c -----------------------------------------------------------------
c
logical function bse_doit(rtdb)
implicit none
#include "errquit.fh"
#include "global.fh"
#include "rtdb.fh"
#include "stdio.fh"
#include "util.fh"
#include "bse.fh"
integer, intent(in) :: rtdb
character(*), parameter :: pname = "bse_doit: "
type(bse_params_t) :: params
character*80 theory
integer me
c
c -----------------------------------------------------------------
c
c Initialize
c
bse_doit = .false.
c
c Determine the theory we are using
c
if (.not. rtdb_cget(rtdb, "dft:theory", 1, theory))
$ call errquit(pname//"theory not specified",0,RTDB_ERR)
c
c Initialize GW calculation
c
params%rtdb = rtdb
params%me = ga_nodeid()
c
call bse_init(params)
c
if (theory.eq."dft") then
if (params%davidson) then
call bse_davidson(params)
else
call bse_analytic(params)
endif
c
elseif (theory.eq.'sodft') then
call errquit(pname//"spin-orbit not implemented",0,0)
else
call errquit(pname//"invalid theory: "//trim(theory),0,0)
endif
c
call bse_finalize(params)
bse_doit = .true.
c
c -----------------------------------------------------------------
c
end function

View file

@ -21,7 +21,7 @@ KNOWN_MODULE_SUBDIRS = NWints atomscf ddscf develop gradients moints nwdft \
property solvation nwpw fft rimp2_grad python argos analyz diana nwmd cafe space \
drdy vscf uccsdt qmmm rism qmd etrans tce geninterface bq mm cons\
perfm dntmc smd nbo dangchang leps ccca lucia dimqm rdmft fcidump\
gwmol xtb nwc_columbus
gwmol xtb nwc_columbus bsemol
# These are the libraries for the high-level modules. They should be
# specified in an order that will link correctly, but that shouldn't
@ -36,7 +36,7 @@ KNOWN_MODULE_LIBS = -ltest -lccsd -llucia -lmcscf -lselci -lmp2 -lmoints \
-lcafe -lspace -lanalyze -ldiana -lqhop \
-lpfft -ldplot -lnwpython -ldrdy -lvscf -lqmmm -lrism -lqmd -letrans\
-lpspw -ltce -lgeninterface -lbq -lmm -lcons -lperfm -lnbo -ldntmc -lccca\
-ldimqm -lrdmft -lfcidump -lgwmol -lnwxtb
-ldimqm -lrdmft -lfcidump -lgwmol -lnwxtb -lbsemol
# This include file handles configuration of the NW_MODULE_SUBDIRS and
# NW_MODULE_LIBS macros for what we actually want to build. It
@ -77,22 +77,22 @@ EXP_alldev = develop $(EXP_all)
EXP_pnnldev = develop $(EXP_pnnl)
EXP_qm = NWints atomscf ddscf gradients moints nwdft rimp2 hessian bq cons \
optim stepper driver ccsd vib mcscf property solvation dplot mp2_grad esp etrans tce vscf dimqm \
gwmol
gwmol bsemol
EXP_smallqm = NWints atomscf ddscf gradients moints nwdft rimp2 hessian bq cons \
optim stepper driver vib qmd mcscf property solvation dplot \
mp2_grad esp etrans vscf nwmd dimqm gwmol
mp2_grad esp etrans vscf nwmd dimqm gwmol bsemol
EXP_smallqmmm = NWints atomscf ddscf gradients moints nwdft rimp2 hessian bq cons \
optim stepper driver vib qmd mcscf property solvation dplot \
mp2_grad esp etrans vscf nwmd qmmm prepar analyz dimqm gwmol
mp2_grad esp etrans vscf nwmd qmmm prepar analyz dimqm gwmol bsemol
EXP_notce = NWints atomscf ddscf gradients moints nwdft rimp2 hessian bq cons \
optim stepper driver ccsd vib mcscf property solvation dplot mp2_grad esp etrans vscf dimqm \
gwmol
gwmol bsemol
EXP_tinyqmpw = nwdft property solvation driver ccsd mp2_grad hessian vib nwpw
EXP_ver2 = prepar gradients mp2_grad stepper vib mcscf \
optim driver analyz diana
#plane_wave ... broken? RJH
EXP_qmandpw = nwpw NWints atomscf ddscf gradients moints nwdft rimp2 hessian bq cons \
optim stepper driver ccsd vib mcscf property solvation dplot mp2_grad esp etrans vscf dimqm gwmol
optim stepper driver ccsd vib mcscf property solvation dplot mp2_grad esp etrans vscf dimqm gwmol bsemol
EXP_tps = prepar nwmd analyz diana
EXP_md = prepar nwmd analyz ddscf nwdft gradients esp
EXP_classical = prepar nwmd analyz
@ -187,6 +187,7 @@ MODULES := $(patsubst nwpw,nwpw fft,$(MODULES) )
MODULES := $(patsubst ddscf,ddscf atomscf nwdft NWints gradients,$(MODULES) )
MODULES := $(patsubst gwmol, gwmol nwdft ddscf gradients hessian atomscf NWints,$(MODULES) )
MODULES := $(patsubst bsemol, bsemol gwmol nwdft ddscf gradients hessian atomscf NWints,$(MODULES) )
MODULES := $(patsubst lucia, lucia mp2_grad ddscf nwdft property cphf gradients hessian atomscf NWints moints,$(MODULES) )
@ -274,6 +275,7 @@ MOD_LIBS := $(patsubst dimqm,-ldimqm,$(MOD_LIBS) )
MOD_LIBS := $(patsubst rdmft,-lrdmft,$(MOD_LIBS) )
MOD_LIBS := $(patsubst fcidump,-lfcidump,$(MOD_LIBS) )
MOD_LIBS := $(patsubst gwmol,-lgwmol,$(MOD_LIBS) )
MOD_LIBS := $(patsubst bsemol,-lbsemol,$(MOD_LIBS) )
MOD_LIBS := $(patsubst xtb,-lnwxtb,$(MOD_LIBS) )
MOD_LIBS := $(patsubst nwc_columbus,-lnwc_columbus,$(MOD_LIBS) )
MOD_LIBS := $(strip $(MOD_LIBS))
@ -339,6 +341,7 @@ STUBROUTINES := $(patsubst solvation, solvation,$(STUBROUTINES))
STUBROUTINES := $(patsubst dimqm,dimqm_addop dimqm_addop_uhf_damp dimqm_addop_uhf dimqm_tidy dimqm_input dimpar_input dimqm_main fock_dim dimqm_prep dim_grid_quadv0 dimqm_lclfld dimqm_polarizability dimqm_broadcast dimqm_Dens2Dipoles dimqm_rohf_calcDIM dimqm_rohf_wrap dimqm_used dimqm_setgo dimqm_getgo dimqm_getenergy dimqm_setlfirst dimqm_getlrsp dimqm_lrsp dimqm_getlclfld dimqm_getgdipeli dimqm_energy dimqm_init,$(STUBROUTINES) )
#STUBROUTINES := $(patsubst cons, cons_end cons_input cons_add_egrad cons_add_energy cons_init cons_print_summary cons_active,$(STUBROUTINES) )
STUBROUTINES := $(patsubst gwmol, gw_input gwmol, $(STUBROUTINES) )
STUBROUTINES := $(patsubst bsemol, bse_input bsemol, $(STUBROUTINES) )
STUBROUTINES := $(patsubst xtb, xtb_input xtb_energy xtb_energy_gradient, $(STUBROUTINES) )
STUBROUTINES := $(patsubst nwc_columbus, sifcfg sifwh sodft_aoints, $(STUBROUTINES) )
STUBROUTINES := $(strip $(STUBROUTINES))

View file

@ -9,7 +9,7 @@
task_save_state.o task.o task_input.o task_shell_input.o task_num_grad.o \
task_vscf.o task_property.o task_dntmc.o task_bsse.o \
task_jefftce.o task_ncc.o task_rfcf.o \
task_gwmol.o \
task_gwmol.o task_bsemol.o
task_aoints.o task_nadtrc.o
USES_BLAS = task_num_grad.F task_bsse.F

View file

@ -163,6 +163,10 @@ c
logical task_nadtrc
external task_nadtrc
logical task_bsemol
external task_bsemol
c logical embed_optimize
c external embed_optimize
c
@ -468,6 +472,8 @@ c status=embed_prepare(rtdb)
status=task_dntmc(rtdb)
else if (inp_compare(.false., 'gw', operation)) then
status = task_gwmol(rtdb)
else if (inp_compare(.false., 'bse', operation)) then
status = task_bsemol(rtdb)
c
c **** check for theories that have specific operations ****
c

140
src/task/task_bsemol.F Normal file
View file

@ -0,0 +1,140 @@
c
c Task handler for Molecular BSE calcualtions
c Daniel Mejia-Rodriguez, 2022
c
logical function task_bsemol(rtdb)
c
implicit none
c
#include "errquit.fh"
#include "rtdb.fh"
#include "mafdecls.fh"
#include "inp.fh"
#include "stdio.fh"
#include "util.fh"
#include "global.fh"
c
integer rtdb
c
character*32 theory
character*255 title
integer iproc, ipol
logical oprint, rodft, dftmp2, status, converged, ignore
double precision cputime, walltime
c
logical bsemol, task_energy
external bsemol, task_energy
c
c call ecce_print_module_entry('task gwmol')
c
c Initialization
c
task_bsemol = .false.
c
c Time at the start
c
cputime = util_cpusec()
walltime = util_wallsec()
c
c Invalidate status
c
if (.not. rtdb_put( rtdb, 'task:status', mt_log, 1, .false.))
& call errquit('task_bse: failed to invalidate status',0,RTDB_ERR)
c
c Print BSE header
c
oprint = util_print('information', print_low)
if (oprint) then
if (.not. rtdb_cget(rtdb, 'title', 1, title)) title = ' '
if (ga_nodeid().eq.0) then
write(LuOut,*)
call util_print_centered(
& LuOut, 'NWChem BSE Module', 40, .true.)
write(LuOut,*)
write(LuOut,*)
if (title.ne.' ') then
call util_print_centered(LuOut, title, 40, .false.)
write(LuOut,*)
write(LuOut,*)
endif
endif
endif
c
c Check if theory has been specified
c
if (.not. rtdb_cget(rtdb, 'task:theory', 1, theory))
& call errquit('task_bsemol: theory not specified',0, RTDB_ERR)
c
c Restrict to DFT for now (SCF could also be used)
c
if (.not.inp_compare(.false.,'dft',theory))
& call errquit('task_bsemol: BSE is only implemented for DFT', 0,
& INPUT_ERR)
c
if (rtdb_get(rtdb, 'dft:rodft', mt_log, 1, rodft))
& call errquit('task_bsemol: BSE is not implemented for RODFT', 0,
& INPUT_ERR)
c
if (rtdb_get(rtdb, 'dft:dftmp2', mt_log, 1, dftmp2))
& call errquit('task_bsemol: BSE is not implemented for double'//
& ' hybrid functionals',0,INPUT_ERR)
c
if (.not.rtdb_get(rtdb, 'dft:ipol', mt_int, 1, ipol))
& call errquit('task_bsemol: failed to read ipol', 0, RTDB_ERR)
if (.not. rtdb_get(rtdb, 'task:ignore', mt_log, 1, ignore))
& ignore = .true.
c
c Make sure a DFT solution is available
c
status = rtdb_get(rtdb, 'dft:converged', mt_log, 1, converged)
if (.not.status) then
status = task_energy(rtdb)
if (.not.rtdb_get(rtdb, 'dft:converged', mt_log, 1, converged))
& call errquit('task_gwmol: failed to read convergence info',0,
& RTDB_ERR)
endif
status = converged
if (ga_nodeid().eq.0) then
write(luout,*) " "
write(luout,'(29X,A)') "***********************"
write(luout,'(29X,A)') "*** BSE ***"
write(luout,'(29X,A)') "***********************"
write(luout,*) " "
if (.not.converged) then
write(LuOut,*)
write(LuOut,'(27X,A)') '!!! WARNING: !!!'
write(LuOut,'(27X,A)') '!!! DFT DID NOT CONVERGE !!!'
write(LuOut,*)
endif
call util_flush(luout)
endif
if ((.not.status) .and. (.not.ignore))
& call errquit('task_bsemol: Energy failure',0,CALC_ERR)
c
c Actually compute GW
c
status = bsemol(rtdb)
c
cputime = util_cpusec() - cputime
walltime = util_wallsec() - walltime
c
if (.not. rtdb_put(rtdb, 'task:status', mt_log, 1, .true.))
& call errquit('task_bse: failed to set status',0,
& RTDB_ERR)
if (.not. rtdb_put(rtdb, 'task:cputime', mt_dbl, 1, cputime))
& call errquit('task_bse: failed storing cputime',0,
& RTDB_ERR)
if (.not. rtdb_put(rtdb, 'task:walltime', mt_dbl, 1, walltime))
& call errquit('task_bse: failed storing walltime',0,
& RTDB_ERR)
c
task_bsemol = status
c
end function