LCPQ
/
quantum_package
mirror of https://github.com/LCPQ/quantum_package
10
0
Fork 0
Code Issues Releases Wiki Activity

Fix travis

This commit is contained in:
TApplencourt 2015-12-14 08:14:06 +01:00
parent 131805da27 c02a7ca735
commit 229b7c5bb2
153 changed files with 5933 additions and 3405 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
.travis.yml
View File

@ -27,7 +27,7 @@ python:
script:
- ./configure --production ./config/gfortran.cfg
- source ./quantum_package.rc
- qp_module.py install Full_CI Hartree_Fock CAS_SD MRCC_CASSD
- qp_install_module.py install Full_CI Hartree_Fock CAS_SD MRCC_CASSD
- ninja
- cd ocaml ; make ; cd -
- cd test ; bats bats/qp.bats

4
config/gfortran.cfg
View File

@ -10,7 +10,7 @@
#
#
[COMMON]
FC : gfortran -g -ffree-line-length-none -I .
FC : gfortran -g -ffree-line-length-none -I . -static-libgcc
LAPACK_LIB : -llapack -lblas
IRPF90 : irpf90
IRPF90_FLAGS : --ninja --align=32
@ -22,7 +22,7 @@ IRPF90_FLAGS : --ninja --align=32
# 0 : Deactivate
#
[OPTION]
MODE : OPT ; [ OPT | PROFILE | DEBUG ] : Chooses the section below
MODE : DEBUG ; [ OPT | PROFILE | DEBUG ] : Chooses the section below
CACHE : 1 ; Enable cache_compile.py
OPENMP : 1 ; Append OpenMP flags

71
configure vendored
View File

@ -33,8 +33,8 @@ import sys
from os.path import join
if not any(i in ["--production", "--development"] for i in sys.argv):
print __doc__
sys.exit()
sys.argv += ["--development"]
if len(sys.argv) != 3:
print __doc__
sys.exit()
@ -46,6 +46,7 @@ if len(sys.argv) != 3:
QP_ROOT = os.getcwd()
QP_ROOT_BIN = join(QP_ROOT, "bin")
QP_ROOT_LIB = join(QP_ROOT, "lib")
QP_ROOT_INSTALL = join(QP_ROOT, "install")
os.environ["PATH"] = os.environ["PATH"] + ":" + QP_ROOT_BIN
@ -62,8 +63,11 @@ d_dependency = {
"resultsFile": ["python"],
"emsl": ["python"],
"gcc": [],
"g++": [],
"zeromq" : [ "g++" ],
"f77zmq" : [ "zeromq", "python" ],
"python": [],
"ninja": ["gcc", "python"],
"ninja": ["g++", "python"],
"make": [],
"p_graphviz": ["python"],
"bats": []
@ -93,12 +97,12 @@ curl = Info(
zlib = Info(
url='http://zlib.net/zlib-1.2.8.tar.gz',
description=' zlib',
default_path=join(QP_ROOT_INSTALL, "zlib"))
default_path=join(QP_ROOT_LIB, "libz.a"))
path = Info(
patch = Info(
url='ftp://ftp.gnu.org/gnu/patch/patch-2.7.5.tar.gz',
description=' path',
default_path=join(QP_ROOT, "lib", "libz.a"))
description=' patch',
default_path=join(QP_ROOT_BIN, "patch"))
irpf90 = Info(
url='{head}/LCPQ/irpf90/{tail}'.format(**path_github),
@ -117,12 +121,11 @@ resultsFile = Info(
ninja = Info(
url='{head}/martine/ninja/{tail}'.format(**path_github),
description=' nina',
description=' ninja',
default_path=join(QP_ROOT_BIN, "ninja"))
emsl = Info(
url='{head}/LCPQ/EMSL_Basis_Set_Exchange_Local/{tail}'.format(**
path_github),
url='{head}/LCPQ/EMSL_Basis_Set_Exchange_Local/{tail}'.format(**path_github),
description=' EMSL basis set library',
default_path=join(QP_ROOT_INSTALL, "emsl"))
@ -131,6 +134,16 @@ ezfio = Info(
description=' EZFIO',
default_path=join(QP_ROOT_INSTALL, "EZFIO"))
zeromq = Info(
url='http://download.zeromq.org/zeromq-4.1.3.tar.gz',
description=' ZeroMQ',
default_path=join(QP_ROOT_LIB, "libzmq.a"))
f77zmq = Info(
url='{head}/zeromq/f77_zmq/{tail}'.format(**path_github),
description=' F77-ZeroMQ',
default_path=join(QP_ROOT_LIB, "libf77zmq.a"))
p_graphviz = Info(
url='https://github.com/xflr6/graphviz/archive/master.tar.gz',
description=' Python library for graphviz',
@ -143,8 +156,14 @@ bats = Info(
d_info = dict()
<<<<<<< HEAD
for m in ["ocaml", "m4", "curl", "zlib", "path", "irpf90", "docopt",
"resultsFile", "ninja", "emsl", "ezfio", "p_graphviz",'bats']:
=======
for m in ["ocaml", "m4", "curl", "zlib", "patch", "irpf90", "docopt",
"resultsFile", "ninja", "emsl", "ezfio", "p_graphviz",
"zeromq", "f77zmq" ]:
>>>>>>> LCPQ-master
exec ("d_info['{0}']={0}".format(m))
@ -195,8 +214,7 @@ def check_output(*popenargs, **kwargs):
def checking(d_dependency):
"""
For each key in d_dependency check if it
is avalabie or not
For each key in d_dependency check if it is avalabie
"""
def check_python():
@ -266,7 +284,7 @@ def checking(d_dependency):
l_installed = dict()
l_needed = []
# Check all the other
# Check all the others
length = max(map(len, d_dependency))
for i in d_dependency.keys():
@ -281,7 +299,7 @@ def checking(d_dependency):
l_needed.append(i)
print ""
# Expend the need_stuff for all the genealogy
# Expand the needed stuff for all the genealogy
l_install_descendant = get_list_descendant(d_dependency, l_installed,
l_needed)
@ -334,7 +352,7 @@ _|_ | | _> |_ (_| | | (_| |_ | (_) | |
d_print = {
"install_ninja": "Install ninja...",
"build": "Creating build.ninja...",
"install": "Installing the dependencies with Ninja..."
"install": "Installing the dependencies using Ninja..."
}
length = max(map(len, d_print.values()))
@ -378,7 +396,7 @@ _|_ | | _> |_ (_| | | (_| |_ | (_) | |
descr = d_info[need].description
default_path = d_info[need].default_path
# Build to dowload
# Build to download
l_build += ["build {0}: download".format(archive_path),
" url = {0}".format(url), " descr = {0}".format(descr),
""]
@ -410,7 +428,16 @@ _|_ | | _> |_ (_| | | (_| |_ | (_) | |
path_ninja = find_path("ninja", l_installed)
subprocess.check_call("cd install ;{0}".format(path_ninja), shell=True)
except:
raise
prefix = os.path.join('install', '_build')
for filename in os.listdir(prefix):
if filename.endswith(".log"):
with open( os.path.join(prefix,filename) ,'r') as f:
print "\n\n"
print "=-=-=-=-=-=- %s =-=-=-=-=-=-" %(filename)
print f.read()
print "=-=-=-=-=-=-=-=-=-=-=-=-=-=-\n\n"
print "Error in installation of dependencies"
sys.exit(1)
else:
print r"""
_________
@ -443,11 +470,12 @@ def create_ninja_and_rc(l_installed):
print str_info("qp_root"),
python_path = [join(QP_ROOT, "scripts"), join(QP_ROOT, "install")]
l_python = [join(QP_ROOT, "scripts")]
l_python = [join("${QP_ROOT}", "scripts")]
for dir_ in python_path:
for folder in os.listdir(dir_):
path = join(dir_, folder)
if os.path.isdir(path):
path = path.replace(QP_ROOT,"${QP_ROOT}")
l_python.append(path)
path_ezfio = find_path('ezfio', l_installed, var_for_qp_root=True)
@ -456,9 +484,9 @@ def create_ninja_and_rc(l_installed):
l_rc = [
'export QP_ROOT={0}'.format(QP_ROOT),
'export QP_EZFIO={0}'.format(path_ezfio),
'export IRPF90={0}'.format(path_irpf90),
'export NINJA={0}'.format(path_ninja),
'export QP_EZFIO={0}'.format(path_ezfio.replace(QP_ROOT,"${QP_ROOT}")),
'export IRPF90={0}'.format(path_irpf90.replace(QP_ROOT,"${QP_ROOT}")),
'export NINJA={0}'.format(path_ninja.replace(QP_ROOT,"${QP_ROOT}")),
'export QP_PYTHON={0}'.format(":".join(l_python)), "",
'export PYTHONPATH="${QP_EZFIO}/Python":"${QP_PYTHON}":"${PYTHONPATH}"',
'export PATH="${QP_PYTHON}":"${QP_ROOT}"/bin:"${QP_ROOT}"/ocaml:"${PATH}"',
@ -532,3 +560,4 @@ if __name__ == '__main__':
create_ninja_and_rc(l_installed)
recommendation()

2
install/scripts/build.sh
View File

@ -7,4 +7,4 @@ mkdir ${BUILD} || exit 1
tar -zxf Downloads/${TARGET}.tar.gz --strip-components=1 --directory=${BUILD} || exit 1
_install || exit 1
rm -rf -- ${BUILD} _build/${TARGET}.log
exit 0
exit 0

23
install/scripts/install_f77zmq.sh Executable file
View File

@ -0,0 +1,23 @@
#!/bin/bash -x
TARGET=f77zmq
function _install()
{
cd ..
QP_ROOT=$PWD
cd -
export C_INCLUDE_PATH="${C_INCLUDE_PATH}":"${QP_ROOT}"/lib
set -e
set -u
export ZMQ_H="${QP_ROOT}"/lib/zmq.h
cd "${BUILD}"
make -j 8 || exit 1
mv libf77zmq.a "${QP_ROOT}"/lib || exit 1
mv libf77zmq.so "${QP_ROOT}"/lib || exit 1
cp f77_zmq.h "${QP_ROOT}"/src/ZMQ/
cd -
return 0
}
source scripts/build.sh

27
install/scripts/install_zeromq.sh Executable file
View File

@ -0,0 +1,27 @@
#!/bin/bash -x
TARGET=zeromq
function _install()
{
cd ..
QP_ROOT=$PWD
cd -
export C_INCLUDE_PATH="${C_INCLUDE_PATH}":./
set -e
set -u
ORIG=$(pwd)
cd "${BUILD}"
./configure --without-libsodium || exit 1
make -j 8 || exit 1
rm -f -- "${QP_ROOT}"/lib/libzmq.a "${QP_ROOT}"/lib/libzmq.so "${QP_ROOT}"/lib/libzmq.so.5
cp .libs/libzmq.a "${QP_ROOT}"/lib
cp .libs/libzmq.so "${QP_ROOT}"/lib/libzmq.so.5
cp include/{zmq.h,zmq_utils.h} "${QP_ROOT}"/lib
cd "${QP_ROOT}"/lib
ln -s libzmq.so.5 libzmq.so
cd ${ORIG}
return 0
}
source scripts/build.sh

12
ocaml/Atom.ml
View File

@ -8,8 +8,8 @@ type t =
coord : Point3d.t ;
} with sexp
(** Read xyz coordinates of the atom with unit u *)
let of_string u s =
(** Read xyz coordinates of the atom *)
let of_string ~units s =
let buffer = s
|> String.split ~on:' '
|> List.filter ~f:(fun x -> x <> "")
@ -18,21 +18,21 @@ let of_string u s =
| [ name; charge; x; y; z ] ->
{ element = Element.of_string name ;
charge = Charge.of_string charge ;
coord = Point3d.of_string u (String.concat [x; y; z] ~sep:" ")
coord = Point3d.of_string ~units (String.concat [x; y; z] ~sep:" ")
}
| [ name; x; y; z ] ->
let e = Element.of_string name in
{ element = e ;
charge = Element.to_charge e;
coord = Point3d.of_string u (String.concat [x; y; z] ~sep:" ")
coord = Point3d.of_string ~units (String.concat [x; y; z] ~sep:" ")
}
| _ -> raise (AtomError s)
;;
let to_string u a =
let to_string ~units a =
[ Element.to_string a.element ;
Charge.to_string a.charge ;
Point3d.to_string u a.coord ]
Point3d.to_string ~units a.coord ]
|> String.concat ~sep:" "
;;

4
ocaml/Atom.mli
View File

@ -5,5 +5,5 @@ type t = { element : Element.t; charge : Charge.t; coord : Point3d.t; }
val t_of_sexp : Sexplib.Sexp.t -> t
val sexp_of_t : t -> Sexplib.Sexp.t
val of_string : Units.units -> string -> t
val to_string : Units.units -> t -> string
val of_string : units:Units.units -> string -> t
val to_string : units:Units.units -> t -> string

218
ocaml/Element.ml
View File

@ -1,4 +1,5 @@
open Core.Std;;
open Core.Std
open Qptypes
exception ElementError of string
@ -8,49 +9,49 @@ type t =
|Li|Be |B |C |N |O |F |Ne
|Na|Mg |Al|Si|P |S |Cl|Ar
|K |Ca|Sc|Ti|V |Cr|Mn|Fe|Co|Ni|Cu|Zn|Ga|Ge|As|Se|Br|Kr
with sexp;;
with sexp
let of_string x =
match (String.capitalize (String.lowercase x)) with
| "X" | "Dummy" -> X
| "H" | "Hydrogen" -> H
| "He" | "Helium" -> He
| "Li" | "Lithium" -> Li
| "Be" | "Beryllium" -> Be
| "B" | "Boron" -> B
| "C" | "Carbon" -> C
| "N" | "Nitrogen" -> N
| "O" | "Oxygen" -> O
| "F" | "Fluorine" -> F
| "Ne" | "Neon" -> Ne
| "Na" | "Sodium" -> Na
| "Mg" | "Magnesium" -> Mg
| "Al" | "Aluminum" -> Al
| "Si" | "Silicon" -> Si
| "P" | "Phosphorus" -> P
| "S" | "Sulfur" -> S
| "Cl" | "Chlorine" -> Cl
| "Ar" | "Argon" -> Ar
| "K" | "Potassium" -> K
| "Ca" | "Calcium" -> Ca
| "Sc" | "Scandium" -> Sc
| "Ti" | "Titanium" -> Ti
| "V" | "Vanadium" -> V
| "Cr" | "Chromium" -> Cr
| "Mn" | "Manganese" -> Mn
| "Fe" | "Iron" -> Fe
| "Co" | "Cobalt" -> Co
| "Ni" | "Nickel" -> Ni
| "Cu" | "Copper" -> Cu
| "Zn" | "Zinc" -> Zn
| "Ga" | "Gallium" -> Ga
| "Ge" | "Germanium" -> Ge
| "As" | "Arsenic" -> As
| "Se" | "Selenium" -> Se
| "Br" | "Bromine" -> Br
| "Kr" | "Krypton" -> Kr
| "X" | "Dummy" -> X
| "H" | "Hydrogen" -> H
| "He" | "Helium" -> He
| "Li" | "Lithium" -> Li
| "Be" | "Beryllium" -> Be
| "B" | "Boron" -> B
| "C" | "Carbon" -> C
| "N" | "Nitrogen" -> N
| "O" | "Oxygen" -> O
| "F" | "Fluorine" -> F
| "Ne" | "Neon" -> Ne
| "Na" | "Sodium" -> Na
| "Mg" | "Magnesium" -> Mg
| "Al" | "Aluminum" -> Al
| "Si" | "Silicon" -> Si
| "P" | "Phosphorus" -> P
| "S" | "Sulfur" -> S
| "Cl" | "Chlorine" -> Cl
| "Ar" | "Argon" -> Ar
| "K" | "Potassium" -> K
| "Ca" | "Calcium" -> Ca
| "Sc" | "Scandium" -> Sc
| "Ti" | "Titanium" -> Ti
| "V" | "Vanadium" -> V
| "Cr" | "Chromium" -> Cr
| "Mn" | "Manganese" -> Mn
| "Fe" | "Iron" -> Fe
| "Co" | "Cobalt" -> Co
| "Ni" | "Nickel" -> Ni
| "Cu" | "Copper" -> Cu
| "Zn" | "Zinc" -> Zn
| "Ga" | "Gallium" -> Ga
| "Ge" | "Germanium" -> Ge
| "As" | "Arsenic" -> As
| "Se" | "Selenium" -> Se
| "Br" | "Bromine" -> Br
| "Kr" | "Krypton" -> Kr
| x -> raise (ElementError ("Element "^x^" unknown"))
;;
let to_string = function
| X -> "X"
@ -90,7 +91,7 @@ let to_string = function
| Se -> "Se"
| Br -> "Br"
| Kr -> "Kr"
;;
let to_long_string = function
| X -> "Dummy"
@ -130,7 +131,7 @@ let to_long_string = function
| Se -> "Selenium"
| Br -> "Bromine"
| Kr -> "Krypton"
;;
let to_charge c =
let result = match c with
@ -172,7 +173,7 @@ let to_charge c =
| Br -> 35
| Kr -> 36
in Charge.of_int result
;;
let of_charge c = match (Charge.to_int c) with
| 0 -> X
@ -213,5 +214,134 @@ let of_charge c = match (Charge.to_int c) with
| 35 -> Br
| 36 -> Kr
| x -> raise (ElementError ("Element of charge "^(string_of_int x)^" unknown"))
;;
let covalent_radius x =
let result = function
| X -> 0.
| H -> 0.37
| He -> 0.70
| Li -> 1.23
| Be -> 0.89
| B -> 0.90
| C -> 0.85
| N -> 0.74
| O -> 0.74
| F -> 0.72
| Ne -> 0.70
| Na -> 1.00
| Mg -> 1.36
| Al -> 1.25
| Si -> 1.17
| P -> 1.10
| S -> 1.10
| Cl -> 0.99
| Ar -> 0.70
| K -> 2.03
| Ca -> 1.74
| Sc -> 1.44
| Ti -> 1.32
| V -> 1.22
| Cr -> 0.00
| Mn -> 1.16
| Fe -> 0.00
| Co -> 1.15
| Ni -> 1.17
| Cu -> 1.25
| Zn -> 1.25
| Ga -> 1.20
| Ge -> 1.21
| As -> 1.16
| Se -> 0.70
| Br -> 1.24
| Kr -> 1.91
in
Units.angstrom_to_bohr *. (result x)
|> Positive_float.of_float
let vdw_radius x =
let result = function
| X -> 0.
| H -> 1.20
| He -> 1.70
| Li -> 1.70
| Be -> 1.70
| B -> 1.70
| C -> 1.70
| N -> 1.55
| O -> 1.52
| F -> 1.47
| Ne -> 1.70
| Na -> 1.70
| Mg -> 1.70
| Al -> 1.94
| Si -> 2.10
| P -> 1.80
| S -> 1.80
| Cl -> 1.75
| Ar -> 1.70
| K -> 1.70
| Ca -> 1.70
| Sc -> 1.70
| Ti -> 1.70
| V -> 1.98
| Cr -> 1.94
| Mn -> 1.93
| Fe -> 1.93
| Co -> 1.92
| Ni -> 1.70
| Cu -> 1.70
| Zn -> 1.70
| Ga -> 2.02
| Ge -> 1.70
| As -> 1.96
| Se -> 1.70
| Br -> 2.10
| Kr -> 1.70
in
Units.angstrom_to_bohr *. (result x)
|> Positive_float.of_float
let mass x =
let result = function
| X -> 0.
| H -> 1.0079
| He -> 4.00260
| Li -> 6.941
| Be -> 9.01218
| B -> 10.81
| C -> 12.011
| N -> 14.0067
| O -> 15.9994
| F -> 18.998403
| Ne -> 20.179
| Na -> 22.98977
| Mg -> 24.305
| Al -> 26.98154
| Si -> 28.0855
| P -> 30.97376
| S -> 32.06
| Cl -> 35.453
| Ar -> 39.948
| K -> 39.0983
| Ca -> 40.08
| Sc -> 44.9559
| Ti -> 47.90
| V -> 50.9415
| Cr -> 51.996
| Mn -> 54.9380
| Fe -> 55.9332
| Co -> 58.9332
| Ni -> 58.70
| Cu -> 63.546
| Zn -> 65.38
| Ga -> 69.72
| Ge -> 72.59
| As -> 74.9216
| Se -> 78.96
| Br -> 79.904
| Kr -> 83.80
in
result x
|> Positive_float.of_float

4
ocaml/Element.mli
View File

@ -13,6 +13,8 @@ val of_string : string -> t
val to_string : t -> string
val to_long_string : t -> string
(** get the positive charge *)
(** Properties *)
val to_charge : t -> Charge.t
val of_charge : Charge.t -> t
val covalent_radius : t -> Qptypes.Positive_float.t
val vdw_radius : t -> Qptypes.Positive_float.t

31
ocaml/Input_determinants_by_hand.ml
View File

@ -7,8 +7,6 @@ module Determinants_by_hand : sig
{ n_int : N_int_number.t;
bit_kind : Bit_kind.t;
n_det : Det_number.t;
n_states : States_number.t;
n_states_diag : States_number.t;
expected_s2 : Positive_float.t;
psi_coef : Det_coef.t array;
psi_det : Determinant.t array;
@ -23,8 +21,6 @@ end = struct
{ n_int : N_int_number.t;
bit_kind : Bit_kind.t;
n_det : Det_number.t;
n_states : States_number.t;
n_states_diag : States_number.t;
expected_s2 : Positive_float.t;
psi_coef : Det_coef.t array;
psi_det : Determinant.t array;
@ -146,11 +142,12 @@ end = struct
|> Array.map ~f:Det_coef.of_float
;;
let write_psi_coef ~n_det ~n_states c =
let write_psi_coef ~n_det c =
let n_det = Det_number.to_int n_det
and c = Array.to_list c
|> List.map ~f:Det_coef.to_float
and n_states = States_number.to_int n_states
and n_states =
read_n_states () |> States_number.to_int
in
Ezfio.ezfio_array_of_list ~rank:2 ~dim:[| n_det ; n_states |] ~data:c
|> Ezfio.set_determinants_psi_coef
@ -214,8 +211,6 @@ end = struct
{ n_int = read_n_int () ;
bit_kind = read_bit_kind () ;
n_det = read_n_det () ;
n_states = read_n_states () ;
n_states_diag = read_n_states_diag () ;
expected_s2 = read_expected_s2 () ;
psi_coef = read_psi_coef () ;
psi_det = read_psi_det () ;
@ -227,8 +222,6 @@ end = struct
let write { n_int ;
bit_kind ;
n_det ;
n_states ;
n_states_diag ;
expected_s2 ;
psi_coef ;
psi_det ;
@ -236,10 +229,8 @@ end = struct
write_n_int n_int ;
write_bit_kind bit_kind;
write_n_det n_det;
write_n_states n_states;
write_n_states_diag ~n_states:n_states n_states_diag;
write_expected_s2 expected_s2;
write_psi_coef ~n_det:n_det psi_coef ~n_states:n_states;
write_psi_coef ~n_det:n_det psi_coef ;
write_psi_det ~n_int:n_int ~n_det:n_det psi_det;
;;
@ -249,7 +240,7 @@ end = struct
let mo_tot_num = MO_number.of_int mo_tot_num ~max:mo_tot_num in
let det_text =
let nstates =
States_number.to_int b.n_states
read_n_states () |> States_number.to_int
and ndet =
Det_number.to_int b.n_det
in
@ -284,12 +275,6 @@ If true, input the expected value of S^2 ::
expected_s2 = %s
Number of requested states, and number of states used for the
Davidson diagonalization ::
n_states = %s
n_states_diag = %s
Number of determinants ::
n_det = %s
@ -299,8 +284,6 @@ Determinants ::
%s
"
(b.expected_s2 |> Positive_float.to_string)
(b.n_states |> States_number.to_string)
(b.n_states_diag |> States_number.to_string)
(b.n_det |> Det_number.to_string)
det_text
|> Rst_string.of_string
@ -313,8 +296,6 @@ Determinants ::
n_int = %s
bit_kind = %s
n_det = %s
n_states = %s
n_states_diag = %s
expected_s2 = %s
psi_coef = %s
psi_det = %s
@ -322,8 +303,6 @@ psi_det = %s
(b.n_int |> N_int_number.to_string)
(b.bit_kind |> Bit_kind.to_string)
(b.n_det |> Det_number.to_string)
(b.n_states |> States_number.to_string)
(b.n_states_diag |> States_number.to_string)
(b.expected_s2 |> Positive_float.to_string)
(b.psi_coef |> Array.to_list |> List.map ~f:Det_coef.to_string
|> String.concat ~sep:", ")

4
ocaml/Input_nuclei.ml
View File

@ -147,7 +147,7 @@ nucl_coord = %s
(b.nucl_charge |> Array.to_list |> List.map
~f:(Charge.to_string) |> String.concat ~sep:", " )
(b.nucl_coord |> Array.to_list |> List.map
~f:(Point3d.to_string Units.Bohr) |> String.concat ~sep:"\n" )
~f:(Point3d.to_string ~units:Units.Bohr) |> String.concat ~sep:"\n" )
;;
@ -161,7 +161,7 @@ nucl_coord = %s
Printf.sprintf " %-3s %d %s"
(b.nucl_label.(i) |> Element.to_string)
(b.nucl_charge.(i) |> Charge.to_int )
(b.nucl_coord.(i) |> Point3d.to_string Units.Angstrom) )
(b.nucl_coord.(i) |> Point3d.to_string ~units:Units.Angstrom) )
) |> String.concat ~sep:"\n"
in
Printf.sprintf "

96
ocaml/Molecule.ml
View File

@ -10,27 +10,36 @@ type t = {
} with sexp
let get_charge { nuclei ; elec_alpha ; elec_beta } =
let result = (Elec_alpha_number.to_int elec_alpha) +
(Elec_beta_number.to_int elec_beta) in
let result =
(Elec_alpha_number.to_int elec_alpha) +
(Elec_beta_number.to_int elec_beta)
in
let rec nucl_charge = function
| a::rest -> (Charge.to_float a.Atom.charge) +. nucl_charge rest
| [] -> 0.
in
Charge.of_float (nucl_charge nuclei -. (Float.of_int result))
;;
let get_multiplicity m =
let elec_alpha = m.elec_alpha in
let elec_alpha =
m.elec_alpha
in
Multiplicity.of_alpha_beta elec_alpha m.elec_beta
;;
let get_nucl_num m =
let nmax = (List.length m.nuclei) in
let nmax =
List.length m.nuclei
in
Nucl_number.of_int nmax ~max:nmax
;;
let name m =
let cm = Charge.to_int (get_charge m) in
let cm =
get_charge m
|> Charge.to_int
in
let c =
match cm with
| 0 -> ""
@ -39,8 +48,12 @@ let name m =
| i when i>1 -> Printf.sprintf " (%d+)" i
| i -> Printf.sprintf " (%d-)" (-i)
in
let mult = Multiplicity.to_string (get_multiplicity m) in
let { nuclei ; elec_alpha ; elec_beta } = m in
let mult =
get_multiplicity m
|> Multiplicity.to_string
in
let { nuclei ; elec_alpha ; elec_beta } = m
in
let rec build_list accu = function
| a::rest ->
begin
@ -53,7 +66,9 @@ let name m =
in
let rec build_name accu = function
| (a, n)::rest ->
let a = Element.to_string a in
let a =
Element.to_string a
in
begin
match n with
| 1 -> build_name (a::accu) rest
@ -64,19 +79,25 @@ let name m =
end
| [] -> accu
in
let result = build_list [] nuclei |> build_name [c ; ", " ; mult]
let result =
build_list [] nuclei |> build_name [c ; ", " ; mult]
in
String.concat (result)
;;
let to_string m =
let { nuclei ; elec_alpha ; elec_beta } = m in
let n = List.length nuclei in
let title = name m in
[ Int.to_string n ; title ] @ (List.map ~f:(fun x -> Atom.to_string
Units.Angstrom x) nuclei)
let { nuclei ; elec_alpha ; elec_beta } = m
in
let n =
List.length nuclei
in
let title =
name m
in
[ Int.to_string n ; title ] @
(List.map ~f:(fun x -> Atom.to_string Units.Angstrom x) nuclei)
|> String.concat ~sep:"\n"
;;
let of_xyz_string
?(charge=(Charge.of_int 0)) ?(multiplicity=(Multiplicity.of_int 1))
@ -94,7 +115,9 @@ let of_xyz_string
) + 1 - (Charge.to_int charge)
|> Elec_number.of_int
in
let (na,nb) = Multiplicity.to_alpha_beta ne multiplicity in
let (na,nb) =
Multiplicity.to_alpha_beta ne multiplicity
in
let result =
{ nuclei = l ;
elec_alpha = na ;
@ -109,7 +132,7 @@ let of_xyz_string
raise (MultiplicityError msg);
else () ;
result
;;
let of_xyz_file
@ -121,8 +144,33 @@ let of_xyz_file
let (_,buffer) = String.lsplit2_exn buffer ~on:'\n' in
of_xyz_string ~charge:charge ~multiplicity:multiplicity
~units:units buffer
;;
include To_md5;;
let distance_matrix molecule =
let coord =
molecule.nuclei
|> List.map ~f:(fun x -> x.Atom.coord)
|> Array.of_list
in
let n =
Array.length coord
in
let result =
Array.make_matrix ~dimx:n ~dimy:n 0.
in
for i = 0 to (n-1)
do
for j = 0 to (n-1)
do
result.(i).(j) <- Point3d.distance coord.(i) coord.(j)
done;
done;
result
include To_md5
let to_md5 = to_md5 sexp_of_t
;;

4
ocaml/Molecule.mli
View File

@ -34,5 +34,9 @@ val of_xyz_string :
?multiplicity:Multiplicity.t ->
?units:Units.units -> string -> t
(** Creates the distance matrix between all the atoms *)
val distance_matrix :
t -> (float array) array
(** Computes the MD5 hash *)
val to_md5 : t -> Qptypes.MD5.t

25
ocaml/Point3d.ml
View File

@ -7,9 +7,16 @@ type t = {
z : float ;
} with sexp
let of_tuple ~units (x,y,z) =
let f = match units with
| Units.Bohr -> 1.
| Units.Angstrom -> Units.angstrom_to_bohr
in
{ x = x *. f ; y = y *. f ; z = z *. f }
(** Read x y z coordinates in string s with units u *)
let of_string u s =
let f = match u with
let of_string ~units s =
let f = match units with
| Units.Bohr -> 1.
| Units.Angstrom -> Units.angstrom_to_bohr
in
@ -22,7 +29,6 @@ let of_string u s =
{ x = l.(0) *. f ;
y = l.(1) *. f ;
z = l.(2) *. f }
;;
let distance2 p1 p2 =
@ -30,17 +36,18 @@ let distance2 p1 p2 =
and { x=x2 ; y=y2 ; z=z2 } = p2 in
(x2-.x1)*.(x2-.x1) +. (y2-.y1)*.(y2-.y1) +. (z2-.z1)*.(z2-.z1)
|> Positive_float.of_float
;;
let distance p1 p2 = sqrt (Positive_float.to_float (distance2 p1 p2))
;;
let to_string u p =
let f = match u with
let distance p1 p2 =
sqrt (Positive_float.to_float (distance2 p1 p2))
let to_string ~units p =
let f = match units with
| Units.Bohr -> 1.
| Units.Angstrom -> Units.bohr_to_angstrom
in
let { x=x ; y=y ; z=z } = p in
Printf.sprintf "%16.8f %16.8f %16.8f" (x*.f) (y*.f) (z*.f)
;;

7
ocaml/Point3d.mli
View File

@ -4,11 +4,14 @@ type t =
z : float;
} with sexp
(** Create from a tuple of floats *)
val of_tuple : units:Units.units -> float*float*float -> t
(** Create from an xyz string *)
val of_string : Units.units -> string -> t
val of_string : units:Units.units -> string -> t
(** Convert to a string for printing *)
val to_string : Units.units -> t -> string
val to_string : units:Units.units -> t -> string
(** Computes the squared distance between 2 points *)
val distance2 : t -> t -> Qptypes.Positive_float.t

1
ocaml/Qputils.ml
View File

@ -12,7 +12,6 @@ let rec transpose = function
;;
*)
let input_to_sexp s =
let result =
String.split_lines s

245
ocaml/qp_create_ezfio_from_xyz.ml
View File

@ -11,76 +11,169 @@ let spec =
~doc:"string Name of basis set."
+> flag "c" (optional_with_default 0 int)
~doc:"int Total charge of the molecule. Default is 0."
+> flag "d" (optional_with_default 0. float)
~doc:"float Add dummy atoms. x * (covalent radii of the atoms)"
+> flag "m" (optional_with_default 1 int)
~doc:"int Spin multiplicity (2S+1) of the molecule. Default is 1."
+> flag "p" no_arg
~doc:" Using pseudopotentials"
+> anon ("xyz_file" %: string)
;;
+> anon ("xyz_file" %: file )
let run ?o b c m p xyz_file =
let dummy_centers ~threshold ~molecule ~nuclei =
let d =
Molecule.distance_matrix molecule
in
let n =
Array.length d
in
let nuclei =
Array.of_list nuclei
in
let rec aux accu = function
| (-1,_) -> accu
| (i,-1) -> aux accu (i-1,i-1)
| (i,j) when (i>j) ->
let new_accu =
let x,y =
Element.covalent_radius (nuclei.(i)).Atom.element |> Positive_float.to_float,
Element.covalent_radius (nuclei.(j)).Atom.element |> Positive_float.to_float
in
let r =
( x +. y ) *. threshold
in
if d.(i).(j) < r then
(i,x,j,y,d.(i).(j)) :: accu
else
accu
in aux new_accu (i,j-1)
| (i,j) when (i=j) -> aux accu (i,j-1)
| _ -> assert false
in
aux [] (n-1,n-1)
|> List.map ~f:(fun (i,x,j,y,r) ->
let f =
x /. (x +. y)
in
let u =
Point3d.of_tuple ~units:Units.Bohr
( nuclei.(i).Atom.coord.Point3d.x +.
(nuclei.(j).Atom.coord.Point3d.x -. nuclei.(i).Atom.coord.Point3d.x) *. f,
nuclei.(i).Atom.coord.Point3d.y +.
(nuclei.(j).Atom.coord.Point3d.y -. nuclei.(i).Atom.coord.Point3d.y) *. f,
nuclei.(i).Atom.coord.Point3d.z +.
(nuclei.(j).Atom.coord.Point3d.z -. nuclei.(i).Atom.coord.Point3d.z) *. f)
in
Atom.{ element = Element.X ; charge = Charge.of_int 0 ; coord = u }
)
let list_basis () =
let basis_list =
Qpackage.root ^ "/install/emsl/EMSL_api.py list_basis"
|> Unix.open_process_in
|> In_channel.input_lines
|> List.map ~f:(fun x ->
match String.split x ~on:'\'' with
| [] -> ""
| a :: []
| _ :: a :: _ -> String.strip a
)
in
List.sort basis_list ~cmp:String.ascending
|> String.concat ~sep:"\t"
let run ?o b c d m p xyz_file =
(* Read molecule *)
let molecule =
(Molecule.of_xyz_file xyz_file ~charge:(Charge.of_int c)
~multiplicity:(Multiplicity.of_int m) )
in
let nuclei = molecule.Molecule.nuclei in
let dummy =
dummy_centers ~threshold:d ~molecule ~nuclei:molecule.Molecule.nuclei
in
(*
List.iter dummy ~f:(fun x ->
Printf.printf "%s\n" (Atom.to_string ~units:Units.Angstrom x)
);
*)
let nuclei =
molecule.Molecule.nuclei @ dummy
in
let basis_table =
Hashtbl.Poly.create ()
in
let basis_table = Hashtbl.Poly.create () in
(* Open basis set channels *)
let basis_channel element =
let key = Element.to_string element in
let key =
Element.to_string element
in
match Hashtbl.find basis_table key with
| Some in_channel ->
in_channel
| None ->
begin
Printf.printf "%s is not defined in basis %s.\nEnter alternate basis : %!"
(Element.to_long_string element) b ;
let bas =
match In_channel.input_line stdin with
| Some line -> String.strip line |> String.lowercase
| None -> failwith "Aborted"
in
let new_channel = In_channel.create
(Qpackage.root ^ "/data/basis/" ^ bas)
in
Hashtbl.add_exn basis_table ~key:key ~data:new_channel;
new_channel
end
| None ->
let msg =
Printf.sprintf "%s is not defined in basis %s.%!"
(Element.to_long_string element) b ;
in
failwith msg
in
let temp_filename =
Filename.temp_file "qp_create_" ".basis"
in
let () =
Sys.remove temp_filename
in
let fetch_channel basis =
let command =
if (p) then
Qpackage.root ^ "/scripts/get_basis.sh \"" ^ temp_filename
^ "." ^ basis ^ "\" \"" ^ basis ^"\" pseudo"
else
Qpackage.root ^ "/scripts/get_basis.sh \"" ^ temp_filename
^ "." ^ basis ^ "\" \"" ^ basis ^"\""
in
match Sys.is_file basis with
| `Yes ->
In_channel.create basis
| _ ->
begin
let filename =
Unix.open_process_in command
|> In_channel.input_all
|> String.strip
in
let new_channel =
In_channel.create filename
in
Unix.unlink filename;
new_channel
end
in
let rec build_basis = function
| [] -> ()
| elem_and_basis_name :: rest ->
begin
match (String.lsplit2 ~on:':' elem_and_basis_name) with
| None -> (* Principal basis *)
let basis = elem_and_basis_name in
let command =
if (p) then
Qpackage.root ^ "/scripts/get_basis.sh \"" ^ temp_filename
^ "\" \"" ^ basis ^"\" pseudo"
else
Qpackage.root ^ "/scripts/get_basis.sh \"" ^ temp_filename
^ "\" \"" ^ basis ^"\""
in
begin
let filename =
Unix.open_process_in command
|> In_channel.input_all
|> String.strip
let basis =
elem_and_basis_name
in
let new_channel =
In_channel.create filename
fetch_channel basis
in
Unix.unlink filename;
List.iter nuclei ~f:(fun elem->
let key = Element.to_string elem.Atom.element
let key =
Element.to_string elem.Atom.element
in
match Hashtbl.add basis_table ~key:key ~data:new_channel with
| `Ok -> ()
@ -89,25 +182,18 @@ let run ?o b c m p xyz_file =
end
| Some (key, basis) -> (*Aux basis *)
begin
let elem = Element.of_string key
and basis = String.lowercase basis
let elem =
Element.of_string key
and basis =
String.lowercase basis
in
let key = Element.to_string elem
let key =
Element.to_string elem
in
let command =
Qpackage.root ^ "/scripts/get_basis.sh \"" ^ temp_filename ^
"\" \"" ^ basis ^ "\" " ^ key
let new_channel =
fetch_channel basis
in
begin
let filename =
Unix.open_process_in command
|> In_channel.input_all
|> String.strip
in
let new_channel =
In_channel.create filename
in
Unix.unlink filename;
match Hashtbl.add basis_table ~key:key ~data:new_channel with
| `Ok -> ()
| `Duplicate -> failwith ("Duplicate definition of basis for "^(Element.to_long_string elem))
@ -164,28 +250,30 @@ let run ?o b c m p xyz_file =
(* Write Basis set *)
let basis =
let nmax = Nucl_number.get_max () in
let nmax =
Nucl_number.get_max ()
in
let rec do_work (accu:(Atom.t*Nucl_number.t) list) (n:int) = function
| [] -> accu
| e::tail ->
let new_accu = (e,(Nucl_number.of_int ~max:nmax n))::accu in
let new_accu =
(e,(Nucl_number.of_int ~max:nmax n))::accu
in
do_work new_accu (n+1) tail
in
let result = do_work [] 1 nuclei
|> List.rev
|> List.map ~f:(fun (x,i) ->
try
Basis.read_element (basis_channel x.Atom.element) i x.Atom.element
let e =
match x.Atom.element with
| Element.X -> Element.H
| e -> e
in
Basis.read_element (basis_channel x.Atom.element) i e
with
| End_of_file ->
begin
let alt_channel = basis_channel x.Atom.element in
try
Basis.read_element alt_channel i x.Atom.element
with
End_of_file -> failwith
("Element "^(Element.to_string x.Atom.element)^" not found")
end
| End_of_file -> failwith
("Element "^(Element.to_string x.Atom.element)^" not found in basis set.")
)
|> List.concat
in
@ -264,28 +352,37 @@ let run ?o b c m p xyz_file =
| None -> failwith "Error in basis"
| Some x -> Input.Ao_basis.write x
;;
let command =
Command.basic
~summary: "Quantum Package command"
~readme:(fun () -> "
Creates an EZFIO directory from a standard xyz file.
The basis set is defined as a single string if all the
atoms are taken from the same basis set, otherwise specific
elements can be defined as follows:
=== Available basis sets ===
" ^ (list_basis ()) ^ "
============================
Creates an EZFIO directory from a standard xyz file. The basis set is defined
as a single string if all the atoms are taken from the same basis set,
otherwise specific elements can be defined as follows:
-b \"cc-pcvdz | H:cc-pvdz | C:6-31g\"
")
If a file with the same name as the basis set exists, this file will be read.
Otherwise, the basis set is obtained from the database.
" )
spec
(fun o b c m p xyz_file () ->
run ?o b c m p xyz_file )
;;
(fun o b c d m p xyz_file () ->
run ?o b c d m p xyz_file )
let () =
Command.run command
;;

308
ocaml/qp_edit.ml Normal file
View File

@ -0,0 +1,308 @@
open Qputils;;
open Qptypes;;
open Core.Std;;
(** Interactive editing of the input.
WARNING
This file is autogenerad by
`${QP_ROOT}/script/ezfio_interface/ei_handler.py`
*)
(** Keywords used to define input sections *)
type keyword =
| Ao_basis
| Determinants_by_hand
| Electrons
| Mo_basis
| Nuclei
| Determinants
| Hartree_fock
| Integrals_bielec
| Perturbation
| Properties
| Pseudo
;;
let keyword_to_string = function
| Ao_basis -> "AO basis"
| Determinants_by_hand -> "Determinants_by_hand"
| Electrons -> "Electrons"
| Mo_basis -> "MO basis"
| Nuclei -> "Molecule"
| Determinants -> "Determinants"
| Hartree_fock -> "Hartree_fock"
| Integrals_bielec -> "Integrals_bielec"
| Perturbation -> "Perturbation"
| Properties -> "Properties"
| Pseudo -> "Pseudo"
;;
(** Create the header of the temporary file *)
let file_header filename =
Printf.sprintf "
==================================================================
Quantum Package
==================================================================
Editing file `%s`
" filename
;;
(** Creates the header of a section *)
let make_header kw =
let s = keyword_to_string kw in
let l = String.length s in
"\n\n"^s^"\n"^(String.init l ~f:(fun _ -> '='))^"\n\n"
;;
(** Returns the rst string of section [s] *)
let get s =
let header = (make_header s) in
let f (read,to_rst) =
match read () with
| Some text -> header ^ (Rst_string.to_string (to_rst text))
| None -> ""
in
let rst =
try
begin
let open Input in
match s with
| Mo_basis ->
f Mo_basis.(read, to_rst)
| Electrons ->
f Electrons.(read, to_rst)
| Nuclei ->
f Nuclei.(read, to_rst)
| Ao_basis ->
f Ao_basis.(read, to_rst)
| Determinants_by_hand ->
f Determinants_by_hand.(read, to_rst)
| Determinants ->
f Determinants.(read, to_rst)
| Hartree_fock ->
f Hartree_fock.(read, to_rst)
| Integrals_bielec ->
f Integrals_bielec.(read, to_rst)
| Perturbation ->
f Perturbation.(read, to_rst)
| Properties ->
f Properties.(read, to_rst)
| Pseudo ->
f Pseudo.(read, to_rst)
end
with
| Sys_error msg -> (Printf.eprintf "Info: %s\n%!" msg ; "")
in
rst
;;
(** Applies the changes from the string [str] corresponding to section [s] *)
let set str s =
let header = (make_header s) in
match String.substr_index ~pos:0 ~pattern:header str with
| None -> ()
| Some idx ->
begin
let index_begin = idx + (String.length header) in
let index_end =
match ( String.substr_index ~pos:(index_begin+(String.length header)+1)
~pattern:"==" str) with
| Some i -> i
| None -> String.length str
in
let l = index_end - index_begin in
let str = String.sub ~pos:index_begin ~len:l str
|> Rst_string.of_string
in
let write (of_rst,w) s =
try
match of_rst str with
| Some data -> w data
| None -> ()
with
| _ -> (Printf.eprintf "Info: Read error in %s\n%!"
(keyword_to_string s); ignore (of_rst str) )
in
let open Input in
match s with
| Determinants -> write Determinants.(of_rst, write) s
| Hartree_fock -> write Hartree_fock.(of_rst, write) s
| Integrals_bielec -> write Integrals_bielec.(of_rst, write) s
| Perturbation -> write Perturbation.(of_rst, write) s
| Properties -> write Properties.(of_rst, write) s
| Pseudo -> write Pseudo.(of_rst, write) s
| Electrons -> write Electrons.(of_rst, write) s
| Determinants_by_hand -> write Determinants_by_hand.(of_rst, write) s
| Nuclei -> write Nuclei.(of_rst, write) s
| Ao_basis -> () (* TODO *)
| Mo_basis -> () (* TODO *)
end
;;
(** Creates the temporary file for interactive editing *)
let create_temp_file ezfio_filename fields =
let temp_filename = Filename.temp_file "qp_edit_" ".rst" in
begin
Out_channel.with_file temp_filename ~f:(fun out_channel ->
(file_header ezfio_filename) :: (List.map ~f:get fields)
|> String.concat ~sep:"\n"
|> Out_channel.output_string out_channel
)
end
; temp_filename
;;
let run check_only ezfio_filename =
(* Open EZFIO *)
if (not (Sys.file_exists_exn ezfio_filename)) then
failwith (ezfio_filename^" does not exists");
Ezfio.set_file ezfio_filename;
(*
let output = (file_header ezfio_filename) :: (
List.map ~f:get [
Ao_basis ;
Mo_basis ;
])
in
String.concat output
|> print_string
*)
let tasks = [
Nuclei ;
Ao_basis;
Electrons ;
Determinants ;
Hartree_fock ;
Integrals_bielec ;
Perturbation ;
Properties ;
Pseudo ;
Mo_basis;
Determinants_by_hand ;
]
in
(* Create the temp file *)
let temp_filename = create_temp_file ezfio_filename tasks in
(* Open the temp file with external editor *)
let editor =
match Sys.getenv "EDITOR" with
| Some editor -> editor
| None -> "vi"
in
match check_only with
| true -> ()
| false ->
Printf.sprintf "%s %s ; tput sgr0 2> /dev/null" editor temp_filename
|> Sys.command_exn
;
(* Re-read the temp file *)
let temp_string =
In_channel.with_file temp_filename ~f:(fun in_channel ->
In_channel.input_all in_channel)
in
List.iter ~f:(fun x -> set temp_string x) tasks;
(* Remove temp_file *)
Sys.remove temp_filename;
;;
(** Create a backup file in case of an exception *)
let create_backup ezfio_filename =
Printf.sprintf "
rm -f %s/backup.tgz ;
tar -zcf .backup.tgz %s && mv .backup.tgz %s/backup.tgz
"
ezfio_filename ezfio_filename ezfio_filename
|> Sys.command_exn
;;
(** Restore the backup file when an exception occuprs *)
let restore_backup ezfio_filename =
Printf.sprintf "tar -zxf %s/backup.tgz"
ezfio_filename
|> Sys.command_exn
;;
let spec =
let open Command.Spec in
empty
+> flag "-c" no_arg
~doc:"Checks the input data"
(*
+> flag "o" (optional string)
~doc:"Prints output data"
*)
+> anon ("ezfio_file" %: string)
;;
let command =
Command.basic
~summary: "Quantum Package command"
~readme:(fun () ->
"
Edit input data
")
spec
(* (fun i o ezfio_file () -> *)
(*fun ezfio_file () ->
try
run ezfio_file
with
| _ msg -> print_string ("\n\nError\n\n"^msg^"\n\n")
*)
(fun c ezfio_file () ->
try
run c ezfio_file ;
(* create_backup ezfio_file; *)
with
| Failure exc
| Invalid_argument exc as e ->
begin
Printf.eprintf "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-\n\n";
Printf.eprintf "%s\n\n" exc;
Printf.eprintf "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-\n\n";
(* restore_backup ezfio_file; *)
raise e
end
| Assert_failure (file, line, ch) as e ->
begin
Printf.eprintf "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-\n\n";
Printf.eprintf "Assert error in file $QP_ROOT/ocaml/%s, line %d, character %d\n\n" file line ch;
Printf.eprintf "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-\n\n";
(* restore_backup ezfio_file; *)
raise e
end
)
;;
let () =
Command.run command;
exit 0
;;

2
ocaml/qp_run.ml
View File

@ -42,7 +42,7 @@ let run exe ezfio_file =
let spec =
let open Command.Spec in
empty
+> anon ("exectuable" %: string)
+> anon ("executable" %: string)
+> anon ("ezfio_file" %: string)
;;

6
ocaml/qptypes_generator.ml
View File

@ -83,6 +83,12 @@ let input_data = "
assert (x >= 0.) ;
assert (x <= 1.) ;
* Energy : float
assert (x <=0.) ;
* S2 : float
assert (x >=0.) ;
* PT2_energy : float
assert (x >=0.) ;

9
ocaml/test_molecule.ml
View File

@ -39,6 +39,15 @@ H 0.54386314 0.00000000 -0.92559535
let m = Molecule.of_xyz_file "c2h6.xyz" in
print_string (Molecule.to_string m);
print_string "\nDistance matrix\n";
print_string "---------------\n";
let d =
Molecule.distance_matrix m
in
Array.iter d ~f:(fun x ->
Array.iter x ~f:(fun y -> Printf.printf "%12.8f " y);
print_newline ();
)
;;
test_molecule ();;

10
plugins/CASSCF/EZFIO.cfg Normal file
View File

@ -0,0 +1,10 @@
[energy]
type: double precision
doc: "Calculated CAS-SCF energy"
interface: ezfio
[energy_pt2]
type: double precision
doc: "Calculated selected CAS-SCF energy with PT2 correction"
interface: ezfio

39
plugins/CASSCF/H_apply.irp.f Normal file
View File

@ -0,0 +1,39 @@
use bitmasks
BEGIN_SHELL [ /usr/bin/env python ]
from generate_h_apply import *
s = H_apply("CAS_SD")
print s
s = H_apply("CAS_SD_selected_no_skip")
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
print s
s = H_apply("CAS_SD_selected")
s.set_selection_pt2("epstein_nesbet_2x2")
print s
s = H_apply("CAS_SD_PT2")
s.set_perturbation("epstein_nesbet_2x2")
print s
s = H_apply("CAS_S",do_double_exc=False)
print s
s = H_apply("CAS_S_selected_no_skip",do_double_exc=False)
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
print s
s = H_apply("CAS_S_selected",do_double_exc=False)
s.set_selection_pt2("epstein_nesbet_2x2")
print s
s = H_apply("CAS_S_PT2",do_double_exc=False)
s.set_perturbation("epstein_nesbet_2x2")
print s
END_SHELL

1
plugins/CASSCF/NEEDED_CHILDREN_MODULES Normal file
View File

@ -0,0 +1 @@
Generators_CAS Perturbation Selectors_full

20
plugins/CASSCF/README.rst Normal file
View File

@ -0,0 +1,20 @@
======
CASSCF
======
This module is not a "real" CAS-SCF. It is an orbital optimization step done by :
1) Doing the CAS+SD
2) Taking one-electron density matrix
3) Cancelling all active-active rotations
4) Finding the order which matches with the input MOs
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.

220
plugins/CASSCF/casscf.irp.f Normal file
View File

@ -0,0 +1,220 @@
program casscf
implicit none
BEGIN_DOC
! Optimize MOs and CI coefficients of the CAS
END_DOC
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer(bit_kind), allocatable :: generators_bitmask_save(:,:,:,:)
integer :: degree, N_generators_bitmask_save, N_det_ci
double precision :: E_old, E_CI
double precision :: selection_criterion_save, selection_criterion_min_save
integer :: N_det_old
integer :: i, j, k, l
integer :: i_bit, j_bit, i_int, j_int
integer(bit_kind), allocatable :: bit_tmp(:), cas_bm(:)
character*(64) :: label
allocate( pt2(N_states), norm_pert(N_states),H_pert_diag(N_states) )
allocate( generators_bitmask_save(N_int,2,6,N_generators_bitmask) )
allocate( bit_tmp(N_int), cas_bm(N_int) )
PROVIDE N_det_cas
N_det_old = 0
pt2 = 1.d0
E_CI = 1.d0
E_old = 0.d0
diag_algorithm = "Lapack"
selection_criterion_save = selection_criterion
selection_criterion_min_save = selection_criterion_min
cas_bm = 0_bit_kind
do i=1,N_cas_bitmask
do j=1,N_int
cas_bm(j) = ior(cas_bm(j), cas_bitmask(j,1,i))
cas_bm(j) = ior(cas_bm(j), cas_bitmask(j,2,i))
enddo
enddo
! Save CAS-SD bitmask
generators_bitmask_save = generators_bitmask
N_generators_bitmask_save = N_generators_bitmask
! Set the CAS bitmask
do i=1,6
generators_bitmask(:,:,i,:) = cas_bitmask
enddo
N_generators_bitmask = N_cas_bitmask
SOFT_TOUCH generators_bitmask N_generators_bitmask
! If the number of dets already in the file is larger than the requested
! number of determinants, truncate the wf
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
! Start MCSCF iteration
! CAS-CI
! ------
E_old = E_CI
! Reset the selection criterion
selection_criterion = selection_criterion_save
selection_criterion_min = selection_criterion_min_save
SOFT_TOUCH selection_criterion_min selection_criterion selection_criterion_factor
! Set the CAS bitmask
do i=1,6
generators_bitmask(:,:,i,:) = cas_bitmask
enddo
N_generators_bitmask = N_cas_bitmask
SOFT_TOUCH generators_bitmask N_generators_bitmask
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_states))) > pt2_max)
N_det_old = N_det
call H_apply_CAS_SD_selected_no_skip(pt2, norm_pert, H_pert_diag, N_states)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det > N_det_max) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
endif
call diagonalize_CI
call save_wavefunction
print *, '======'
print *, 'CAS-CI'
print *, '======'
print *, ''
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E(CAS) = ', CI_energy
print *, 'E(CAS)+PT2 = ', CI_energy+pt2
print *, '-----'
print *, ''
E_CI = sum(CI_energy(1:N_states)+pt2(1:N_states))/dble(N_states)
call ezfio_set_casscf_energy(CI_energy(1))
if (abort_all) then
exit
endif
if (N_det == N_det_old) then
exit
endif
enddo
! Super-CI
! --------
selection_criterion_min = 1.d-12
selection_criterion = 1.d-12
! Set the CAS bitmask
generators_bitmask = generators_bitmask_save
N_generators_bitmask = N_generators_bitmask_save
SOFT_TOUCH generators_bitmask N_generators_bitmask selection_criterion selection_criterion_min selection_criterion_factor
N_det_ci = N_det
call H_apply_CAS_SD_selected(pt2, norm_pert, H_pert_diag, N_states)
do i=1,mo_tot_num
i_int = ishft(i-1,-bit_kind_shift)+1
i_bit = j-ishft(i_int-1,bit_kind_shift)-1
bit_tmp(:) = 0_bit_kind
bit_tmp(i_int) = ibset(0_bit_kind,i_bit)
if (iand(bit_tmp(i_int), cas_bm(i_int)) == 0_bit_kind) then
! Not a CAS MO
cycle
endif
do j=1,mo_tot_num
if (j == i) then
cycle
endif
j_int = ishft(j-1,-bit_kind_shift)+1
j_bit = j-ishft(j_int-1,bit_kind_shift)-1
bit_tmp(:) = 0_bit_kind
bit_tmp(j_int) = ibset(0_bit_kind,j_bit)
if (iand(bit_tmp(j_int), cas_bm(j_int)) == 0_bit_kind) then
! Not a CAS MO
cycle
endif
! Now, both i and j are MOs of the CAS. De-couple them in the DM
one_body_dm_mo(i,j) = 0.d0
enddo
enddo
SOFT_TOUCH one_body_dm_mo
double precision :: mx, ov
double precision, allocatable :: mo_coef_old(:,:)
integer, allocatable :: iorder(:)
logical, allocatable :: selected(:)
allocate( mo_coef_old(size(mo_coef,1), size(mo_coef,2)), iorder(mo_tot_num), selected(mo_tot_num) )
mo_coef_old = mo_coef
label = "Canonical"
call mo_as_eigvectors_of_mo_matrix(one_body_dm_mo,size(one_body_dm_mo,1),size(one_body_dm_mo,2),label,-1)
selected = .False.
do j=1,mo_tot_num
mx = -1.d0
iorder(j) = j
do i=1,mo_tot_num
if (selected(i)) then
cycle
endif
ov = 0.d0
do l=1,ao_num
do k=1,ao_num
ov = ov + mo_coef_old(k,j) * ao_overlap(k,l) * mo_coef(l,i)
enddo
enddo
ov= dabs(ov)
if (ov > mx) then
mx = ov
iorder(j) = i
endif
enddo
selected( iorder(j) ) = .True.
enddo
mo_coef_old = mo_coef
do i=1,mo_tot_num
mo_coef(:,i) = mo_coef_old(:,iorder(i))
enddo
call save_mos
call write_double(6,E_CI,"Energy(CAS)")
deallocate( mo_coef_old )
deallocate( pt2, norm_pert,H_pert_diag )
deallocate( generators_bitmask_save )
deallocate( bit_tmp, cas_bm, iorder )
end

18
plugins/CAS_SD/H_apply.irp.f
View File

@ -5,7 +5,6 @@ from generate_h_apply import *
s = H_apply("CAS_SD")
print s
s = H_apply("CAS_SD_selected_no_skip")
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
@ -19,5 +18,22 @@ s = H_apply("CAS_SD_PT2")
s.set_perturbation("epstein_nesbet_2x2")
print s
s = H_apply("CAS_S",do_double_exc=False)
print s
s = H_apply("CAS_S_selected_no_skip",do_double_exc=False)
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
print s
s = H_apply("CAS_S_selected",do_double_exc=False)
s.set_selection_pt2("epstein_nesbet_2x2")
print s
s = H_apply("CAS_S_PT2",do_double_exc=False)
s.set_perturbation("epstein_nesbet_2x2")
print s
END_SHELL

56
plugins/CAS_SD/README.rst
View File

@ -118,69 +118,101 @@ Documentation
Undocumented
`h_apply_cas_sd <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L414>`_
h_apply_cas_sd
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cas_sd_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L1>`_
h_apply_cas_sd_diexc
Undocumented
h_apply_cas_sd_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L269>`_
h_apply_cas_sd_diexcp
Undocumented
h_apply_cas_sd_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_pt2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L2610>`_
h_apply_cas_sd_pt2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cas_sd_pt2_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L2118>`_
h_apply_cas_sd_pt2_diexc
Undocumented
h_apply_cas_sd_pt2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_pt2_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L2427>`_
h_apply_cas_sd_pt2_diexcp
Undocumented
h_apply_cas_sd_pt2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_selected <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L1872>`_
h_apply_cas_sd_selected
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cas_sd_selected_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L1346>`_
h_apply_cas_sd_selected_diexc
Undocumented
h_apply_cas_sd_selected_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_selected_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L1675>`_
h_apply_cas_sd_selected_diexcp
Undocumented
h_apply_cas_sd_selected_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_selected_no_skip <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L1128>`_
h_apply_cas_sd_selected_no_skip
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cas_sd_selected_no_skip_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L602>`_
h_apply_cas_sd_selected_no_skip_diexc
Undocumented
h_apply_cas_sd_selected_no_skip_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_selected_no_skip_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/H_apply.irp.f_shell_22#L931>`_
h_apply_cas_sd_selected_no_skip_diexcp
Undocumented
h_apply_cas_sd_selected_no_skip_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.

95
plugins/CAS_SD/cas_s.irp.f Normal file
View File

@ -0,0 +1,95 @@
program full_ci
implicit none
integer :: i,k
integer :: N_det_old
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer :: N_st, degree
N_st = N_states
allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
character*(64) :: perturbation
PROVIDE N_det_cas
N_det_old = 0
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_st))) > pt2_max)
N_det_old = N_det
call H_apply_CAS_S(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det > N_det_max) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
endif
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
call ezfio_set_cas_sd_energy(CI_energy(1))
if (abort_all) then
exit
endif
if (N_det == N_det_old) then
exit
endif
enddo
call diagonalize_CI
if(do_pt2_end)then
print*,'Last iteration only to compute the PT2'
threshold_selectors = 1.d0
threshold_generators = 0.999d0
call H_apply_CAS_S_PT2(pt2, norm_pert, H_pert_diag, N_st)
print *, 'Final step'
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
call ezfio_set_cas_sd_energy_pt2(CI_energy(1)+pt2(1))
endif
integer :: exc_max, degree_min
exc_max = 0
print *, 'CAS determinants : ', N_det_cas
do i=1,min(N_det_cas,10)
do k=i,N_det_cas
call get_excitation_degree(psi_cas(1,1,k),psi_cas(1,1,i),degree,N_int)
exc_max = max(exc_max,degree)
enddo
call debug_det(psi_cas(1,1,i),N_int)
print *, ''
enddo
print *, 'Max excitation degree in the CAS :', exc_max
end

89
plugins/CAS_SD/cas_s_selected.irp.f Normal file
View File

@ -0,0 +1,89 @@
program full_ci
implicit none
integer :: i,k
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer :: N_st, degree
N_st = N_states
allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
character*(64) :: perturbation
PROVIDE N_det_cas
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_st))) > pt2_max)
call H_apply_CAS_S_selected(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det > N_det_max) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
endif
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
call ezfio_set_cas_sd_energy(CI_energy(1))
if (abort_all) then
exit
endif
enddo
call diagonalize_CI
if(do_pt2_end)then
print*,'Last iteration only to compute the PT2'
threshold_selectors = 1.d0
threshold_generators = 0.999d0
call H_apply_CAS_S_PT2(pt2, norm_pert, H_pert_diag, N_st)
print *, 'Final step'
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
call ezfio_set_cas_sd_energy_pt2(CI_energy(1)+pt2(1))
endif
integer :: exc_max, degree_min
exc_max = 0
print *, 'CAS determinants : ', N_det_cas
do i=1,min(N_det_cas,10)
do k=i,N_det_cas
call get_excitation_degree(psi_cas(1,1,k),psi_cas(1,1,i),degree,N_int)
exc_max = max(exc_max,degree)
enddo
call debug_det(psi_cas(1,1,i),N_int)
print *, ''
enddo
print *, 'Max excitation degree in the CAS :', exc_max
end

45
plugins/CID/README.rst
View File

@ -203,3 +203,48 @@ Documentation
particles.
Assume N_int is already provided.
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/plugins/Selectors_full>`_
* `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/plugins/SingleRefMethod>`_
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
`cid <http://github.com/LCPQ/quantum_package/tree/master/plugins/CID/cid_lapack.irp.f#L1>`_
Undocumented
h_apply_cid
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_cid_diexc
Undocumented
h_apply_cid_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cid_diexcp
Undocumented
h_apply_cid_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.

1
plugins/CISD/.gitignore vendored
View File

@ -20,6 +20,7 @@ Pseudo
Selectors_full
SingleRefMethod
Utils
cisd
cisd_lapack
ezfio_interface.irp.f
irpf90.make

10
plugins/CISD/EZFIO.cfg Normal file
View File

@ -0,0 +1,10 @@
[energy]
type: double precision
doc: Variational CISD energy
interface: ezfio
[energy_pt2]
type: double precision
doc: Estimated CISD energy (including PT2)
interface: ezfio

18
plugins/CISD/README.rst
View File

@ -59,22 +59,26 @@ Documentation
.. by the `update_README.py` script.
`cisd <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD/cisd_lapack.irp.f#L1>`_
Undocumented
`h_apply_cisd <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD/H_apply.irp.f_shell_8#L414>`_
h_apply_cisd
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD/H_apply.irp.f_shell_8#L1>`_
h_apply_cisd_diexc
Undocumented
h_apply_cisd_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD/H_apply.irp.f_shell_8#L269>`_
h_apply_cisd_diexcp
Undocumented
h_apply_cisd_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.

BIN
plugins/CISD/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 87 KiB

After

Width:  |  Height:  |  Size: 84 KiB

119
plugins/CISD_SC2_selected/README.rst
View File

@ -72,3 +72,122 @@ Needed Modules
* `CISD_selected <http://github.com/LCPQ/quantum_package/tree/master/src/CISD_selected>`_
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `CISD_selected <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected>`_
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
`cisd_sc2_selected <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_SC2_selected/cisd_sc2_selection.irp.f#L1>`_
Undocumented
h_apply_cisd
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_cisd_diexc
Undocumented
h_apply_cisd_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cisd_diexcp
Undocumented
h_apply_cisd_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_pt2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_pt2_diexc
Undocumented
h_apply_pt2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_pt2_diexcp
Undocumented
h_apply_pt2_en_sc2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_pt2_en_sc2_diexc
Undocumented
h_apply_pt2_en_sc2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_pt2_en_sc2_diexcp
Undocumented
h_apply_pt2_en_sc2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_pt2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_sc2_selected
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_sc2_selected_diexc
Undocumented
h_apply_sc2_selected_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_sc2_selected_diexcp
Undocumented
h_apply_sc2_selected_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.

144
plugins/CISD_selected/README.rst
View File

@ -196,22 +196,26 @@ Documentation
.. by the `update_README.py` script.
`cisd <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/cisd_selection.irp.f#L1>`_
Undocumented
`h_apply_cisd <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_8#L414>`_
h_apply_cisd
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_8#L1>`_
h_apply_cisd_diexc
Undocumented
h_apply_cisd_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_8#L269>`_
h_apply_cisd_diexcp
Undocumented
h_apply_cisd_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
@ -221,154 +225,226 @@ Documentation
Undocumented
`h_apply_cisd_selection_delta_rho_one_point <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L5931>`_
h_apply_cisd_selection_delta_rho_one_point
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_selection_delta_rho_one_point_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L5405>`_
h_apply_cisd_selection_delta_rho_one_point_diexc
Undocumented
h_apply_cisd_selection_delta_rho_one_point_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_delta_rho_one_point_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L5734>`_
h_apply_cisd_selection_delta_rho_one_point_diexcp
Undocumented
h_apply_cisd_selection_delta_rho_one_point_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_dipole_moment_z <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L5159>`_
h_apply_cisd_selection_dipole_moment_z
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_selection_dipole_moment_z_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L4633>`_
h_apply_cisd_selection_dipole_moment_z_diexc
Undocumented
h_apply_cisd_selection_dipole_moment_z_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_dipole_moment_z_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L4962>`_
h_apply_cisd_selection_dipole_moment_z_diexcp
Undocumented
h_apply_cisd_selection_dipole_moment_z_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L3615>`_
h_apply_cisd_selection_epstein_nesbet
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_selection_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L4387>`_
h_apply_cisd_selection_epstein_nesbet_2x2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_selection_epstein_nesbet_2x2_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L3861>`_
h_apply_cisd_selection_epstein_nesbet_2x2_diexc
Undocumented
h_apply_cisd_selection_epstein_nesbet_2x2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet_2x2_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L4190>`_
h_apply_cisd_selection_epstein_nesbet_2x2_diexcp
Undocumented
h_apply_cisd_selection_epstein_nesbet_2x2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L3089>`_
h_apply_cisd_selection_epstein_nesbet_diexc
Undocumented
h_apply_cisd_selection_epstein_nesbet_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L3418>`_
h_apply_cisd_selection_epstein_nesbet_diexcp
Undocumented
h_apply_cisd_selection_epstein_nesbet_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L2843>`_
h_apply_cisd_selection_epstein_nesbet_sc2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_selection_epstein_nesbet_sc2_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L2317>`_
h_apply_cisd_selection_epstein_nesbet_sc2_diexc
Undocumented
h_apply_cisd_selection_epstein_nesbet_sc2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet_sc2_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L2646>`_
h_apply_cisd_selection_epstein_nesbet_sc2_diexcp
Undocumented
h_apply_cisd_selection_epstein_nesbet_sc2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L2071>`_
h_apply_cisd_selection_epstein_nesbet_sc2_no_projected
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_selection_epstein_nesbet_sc2_no_projected_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L1545>`_
h_apply_cisd_selection_epstein_nesbet_sc2_no_projected_diexc
Undocumented
h_apply_cisd_selection_epstein_nesbet_sc2_no_projected_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet_sc2_no_projected_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L1874>`_
h_apply_cisd_selection_epstein_nesbet_sc2_no_projected_diexcp
Undocumented
h_apply_cisd_selection_epstein_nesbet_sc2_no_projected_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L1299>`_
h_apply_cisd_selection_epstein_nesbet_sc2_projected
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_selection_epstein_nesbet_sc2_projected_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L773>`_
h_apply_cisd_selection_epstein_nesbet_sc2_projected_diexc
Undocumented
h_apply_cisd_selection_epstein_nesbet_sc2_projected_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_epstein_nesbet_sc2_projected_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L1102>`_
h_apply_cisd_selection_epstein_nesbet_sc2_projected_diexcp
Undocumented
h_apply_cisd_selection_epstein_nesbet_sc2_projected_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_h_core <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L527>`_
h_apply_cisd_selection_h_core
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_selection_h_core_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L1>`_
h_apply_cisd_selection_h_core_diexc
Undocumented
h_apply_cisd_selection_h_core_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_h_core_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L330>`_
h_apply_cisd_selection_h_core_diexcp
Undocumented
h_apply_cisd_selection_h_core_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L6703>`_
h_apply_cisd_selection_moller_plesset
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cisd_selection_moller_plesset_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L6177>`_
h_apply_cisd_selection_moller_plesset_diexc
Undocumented
h_apply_cisd_selection_moller_plesset_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cisd_selection_moller_plesset_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/CISD_selected/H_apply.irp.f_shell_10#L6506>`_
h_apply_cisd_selection_moller_plesset_diexcp
Undocumented
h_apply_cisd_selection_moller_plesset_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.

4
plugins/CISD_selected/cisd_selection.irp.f
View File

@ -41,8 +41,8 @@ program cisd
N_det = min(N_det,N_det_max)
touch N_det psi_det psi_coef
call diagonalize_CI
deallocate(pt2,norm_pert,H_pert_diag)
call save_wavefunction
call save_wavefunction
call ezfio_set_cisd_selected_energy(CI_energy)
call ezfio_set_cisd_selected_energy_pt2(CI_energy+pt2)
deallocate(pt2,norm_pert,H_pert_diag)
end

71
plugins/DDCI_selected/README.rst
View File

@ -57,3 +57,74 @@ Needed Modules
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/src/Generators_CAS>`_
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation>`_
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/plugins/Selectors_full>`_
* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/plugins/Generators_CAS>`_
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
`ddci <http://github.com/LCPQ/quantum_package/tree/master/plugins/DDCI_selected/ddci.irp.f#L1>`_
Undocumented
h_apply_ddci_pt2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_ddci_pt2_diexc
Undocumented
h_apply_ddci_pt2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_ddci_pt2_diexcp
Undocumented
h_apply_ddci_pt2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_ddci_selection
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_ddci_selection_diexc
Undocumented
h_apply_ddci_selection_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_ddci_selection_diexcp
Undocumented
h_apply_ddci_selection_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.

BIN
plugins/FCIdump/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 64 KiB

After

Width:  |  Height:  |  Size: 59 KiB

32
plugins/Full_CI/.gitignore vendored Normal file
View File

@ -0,0 +1,32 @@
# Automatically created by $QP_ROOT/scripts/module/module_handler.py
.ninja_deps
.ninja_log
AO_Basis
Bitmask
Determinants
Electrons
Ezfio_files
Generators_full
Hartree_Fock
IRPF90_man
IRPF90_temp
Integrals_Bielec
Integrals_Monoelec
MOGuess
MO_Basis
Makefile
Makefile.depend
Nuclei
Perturbation
Properties
Pseudo
Selectors_full
Utils
ezfio_interface.irp.f
full_ci
full_ci_no_skip
irpf90.make
irpf90_entities
tags
target_pt2
var_pt2_ratio

112
plugins/Full_CI/README.rst
View File

@ -27,137 +27,201 @@ Documentation
Undocumented
`h_apply_fci <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L527>`_
h_apply_fci
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_fci_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L1>`_
h_apply_fci_diexc
Undocumented
h_apply_fci_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_fci_mono <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L2744>`_
h_apply_fci_diexcp
Undocumented
h_apply_fci_mono
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_fci_mono_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L2216>`_
h_apply_fci_mono_diexc
Undocumented
h_apply_fci_mono_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_fci_mono_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L2545>`_
h_apply_fci_mono_diexcp
Undocumented
h_apply_fci_mono_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_fci_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L330>`_
h_apply_fci_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_fci_no_skip <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L1998>`_
h_apply_fci_no_skip
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_fci_no_skip_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L1472>`_
h_apply_fci_no_skip_diexc
Undocumented
h_apply_fci_no_skip_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_fci_no_skip_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L1801>`_
h_apply_fci_no_skip_diexcp
Undocumented
h_apply_fci_no_skip_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_fci_pt2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L1265>`_
h_apply_fci_pt2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_fci_pt2_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L773>`_
h_apply_fci_pt2_diexc
Undocumented
h_apply_fci_pt2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_fci_pt2_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L1082>`_
h_apply_fci_pt2_diexcp
Undocumented
h_apply_fci_pt2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_pt2_mono_delta_rho <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L4258>`_
h_apply_pt2_mono_delta_rho
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_pt2_mono_delta_rho_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L3764>`_
h_apply_pt2_mono_delta_rho_diexc
Undocumented
h_apply_pt2_mono_delta_rho_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_pt2_mono_delta_rho_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L4073>`_
h_apply_pt2_mono_delta_rho_diexcp
Undocumented
h_apply_pt2_mono_delta_rho_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_pt2_mono_di_delta_rho <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L5729>`_
h_apply_pt2_mono_di_delta_rho
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_pt2_mono_di_delta_rho_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L5237>`_
h_apply_pt2_mono_di_delta_rho_diexc
Undocumented
h_apply_pt2_mono_di_delta_rho_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_pt2_mono_di_delta_rho_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L5546>`_
h_apply_pt2_mono_di_delta_rho_diexcp
Undocumented
h_apply_pt2_mono_di_delta_rho_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_select_mono_delta_rho <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L3518>`_
h_apply_select_mono_delta_rho
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_select_mono_delta_rho_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L2990>`_
h_apply_select_mono_delta_rho_diexc
Undocumented
h_apply_select_mono_delta_rho_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_select_mono_delta_rho_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L3319>`_
h_apply_select_mono_delta_rho_diexcp
Undocumented
h_apply_select_mono_delta_rho_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_select_mono_di_delta_rho <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L4991>`_
h_apply_select_mono_di_delta_rho
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_select_mono_di_delta_rho_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L4465>`_
h_apply_select_mono_di_delta_rho_diexc
Undocumented
h_apply_select_mono_di_delta_rho_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_select_mono_di_delta_rho_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/Full_CI/H_apply.irp.f_shell_43#L4794>`_
h_apply_select_mono_di_delta_rho_diexcp
Undocumented
h_apply_select_mono_di_delta_rho_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.

2
plugins/Full_CI/target_pt2.irp.f
View File

@ -27,6 +27,8 @@ program var_pt2_ratio_run
call diagonalize_CI
ratio = (CI_energy(1) - HF_energy) / (CI_energy(1)+pt2(1) - HF_energy)
if (N_det > 20000) then
N_det = 20000
TOUCH N_det
exit
endif
enddo

BIN
plugins/Full_CI/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 100 KiB

After

Width:  |  Height:  |  Size: 98 KiB

BIN
plugins/Generators_full/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 74 KiB

After

Width:  |  Height:  |  Size: 73 KiB

6
plugins/Hartree_Fock/EZFIO.cfg
View File

@ -10,6 +10,12 @@ doc: Maximum number of SCF iterations
interface: ezfio,provider,ocaml
default: 200
[level_shift]
type: Positive_float
doc: Energy shift on the virtual MOs to improve SCF convergence
interface: ezfio,provider,ocaml
default: 0.5
[mo_guess_type]
type: MO_guess
doc: Initial MO guess. Can be [ Huckel | HCore ]

99
plugins/Hartree_Fock/Fock_matrix.irp.f
View File

@ -73,8 +73,13 @@
enddo
endif
! Introduce level shift here
do i = elec_alpha_num+1, mo_tot_num
Fock_matrix_mo(i,i) += level_shift
enddo
do i = 1, mo_tot_num
Fock_matrix_diag_mo(i) = Fock_matrix_mo(i,i)
Fock_matrix_diag_mo(i) = Fock_matrix_mo(i,i)
enddo
END_PROVIDER
@ -108,9 +113,10 @@ END_PROVIDER
END_DOC
integer :: i,j,k,l,k1,r,s
integer :: i0,j0,k0,l0
integer*8 :: p,q
double precision :: integral
double precision :: ao_bielec_integral
double precision :: integral, c0, c1, c2
double precision :: ao_bielec_integral, local_threshold
double precision, allocatable :: ao_bi_elec_integral_alpha_tmp(:,:)
double precision, allocatable :: ao_bi_elec_integral_beta_tmp(:,:)
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: ao_bi_elec_integral_beta_tmp
@ -121,11 +127,12 @@ END_PROVIDER
if (do_direct_integrals) then
!$OMP PARALLEL DEFAULT(NONE) &
!$OMP PRIVATE(i,j,l,k1,k,integral,ii,jj,kk,ll,i8,keys,values,p,q,r,s, &
!$OMP ao_bi_elec_integral_alpha_tmp,ao_bi_elec_integral_beta_tmp)&
!$OMP PRIVATE(i,j,l,k1,k,integral,ii,jj,kk,ll,i8,keys,values,p,q,r,s,i0,j0,k0,l0, &
!$OMP ao_bi_elec_integral_alpha_tmp,ao_bi_elec_integral_beta_tmp, c0, c1, c2, &
!$OMP local_threshold)&
!$OMP SHARED(ao_num,ao_num_align,HF_density_matrix_ao_alpha,HF_density_matrix_ao_beta,&
!$OMP ao_integrals_map,ao_integrals_threshold, ao_bielec_integral_schwartz, &
!$OMP ao_overlap_abs, ao_bi_elec_integral_alpha, ao_bi_elec_integral_beta)
!$OMP ao_overlap_abs, ao_bi_elec_integral_alpha, ao_bi_elec_integral_beta)
allocate(keys(1), values(1))
allocate(ao_bi_elec_integral_alpha_tmp(ao_num_align,ao_num), &
@ -152,14 +159,16 @@ END_PROVIDER
< ao_integrals_threshold) then
cycle
endif
if (ao_bielec_integral_schwartz(k,l)*ao_bielec_integral_schwartz(i,j) &
< ao_integrals_threshold) then
cycle
endif
values(1) = ao_bielec_integral(k,l,i,j)
if (abs(values(1)) < ao_integrals_threshold) then
local_threshold = ao_bielec_integral_schwartz(k,l)*ao_bielec_integral_schwartz(i,j)
if (local_threshold < ao_integrals_threshold) then
cycle
endif
i0 = i
j0 = j
k0 = k
l0 = l
values(1) = 0.d0
local_threshold = ao_integrals_threshold/local_threshold
do k2=1,8
if (kk(k2)==0) then
cycle
@ -168,12 +177,21 @@ END_PROVIDER
j = jj(k2)
k = kk(k2)
l = ll(k2)
integral = (HF_density_matrix_ao_alpha(k,l)+HF_density_matrix_ao_beta(k,l)) * values(1)
c0 = HF_density_matrix_ao_alpha(k,l)+HF_density_matrix_ao_beta(k,l)
c1 = HF_density_matrix_ao_alpha(k,i)
c2 = HF_density_matrix_ao_beta(k,i)
if ( dabs(c0)+dabs(c1)+dabs(c2) < local_threshold) then
cycle
endif
if (values(1) == 0.d0) then
values(1) = ao_bielec_integral(k0,l0,i0,j0)
endif
integral = c0 * values(1)
ao_bi_elec_integral_alpha_tmp(i,j) += integral
ao_bi_elec_integral_beta_tmp (i,j) += integral
integral = values(1)
ao_bi_elec_integral_alpha_tmp(l,j) -= HF_density_matrix_ao_alpha(k,i) * integral
ao_bi_elec_integral_beta_tmp (l,j) -= HF_density_matrix_ao_beta (k,i) * integral
ao_bi_elec_integral_alpha_tmp(l,j) -= c1 * integral
ao_bi_elec_integral_beta_tmp (l,j) -= c2 * integral
enddo
enddo
!$OMP END DO NOWAIT
@ -315,7 +333,9 @@ BEGIN_PROVIDER [ double precision, Fock_matrix_ao, (ao_num_align, ao_num) ]
! Fock matrix in AO basis set
END_DOC
if (elec_alpha_num == elec_beta_num) then
if ( (elec_alpha_num == elec_beta_num).and. &
(level_shift == 0.) ) &
then
integer :: i,j
do j=1,ao_num
!DIR$ VECTOR ALIGNED
@ -324,30 +344,47 @@ BEGIN_PROVIDER [ double precision, Fock_matrix_ao, (ao_num_align, ao_num) ]
enddo
enddo
else
double precision, allocatable :: T(:,:), M(:,:)
double precision, allocatable :: T(:,:), M(:,:)
integer :: ierr
! F_ao = S C F_mo C^t S
allocate (T(ao_num_align,ao_num),M(ao_num_align,ao_num))
call dgemm('N','N', ao_num,ao_num,ao_num, 1.d0, &
allocate (T(ao_num_align,ao_num),M(ao_num_align,ao_num),stat=ierr)
if (ierr /=0 ) then
print *, irp_here, ' : allocation failed'
endif
! ao_overlap (ao_num,ao_num) . mo_coef (ao_num,mo_tot_num)
! -> M(ao_num,mo_tot_num)
call dgemm('N','N', ao_num,mo_tot_num,ao_num, 1.d0, &
ao_overlap, size(ao_overlap,1), &
mo_coef, size(mo_coef,1), &
0.d0, &
M, size(M,1))
! M(ao_num,mo_tot_num) . Fock_matrix_mo (mo_tot_num,mo_tot_num)
! -> T(ao_num,mo_tot_num)
call dgemm('N','N', ao_num,mo_tot_num,mo_tot_num, 1.d0, &
M, size(M,1), &
Fock_matrix_mo, size(Fock_matrix_mo,1), &
0.d0, &
T, size(T,1))
call dgemm('N','T', mo_tot_num,ao_num,mo_tot_num, 1.d0, &
! T(ao_num,mo_tot_num) . mo_coef^T (mo_tot_num,ao_num)
! -> M(ao_num,ao_num)
call dgemm('N','T', ao_num,ao_num,mo_tot_num, 1.d0, &
T, size(T,1), &
mo_coef, size(mo_coef,1), &
0.d0, &
M, size(M,1))
! M(ao_num,ao_num) . ao_overlap (ao_num,ao_num)
! -> Fock_matrix_ao(ao_num,ao_num)
call dgemm('N','N', ao_num,ao_num,ao_num, 1.d0, &
M, size(M,1), &
ao_overlap, size(ao_overlap,1), &
0.d0, &
Fock_matrix_ao, size(Fock_matrix_ao,1))
deallocate(T)
endif
END_PROVIDER
@ -360,23 +397,39 @@ subroutine Fock_mo_to_ao(FMO,LDFMO,FAO,LDFAO)
double precision, intent(out) :: FAO(LDFAO,*)
double precision, allocatable :: T(:,:), M(:,:)
integer :: ierr
! F_ao = S C F_mo C^t S
allocate (T(ao_num_align,ao_num),M(ao_num_align,ao_num))
call dgemm('N','N', ao_num,ao_num,ao_num, 1.d0, &
allocate (T(ao_num_align,ao_num),M(ao_num_align,ao_num),stat=ierr)
if (ierr /=0 ) then
print *, irp_here, ' : allocation failed'
endif
! ao_overlap (ao_num,ao_num) . mo_coef (ao_num,mo_tot_num)
! -> M(ao_num,mo_tot_num)
call dgemm('N','N', ao_num,mo_tot_num,ao_num, 1.d0, &
ao_overlap, size(ao_overlap,1), &
mo_coef, size(mo_coef,1), &
0.d0, &
M, size(M,1))
! M(ao_num,mo_tot_num) . FMO (mo_tot_num,mo_tot_num)
! -> T(ao_num,mo_tot_num)
call dgemm('N','N', ao_num,mo_tot_num,mo_tot_num, 1.d0, &
M, size(M,1), &
FMO, size(FMO,1), &
0.d0, &
T, size(T,1))
call dgemm('N','T', mo_tot_num,ao_num,mo_tot_num, 1.d0, &
! T(ao_num,mo_tot_num) . mo_coef^T (mo_tot_num,ao_num)
! -> M(ao_num,ao_num)
call dgemm('N','T', ao_num,ao_num,mo_tot_num, 1.d0, &
T, size(T,1), &
mo_coef, size(mo_coef,1), &
0.d0, &
M, size(M,1))
! M(ao_num,ao_num) . ao_overlap (ao_num,ao_num)
! -> Fock_matrix_ao(ao_num,ao_num)
call dgemm('N','N', ao_num,ao_num,ao_num, 1.d0, &
M, size(M,1), &
ao_overlap, size(ao_overlap,1), &

2
plugins/Hartree_Fock/NEEDED_CHILDREN_MODULES
View File

@ -1 +1 @@
Integrals_Bielec MOGuess
Integrals_Bielec MOGuess

26
plugins/Hartree_Fock/README.rst
View File

@ -32,11 +32,11 @@ Documentation
.. by the `update_README.py` script.
`ao_bi_elec_integral_alpha <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L102>`_
`ao_bi_elec_integral_alpha <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L107>`_
Alpha Fock matrix in AO basis set
`ao_bi_elec_integral_beta <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L103>`_
`ao_bi_elec_integral_beta <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L108>`_
Alpha Fock matrix in AO basis set
@ -62,23 +62,23 @@ Documentation
Diagonal Fock matrix in the MO basis
`fock_matrix_alpha_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L83>`_
`fock_matrix_alpha_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L88>`_
Alpha Fock matrix in AO basis set
`fock_matrix_alpha_mo <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L254>`_
`fock_matrix_alpha_mo <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L272>`_
Fock matrix on the MO basis
`fock_matrix_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L312>`_
`fock_matrix_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L330>`_
Fock matrix in AO basis set
`fock_matrix_beta_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L84>`_
`fock_matrix_beta_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L89>`_
Alpha Fock matrix in AO basis set
`fock_matrix_beta_mo <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L274>`_
`fock_matrix_beta_mo <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L292>`_
Fock matrix on the MO basis
@ -114,7 +114,7 @@ Documentation
.br
`fock_mo_to_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L355>`_
`fock_mo_to_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L375>`_
Undocumented
@ -134,7 +134,7 @@ Documentation
S^-1 Beta density matrix in the AO basis x S^-1
`hf_energy <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L293>`_
`hf_energy <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L311>`_
Hartree-Fock energy
@ -142,7 +142,11 @@ Documentation
Build the MOs using the extended Huckel model
`mo_guess_type <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/ezfio_interface.irp.f#L28>`_
`level_shift <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/ezfio_interface.irp.f#L28>`_
Energy shift on the virtual MOs to improve SCF convergence
`mo_guess_type <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/ezfio_interface.irp.f#L50>`_
Initial MO guess. Can be [ Huckel | HCore ]
@ -161,6 +165,6 @@ Documentation
optional: mo_basis.mo_coef
`thresh_scf <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/ezfio_interface.irp.f#L46>`_
`thresh_scf <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/ezfio_interface.irp.f#L68>`_
Threshold on the convergence of the Hartree Fock energy

7
plugins/Hartree_Fock/SCF.irp.f
View File

@ -42,16 +42,13 @@ subroutine run
BEGIN_DOC
! Run SCF calculation
END_DOC
double precision :: SCF_energy_before,SCF_energy_after,diag_H_mat_elem,get_mo_bielec_integral
double precision :: SCF_energy_before,SCF_energy_after,diag_H_mat_elem
double precision :: E0
integer :: i_it, i, j, k
E0 = HF_energy
thresh_SCF = 1.d-10
call damping_SCF
mo_label = "Canonical"
TOUCH mo_label mo_coef
call save_mos
call damping_SCF
end

4
plugins/Hartree_Fock/damping_SCF.irp.f
View File

@ -86,7 +86,7 @@ subroutine damping_SCF
if ((E_half > E).and.(E_new < E)) then
lambda = 1.d0
exit
else if ((E_half > E).and.(lambda > 5.d-2)) then
else if ((E_half > E).and.(lambda > 5.d-4)) then
lambda = 0.5d0 * lambda
E_new = E_half
else
@ -119,7 +119,7 @@ subroutine damping_SCF
write(output_hartree_fock,'(A4,X,A16, X, A16, X, A16, X, A4 )'), '====','================','================','================', '===='
write(output_hartree_fock,*)
call mo_as_eigvectors_of_mo_matrix(Fock_matrix_mo,size(Fock_matrix_mo,1),size(Fock_matrix_mo,2),mo_label)
call mo_as_eigvectors_of_mo_matrix(Fock_matrix_mo,size(Fock_matrix_mo,1),size(Fock_matrix_mo,2),mo_label,1)
call write_double(output_hartree_fock, E_min, 'Hartree-Fock energy')
call ezfio_set_hartree_fock_energy(E_min)

141
plugins/Hartree_Fock/diagonalize_fock.irp.f
View File

@ -10,58 +10,103 @@
integer, allocatable :: iwork(:)
double precision, allocatable :: work(:), F(:,:), S(:,:)
allocate(F(ao_num_align,ao_num), S(ao_num_align,ao_num) )
do j=1,ao_num
do i=1,ao_num
S(i,j) = ao_overlap(i,j)
F(i,j) = Fock_matrix_ao(i,j)
enddo
enddo
n = ao_num
lwork = 1+6*n + 2*n*n
liwork = 3 + 5*n
allocate(work(lwork), iwork(liwork) )
if (mo_tot_num == ao_num) then
! Solve H.C = E.S.C in AO basis set
lwork = -1
liwork = -1
allocate(F(ao_num_align,ao_num), S(ao_num_align,ao_num) )
do j=1,ao_num
do i=1,ao_num
S(i,j) = ao_overlap(i,j)
F(i,j) = Fock_matrix_ao(i,j)
enddo
enddo
call dsygvd(1,'v','u',ao_num,F,size(F,1),S,size(S,1),&
diagonal_Fock_matrix_mo, work, lwork, iwork, liwork, info)
! call dsygv(1, 'v', 'u',ao_num,F,size(F,1),S,size(S,1),&
! diagonal_Fock_matrix_mo, work, lwork, info)
n = ao_num
lwork = 1+6*n + 2*n*n
liwork = 3 + 5*n
allocate(work(lwork), iwork(liwork) )
lwork = -1
liwork = -1
call dsygvd(1,'v','u',ao_num,F,size(F,1),S,size(S,1),&
diagonal_Fock_matrix_mo, work, lwork, iwork, liwork, info)
if (info /= 0) then
print *, irp_here//' failed : ', info
stop 1
endif
lwork = int(work(1))
liwork = iwork(1)
deallocate(work,iwork)
allocate(work(lwork), iwork(liwork) )
! deallocate(work)
! allocate(work(lwork))
if (info /= 0) then
print *, irp_here//' failed : ', info
stop 1
endif
lwork = int(work(1))
liwork = iwork(1)
deallocate(work,iwork)
allocate(work(lwork), iwork(liwork) )
call dsygvd(1,'v','u',ao_num,F,size(F,1),S,size(S,1),&
diagonal_Fock_matrix_mo, work, lwork, iwork, liwork, info)
call dsygvd(1,'v','u',ao_num,F,size(F,1),S,size(S,1),&
diagonal_Fock_matrix_mo, work, lwork, iwork, liwork, info)
! call dsygv(1, 'v', 'u',ao_num,F,size(F,1),S,size(S,1),&
! diagonal_Fock_matrix_mo, work, lwork, info)
if (info /= 0) then
print *, irp_here//' failed : ', info
stop 1
endif
do j=1,mo_tot_num
do i=1,ao_num
eigenvectors_Fock_matrix_mo(i,j) = F(i,j)
enddo
enddo
if (info /= 0) then
print *, irp_here//' failed : ', info
stop 1
endif
do j=1,mo_tot_num
do i=1,ao_num
eigenvectors_Fock_matrix_mo(i,j) = F(i,j)
enddo
enddo
deallocate(work, iwork, F, S)
else
! Solve H.C = E.C in MO basis set
allocate( F(mo_tot_num_align,mo_tot_num) )
do j=1,mo_tot_num
do i=1,mo_tot_num
F(i,j) = Fock_matrix_mo(i,j)
enddo
enddo
n = mo_tot_num
lwork = 1+6*n + 2*n*n
liwork = 3 + 5*n
allocate(work(lwork), iwork(liwork) )
lwork = -1
liwork = -1
call dsyevd( 'V', 'U', mo_tot_num, F, &
size(F,1), diagonal_Fock_matrix_mo, &
work, lwork, iwork, liwork, info)
if (info /= 0) then
print *, irp_here//' failed : ', info
stop 1
endif
lwork = int(work(1))
liwork = iwork(1)
deallocate(work,iwork)
allocate(work(lwork), iwork(liwork) )
call dsyevd( 'V', 'U', mo_tot_num, F, &
size(F,1), diagonal_Fock_matrix_mo, &
work, lwork, iwork, liwork, info)
if (info /= 0) then
print *, irp_here//' failed : ', info
stop 1
endif
call dgemm('N','N',ao_num,mo_tot_num,mo_tot_num, 1.d0, &
mo_coef, size(mo_coef,1), F, size(F,1), &
0.d0, eigenvectors_Fock_matrix_mo, size(eigenvectors_Fock_matrix_mo,1))
deallocate(work, iwork, F)
endif
deallocate(work, iwork, F, S)
END_PROVIDER
BEGIN_PROVIDER [double precision, diagonal_Fock_matrix_mo_sum, (mo_tot_num)]
@ -73,19 +118,19 @@ BEGIN_PROVIDER [double precision, diagonal_Fock_matrix_mo_sum, (mo_tot_num)]
END_DOC
integer :: i,j
double precision :: accu
do i = 1,elec_alpha_num
do j = 1,elec_alpha_num
accu = 0.d0
do j = 1, elec_alpha_num
do i = 1, elec_alpha_num
accu += 2.d0 * mo_bielec_integral_jj_from_ao(i,j) - mo_bielec_integral_jj_exchange_from_ao(i,j)
enddo
diagonal_Fock_matrix_mo_sum(i) = accu + mo_mono_elec_integral(i,i)
diagonal_Fock_matrix_mo_sum(j) = accu + mo_mono_elec_integral(j,j)
enddo
do i = elec_alpha_num+1,mo_tot_num
do j = elec_alpha_num+1,mo_tot_num
accu = 0.d0
do j = 1, elec_alpha_num
do i = 1, elec_alpha_num
accu += 2.d0 * mo_bielec_integral_jj_from_ao(i,j) - mo_bielec_integral_jj_exchange_from_ao(i,j)
enddo
diagonal_Fock_matrix_mo_sum(i) = accu + mo_mono_elec_integral(i,i)
diagonal_Fock_matrix_mo_sum(j) = accu + mo_mono_elec_integral(j,j)
enddo
END_PROVIDER

14
plugins/Hartree_Fock/huckel.irp.f
View File

@ -4,7 +4,7 @@ subroutine huckel_guess
! Build the MOs using the extended Huckel model
END_DOC
integer :: i,j
double precision :: tmp_matrix(ao_num_align,ao_num),accu
double precision :: accu
double precision :: c
character*(64) :: label
@ -13,20 +13,18 @@ subroutine huckel_guess
label = "Guess"
call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral, &
size(mo_mono_elec_integral,1), &
size(mo_mono_elec_integral,2),label)
size(mo_mono_elec_integral,2),label,1)
TOUCH mo_coef
c = 0.5d0 * 1.75d0
do j=1,ao_num
!DIR$ VECTOR ALIGNED
do i=1,ao_num
if (i.ne.j) then
Fock_matrix_ao(i,j) = c*ao_overlap(i,j)*(ao_mono_elec_integral(i,i) + &
ao_mono_elec_integral(j,j))
else
Fock_matrix_ao(i,j) = Fock_matrix_alpha_ao(i,j)
endif
Fock_matrix_ao(i,j) = c*ao_overlap(i,j)*(ao_mono_elec_integral_diag(i) + &
ao_mono_elec_integral_diag(j))
enddo
Fock_matrix_ao(j,j) = Fock_matrix_alpha_ao(j,j)
enddo
TOUCH Fock_matrix_ao
mo_coef = eigenvectors_fock_matrix_mo

BIN
plugins/Hartree_Fock/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 59 KiB

After

Width:  |  Height:  |  Size: 57 KiB

5
plugins/MP2/EZFIO.cfg Normal file
View File

@ -0,0 +1,5 @@
[energy]
type: double precision
doc: MP2 energy
interface: ezfio

4
plugins/MP2/H_apply.irp.f
View File

@ -6,5 +6,9 @@ from perturbation import perturbations
s = H_apply("mp2")
s.set_perturbation("Moller_plesset")
print s
s = H_apply("mp2_selection")
s.set_selection_pt2("Moller_plesset")
print s
END_SHELL

46
plugins/MP2/README.rst
View File

@ -40,3 +40,49 @@ Needed Modules
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
* `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/src/SingleRefMethod>`_
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation>`_
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/plugins/Selectors_full>`_
* `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/plugins/SingleRefMethod>`_
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
h_apply_mp2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_mp2_diexc
Undocumented
h_apply_mp2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_mp2_diexcp
Undocumented
h_apply_mp2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`mp2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/MP2/mp2.irp.f#L1>`_
Undocumented

31
plugins/MP2/mp2_wf.irp.f Normal file
View File

@ -0,0 +1,31 @@
program mp2_wf
implicit none
BEGIN_DOC
! Save the MP2 wave function
END_DOC
integer :: i,k
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer :: N_st, iter
N_st = N_states
allocate (pt2(N_st), norm_pert(N_st), H_pert_diag(N_st))
pt2 = 1.d0
selection_criterion = 1.e-12
selection_criterion_min = 1.e-12
TOUCH selection_criterion_min selection_criterion selection_criterion_factor
call H_apply_mp2_selection(pt2, norm_pert, H_pert_diag, N_st)
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
touch N_det psi_det psi_coef
print*,'N_det = ',N_det
print*,'-----'
print *, 'PT2 = ', pt2(1)
print *, 'E = ', HF_energy
print *, 'E_before +PT2 = ', HF_energy+pt2(1)
N_det = min(N_det,N_det_max)
call save_wavefunction
call ezfio_set_mp2_energy(HF_energy+pt2(1))
deallocate(pt2,norm_pert,H_pert_diag)
end

BIN
plugins/MRCC_CASSD/tree_dependency.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 127 KiB

214
plugins/MRCC_Utils/README.rst
View File

@ -1,171 +1,3 @@
===========
MRCC Module
===========
Needed Modules
==============
.. Do not edit this section. It was auto-generated from the
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/src/Perturbation>`_
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
* `Generators_full <http://github.com/LCPQ/quantum_package/tree/master/src/Generators_full>`_
* `Psiref_Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Psiref_Utils>`_
Documentation
=============
.. Do not edit this section. It was auto-generated from the
.. by the `update_README.py` script.
`ci_eigenvectors_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L166>`_
Eigenvectors/values of the CI matrix
`ci_eigenvectors_s2_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L167>`_
Eigenvectors/values of the CI matrix
`ci_electronic_energy_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L165>`_
Eigenvectors/values of the CI matrix
`ci_energy_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L232>`_
N_states lowest eigenvalues of the dressed CI matrix
`davidson_diag_hjj_mrcc <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/davidson.irp.f#L56>`_
Davidson diagonalization with specific diagonal elements of the H matrix
.br
H_jj : specific diagonal H matrix elements to diagonalize de Davidson
.br
dets_in : bitmasks corresponding to determinants
.br
u_in : guess coefficients on the various states. Overwritten
on exit
.br
dim_in : leftmost dimension of u_in
.br
sze : Number of determinants
.br
N_st : Number of eigenstates
.br
iunit : Unit for the I/O
.br
Initial guess vectors are not necessarily orthonormal
`davidson_diag_mrcc <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/davidson.irp.f#L1>`_
Davidson diagonalization.
.br
dets_in : bitmasks corresponding to determinants
.br
u_in : guess coefficients on the various states. Overwritten
on exit
.br
dim_in : leftmost dimension of u_in
.br
sze : Number of determinants
.br
N_st : Number of eigenstates
.br
iunit : Unit number for the I/O
.br
Initial guess vectors are not necessarily orthonormal
`delta_ii <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L104>`_
Dressing matrix in N_det basis
`delta_ij <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L103>`_
Dressing matrix in N_det basis
`diagonalize_ci_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L247>`_
Replace the coefficients of the CI states by the coefficients of the
eigenstates of the CI matrix
`find_triples_and_quadruples <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_dress.irp.f#L206>`_
Undocumented
`h_apply_mrcc <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/H_apply.irp.f_shell_27#L422>`_
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_mrcc_diexc <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/H_apply.irp.f_shell_27#L1>`_
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_mrcc_monoexc <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/H_apply.irp.f_shell_27#L273>`_
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_matrix_dressed <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L140>`_
Dressed H with Delta_ij
`h_u_0_mrcc <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/davidson.irp.f#L360>`_
Computes v_0 = H|u_0>
.br
n : number of determinants
.br
H_jj : array of <j|H|j>
`lambda_mrcc <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L5>`_
cm/<Psi_0|H|D_m> or perturbative 1/Delta_E(m)
`lambda_mrcc_tmp <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L81>`_
Undocumented
`lambda_pert <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L6>`_
cm/<Psi_0|H|D_m> or perturbative 1/Delta_E(m)
`mrcc_dress <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_dress.irp.f#L15>`_
Undocumented
`mrcc_dress_simple <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_dress.irp.f#L160>`_
Undocumented
`mrcc_iterations <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_general.irp.f#L7>`_
Undocumented
`oscillations <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L86>`_
Undocumented
`pert_determinants <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_utils.irp.f#L1>`_
Undocumented
`psi_ref_lock <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_dress.irp.f#L3>`_
Locks on ref determinants to fill delta_ij
`run_mrcc <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_general.irp.f#L1>`_
Undocumented
`set_generators_bitmasks_as_holes_and_particles <http://github.com/LCPQ/quantum_package/tree/master/src/MRCC_Utils/mrcc_general.irp.f#L69>`_
Undocumented
Needed Modules
==============
.. Do not edit this section It was auto-generated
@ -258,10 +90,6 @@ Documentation
N_states lowest eigenvalues of the dressed CI matrix
`create_minilist <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/mrcc_dress.irp.f#L17>`_
Undocumented
`davidson_diag_hjj_mrcc <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/davidson.irp.f#L59>`_
Davidson diagonalization with specific diagonal elements of the H matrix
.br
@ -374,7 +202,7 @@ Documentation
Find A.C = B
`find_triples_and_quadruples <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/mrcc_dress.irp.f#L271>`_
`find_triples_and_quadruples <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/mrcc_dress.irp.f#L245>`_
Undocumented
@ -400,6 +228,22 @@ Documentation
Undocumented
`gen_det_idx <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/mrcc_dress.irp.f#L237>`_
Undocumented
`gen_det_shortcut <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/mrcc_dress.irp.f#L235>`_
Undocumented
`gen_det_sorted <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/mrcc_dress.irp.f#L234>`_
Undocumented
`gen_det_version <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/mrcc_dress.irp.f#L236>`_
Undocumented
`get_pseudo_inverse <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/LinearAlgebra.irp.f#L95>`_
Find C = A^-1
@ -435,26 +279,26 @@ Documentation
Undocumented
`h_apply_mrcc <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/H_apply.irp.f_shell_27#L553>`_
h_apply_mrcc
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_mrcc_diexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/H_apply.irp.f_shell_27#L3>`_
h_apply_mrcc_diexc
Undocumented
`h_apply_mrcc_diexcorg <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/H_apply.irp.f_shell_27#L127>`_
h_apply_mrcc_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_mrcc_diexcp <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/H_apply.irp.f_shell_27#L100>`_
h_apply_mrcc_diexcp
Undocumented
`h_apply_mrcc_monoexc <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/H_apply.irp.f_shell_27#L399>`_
h_apply_mrcc_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
@ -464,15 +308,7 @@ Documentation
Dressed H with Delta_ij
`h_u_0_mrcc <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/davidson.irp.f#L371>`_
Computes v_0 = H|u_0>
.br
n : number of determinants
.br
H_jj : array of <j|H|j>
`h_u_0_mrcc_org <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/davidson.irp.f#L475>`_
`h_u_0_mrcc <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/davidson.irp.f#L367>`_
Computes v_0 = H|u_0>
.br
n : number of determinants
@ -781,10 +617,6 @@ Documentation
Undocumented
`mrcc_dress_simple <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/mrcc_dress.irp.f#L225>`_
Undocumented
`mrcc_iterations <http://github.com/LCPQ/quantum_package/tree/master/plugins/MRCC_Utils/mrcc_general.irp.f#L7>`_
Undocumented

198
plugins/MRCC_Utils/davidson.irp.f
View File

@ -105,10 +105,7 @@ subroutine davidson_diag_hjj_mrcc(dets_in,u_in,H_jj,energies,dim_in,sze,N_st,Nin
double precision :: to_print(2,N_st)
double precision :: cpu, wall
integer(bit_kind) :: dets_in_sorted(Nint,2,sze)
integer :: idx(sze), shortcut(0:sze+1),sh,ii,tmp
PROVIDE det_connections
!PROVIDE det_connections
call write_time(iunit)
call wall_time(wall)
@ -154,8 +151,6 @@ subroutine davidson_diag_hjj_mrcc(dets_in,u_in,H_jj,energies,dim_in,sze,N_st,Nin
! ==============
dets_in_sorted(:,:,:) = dets_in(:,:,:)
call sort_dets_ab(dets_in_sorted, idx, shortcut, sze, Nint)
k_pairs=0
do l=1,N_st
@ -215,8 +210,7 @@ subroutine davidson_diag_hjj_mrcc(dets_in,u_in,H_jj,energies,dim_in,sze,N_st,Nin
! ----------------------
do k=1,N_st
call H_u_0_mrcc(W(1,k,iter),U(1,k,iter),H_jj,sze,dets_in_sorted,shortcut,idx,Nint,istate)
!call H_u_0_mrcc_org(W(1,k,iter),U(1,k,iter),H_jj,sze,dets_in,Nint,istate)
call H_u_0_mrcc(W(1,k,iter),U(1,k,iter),H_jj,sze,dets_in,Nint,istate)
enddo
! Compute h_kl = <u_k | W_l> = <u_k| H |u_l>
@ -368,7 +362,9 @@ subroutine davidson_diag_hjj_mrcc(dets_in,u_in,H_jj,energies,dim_in,sze,N_st,Nin
abort_here = abort_all
end
subroutine H_u_0_mrcc(v_0,u_0,H_jj,n,keys_tmp,shortcut,sort_idx,Nint,istate)
subroutine H_u_0_mrcc(v_0,u_0,H_jj,n,keys_tmp,Nint,istate)
use bitmasks
implicit none
BEGIN_DOC
@ -389,8 +385,11 @@ subroutine H_u_0_mrcc(v_0,u_0,H_jj,n,keys_tmp,shortcut,sort_idx,Nint,istate)
integer :: i,j,k,l, jj,ii
integer :: i0, j0
integer,intent(in) :: shortcut(0:n+1), sort_idx(n)
integer :: tmp, warp(2,0:n+1), sh, ni
integer :: shortcut(0:n+1), sort_idx(n)
integer(bit_kind) :: sorted(Nint,n), version(Nint,n)
integer :: sh, sh2, ni, exa, ext, org_i, org_j, endi, pass
!
ASSERT (Nint > 0)
@ -399,12 +398,14 @@ subroutine H_u_0_mrcc(v_0,u_0,H_jj,n,keys_tmp,shortcut,sort_idx,Nint,istate)
PROVIDE ref_bitmask_energy delta_ij
integer, parameter :: block_size = 157
!$OMP PARALLEL DEFAULT(NONE) &
!$OMP PRIVATE(i,hij,j,k,idx,jj,vt,ii,warp,tmp,sh) &
!$OMP SHARED(n_det_ref,n_det_non_ref,idx_ref,idx_non_ref,n,H_jj,u_0,keys_tmp,Nint,v_0,istate,delta_ij,shortcut,sort_idx)
!$OMP PARALLEL DEFAULT(NONE) &
!$OMP PRIVATE(i,hij,j,k,idx,jj,vt,ii,sh, sh2, ni, exa, ext, org_i, org_j, endi, pass) &
!$OMP SHARED(n_det_ref,n_det_non_ref,idx_ref,idx_non_ref,n,H_jj,u_0,keys_tmp,Nint,v_0,istate,delta_ij,sorted,shortcut,sort_idx,version)
!$OMP DO SCHEDULE(static)
!$OMP DO SCHEDULE(static)
do i=1,n
v_0(i) = H_jj(i) * u_0(i)
enddo
@ -412,43 +413,73 @@ subroutine H_u_0_mrcc(v_0,u_0,H_jj,n,keys_tmp,shortcut,sort_idx,Nint,istate)
allocate(idx(0:n), vt(n))
Vt = 0.d0
!$OMP SINGLE
call sort_dets_ab_v(keys_tmp, sorted, sort_idx, shortcut, version, n, Nint)
!$OMP END SINGLE
!$OMP DO SCHEDULE(dynamic)
do sh=1,shortcut(0)
warp(1,0) = 0
do ii=1,sh!shortcut(0)
tmp = 0
do ni=1,Nint
tmp = popcnt(xor(keys_tmp(ni,1, shortcut(ii)), keys_tmp(ni,1,shortcut(sh))))
end do
if(tmp <= 4) then
tmp = warp(1,0) + 1
warp(1,0) = tmp
warp(1,tmp) = shortcut(ii)
warp(2,tmp) = shortcut(ii+1)-1
end if
do sh2=1,sh
exa = 0
do ni=1,Nint
exa += popcnt(xor(version(ni,sh), version(ni,sh2)))
end do
if(exa > 2) then
cycle
end if
do ii=shortcut(sh),shortcut(sh+1)-1
idx(0) = ii
do i=shortcut(sh),shortcut(sh+1)-1
if(sh==sh2) then
endi = i-1
else
endi = shortcut(sh2+1)-1
end if
!call filter_connected_davidson_mwen(keys_tmp,shortcut,keys_tmp(1,1,ii),Nint,ii-1,idx)
call filter_connected_davidson_warp(keys_tmp,warp,keys_tmp(1,1,ii),Nint,ii-1,idx)
i = sort_idx(ii)
do jj=1,idx(0)
j = sort_idx(idx(jj))
!j = idx(jj)
if ( (dabs(u_0(j)) > 1.d-7).or.((dabs(u_0(i)) > 1.d-7)) ) then
call i_H_j(keys_tmp(1,1,idx(jj)),keys_tmp(1,1,ii),Nint,hij)
vt (i) = vt (i) + hij*u_0(j)
vt (j) = vt (j) + hij*u_0(i)
endif
enddo
enddo
do j=shortcut(sh2),endi
ext = exa
do ni=1,Nint
ext += popcnt(xor(sorted(ni,i), sorted(ni,j)))
end do
if(ext <= 4) then
org_i = sort_idx(i)
org_j = sort_idx(j)
call i_H_j(keys_tmp(1,1,org_j),keys_tmp(1,1,org_i),Nint,hij)
vt (org_i) = vt (org_i) + hij*u_0(org_j)
vt (org_j) = vt (org_j) + hij*u_0(org_i)
end if
end do
end do
end do
enddo
!$OMP END DO
!$OMP SINGLE
call sort_dets_ba_v(keys_tmp, sorted, sort_idx, shortcut, version, n, Nint)
!$OMP END SINGLE
!$OMP DO SCHEDULE(dynamic)
do sh=1,shortcut(0)
do i=shortcut(sh),shortcut(sh+1)-1
do j=shortcut(sh),i-1
ext = 0
do ni=1,Nint
ext += popcnt(xor(sorted(ni,i), sorted(ni,j)))
end do
if(ext == 4) then
org_i = sort_idx(i)
org_j = sort_idx(j)
call i_H_j(keys_tmp(1,1,org_j),keys_tmp(1,1,org_i),Nint,hij)
vt (org_i) = vt (org_i) + hij*u_0(org_j)
vt (org_j) = vt (org_j) + hij*u_0(org_i)
end if
end do
end do
enddo
!$OMP END DO
!$OMP DO SCHEDULE(guided)
@ -470,80 +501,3 @@ subroutine H_u_0_mrcc(v_0,u_0,H_jj,n,keys_tmp,shortcut,sort_idx,Nint,istate)
!$OMP END PARALLEL
end
subroutine H_u_0_mrcc_org(v_0,u_0,H_jj,n,keys_tmp,Nint,istate)
use bitmasks
implicit none
BEGIN_DOC
! Computes v_0 = H|u_0>
!
! n : number of determinants
!
! H_jj : array of <j|H|j>
END_DOC
integer, intent(in) :: n,Nint,istate
double precision, intent(out) :: v_0(n)
double precision, intent(in) :: u_0(n)
double precision, intent(in) :: H_jj(n)
integer(bit_kind),intent(in) :: keys_tmp(Nint,2,n)
integer, allocatable :: idx(:)
double precision :: hij
double precision, allocatable :: vt(:)
integer :: i,j,k,l, jj,ii
integer :: i0, j0
ASSERT (Nint > 0)
ASSERT (Nint == N_int)
ASSERT (n>0)
PROVIDE ref_bitmask_energy delta_ij
integer, parameter :: block_size = 157
!$OMP PARALLEL DEFAULT(NONE) &
!$OMP PRIVATE(i,hij,j,k,idx,jj,ii,vt) &
!$OMP SHARED(n_det_ref,n_det_non_ref,idx_ref,idx_non_ref,n,H_jj,u_0,keys_tmp,Nint,v_0,istate,delta_ij)
!$OMP DO SCHEDULE(static)
do i=1,n
v_0(i) = H_jj(i) * u_0(i)
enddo
!$OMP END DO
allocate(idx(0:n), vt(n))
Vt = 0.d0
!$OMP DO SCHEDULE(guided)
do i=1,n
idx(0) = i
call filter_connected_davidson(keys_tmp,keys_tmp(1,1,i),Nint,i-1,idx)
do jj=1,idx(0)
j = idx(jj)
! if ( (dabs(u_0(j)) > 1.d-7).or.((dabs(u_0(i)) > 1.d-7)) ) then
call i_H_j(keys_tmp(1,1,j),keys_tmp(1,1,i),Nint,hij)
hij = hij
vt (i) = vt (i) + hij*u_0(j)
vt (j) = vt (j) + hij*u_0(i)
! endif
enddo
enddo
!$OMP END DO
!$OMP DO SCHEDULE(guided)
do ii=1,n_det_ref
i = idx_ref(ii)
do jj = 1, n_det_non_ref
j = idx_non_ref(jj)
vt (i) = vt (i) + delta_ij(ii,jj,istate)*u_0(j)
vt (j) = vt (j) + delta_ij(ii,jj,istate)*u_0(i)
enddo
enddo
!$OMP END DO
!$OMP CRITICAL
do i=1,n
v_0(i) = v_0(i) + vt(i)
enddo
!$OMP END CRITICAL
deallocate(idx,vt)
!$OMP END PARALLEL
end

154
plugins/MRCC_Utils/mrcc_dress.irp.f
View File

@ -14,54 +14,6 @@ BEGIN_PROVIDER [ integer(omp_lock_kind), psi_ref_lock, (psi_det_size) ]
END_PROVIDER
subroutine create_minilist(key_mask, fullList, miniList, idx_miniList, N_fullList, N_miniList, Nint)
use bitmasks
implicit none
integer(bit_kind), intent(in) :: fullList(Nint, 2, N_fullList)
integer, intent(in) :: N_fullList
integer(bit_kind),intent(out) :: miniList(Nint, 2, N_fullList)
integer,intent(out) :: idx_miniList(N_fullList), N_miniList
integer, intent(in) :: Nint
integer(bit_kind) :: key_mask(Nint, 2)
integer :: ni, i, n_a, n_b, e_a, e_b
n_a = 0
n_b = 0
do ni=1,nint
n_a = n_a + popcnt(key_mask(ni,1))
n_b = n_b + popcnt(key_mask(ni,2))
end do
if(n_a == 0) then
N_miniList = N_fullList
miniList(:,:,:) = fullList(:,:,:)
do i=1,N_fullList
idx_miniList(i) = i
end do
return
end if
N_miniList = 0
do i=1,N_fullList
e_a = n_a
e_b = n_b
do ni=1,nint
e_a -= popcnt(iand(fullList(ni, 1, i), key_mask(ni, 1)))
e_b -= popcnt(iand(fullList(ni, 2, i), key_mask(ni, 2)))
end do
if(e_a + e_b <= 2) then
N_miniList = N_miniList + 1
miniList(:,:,N_miniList) = fullList(:,:,i)
idx_miniList(N_miniList) = i
end if
end do
end subroutine
subroutine mrcc_dress(delta_ij_, delta_ii_, Ndet_ref, Ndet_non_ref,i_generator,n_selected,det_buffer,Nint,iproc,key_mask)
use bitmasks
implicit none
@ -75,11 +27,10 @@ subroutine mrcc_dress(delta_ij_, delta_ii_, Ndet_ref, Ndet_non_ref,i_generator,n
integer :: i,j,k,l
integer :: degree_alpha(psi_det_size)
integer :: idx_alpha(0:psi_det_size)
logical :: good
logical :: good, fullMatch
integer(bit_kind) :: tq(Nint,2,n_selected)
integer :: N_tq, c_ref ,degree
integer :: connected_to_ref
double precision :: hIk, hla, hIl, dIk(N_states), dka(N_states), dIa(N_states)
double precision, allocatable :: dIa_hla(:,:)
@ -91,11 +42,24 @@ subroutine mrcc_dress(delta_ij_, delta_ii_, Ndet_ref, Ndet_non_ref,i_generator,n
integer :: iint, ipos
integer :: i_state, k_sd, l_sd, i_I, i_alpha
integer(bit_kind) :: miniList(Nint, 2, N_det_non_ref), key_mask(Nint, 2)
integer :: idx_miniList(N_det_non_ref), N_miniList
integer(bit_kind),allocatable :: miniList(:,:,:)
integer(bit_kind),intent(in) :: key_mask(Nint, 2)
integer,allocatable :: idx_miniList(:)
integer :: N_miniList, ni, leng
call find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq)
leng = max(N_det_generators, N_det_non_ref)
allocate(miniList(Nint, 2, leng), idx_miniList(leng))
!create_minilist_find_previous(key_mask, fullList, miniList, N_fullList, N_miniList, fullMatch, Nint)
call create_minilist_find_previous(key_mask, psi_det_generators, miniList, i_generator-1, N_miniList, fullMatch, Nint)
if(fullMatch) then
return
end if
call find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq,miniList,N_minilist)
allocate (dIa_hla(N_states,Ndet_non_ref))
@ -214,61 +178,24 @@ subroutine mrcc_dress(delta_ij_, delta_ii_, Ndet_ref, Ndet_non_ref,i_generator,n
enddo
enddo
deallocate (dIa_hla)
deallocate(miniList, idx_miniList)
end
BEGIN_PROVIDER [ integer(bit_kind), gen_det_sorted, (N_int,2,N_det_generators,2) ]
&BEGIN_PROVIDER [ integer, gen_det_shortcut, (0:N_det_generators,2) ]
&BEGIN_PROVIDER [ integer, gen_det_version, (N_int, N_det_generators,2) ]
&BEGIN_PROVIDER [ integer, gen_det_idx, (N_det_generators,2) ]
gen_det_sorted(:,:,:,1) = psi_det_generators(:,:,:N_det_generators)
gen_det_sorted(:,:,:,2) = psi_det_generators(:,:,:N_det_generators)
call sort_dets_ab_v(gen_det_sorted(:,:,:,1), gen_det_idx(:,1), gen_det_shortcut(0:,1), gen_det_version(:,:,1), N_det_generators, N_int)
call sort_dets_ba_v(gen_det_sorted(:,:,:,2), gen_det_idx(:,2), gen_det_shortcut(0:,2), gen_det_version(:,:,2), N_det_generators, N_int)
END_PROVIDER
subroutine find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq,miniList,N_miniList)
subroutine mrcc_dress_simple(delta_ij_non_ref_,Ndet_non_ref,i_generator,n_selected,det_buffer,Nint,iproc)
use bitmasks
implicit none
integer, intent(in) :: i_generator,n_selected, Nint, iproc
integer, intent(in) :: Ndet_non_ref
double precision, intent(inout) :: delta_ij_non_ref_(Ndet_non_ref,Ndet_non_ref,*)
integer(bit_kind), intent(in) :: det_buffer(Nint,2,n_selected)
integer :: i,j,k,m
integer :: new_size
integer :: degree(psi_det_size)
integer :: idx(0:psi_det_size)
logical :: good
integer(bit_kind) :: tq(Nint,2,n_selected)
integer :: N_tq, c_ref
integer :: connected_to_ref
call find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq)
! Compute <k|H|a><a|H|j> / (E0 - Haa)
double precision :: hka, haa
double precision :: haj
double precision :: f(N_states)
do i=1,N_tq
call get_excitation_degree_vector(psi_non_ref,tq(1,1,i),degree,Nint,Ndet_non_ref,idx)
call i_h_j(tq(1,1,i),tq(1,1,i),Nint,haa)
do m=1,N_states
f(m) = 1.d0/(ci_electronic_energy(m)-haa)
enddo
do k=1,idx(0)
call i_h_j(tq(1,1,i),psi_non_ref(1,1,idx(k)),Nint,hka)
do j=k,idx(0)
call i_h_j(tq(1,1,i),psi_non_ref(1,1,idx(j)),Nint,haj)
do m=1,N_states
delta_ij_non_ref_(idx(k), idx(j),m) += haj*hka* f(m)
delta_ij_non_ref_(idx(j), idx(k),m) += haj*hka* f(m)
enddo
enddo
enddo
enddo
end
subroutine find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq)
use bitmasks
implicit none
@ -283,18 +210,24 @@ subroutine find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq
integer(bit_kind), intent(out) :: tq(Nint,2,n_selected)
integer, intent(out) :: N_tq
integer :: c_ref
integer :: connected_to_ref
integer :: nt,ni
logical, external :: is_connected_to
integer(bit_kind),intent(in) :: miniList(Nint,2,N_det_generators)
integer,intent(in) :: N_miniList
N_tq = 0
do i=1,N_selected
c_ref = connected_to_ref(det_buffer(1,1,i),psi_det_generators,Nint, &
i_generator,N_det_generators)
if (c_ref /= 0) then
i_loop : do i=1,N_selected
if(is_connected_to(det_buffer(1,1,i), miniList, Nint, N_miniList)) then
cycle
endif
end if
! Select determinants that are triple or quadruple excitations
! from the ref
@ -316,8 +249,7 @@ subroutine find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq
enddo
endif
endif
enddo
enddo i_loop
end

BIN
plugins/MRCC_Utils/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 110 KiB

After

Width:  |  Height:  |  Size: 106 KiB

36
plugins/Molden/README.rst
View File

@ -38,3 +38,39 @@ Needed Modules
* `MO_Basis <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis>`_
* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `MO_Basis <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis>`_
* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
`print_mos <http://github.com/LCPQ/quantum_package/tree/master/plugins/Molden/print_mo.irp.f#L1>`_
Undocumented
`write_ao_basis <http://github.com/LCPQ/quantum_package/tree/master/plugins/Molden/print_mo.irp.f#L63>`_
Undocumented
`write_geometry <http://github.com/LCPQ/quantum_package/tree/master/plugins/Molden/print_mo.irp.f#L45>`_
Undocumented
`write_intro_gamess <http://github.com/LCPQ/quantum_package/tree/master/plugins/Molden/print_mo.irp.f#L26>`_
Undocumented
`write_mo_basis <http://github.com/LCPQ/quantum_package/tree/master/plugins/Molden/print_mo.irp.f#L112>`_
Undocumented

42
plugins/Perturbation/Moller_plesset.irp.f
View File

@ -1,42 +0,0 @@
subroutine pt2_moller_plesset(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st)
use bitmasks
implicit none
integer, intent(in) :: Nint,ndet,n_st
integer(bit_kind), intent(in) :: det_pert(Nint,2)
double precision , intent(out) :: c_pert(n_st),e_2_pert(n_st),H_pert_diag(N_st)
double precision :: i_H_psi_array(N_st)
BEGIN_DOC
! compute the standard Moller-Plesset perturbative first order coefficient and second order energetic contribution
!
! for the various n_st states.
!
! c_pert(i) = <psi(i)|H|det_pert>/(difference of orbital energies)
!
! e_2_pert(i) = <psi(i)|H|det_pert>^2/(difference of orbital energies)
!
END_DOC
integer :: i,j
double precision :: diag_H_mat_elem
integer :: exc(0:2,2,2)
integer :: degree
double precision :: phase,delta_e,h
integer :: h1,h2,p1,p2,s1,s2
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
call get_excitation(ref_bitmask,det_pert,exc,degree,phase,Nint)
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
delta_e = Fock_matrix_diag_mo(h1) + Fock_matrix_diag_mo(h2) - &
(Fock_matrix_diag_mo(p1) + Fock_matrix_diag_mo(p2))
delta_e = 1.d0/delta_e
call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det,psi_selectors_size,n_st,i_H_psi_array)
h = diag_H_mat_elem(det_pert,Nint)
do i =1,n_st
H_pert_diag(i) = h
c_pert(i) = i_H_psi_array(i) *delta_e
e_2_pert(i) = c_pert(i) * i_H_psi_array(i)
enddo
end

50
plugins/Perturbation/README.rst
View File

@ -107,92 +107,92 @@ Documentation
Undocumented
`perturb_buffer_by_mono_delta_rho_one_point <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L161>`_
perturb_buffer_by_mono_delta_rho_one_point
Applly pertubration ``delta_rho_one_point`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_by_mono_dipole_moment_z <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L266>`_
perturb_buffer_by_mono_dipole_moment_z
Applly pertubration ``dipole_moment_z`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_by_mono_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L371>`_
perturb_buffer_by_mono_epstein_nesbet
Applly pertubration ``epstein_nesbet`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_by_mono_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L476>`_
perturb_buffer_by_mono_epstein_nesbet_2x2
Applly pertubration ``epstein_nesbet_2x2`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_by_mono_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L791>`_
perturb_buffer_by_mono_epstein_nesbet_sc2
Applly pertubration ``epstein_nesbet_sc2`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_by_mono_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L686>`_
perturb_buffer_by_mono_epstein_nesbet_sc2_no_projected
Applly pertubration ``epstein_nesbet_sc2_no_projected`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_by_mono_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L581>`_
perturb_buffer_by_mono_epstein_nesbet_sc2_projected
Applly pertubration ``epstein_nesbet_sc2_projected`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_by_mono_h_core <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L896>`_
perturb_buffer_by_mono_h_core
Applly pertubration ``h_core`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_by_mono_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L56>`_
perturb_buffer_by_mono_moller_plesset
Applly pertubration ``moller_plesset`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_delta_rho_one_point <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L110>`_
perturb_buffer_delta_rho_one_point
Applly pertubration ``delta_rho_one_point`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_dipole_moment_z <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L215>`_
perturb_buffer_dipole_moment_z
Applly pertubration ``dipole_moment_z`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L320>`_
perturb_buffer_epstein_nesbet
Applly pertubration ``epstein_nesbet`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L425>`_
perturb_buffer_epstein_nesbet_2x2
Applly pertubration ``epstein_nesbet_2x2`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L740>`_
perturb_buffer_epstein_nesbet_sc2
Applly pertubration ``epstein_nesbet_sc2`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L635>`_
perturb_buffer_epstein_nesbet_sc2_no_projected
Applly pertubration ``epstein_nesbet_sc2_no_projected`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L530>`_
perturb_buffer_epstein_nesbet_sc2_projected
Applly pertubration ``epstein_nesbet_sc2_projected`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_h_core <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L845>`_
perturb_buffer_h_core
Applly pertubration ``h_core`` to the buffer of determinants generated in the H_apply
routine.
`perturb_buffer_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/perturbation.irp.f_shell_13#L5>`_
perturb_buffer_moller_plesset
Applly pertubration ``moller_plesset`` to the buffer of determinants generated in the H_apply
routine.
@ -239,7 +239,7 @@ Documentation
.br
`pt2_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/epstein_nesbet.irp.f#L1>`_
`pt2_epstein_nesbet <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pt2_equations.irp.f_template_363#L3>`_
compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
.br
for the various N_st states.
@ -250,7 +250,7 @@ Documentation
.br
`pt2_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/epstein_nesbet.irp.f#L45>`_
`pt2_epstein_nesbet_2x2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pt2_equations.irp.f_template_363#L59>`_
compute the Epstein-Nesbet 2x2 diagonalization coefficient and energetic contribution
.br
for the various N_st states.
@ -261,7 +261,7 @@ Documentation
.br
`pt2_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pert_sc2.irp.f#L186>`_
`pt2_epstein_nesbet_sc2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pt2_equations.irp.f_template_363#L353>`_
compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
.br
for the various N_st states, but with the CISD_SC2 energies and coefficients
@ -272,7 +272,7 @@ Documentation
.br
`pt2_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pert_sc2.irp.f#L87>`_
`pt2_epstein_nesbet_sc2_no_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pt2_equations.irp.f_template_363#L275>`_
compute the Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
.br
for the various N_st states,
@ -296,7 +296,7 @@ Documentation
H_pert_diag = <HF|H|det_pert> c_pert
`pt2_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pert_sc2.irp.f#L2>`_
`pt2_epstein_nesbet_sc2_projected <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pt2_equations.irp.f_template_363#L181>`_
compute the Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
.br
for the various N_st states,
@ -336,7 +336,7 @@ Documentation
than pt2_max in absolute value
`pt2_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/Moller_plesset.irp.f#L1>`_
`pt2_moller_plesset <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pt2_equations.irp.f_template_363#L120>`_
compute the standard Moller-Plesset perturbative first order coefficient and second order energetic contribution
.br
for the various n_st states.
@ -352,7 +352,7 @@ Documentation
provided.
`repeat_all_e_corr <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pert_sc2.irp.f#L156>`_
`repeat_all_e_corr <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation/pert_sc2.irp.f#L1>`_
Undocumented

9
plugins/Perturbation/delta_rho_perturbation.irp.f
View File

@ -1,4 +1,4 @@
subroutine pt2_delta_rho_one_point(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st)
subroutine pt2_delta_rho_one_point(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st,minilist,idx_minilist,N_minilist)
use bitmasks
implicit none
integer, intent(in) :: Nint,ndet,n_st
@ -7,6 +7,10 @@ subroutine pt2_delta_rho_one_point(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,nde
double precision :: i_O1_psi_array(N_st)
double precision :: i_H_psi_array(N_st)
integer, intent(in) :: N_minilist
integer, intent(in) :: idx_minilist(0:N_det_selectors)
integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors)
BEGIN_DOC
! compute the perturbatibe contribution to the Integrated Spin density at z = z_one point of one determinant
!
@ -46,7 +50,8 @@ subroutine pt2_delta_rho_one_point(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,nde
! endif
call i_O1_psi_alpha_beta(mo_integrated_delta_rho_one_point,det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_O1_psi_array)
call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
!call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
call i_H_psi_minilist(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem(det_pert,Nint)
oii = diag_O1_mat_elem_alpha_beta(mo_integrated_delta_rho_one_point,det_pert,N_int)

10
plugins/Perturbation/dipole_moment.irp.f
View File

@ -1,4 +1,4 @@
subroutine pt2_dipole_moment_z(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st)
subroutine pt2_dipole_moment_z(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_st,minilist,idx_minilist,N_minilist)
use bitmasks
implicit none
integer, intent(in) :: Nint,ndet,n_st
@ -7,6 +7,10 @@ subroutine pt2_dipole_moment_z(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_
double precision :: i_O1_psi_array(N_st)
double precision :: i_H_psi_array(N_st)
integer, intent(in) :: N_minilist
integer, intent(in) :: idx_minilist(0:N_det_selectors)
integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors)
BEGIN_DOC
! compute the perturbatibe contribution to the dipole moment of one determinant
!
@ -46,7 +50,9 @@ subroutine pt2_dipole_moment_z(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,n_
! endif
call i_O1_psi(mo_dipole_z,det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_O1_psi_array)
call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
!call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
call i_H_psi_minilist(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem(det_pert,Nint)
oii = diag_O1_mat_elem(mo_dipole_z,det_pert,N_int)

93
plugins/Perturbation/epstein_nesbet.irp.f
View File

@ -1,93 +0,0 @@
subroutine pt2_epstein_nesbet(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st)
use bitmasks
implicit none
integer, intent(in) :: Nint,ndet,N_st
integer(bit_kind), intent(in) :: det_pert(Nint,2)
double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st)
double precision :: i_H_psi_array(N_st)
BEGIN_DOC
! compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
!
! for the various N_st states.
!
! c_pert(i) = <psi(i)|H|det_pert>/( E(i) - <det_pert|H|det_pert> )
!
! e_2_pert(i) = <psi(i)|H|det_pert>^2/( E(i) - <det_pert|H|det_pert> )
!
END_DOC
integer :: i,j
double precision :: diag_H_mat_elem, h
PROVIDE selection_criterion
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem(det_pert,Nint)
do i =1,N_st
if(CI_electronic_energy(i)>h.and.CI_electronic_energy(i).ne.0.d0)then
c_pert(i) = -1.d0
e_2_pert(i) = selection_criterion*selection_criterion_factor*2.d0
else if (dabs(CI_electronic_energy(i) - h) > 1.d-6) then
c_pert(i) = i_H_psi_array(i) / (CI_electronic_energy(i) - h)
H_pert_diag(i) = h*c_pert(i)*c_pert(i)
e_2_pert(i) = c_pert(i) * i_H_psi_array(i)
else
c_pert(i) = -1.d0
e_2_pert(i) = -dabs(i_H_psi_array(i))
H_pert_diag(i) = h
endif
enddo
end
subroutine pt2_epstein_nesbet_2x2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st)
use bitmasks
implicit none
integer, intent(in) :: Nint,ndet,N_st
integer(bit_kind), intent(in) :: det_pert(Nint,2)
double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st)
double precision :: i_H_psi_array(N_st)
BEGIN_DOC
! compute the Epstein-Nesbet 2x2 diagonalization coefficient and energetic contribution
!
! for the various N_st states.
!
! e_2_pert(i) = 0.5 * (( <det_pert|H|det_pert> - E(i) ) - sqrt( ( <det_pert|H|det_pert> - E(i)) ^2 + 4 <psi(i)|H|det_pert>^2 )
!
! c_pert(i) = e_2_pert(i)/ <psi(i)|H|det_pert>
!
END_DOC
integer :: i,j
double precision :: diag_H_mat_elem,delta_e, h
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
PROVIDE CI_electronic_energy
call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem(det_pert,Nint)
do i =1,N_st
if (i_H_psi_array(i) /= 0.d0) then
delta_e = h - CI_electronic_energy(i)
if (delta_e > 0.d0) then
e_2_pert(i) = 0.5d0 * (delta_e - dsqrt(delta_e * delta_e + 4.d0 * i_H_psi_array(i) * i_H_psi_array(i)))
else
e_2_pert(i) = 0.5d0 * (delta_e + dsqrt(delta_e * delta_e + 4.d0 * i_H_psi_array(i) * i_H_psi_array(i)))
endif
if (dabs(i_H_psi_array(i)) > 1.d-6) then
c_pert(i) = e_2_pert(i)/i_H_psi_array(i)
else
c_pert(i) = 0.d0
endif
H_pert_diag(i) = h*c_pert(i)*c_pert(i)
else
e_2_pert(i) = 0.d0
c_pert(i) = 0.d0
H_pert_diag(i) = 0.d0
endif
enddo
end

199
plugins/Perturbation/pert_sc2.irp.f
View File

@ -1,158 +1,3 @@
subroutine pt2_epstein_nesbet_SC2_projected(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st)
use bitmasks
implicit none
integer, intent(in) :: Nint,ndet,N_st
integer(bit_kind), intent(in) :: det_pert(Nint,2)
double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st)
double precision :: i_H_psi_array(N_st)
integer :: idx_repeat(0:ndet)
BEGIN_DOC
! compute the Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
!
! for the various N_st states,
!
! but with the correction in the denominator
!
! comming from the interaction of that determinant with all the others determinants
!
! that can be repeated by repeating all the double excitations
!
! : you repeat all the correlation energy already taken into account in CI_electronic_energy(1)
!
! that could be repeated to this determinant.
!
! In addition, for the perturbative energetic contribution you have the standard second order
!
! e_2_pert = <psi_i|H|det_pert>^2/(Delta_E)
!
! and also the purely projected contribution
!
! H_pert_diag = <HF|H|det_pert> c_pert
END_DOC
integer :: i,j,degree,l
double precision :: diag_H_mat_elem,accu_e_corr,hij,h0j,h,delta_E
double precision :: repeat_all_e_corr,accu_e_corr_tmp,e_2_pert_fonda
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
call i_H_psi_SC2(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array,idx_repeat)
accu_e_corr = 0.d0
!$IVDEP
do i = 1, idx_repeat(0)
accu_e_corr = accu_e_corr + E_corr_per_selectors(idx_repeat(i))
enddo
h = diag_H_mat_elem(det_pert,Nint) + accu_e_corr
delta_E = 1.d0/(CI_SC2_electronic_energy(1) - h)
c_pert(1) = i_H_psi_array(1) /(CI_SC2_electronic_energy(1) - h)
e_2_pert(1) = i_H_psi_array(1) * c_pert(1)
do i =2,N_st
H_pert_diag(i) = h
if (dabs(CI_SC2_electronic_energy(i) - h) > 1.d-6) then
c_pert(i) = i_H_psi_array(i) / (-dabs(CI_SC2_electronic_energy(i) - h))
e_2_pert(i) = (c_pert(i) * i_H_psi_array(i))
else
c_pert(i) = i_H_psi_array(i)
e_2_pert(i) = -dabs(i_H_psi_array(i))
endif
enddo
degree = popcnt(xor( ref_bitmask(1,1), det_pert(1,1))) + &
popcnt(xor( ref_bitmask(1,2), det_pert(1,2)))
!DEC$ NOUNROLL
do l=2,Nint
degree = degree+ popcnt(xor( ref_bitmask(l,1), det_pert(l,1))) + &
popcnt(xor( ref_bitmask(l,2), det_pert(l,2)))
enddo
if(degree==4)then
! <psi|delta_H|psi>
e_2_pert_fonda = e_2_pert(1)
H_pert_diag(1) = e_2_pert(1) * c_pert(1) * c_pert(1)
do i = 1, N_st
do j = 1, idx_repeat(0)
e_2_pert(i) += e_2_pert_fonda * psi_selectors_coef(idx_repeat(j),i) * psi_selectors_coef(idx_repeat(j),i)
enddo
enddo
endif
end
subroutine pt2_epstein_nesbet_SC2_no_projected(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st)
use bitmasks
implicit none
integer, intent(in) :: Nint,ndet,N_st
integer(bit_kind), intent(in) :: det_pert(Nint,2)
double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st)
double precision :: i_H_psi_array(N_st)
integer :: idx_repeat(0:ndet)
BEGIN_DOC
! compute the Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
!
! for the various N_st states,
!
! but with the correction in the denominator
!
! comming from the interaction of that determinant with all the others determinants
!
! that can be repeated by repeating all the double excitations
!
! : you repeat all the correlation energy already taken into account in CI_electronic_energy(1)
!
! that could be repeated to this determinant.
!
! In addition, for the perturbative energetic contribution you have the standard second order
!
! e_2_pert = <psi_i|H|det_pert>^2/(Delta_E)
!
! and also the purely projected contribution
!
! H_pert_diag = <HF|H|det_pert> c_pert
END_DOC
integer :: i,j,degree,l
double precision :: diag_H_mat_elem,accu_e_corr,hij,h0j,h,delta_E
double precision :: repeat_all_e_corr,accu_e_corr_tmp,e_2_pert_fonda
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
call i_H_psi_SC2(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array,idx_repeat)
accu_e_corr = 0.d0
!$IVDEP
do i = 1, idx_repeat(0)
accu_e_corr = accu_e_corr + E_corr_per_selectors(idx_repeat(i))
enddo
h = diag_H_mat_elem(det_pert,Nint) + accu_e_corr
delta_E = 1.d0/(CI_SC2_electronic_energy(1) - h)
c_pert(1) = i_H_psi_array(1) /(CI_SC2_electronic_energy(1) - h)
e_2_pert(1) = i_H_psi_array(1) * c_pert(1)
do i =2,N_st
H_pert_diag(i) = h
if (dabs(CI_SC2_electronic_energy(i) - h) > 1.d-6) then
c_pert(i) = i_H_psi_array(i) / (-dabs(CI_SC2_electronic_energy(i) - h))
e_2_pert(i) = (c_pert(i) * i_H_psi_array(i))
else
c_pert(i) = i_H_psi_array(i)
e_2_pert(i) = -dabs(i_H_psi_array(i))
endif
enddo
end
double precision function repeat_all_e_corr(key_in)
implicit none
integer(bit_kind), intent(in) :: key_in(N_int,2)
@ -182,47 +27,3 @@ double precision function repeat_all_e_corr(key_in)
end
subroutine pt2_epstein_nesbet_sc2(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st)
use bitmasks
implicit none
integer, intent(in) :: Nint,ndet,N_st
integer(bit_kind), intent(in) :: det_pert(Nint,2)
double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st)
double precision :: i_H_psi_array(N_st)
BEGIN_DOC
! compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
!
! for the various N_st states, but with the CISD_SC2 energies and coefficients
!
! c_pert(i) = <psi(i)|H|det_pert>/( E(i) - <det_pert|H|det_pert> )
!
! e_2_pert(i) = <psi(i)|H|det_pert>^2/( E(i) - <det_pert|H|det_pert> )
!
END_DOC
integer :: i,j
double precision :: diag_H_mat_elem, h
PROVIDE selection_criterion
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem(det_pert,Nint)
do i =1,N_st
if(CI_SC2_electronic_energy(i)>h.and.CI_SC2_electronic_energy(i).ne.0.d0)then
c_pert(i) = -1.d0
e_2_pert(i) = selection_criterion*selection_criterion_factor*2.d0
else if (dabs(CI_SC2_electronic_energy(i) - h) > 1.d-6) then
c_pert(i) = i_H_psi_array(i) / (CI_SC2_electronic_energy(i) - h)
H_pert_diag(i) = h*c_pert(i)*c_pert(i)
e_2_pert(i) = c_pert(i) * i_H_psi_array(i)
else
c_pert(i) = -1.d0
e_2_pert(i) = -dabs(i_H_psi_array(i))
H_pert_diag(i) = h
endif
enddo
end

6
plugins/Perturbation/pert_single.irp.f
View File

@ -1,4 +1,4 @@
subroutine pt2_h_core(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st)
subroutine pt2_h_core(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist)
use bitmasks
implicit none
integer, intent(in) :: Nint,ndet,N_st
@ -6,6 +6,10 @@ subroutine pt2_h_core(det_pert,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st)
double precision , intent(out) :: c_pert(N_st),e_2_pert(N_st),H_pert_diag(N_st)
double precision :: i_H_psi_array(N_st)
integer, intent(in) :: N_minilist
integer, intent(in) :: idx_minilist(0:N_det_selectors)
integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors)
BEGIN_DOC
! compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
!

112
plugins/Perturbation/perturbation.template.f
View File

@ -2,7 +2,7 @@ BEGIN_SHELL [ /usr/bin/env python ]
import perturbation
END_SHELL
subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,coef_pert_buffer,sum_e_2_pert,sum_norm_pert,sum_H_pert_diag,N_st,Nint)
subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,coef_pert_buffer,sum_e_2_pert,sum_norm_pert,sum_H_pert_diag,N_st,Nint,key_mask,fock_diag_tmp)
implicit none
BEGIN_DOC
! Applly pertubration ``$PERT`` to the buffer of determinants generated in the H_apply
@ -11,49 +11,73 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c
integer, intent(in) :: Nint, N_st, buffer_size, i_generator
integer(bit_kind), intent(in) :: buffer(Nint,2,buffer_size)
integer(bit_kind),intent(in) :: key_mask(Nint,2)
double precision, intent(in) :: fock_diag_tmp(2,0:mo_tot_num)
double precision, intent(inout) :: sum_norm_pert(N_st),sum_e_2_pert(N_st)
double precision, intent(inout) :: coef_pert_buffer(N_st,buffer_size),e_2_pert_buffer(N_st,buffer_size),sum_H_pert_diag(N_st)
double precision :: c_pert(N_st), e_2_pert(N_st), H_pert_diag(N_st)
integer :: i,k, c_ref
integer :: i,k, c_ref, ni, ex
integer, external :: connected_to_ref
logical, external :: is_in_wavefunction
integer(bit_kind), allocatable :: minilist(:,:,:)
integer, allocatable :: idx_minilist(:)
integer :: N_minilist
integer(bit_kind), allocatable :: minilist_gen(:,:,:)
integer :: N_minilist_gen
logical :: fullMatch
logical, external :: is_connected_to
allocate( minilist(Nint,2,N_det_selectors), &
minilist_gen(Nint,2,N_det_generators), &
idx_minilist(N_det_selectors) )
ASSERT (Nint > 0)
ASSERT (Nint == N_int)
ASSERT (buffer_size >= 0)
ASSERT (minval(sum_norm_pert) >= 0.d0)
ASSERT (N_st > 0)
do i = 1,buffer_size
c_ref = connected_to_ref(buffer(1,1,i),psi_det_generators,Nint,i_generator,N_det_generators)
if (c_ref /= 0) then
call create_minilist(key_mask, psi_selectors, miniList, idx_miniList, N_det_selectors, N_minilist, Nint)
call create_minilist_find_previous(key_mask, psi_det_generators, miniList_gen, i_generator-1, N_minilist_gen, fullMatch, Nint)
if(fullMatch) then
deallocate( minilist, minilist_gen, idx_minilist )
return
end if
do i=1,buffer_size
if(is_connected_to(buffer(1,1,i), miniList_gen, Nint, N_minilist_gen)) then
cycle
endif
end if
if (is_in_wavefunction(buffer(1,1,i),Nint)) then
cycle
endif
integer :: degree
call get_excitation_degree(HF_bitmask,buffer(1,1,i),degree,N_int)
call pt2_$PERT(buffer(1,1,i), &
c_pert,e_2_pert,H_pert_diag,Nint,N_det_selectors,n_st)
call pt2_$PERT(psi_det_generators(1,1,i_generator),buffer(1,1,i), fock_diag_tmp, &
c_pert,e_2_pert,H_pert_diag,Nint,N_minilist,n_st,minilist,idx_minilist,N_minilist)
do k = 1,N_st
e_2_pert_buffer(k,i) = e_2_pert(k)
coef_pert_buffer(k,i) = c_pert(k)
sum_norm_pert(k) += c_pert(k) * c_pert(k)
sum_e_2_pert(k) += e_2_pert(k)
sum_H_pert_diag(k) += H_pert_diag(k)
e_2_pert_buffer(k,i) = e_2_pert(k)
coef_pert_buffer(k,i) = c_pert(k)
sum_norm_pert(k) = sum_norm_pert(k) + c_pert(k) * c_pert(k)
sum_e_2_pert(k) = sum_e_2_pert(k) + e_2_pert(k)
sum_H_pert_diag(k) = sum_H_pert_diag(k) + H_pert_diag(k)
enddo
enddo
enddo
deallocate( minilist, minilist_gen, idx_minilist )
end
subroutine perturb_buffer_by_mono_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,coef_pert_buffer,sum_e_2_pert,sum_norm_pert,sum_H_pert_diag,N_st,Nint)
subroutine perturb_buffer_by_mono_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,coef_pert_buffer,sum_e_2_pert,sum_norm_pert,sum_H_pert_diag,N_st,Nint,key_mask,fock_diag_tmp)
implicit none
BEGIN_DOC
! Applly pertubration ``$PERT`` to the buffer of determinants generated in the H_apply
@ -62,20 +86,46 @@ subroutine perturb_buffer_by_mono_$PERT(i_generator,buffer,buffer_size,e_2_pert_
integer, intent(in) :: Nint, N_st, buffer_size, i_generator
integer(bit_kind), intent(in) :: buffer(Nint,2,buffer_size)
integer(bit_kind),intent(in) :: key_mask(Nint,2)
double precision, intent(in) :: fock_diag_tmp(2,0:mo_tot_num)
double precision, intent(inout) :: sum_norm_pert(N_st),sum_e_2_pert(N_st)
double precision, intent(inout) :: coef_pert_buffer(N_st,buffer_size),e_2_pert_buffer(N_st,buffer_size),sum_H_pert_diag(N_st)
double precision :: c_pert(N_st), e_2_pert(N_st), H_pert_diag(N_st)
integer :: i,k, c_ref
integer :: i,k, c_ref, ni, ex
integer, external :: connected_to_ref_by_mono
logical, external :: is_in_wavefunction
integer(bit_kind), allocatable :: minilist(:,:,:)
integer, allocatable :: idx_minilist(:)
integer :: N_minilist
integer(bit_kind), allocatable :: minilist_gen(:,:,:)
integer :: N_minilist_gen
logical :: fullMatch
logical, external :: is_connected_to
allocate( minilist(Nint,2,N_det_selectors), &
minilist_gen(Nint,2,N_det_generators), &
idx_minilist(N_det_selectors) )
ASSERT (Nint > 0)
ASSERT (Nint == N_int)
ASSERT (buffer_size >= 0)
ASSERT (minval(sum_norm_pert) >= 0.d0)
ASSERT (N_st > 0)
do i = 1,buffer_size
call create_minilist(key_mask, psi_selectors, miniList, idx_miniList, N_det_selectors, N_minilist, Nint)
call create_minilist_find_previous(key_mask, psi_det_generators, miniList_gen, i_generator-1, N_minilist_gen, fullMatch, Nint)
if(fullMatch) then
deallocate( minilist, minilist_gen, idx_minilist )
return
end if
do i=1,buffer_size
c_ref = connected_to_ref_by_mono(buffer(1,1,i),psi_det_generators,Nint,i_generator,N_det)
if (c_ref /= 0) then
@ -86,19 +136,19 @@ subroutine perturb_buffer_by_mono_$PERT(i_generator,buffer,buffer_size,e_2_pert_
cycle
endif
integer :: degree
call pt2_$PERT(buffer(1,1,i), &
c_pert,e_2_pert,H_pert_diag,Nint,N_det_selectors,n_st)
call pt2_$PERT(psi_det_generators(1,1,i_generator),buffer(1,1,i), fock_diag_tmp, &
c_pert,e_2_pert,H_pert_diag,Nint,N_minilist,n_st,minilist,idx_minilist,N_minilist)
do k = 1,N_st
e_2_pert_buffer(k,i) = e_2_pert(k)
coef_pert_buffer(k,i) = c_pert(k)
sum_norm_pert(k) += c_pert(k) * c_pert(k)
sum_e_2_pert(k) += e_2_pert(k)
sum_H_pert_diag(k) += H_pert_diag(k)
e_2_pert_buffer(k,i) = e_2_pert(k)
coef_pert_buffer(k,i) = c_pert(k)
sum_norm_pert(k) = sum_norm_pert(k) + c_pert(k) * c_pert(k)
sum_e_2_pert(k) = sum_e_2_pert(k) + e_2_pert(k)
sum_H_pert_diag(k) = sum_H_pert_diag(k) + H_pert_diag(k)
enddo
enddo
enddo
deallocate( minilist, minilist_gen, idx_minilist )
end

367
plugins/Perturbation/pt2_equations.irp.f Normal file
View File

@ -0,0 +1,367 @@
BEGIN_TEMPLATE
subroutine pt2_epstein_nesbet ($arguments)
use bitmasks
implicit none
$declarations
BEGIN_DOC
! compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
!
! for the various N_st states.
!
! c_pert(i) = <psi(i)|H|det_pert>/( E(i) - <det_pert|H|det_pert> )
!
! e_2_pert(i) = <psi(i)|H|det_pert>^2/( E(i) - <det_pert|H|det_pert> )
!
END_DOC
integer :: i,j
double precision :: diag_H_mat_elem_fock, h
double precision :: i_H_psi_array(N_st)
PROVIDE selection_criterion
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
!call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
call i_H_psi_minilist(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
do i =1,N_st
if(CI_electronic_energy(i)>h.and.CI_electronic_energy(i).ne.0.d0)then
c_pert(i) = -1.d0
e_2_pert(i) = selection_criterion*selection_criterion_factor*2.d0
else if (dabs(CI_electronic_energy(i) - h) > 1.d-6) then
c_pert(i) = i_H_psi_array(i) / (CI_electronic_energy(i) - h)
H_pert_diag(i) = h*c_pert(i)*c_pert(i)
e_2_pert(i) = c_pert(i) * i_H_psi_array(i)
else
c_pert(i) = -1.d0
e_2_pert(i) = -dabs(i_H_psi_array(i))
H_pert_diag(i) = h
endif
enddo
end
subroutine pt2_epstein_nesbet_2x2 ($arguments)
use bitmasks
implicit none
$declarations
BEGIN_DOC
! compute the Epstein-Nesbet 2x2 diagonalization coefficient and energetic contribution
!
! for the various N_st states.
!
! e_2_pert(i) = 0.5 * (( <det_pert|H|det_pert> - E(i) ) - sqrt( ( <det_pert|H|det_pert> - E(i)) ^2 + 4 <psi(i)|H|det_pert>^2 )
!
! c_pert(i) = e_2_pert(i)/ <psi(i)|H|det_pert>
!
END_DOC
integer :: i,j
double precision :: diag_H_mat_elem_fock,delta_e, h
double precision :: i_H_psi_array(N_st)
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
PROVIDE CI_electronic_energy
!call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
call i_H_psi_minilist(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
do i =1,N_st
if (i_H_psi_array(i) /= 0.d0) then
delta_e = h - CI_electronic_energy(i)
if (delta_e > 0.d0) then
e_2_pert(i) = 0.5d0 * (delta_e - dsqrt(delta_e * delta_e + 4.d0 * i_H_psi_array(i) * i_H_psi_array(i)))
else
e_2_pert(i) = 0.5d0 * (delta_e + dsqrt(delta_e * delta_e + 4.d0 * i_H_psi_array(i) * i_H_psi_array(i)))
endif
if (dabs(i_H_psi_array(i)) > 1.d-6) then
c_pert(i) = e_2_pert(i)/i_H_psi_array(i)
else
c_pert(i) = 0.d0
endif
H_pert_diag(i) = h*c_pert(i)*c_pert(i)
else
e_2_pert(i) = 0.d0
c_pert(i) = 0.d0
H_pert_diag(i) = 0.d0
endif
enddo
end
subroutine pt2_moller_plesset ($arguments)
use bitmasks
implicit none
$declarations
BEGIN_DOC
! compute the standard Moller-Plesset perturbative first order coefficient and second order energetic contribution
!
! for the various n_st states.
!
! c_pert(i) = <psi(i)|H|det_pert>/(difference of orbital energies)
!
! e_2_pert(i) = <psi(i)|H|det_pert>^2/(difference of orbital energies)
!
END_DOC
integer :: i,j
double precision :: diag_H_mat_elem_fock
integer :: exc(0:2,2,2)
integer :: degree
double precision :: phase,delta_e,h
double precision :: i_H_psi_array(N_st)
integer :: h1,h2,p1,p2,s1,s2
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
call get_excitation(ref_bitmask,det_pert,exc,degree,phase,Nint)
if (degree == 2) then
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
delta_e = Fock_matrix_diag_mo(h1) + Fock_matrix_diag_mo(h2) - &
(Fock_matrix_diag_mo(p1) + Fock_matrix_diag_mo(p2))
delta_e = 1.d0/delta_e
else if (degree == 1) then
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
delta_e = Fock_matrix_diag_mo(h1) - Fock_matrix_diag_mo(p1)
delta_e = 1.d0/delta_e
else
delta_e = 0.d0
endif
call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det,psi_selectors_size,n_st,i_H_psi_array)
h = diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
do i =1,n_st
H_pert_diag(i) = h
c_pert(i) = i_H_psi_array(i) *delta_e
e_2_pert(i) = c_pert(i) * i_H_psi_array(i)
enddo
end
subroutine pt2_epstein_nesbet_SC2_projected ($arguments)
use bitmasks
implicit none
$declarations
BEGIN_DOC
! compute the Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
!
! for the various N_st states,
!
! but with the correction in the denominator
!
! comming from the interaction of that determinant with all the others determinants
!
! that can be repeated by repeating all the double excitations
!
! : you repeat all the correlation energy already taken into account in CI_electronic_energy(1)
!
! that could be repeated to this determinant.
!
! In addition, for the perturbative energetic contribution you have the standard second order
!
! e_2_pert = <psi_i|H|det_pert>^2/(Delta_E)
!
! and also the purely projected contribution
!
! H_pert_diag = <HF|H|det_pert> c_pert
END_DOC
double precision :: i_H_psi_array(N_st)
integer :: idx_repeat(0:ndet)
integer :: i,j,degree,l
double precision :: diag_H_mat_elem_fock,accu_e_corr,hij,h0j,h,delta_E
double precision :: repeat_all_e_corr,accu_e_corr_tmp,e_2_pert_fonda
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
call i_H_psi_SC2(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array,idx_repeat)
accu_e_corr = 0.d0
!$IVDEP
do i = 1, idx_repeat(0)
accu_e_corr = accu_e_corr + E_corr_per_selectors(idx_repeat(i))
enddo
h = diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
h = h + accu_e_corr
delta_E = 1.d0/(CI_SC2_electronic_energy(1) - h)
c_pert(1) = i_H_psi_array(1) /(CI_SC2_electronic_energy(1) - h)
e_2_pert(1) = i_H_psi_array(1) * c_pert(1)
do i =2,N_st
H_pert_diag(i) = h
if (dabs(CI_SC2_electronic_energy(i) - h) > 1.d-6) then
c_pert(i) = i_H_psi_array(i) / (-dabs(CI_SC2_electronic_energy(i) - h))
e_2_pert(i) = (c_pert(i) * i_H_psi_array(i))
else
c_pert(i) = i_H_psi_array(i)
e_2_pert(i) = -dabs(i_H_psi_array(i))
endif
enddo
degree = popcnt(xor( ref_bitmask(1,1), det_pert(1,1))) + &
popcnt(xor( ref_bitmask(1,2), det_pert(1,2)))
!DEC$ NOUNROLL
do l=2,Nint
degree = degree+ popcnt(xor( ref_bitmask(l,1), det_pert(l,1))) + &
popcnt(xor( ref_bitmask(l,2), det_pert(l,2)))
enddo
if(degree==4)then
! <psi|delta_H|psi>
e_2_pert_fonda = e_2_pert(1)
H_pert_diag(1) = e_2_pert(1) * c_pert(1) * c_pert(1)
do i = 1, N_st
do j = 1, idx_repeat(0)
e_2_pert(i) += e_2_pert_fonda * psi_selectors_coef(idx_repeat(j),i) * psi_selectors_coef(idx_repeat(j),i)
enddo
enddo
endif
end
subroutine pt2_epstein_nesbet_SC2_no_projected ($arguments)
use bitmasks
implicit none
$declarations
BEGIN_DOC
! compute the Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
!
! for the various N_st states,
!
! but with the correction in the denominator
!
! comming from the interaction of that determinant with all the others determinants
!
! that can be repeated by repeating all the double excitations
!
! : you repeat all the correlation energy already taken into account in CI_electronic_energy(1)
!
! that could be repeated to this determinant.
!
! In addition, for the perturbative energetic contribution you have the standard second order
!
! e_2_pert = <psi_i|H|det_pert>^2/(Delta_E)
!
! and also the purely projected contribution
!
! H_pert_diag = <HF|H|det_pert> c_pert
END_DOC
double precision :: i_H_psi_array(N_st)
integer :: idx_repeat(0:ndet)
integer :: i,j,degree,l
double precision :: diag_H_mat_elem_fock,accu_e_corr,hij,h0j,h,delta_E
double precision :: repeat_all_e_corr,accu_e_corr_tmp,e_2_pert_fonda
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
call i_H_psi_SC2(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array,idx_repeat)
accu_e_corr = 0.d0
!$IVDEP
do i = 1, idx_repeat(0)
accu_e_corr = accu_e_corr + E_corr_per_selectors(idx_repeat(i))
enddo
h = diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
h = h + accu_e_corr
delta_E = 1.d0/(CI_SC2_electronic_energy(1) - h)
c_pert(1) = i_H_psi_array(1) /(CI_SC2_electronic_energy(1) - h)
e_2_pert(1) = i_H_psi_array(1) * c_pert(1)
do i =2,N_st
H_pert_diag(i) = h
if (dabs(CI_SC2_electronic_energy(i) - h) > 1.d-6) then
c_pert(i) = i_H_psi_array(i) / (-dabs(CI_SC2_electronic_energy(i) - h))
e_2_pert(i) = (c_pert(i) * i_H_psi_array(i))
else
c_pert(i) = i_H_psi_array(i)
e_2_pert(i) = -dabs(i_H_psi_array(i))
endif
enddo
end
subroutine pt2_epstein_nesbet_sc2 ($arguments)
use bitmasks
implicit none
$declarations
BEGIN_DOC
! compute the standard Epstein-Nesbet perturbative first order coefficient and second order energetic contribution
!
! for the various N_st states, but with the CISD_SC2 energies and coefficients
!
! c_pert(i) = <psi(i)|H|det_pert>/( E(i) - <det_pert|H|det_pert> )
!
! e_2_pert(i) = <psi(i)|H|det_pert>^2/( E(i) - <det_pert|H|det_pert> )
!
END_DOC
integer :: i,j
double precision :: i_H_psi_array(N_st)
double precision :: diag_H_mat_elem_fock, h
PROVIDE selection_criterion
ASSERT (Nint == N_int)
ASSERT (Nint > 0)
!call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det_selectors,psi_selectors_size,N_st,i_H_psi_array)
call i_H_psi_minilist(det_pert,minilist,idx_minilist,N_minilist,psi_selectors_coef,Nint,N_minilist,psi_selectors_size,N_st,i_H_psi_array)
h = diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
do i =1,N_st
if(CI_SC2_electronic_energy(i)>h.and.CI_SC2_electronic_energy(i).ne.0.d0)then
c_pert(i) = -1.d0
e_2_pert(i) = selection_criterion*selection_criterion_factor*2.d0
else if (dabs(CI_SC2_electronic_energy(i) - h) > 1.d-6) then
c_pert(i) = i_H_psi_array(i) / (CI_SC2_electronic_energy(i) - h)
H_pert_diag(i) = h*c_pert(i)*c_pert(i)
e_2_pert(i) = c_pert(i) * i_H_psi_array(i)
else
c_pert(i) = -1.d0
e_2_pert(i) = -dabs(i_H_psi_array(i))
H_pert_diag(i) = h
endif
enddo
end
SUBST [ arguments, declarations ]
det_ref,det_pert,fock_diag_tmp,c_pert,e_2_pert,H_pert_diag,Nint,ndet,N_st,minilist,idx_minilist,N_minilist ;
integer, intent(in) :: Nint
integer, intent(in) :: ndet
integer, intent(in) :: N_st
integer, intent(in) :: N_minilist
integer(bit_kind), intent(in) :: det_ref (Nint,2)
integer(bit_kind), intent(in) :: det_pert(Nint,2)
double precision , intent(in) :: fock_diag_tmp(2,mo_tot_num+1)
double precision , intent(out) :: c_pert(N_st)
double precision , intent(out) :: e_2_pert(N_st)
double precision, intent(out) :: H_pert_diag(N_st)
integer, intent(in) :: idx_minilist(0:N_det_selectors)
integer(bit_kind), intent(in) :: minilist(Nint,2,N_det_selectors)
;;
END_TEMPLATE
! Note : If the arguments are changed here, they should also be changed accordingly in
! the perturbation.template.f file.

BIN
plugins/Perturbation/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 78 KiB

After

Width:  |  Height:  |  Size: 76 KiB

BIN
plugins/Properties/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 62 KiB

After

Width:  |  Height:  |  Size: 60 KiB

BIN
plugins/Psiref_CAS/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 7.6 KiB

After

Width:  |  Height:  |  Size: 7.4 KiB

BIN
plugins/Psiref_Utils/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 3.5 KiB

After

Width:  |  Height:  |  Size: 3.4 KiB

3
plugins/QmcChem/.gitignore vendored
View File

@ -20,4 +20,5 @@ ezfio_interface.irp.f
irpf90.make
irpf90_entities
save_for_qmcchem
tags
tags
target_pt2_qmc

8
plugins/QmcChem/README.rst
View File

@ -66,6 +66,14 @@ Documentation
title="f(|r-r_A|) = \int Y_{lm}^{C} (|r-r_C|, \Omega_C) \chi_i^{A} (r-r_A) d\Omega_C" />
`compute_energy <http://github.com/LCPQ/quantum_package/tree/master/plugins/QmcChem/target_pt2_qmc.irp.f#L80>`_
Compute an energy when a threshold is applied
`e_curve <http://github.com/LCPQ/quantum_package/tree/master/plugins/QmcChem/target_pt2_qmc.irp.f#L1>`_
Undocumented
`mo_pseudo_grid <http://github.com/LCPQ/quantum_package/tree/master/plugins/QmcChem/pot_ao_pseudo_ints.irp.f#L56>`_
Grid points for f(|r-r_A|) = \int Y_{lm}^{C} (|r-r_C|, \Omega_C) \phi_i^{A} (r-r_A) d\Omega_C
.br

1
plugins/QmcChem/save_for_qmcchem.irp.f
View File

@ -1,6 +1,7 @@
program save_for_qmc
read_wf = .True.
TOUCH read_wf
print *, "N_det = ", N_det
call write_spindeterminants
if (do_pseudo) then
call write_pseudopotential

121
plugins/QmcChem/target_pt2_qmc.irp.f Normal file
View File

@ -0,0 +1,121 @@
program e_curve
use bitmasks
implicit none
integer :: i,j,k, nab, m, l, n_up, n_dn, n
double precision :: norm, E, hij, num, ci, cj
integer, allocatable :: iorder(:)
double precision , allocatable :: norm_sort(:), psi_bilinear_matrix_values_save(:)
nab = n_det_alpha_unique+n_det_beta_unique
allocate ( norm_sort(0:nab), iorder(0:nab), psi_bilinear_matrix_values_save(N_det) )
norm_sort(0) = 0.d0
iorder(0) = 0
do i=1,n_det_alpha_unique
norm_sort(i) = det_alpha_norm(i)
iorder(i) = i
enddo
do i=1,n_det_beta_unique
norm_sort(i+n_det_alpha_unique) = det_beta_norm(i)
iorder(i+n_det_alpha_unique) = -i
enddo
call dsort(norm_sort(1),iorder(1),nab)
if (.not.read_wf) then
stop 'Please set read_wf to true'
endif
psi_bilinear_matrix_values_save = psi_bilinear_matrix_values(:,1)
print *, '=========================================================='
print '(A8,2X,A8,2X,A12,2X,A10,2X,A12)', 'Thresh.', 'Ndet', 'Cost', 'Norm', 'E'
print *, '=========================================================='
integer(bit_kind), allocatable :: det_i(:,:), det_j(:,:)
double precision :: thresh, E_min, E_max, E_prev
thresh = 0.d0
call compute_energy(psi_bilinear_matrix_values_save,E_max,m,norm)
call i_h_j(psi_det_sorted(1,1,1), psi_det_sorted(1,1,1), N_int, E_min)
print *, E_min, E_max
n_up = nab
n_dn = 0
do while (n_up > n_dn)
n = n_dn + (n_up-n_dn)/2
psi_bilinear_matrix_values_save = psi_bilinear_matrix_values(:,1)
do j=1,n
i = iorder(j)
if (i<0) then
do k=1,n_det
if (psi_bilinear_matrix_columns(k) == -i) then
psi_bilinear_matrix_values_save(k) = 0.d0
endif
enddo
else
do k=1,n_det
if (psi_bilinear_matrix_rows(k) == i) then
psi_bilinear_matrix_values_save(k) = 0.d0
endif
enddo
endif
enddo
call compute_energy(psi_bilinear_matrix_values_save,E,m,norm)
print '(E9.1,2X,I8,2X,F10.2,2X,F10.6,2X,F12.6)', norm_sort(n), m, &
dble( elec_alpha_num**3 + elec_alpha_num**2 * m ) / &
dble( elec_alpha_num**3 + elec_alpha_num**2 * n ), norm, E
if (E < target_energy) then
n_dn = n+1
else
n_up = n
endif
enddo
print *, '=========================================================='
print *, norm_sort(n), target_energy
deallocate (iorder, norm_sort, psi_bilinear_matrix_values_save)
end
subroutine compute_energy(psi_bilinear_matrix_values_save, E, m, norm)
implicit none
BEGIN_DOC
! Compute an energy when a threshold is applied
END_DOC
double precision, intent(in) :: psi_bilinear_matrix_values_save(n_det)
integer(bit_kind), allocatable :: det_i(:,:), det_j(:,:)
integer :: i,j, k, l, m
double precision :: num, norm, ci, cj, hij, E
num = 0.d0
norm = 0.d0
m = 0
!$OMP PARALLEL DEFAULT(SHARED) PRIVATE(k,l,det_i,det_j,ci,cj,hij) REDUCTION(+:norm,m,num)
allocate( det_i(N_int,2), det_j(N_int,2))
!$OMP DO schedule(guided)
do k=1,n_det
if (psi_bilinear_matrix_values_save(k) == 0.d0) then
cycle
endif
ci = psi_bilinear_matrix_values_save(k)
det_i(:,1) = psi_det_alpha_unique(:,psi_bilinear_matrix_rows(k))
det_i(:,2) = psi_det_beta_unique(:,psi_bilinear_matrix_columns(k))
do l=1,n_det
if (psi_bilinear_matrix_values_save(l) == 0.d0) then
cycle
endif
cj = psi_bilinear_matrix_values_save(l)
det_j(:,1) = psi_det_alpha_unique(:,psi_bilinear_matrix_rows(l))
det_j(:,2) = psi_det_beta_unique(:,psi_bilinear_matrix_columns(l))
call i_h_j(det_i, det_j, N_int, hij)
num = num + ci*cj*hij
enddo
norm = norm + ci*ci
m = m+1
enddo
!$OMP END DO
deallocate (det_i,det_j)
!$OMP END PARALLEL
E = num / norm + nuclear_repulsion
end

BIN
plugins/Selectors_full/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 74 KiB

After

Width:  |  Height:  |  Size: 72 KiB

BIN
plugins/SingleRefMethod/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 33 KiB

After

Width:  |  Height:  |  Size: 33 KiB

22
plugins/loc_cele/README.rst
View File

@ -20,3 +20,25 @@ Needed Modules
* `MO_Basis <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis>`_
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `MO_Basis <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis>`_
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
`loc_rasorb <http://github.com/LCPQ/quantum_package/tree/master/plugins/loc_cele/loc_cele.irp.f#L1>`_
This program performs a localization of the active orbitals
of a CASSCF wavefunction, reading the orbitals from a RASORB
file of molcas.
id1=max is the number of MO in a given symmetry.

2
plugins/loc_cele/loc_cele.irp.f
View File

@ -280,7 +280,7 @@
enddo !big loop over symmetry
10 format (4E18.12)
10 format (4E20.12)
! Now we copyt the newcmo into the mo_coef

10
scripts/compilation/qp_create_ninja.py
View File

@ -35,7 +35,9 @@ except ImportError:
from qp_path import QP_ROOT, QP_SRC, QP_EZFIO
EZFIO_LIB = join(QP_ROOT, "lib", "libezfio.a")
LIB = "" # join(QP_ROOT, "lib", "rdtsc.o")
EZFIO_LIB = join(QP_ROOT, "lib", "libezfio_irp.a")
ZMQ_LIB = join(QP_ROOT, "lib", "libf77zmq.a") + " " + join(QP_ROOT, "lib", "libzmq.a") + " -lstdc++ -lrt"
ROOT_BUILD_NINJA = join(QP_ROOT, "config", "build.ninja")
header = r"""#
@ -94,7 +96,7 @@ def ninja_create_env_variable(pwd_config_file):
l_string.append(str_)
lib_lapack = get_compilation_option(pwd_config_file, "LAPACK_LIB")
l_string.append("{0} = {1} {2}".format("LIB", lib_lapack, EZFIO_LIB))
l_string.append("LIB = {0} {1} {2} {3}".format(LIB, lib_lapack, EZFIO_LIB, ZMQ_LIB))
l_string.append("")
@ -260,7 +262,7 @@ def ninja_ezfio_rule():
l_flag = ["export {0}='${0}'".format(flag)
for flag in ["FC", "FCFLAGS", "IRPF90"]]
install_lib_ezfio = join(QP_ROOT, 'install', 'EZFIO', "lib", "libezfio.a")
install_lib_ezfio = join(QP_ROOT, 'install', 'EZFIO', "lib", "libezfio_irp.a")
l_cmd = ["cd {0}".format(QP_EZFIO)] + l_flag
l_cmd += ["rm -f make.config ; ninja && ln -sf {0} {1}".format(install_lib_ezfio, EZFIO_LIB)]
@ -707,7 +709,7 @@ def ninja_dot_tree_rule():
l_string = [
"rule build_dot_tree", " command = {0}".format(" ; ".join(l_cmd)),
" generator = 1",
" description = Generate Png representtion of the Tree Dependencies of $module_rel",
" description = Generating Png representation of the Tree Dependencies of $module_rel",
""
]

21
scripts/generate_h_apply.py
View File

@ -60,7 +60,7 @@ class H_apply(object):
s["omp_master"] = "!$OMP MASTER"
s["omp_end_master"] = "!$OMP END MASTER"
s["omp_barrier"] = "!$OMP BARRIER"
s["omp_do"] = "!$OMP DO SCHEDULE (static)"
s["omp_do"] = "!$OMP DO SCHEDULE (static,1)"
s["omp_enddo"] = "!$OMP ENDDO NOWAIT"
d = { True : '.True.', False : '.False.'}
@ -99,7 +99,7 @@ class H_apply(object):
deallocate(H_jj,iorder)
"""
s["size_max"] = str(1024*128)
s["size_max"] = "256"
s["copy_buffer"] = """call copy_H_apply_buffer_to_wf
if (s2_eig) then
call make_s2_eigenfunction
@ -131,10 +131,10 @@ class H_apply(object):
def filter_vvvv_excitation(self):
self["filter_vvvv_excitation"] = """
key_union_hole_part = 0_bit_kind
call set_bite_to_integer(i_a,key_union_hole_part,N_int)
call set_bite_to_integer(j_a,key_union_hole_part,N_int)
call set_bite_to_integer(i_b,key_union_hole_part,N_int)
call set_bite_to_integer(j_b,key_union_hole_part,N_int)
call set_bit_to_integer(i_a,key_union_hole_part,N_int)
call set_bit_to_integer(j_a,key_union_hole_part,N_int)
call set_bit_to_integer(i_b,key_union_hole_part,N_int)
call set_bit_to_integer(j_b,key_union_hole_part,N_int)
do jtest_vvvv = 1, N_int
if(iand(key_union_hole_part(jtest_vvvv),virt_bitmask(jtest_vvvv,1).ne.key_union_hole_part(jtest_vvvv)))then
b_cycle = .False.
@ -157,7 +157,6 @@ class H_apply(object):
def set_filter_2h_2p(self):
self["filter2h2p"] = """
! ! DIR$ FORCEINLINE
if (is_a_two_holes_two_particles(key)) cycle
"""
@ -201,11 +200,11 @@ class H_apply(object):
"""
self.data["size_max"] = "256"
self.data["initialization"] = """
PROVIDE CI_electronic_energy psi_selectors_coef psi_selectors E_corr_per_selectors psi_det_sorted_bit
PROVIDE psi_selectors_coef psi_selectors E_corr_per_selectors psi_det_sorted_bit
"""
self.data["keys_work"] = """
call perturb_buffer_%s(i_generator,keys_out,key_idx,e_2_pert_buffer,coef_pert_buffer,sum_e_2_pert, &
sum_norm_pert,sum_H_pert_diag,N_st,N_int)
sum_norm_pert,sum_H_pert_diag,N_st,N_int,key_mask,fock_diag_tmp)
"""%(pert,)
self.data["finalization"] = """
"""
@ -219,7 +218,7 @@ class H_apply(object):
double precision, intent(inout):: norm_pert(N_st)
double precision, intent(inout):: H_pert_diag(N_st)
double precision :: delta_pt2(N_st), norm_psi(N_st), pt2_old(N_st)
PROVIDE CI_electronic_energy N_det_generators
PROVIDE N_det_generators
do k=1,N_st
pt2(k) = 0.d0
norm_pert(k) = 0.d0
@ -266,7 +265,7 @@ class H_apply(object):
double precision, intent(inout) :: select_max_out"""
self.data["params_post"] += ", select_max(min(i_generator,size(select_max,1)))"
self.data["size_max"] = str(1024*128)
self.data["size_max"] = "256"
self.data["copy_buffer"] = """
call copy_H_apply_buffer_to_wf
if (s2_eig) then

3
scripts/get_basis.sh
View File

@ -8,7 +8,8 @@
# Prints in stdout the name of a temporary file containing the basis set.
#
#DEBUG:
#echo $0 $@ 1>&2
if [[ -z ${QP_ROOT} ]]
then

1
scripts/module/create_gitignore.sh
View File

@ -13,7 +13,6 @@ then
fi
source ${QP_ROOT}/scripts/qp_include.sh
function do_gitingore()
{
cat << EOF > .gitignore

13
scripts/module/module_handler.py
View File

@ -88,7 +88,7 @@ def get_l_module_descendant(d_child, l_module):
except KeyError:
print >> sys.stderr, "Error: "
print >> sys.stderr, "`{0}` is not a submodule".format(module)
print >> sys.stderr, "Check the typo (orthograph, case, '/', etc.) "
print >> sys.stderr, "Check the typo (spelling, case, '/', etc.) "
sys.exit(1)
return list(set(l))
@ -180,6 +180,11 @@ class ModuleHandler():
def create_png(self, l_module):
"""Create the png of the dependency tree for a l_module"""
# Don't update if we are not in the main repository
from is_master_repository import is_master_repository
if not is_master_repository:
return
basename = "tree_dependency"
path = '{0}.png'.format(basename)
@ -289,6 +294,12 @@ if __name__ == '__main__':
pass
if arguments["create_git_ignore"]:
# Don't update if we are not in the main repository
from is_master_repository import is_master_repository
if not is_master_repository:
sys.exit()
path = os.path.join(module_abs, ".gitignore")
with open(path, "w+") as f:

55
scripts/module/qp_module.py
View File

@ -57,10 +57,48 @@ def save_new_module(path, l_child):
f.write(D_KEY["needed_module"])
f.write(D_KEY["documentation"])
with open(os.path.join(path, "%s.main.irp.f"%(module_name) ), "w") as f:
f.write("program {0}".format(module_name) )
f.write(""" implicit none
BEGIN_DOC
! TODO
END_DOC
print *, ' _/ '
print *, ' -:\_?, _Jm####La '
print *, 'J"(:" > _]#AZ#Z#UUZ##, '
print *, '_,::./ %(|i%12XmX1*1XL _?, '
print *, ' \..\ _\(vmWQwodY+ia%lnL _",/ ( '
print *, ' .:< ]J=mQD?WXn<uQWmmvd, -.-:=!\'
print *, ' "{Z jC]QW|=3Zv)Bi3BmXv3 = _7'
print *, ' ]h[Z6)WQ;)jZs]C;|$BZv+, : ./ '
print *, ' -#sJX%$Wmm#ev]hinW#Xi:` c ; '
print *, ' #X#X23###1}vI$WWmX1>|,)nr" '
print *, ' 4XZ#Xov1v}=)vnXAX1nnv;1n" '
print *, ' ]XX#ZXoovvvivnnnlvvo2*i7 '
print *, ' "23Z#1S2oo2XXSnnnoSo2>v" '
print *, ' miX#L -~`""!!1}oSoe|i7 '
print *, ' 4cn#m, v221=|v[ '
print *, ' ]hI3Zma,;..__wXSe=+vo '
print *, ' ]Zov*XSUXXZXZXSe||vo2 '
print *, ' ]Z#><iiii|i||||==vn2( '
print *, ' ]Z#i<ii||+|=||=:{no2[ '
print *, ' ]ZUsiiiiivi|=||=vo22[ '
print *, ' ]XZvlliiIi|i=|+|vooo '
print *, ' =v1llli||||=|||||lii( '
print *, ' ]iillii||||||||=>=|< '
print *, ' -ziiiii||||||+||==+> '
print *, ' -%|+++||=|=+|=|==/ '
print *, ' -a>====+|====-:- '
print *, ' "~,- -- /- '
print *, ' -. )> '
print *, ' .~ +- '
print *, ' . .... : . '
print *, ' -------~ '
print *, ''
end
""")
if __name__ == '__main__':
arguments = docopt(__doc__)
def main(arguments):
if arguments["list"]:
if arguments["--installed"]:
@ -107,12 +145,14 @@ if __name__ == '__main__':
save_new_module(path, l_child_reduce)
print " [ OK ]"
print "Your module is created in the `plugins` directory."
print "You need to create some `.irp.f` to be able to install it."
# print "` {0} install {1} `".format(os.path.basename(__file__), name)
print ""
arguments["create"]=False
arguments["install"]=True
main(arguments)
elif arguments["download"]:
print "Not yet implemented"
pass
# d_local = get_dict_child([QP_SRC])
# d_remote = get_dict_child(arguments["<path_folder>"])
@ -205,3 +245,8 @@ if __name__ == '__main__':
except OSError:
print "%s is a core module which can't be removed" % module
if __name__ == '__main__':
arguments = docopt(__doc__)
main(arguments)

26
scripts/module/qp_update_readme.py
View File

@ -44,7 +44,7 @@ def get_url(path_module_rel):
elif is_module(path_module_rel):
url = "http://github.com/LCPQ/quantum_package/tree/master/src"
else:
print "{0} Is not a valide module nor plugin".format(path_module_rel)
print "{0} Is not a valid module nor plugin".format(path_module_rel)
sys.exit(1)
return os.path.join(url, path_module_rel)
@ -155,20 +155,32 @@ def update_documentation(d_readmen, root_module):
l_doc = []
for irp in d_info[path]:
url = os.path.join(get_url(os.path.basename(path)), irp.file)
doc = extract_doc(root_module, irp.provider)
l_doc += ["`{0} <{1}#L{2}>`_".format(irp.provider, url, irp.line),
doc,
""]
if ".irp.f_shell_" in irp.file:
l_doc += ["{0}".format(irp.provider),
doc,
""]
else:
l_doc += ["`{0} <{1}#L{2}>`_".format(irp.provider, url, irp.line),
doc,
""]
l_doc_section = [D_KEY["documentation"], '',
"\n".join(l_doc)]
d_readme[path]["documentation"] = "\n".join(l_doc_section)
if __name__ == '__main__':
# Update documentation only if the remote repository is
# the main repository
from is_master_repository import is_master_repository
if not is_master_repository:
sys.exit(0)
arguments = docopt(__doc__)
if arguments["--root_module"]:
@ -188,8 +200,8 @@ if __name__ == '__main__':
fetch_splitted_data(d_readme, l_module_readme)
except IOError:
print l_module_readme, "is not a module and/or",
print "have not a `README.rst` file inside"
print "Abort..."
print "has not a `README.rst` file inside"
print "Aborting..."
sys.exit(1)
update_needed(d_readme)

3
scripts/pseudo/elts_num_ele.py
View File

@ -1,4 +1,5 @@
name_to_elec = {"H": 1,
name_to_elec = {"X": 0,
"H": 1,
"He": 2,
"Li": 3,
"Be": 4,

32
scripts/pseudo/put_pseudo_in_ezfio.py
View File

@ -58,17 +58,35 @@ def get_pseudo_str(l_atom):
str_ = ""
for a in l_atom:
l_cmd_atom = ["--atom", a]
l_cmd_head = [EMSL_path, "get_basis_data",
"--db_path", db_path,
"--basis", "BFD-Pseudo"]
if a is not 'X':
l_cmd_atom = ["--atom", a]
process = Popen(l_cmd_head + l_cmd_atom, stdout=PIPE, stderr=PIPE)
l_cmd_head = [EMSL_path, "get_basis_data",
"--db_path", db_path,
"--basis", "BFD-Pseudo"]
stdout, _ = process.communicate()
str_ += stdout.strip() + "\n"
process = Popen(l_cmd_head + l_cmd_atom, stdout=PIPE, stderr=PIPE)
stdout, _ = process.communicate()
str_ += stdout.strip() + "\n"
else: # Dummy atoms
str_ += """Element Symbol: X
Number of replaced protons: 0
Number of projectors: 0
Pseudopotential data:
Local component:
Coeff. r^n Exp.
0.0 -1 0.
0.0 1 0.
0.0 0 0.
Non-local component:
Coeff. r^n Exp. Proj.
"""
return str_

14
scripts/utility/is_master_repository.py Executable file
View File

@ -0,0 +1,14 @@
#!/usr/bin/env python
import subprocess
pipe = subprocess.Popen("git config --get remote.origin.url", \
shell=True, stdout=subprocess.PIPE)
result = pipe.stdout.read()
is_master_repository = "LCPQ/quantum_package" in result
if __name__ == "__main__":
import sys
if is_master_repository:
sys.exit(0)
else:
sys.exit(-1)

1
src/AO_Basis/ao_overlap.irp.f
View File

@ -59,6 +59,7 @@
enddo
enddo
!$OMP END PARALLEL DO
END_PROVIDER

BIN
src/AO_Basis/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 14 KiB

After

Width:  |  Height:  |  Size: 14 KiB

8
src/Bitmask/bitmasks.irp.f
View File

@ -262,13 +262,7 @@ END_PROVIDER
logical :: exists
integer :: j,i
integer :: i_hole,i_part,i_gen
PROVIDE ezfio_filename
!do j = 1, N_int
! inact_bitmask(j,1) = xor(generators_bitmask(j,1,1,1),cas_bitmask(j,1,1))
! inact_bitmask(j,2) = xor(generators_bitmask(j,2,1,1),cas_bitmask(j,2,1))
! virt_bitmask(j,1) = xor(generators_bitmask(j,1,2,1),cas_bitmask(j,1,1))
! virt_bitmask(j,2) = xor(generators_bitmask(j,2,2,1),cas_bitmask(j,2,1))
!enddo
n_inact_orb = 0
n_virt_orb = 0
if(N_generators_bitmask == 1)then

BIN
src/Bitmask/tree_dependency.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 28 KiB

After

Width:  |  Height:  |  Size: 28 KiB

Some files were not shown because too many files have changed in this diff Show More