LCPQ/quantum_package
10
0
Fork 0
mirror of https://github.com/LCPQ/quantum_package synced 2025-04-22 16:30:27 +02:00
Code Issues Releases Wiki Activity

merged

Browse Source
This commit is contained in:
Yann Garniron 2016-01-25 16:13:43 +01:00
commit f5a32e698f
109 changed files with 750006 additions and 1024 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

11
.travis.yml
View File

@ -22,12 +22,9 @@ language: python
python: python:
- "2.6" - "2.6"
script: script:
- ./configure --production ./config/gfortran.cfg - ./configure --production ./config/gfortran.cfg
- source ./quantum_package.rc - 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 - source ./quantum_package.rc ; ninja
- ninja - source ./quantum_package.rc ; cd ocaml ; make ; cd -
- cd ocaml ; make ; cd - - source ./quantum_package.rc ; cd tests ; bats bats/qp.bats
- cd testing_no_regression ; ./unit_test.py

17
README.md
View File

@ -18,6 +18,9 @@ For more information, you can visit the [wiki of the project](http://github.com/
* Python >= 2.6 * Python >= 2.6
* GNU make * GNU make
* Bash * Bash
* Blast/Lapack
* unzip
* g++ (For ninja)
## Standard installation ## Standard installation
@ -51,19 +54,17 @@ This file contains all the environment variables needed by the quantum package b
### Optional) Add some new module ### Optional) Add some new module
Usage: qp_install_module.py list (--installed|--avalaible-local|--avalaible-remote) Usage: qp_module.py list (--installed|--avalaible-local|--avalaible-remote)
qp_install_module.py install <name>... qp_module.py install <name>...
qp_install_module.py create -n <name> [<children_module>...] qp_module.py create -n <name> [<children_module>...]
qp_install_module.py download -n <name> [<path_folder>...] qp_module.py download -n <name> [<path_folder>...]
For exemple you can type : For exemple you can type :
`qp_install_module.py install Full_CI` `qp_module.py install Full_CI`
### 3) Compiling the fortran ### 3) Compiling the fortran
ninja Just type `ninja` if you are in `$QP_ROOT` (or `ninja -f $QP_ROOT/build.ninja` elsewhere). The compilation will take approximately 3 min.
Just type `ninja` if you are in `$QP_ROOT` (or `ninja -f $QP_ROOT/build.ninja`
elsewhere). The compilation will take approximately 3 min.
If you have set the `--developement` flag in a specific module you can go in If you have set the `--developement` flag in a specific module you can go in
the corresponding module directory and run `ninja` to build only this module. the corresponding module directory and run `ninja` to build only this module.

2
config/gfortran.cfg
View File

@ -35,7 +35,7 @@ OPENMP : 1 ; Append OpenMP flags
# -ffast-math and the Fortran-specific # -ffast-math and the Fortran-specific
# -fno-protect-parens and -fstack-arrays. # -fno-protect-parens and -fstack-arrays.
[OPT] [OPT]
FCFLAGS : -Ofast FCFLAGS : -Ofast -march=native
# Profiling flags # Profiling flags
################# #################

5
config/ifort.cfg
View File

@ -31,14 +31,14 @@ OPENMP : 1 ; Append OpenMP flags
# -ftz : Flushes denormal results to zero # -ftz : Flushes denormal results to zero
# #
[OPT] [OPT]
FCFLAGS : -xSSE4.2 -O2 -ip -opt-prefetch -ftz -g FCFLAGS : -xHost -O2 -ip -ftz -g
# Profiling flags # Profiling flags
################# #################
# #
[PROFILE] [PROFILE]
FC : -p -g FC : -p -g
FCFLAGS : -xSSE4.2 -O2 -ip -opt-prefetch -ftz FCFLAGS : -xSSE4.2 -O2 -ip -ftz
# Debugging flags # Debugging flags
################# #################
@ -52,6 +52,7 @@ FCFLAGS : -xSSE4.2 -O2 -ip -opt-prefetch -ftz
[DEBUG] [DEBUG]
FC : -g -traceback FC : -g -traceback
FCFLAGS : -xSSE2 -C -fpe0 FCFLAGS : -xSSE2 -C -fpe0
IRPF90_FLAGS : --openmp
# OpenMP flags # OpenMP flags
################# #################

38
configure vendored
View File

@ -26,6 +26,8 @@ Examples:
""" """
OK="✓"
FAIL="✗"
import subprocess import subprocess
import os import os
import sys import sys
@ -52,7 +54,7 @@ QP_ROOT_INSTALL = join(QP_ROOT, "install")
os.environ["PATH"] = os.environ["PATH"] + ":" + QP_ROOT_BIN os.environ["PATH"] = os.environ["PATH"] + ":" + QP_ROOT_BIN
d_dependency = { d_dependency = {
"ocaml": ["m4", "curl", "zlib", "patch", "gcc"], "ocaml": ["m4", "curl", "zlib", "patch", "gcc", "zeromq"],
"m4": ["make"], "m4": ["make"],
"curl": ["make"], "curl": ["make"],
"zlib": ["gcc", "make"], "zlib": ["gcc", "make"],
@ -69,7 +71,8 @@ d_dependency = {
"python": [], "python": [],
"ninja": ["g++", "python"], "ninja": ["g++", "python"],
"make": [], "make": [],
"p_graphviz": ["python"] "p_graphviz": ["python"],
"bats": []
} }
from collections import namedtuple from collections import namedtuple
@ -134,25 +137,31 @@ ezfio = Info(
default_path=join(QP_ROOT_INSTALL, "EZFIO")) default_path=join(QP_ROOT_INSTALL, "EZFIO"))
zeromq = Info( zeromq = Info(
url='http://download.zeromq.org/zeromq-4.1.3.tar.gz', url='http://download.zeromq.org/zeromq-4.0.7.tar.gz',
description=' ZeroMQ', description=' ZeroMQ',
default_path=join(QP_ROOT_LIB, "libzmq.a")) default_path=join(QP_ROOT_LIB, "libzmq.a"))
f77zmq = Info( f77zmq = Info(
url='{head}/zeromq/f77_zmq/{tail}'.format(**path_github), url='{head}/zeromq/f77_zmq/{tail}'.format(**path_github),
description=' F77-ZeroMQ', description=' F77-ZeroMQ',
default_path=join(QP_ROOT_LIB, "libf77zmq.a")) default_path=join(QP_ROOT_LIB, "libf77zmq.a") + " " + \
join(QP_ROOT, "src", "ZMQ", "f77zmq.h") )
p_graphviz = Info( p_graphviz = Info(
url='https://github.com/xflr6/graphviz/archive/master.tar.gz', url='https://github.com/xflr6/graphviz/archive/master.tar.gz',
description=' Python library for graphviz', description=' Python library for graphviz',
default_path=join(QP_ROOT_INSTALL, "p_graphviz")) default_path=join(QP_ROOT_INSTALL, "p_graphviz"))
bats = Info(
url='https://github.com/sstephenson/bats/archive/master.tar.gz',
description=' Bash Automated Testing System',
default_path=join(QP_ROOT_INSTALL, "bats"))
d_info = dict() d_info = dict()
for m in ["ocaml", "m4", "curl", "zlib", "patch", "irpf90", "docopt", for m in ["ocaml", "m4", "curl", "zlib", "patch", "irpf90", "docopt",
"resultsFile", "ninja", "emsl", "ezfio", "p_graphviz", "resultsFile", "ninja", "emsl", "ezfio", "p_graphviz",
"zeromq", "f77zmq" ]: "zeromq", "f77zmq","bats" ]:
exec ("d_info['{0}']={0}".format(m)) exec ("d_info['{0}']={0}".format(m))
@ -281,10 +290,10 @@ def checking(d_dependency):
r = check_availability(i) r = check_availability(i)
if r: if r:
print "[ OK ] ( {0} )".format(r.strip()) print OK+" ( {0} )".format(r.strip())
l_installed[i] = r.strip() l_installed[i] = r.strip()
else: else:
print "[ FAIL ]" print FAIL
l_needed.append(i) l_needed.append(i)
print "" print ""
@ -366,7 +375,7 @@ _|_ | | _> |_ (_| | | (_| |_ | (_) | |
except: except:
raise raise
else: else:
print "[ OK ]" print OK
l_install_descendant.remove("ninja") l_install_descendant.remove("ninja")
@ -409,10 +418,10 @@ _|_ | | _> |_ (_| | | (_| |_ | (_) | |
with open(path, "w+") as f: with open(path, "w+") as f:
f.write("\n".join(l_string)) f.write("\n".join(l_string))
print " [ OK ] ({0})".format(path) print OK+" ({0})".format(path)
print str_info("install"), print str_info("install"),
print " [ Running ]" print "Running"
try: try:
path_ninja = find_path("ninja", l_installed) path_ninja = find_path("ninja", l_installed)
subprocess.check_call("cd install ;{0}".format(path_ninja), shell=True) subprocess.check_call("cd install ;{0}".format(path_ninja), shell=True)
@ -477,7 +486,7 @@ def create_ninja_and_rc(l_installed):
'export IRPF90={0}'.format(path_irpf90.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 NINJA={0}'.format(path_ninja.replace(QP_ROOT,"${QP_ROOT}")),
'export QP_PYTHON={0}'.format(":".join(l_python)), "", 'export QP_PYTHON={0}'.format(":".join(l_python)), "",
'export PYTHONPATH="${QP_EZFIO}":"${QP_PYTHON}":"${PYTHONPATH}"', 'export PYTHONPATH="${QP_EZFIO}/Python":"${QP_PYTHON}":"${PYTHONPATH}"',
'export PATH="${QP_PYTHON}":"${QP_ROOT}"/bin:"${QP_ROOT}"/ocaml:"${PATH}"', 'export PATH="${QP_PYTHON}":"${QP_ROOT}"/bin:"${QP_ROOT}"/ocaml:"${PATH}"',
'export LD_LIBRARY_PATH="${QP_ROOT}"/lib:"${LD_LIBRARY_PATH}"', 'export LD_LIBRARY_PATH="${QP_ROOT}"/lib:"${LD_LIBRARY_PATH}"',
'export LIBRARY_PATH="${QP_ROOT}"/lib:"${LIBRARY_PATH}"', "", 'export LIBRARY_PATH="${QP_ROOT}"/lib:"${LIBRARY_PATH}"', "",
@ -490,7 +499,7 @@ def create_ninja_and_rc(l_installed):
with open(path, "w+") as f: with open(path, "w+") as f:
f.write("\n".join(l_rc)) f.write("\n".join(l_rc))
print "[ OK ] ({0})".format(path) print OK+" ({0})".format(path)
command = ['bash', '-c', 'source {0} && env'.format(path)] command = ['bash', '-c', 'source {0} && env'.format(path)]
proc = subprocess.Popen(command, stdout=subprocess.PIPE) proc = subprocess.Popen(command, stdout=subprocess.PIPE)
@ -515,7 +524,7 @@ def create_ninja_and_rc(l_installed):
sys.exit(1) sys.exit(1)
else: else:
print "[ OK ]" print OK
def recommendation(): def recommendation():
@ -530,7 +539,8 @@ def recommendation():
print " ninja" print " ninja"
print " make -C ocaml" print " make -C ocaml"
print "" print ""
print "PS : For more info on compiling the code, read the COMPILE_RUN.md file." print "You can install more plugin with the qp_install command"
print "PS : For more info on compiling the code, read the README.md"
if __name__ == '__main__': if __name__ == '__main__':

28245
data/BFD-Orbital.dump Normal file
View File

File diff suppressed because it is too large Load Diff

895
data/BFD-Pseudo.dump Normal file
View File

@ -0,0 +1,895 @@
PRAGMA foreign_keys=OFF;
BEGIN TRANSACTION;
CREATE TABLE basis_tab(
basis_id INTEGER PRIMARY KEY AUTOINCREMENT,
name text,
description text,
UNIQUE(name)
);
INSERT INTO "basis_tab" VALUES(1,'BFD-Pseudo','http://ftp.aip.org/epaps/journ_chem_phys/E-JCPSA6-126-315722/epaps_material.html');
CREATE TABLE data_tab(
basis_id INTEGER,
elt TEXT,
data TEXT,
FOREIGN KEY(basis_id)
REFERENCES basis_tab(basis_id)
);
INSERT INTO "data_tab" VALUES(1,'H','Element Symbol: H
Number of replaced protons: 0
Number of projectors: 0
Pseudopotential data:
Local component:
Coeff. r^n Exp.
1.00000000 -1 4.47692410
4.47692410 1 2.97636451
-4.32112340 0 3.01841596
Non-local component:
Coeff. r^n Exp. Proj.');
INSERT INTO "data_tab" VALUES(1,'He','Element Symbol: He
Number of replaced protons: 0
Number of projectors: 0
Pseudopotential data:
Local component:
Coeff. r^n Exp.
2.00000000 -1 9.84368811
19.68737622 1 10.94516233
-17.20570338 0 12.03715264
Non-local component:
Coeff. r^n Exp. Proj.');
INSERT INTO "data_tab" VALUES(1,'Li','Element Symbol: Li
Number of replaced protons: 2
Number of projectors: 1
Pseudopotential data:
Local component:
Coeff. r^n Exp.
1.00000000 -1 5.41040609
5.41040609 1 2.70520138
-4.60151975 0 2.07005488
Non-local component:
Coeff. r^n Exp. Proj.
7.09172172 0 1.34319829 |0><0|');
INSERT INTO "data_tab" VALUES(1,'Be','Element Symbol: Be
Number of replaced protons: 2
Number of projectors: 1
Pseudopotential data:
Local component:
Coeff. r^n Exp.
2.00000000 -1 4.58686001
9.17372003 1 2.29371778
-8.12599146 0 2.10401964
Non-local component:
Coeff. r^n Exp. Proj.
14.90499810 0 2.71723988 |0><0|');
INSERT INTO "data_tab" VALUES(1,'B','Element Symbol: B
Number of replaced protons: 2
Number of projectors: 1
Pseudopotential data:
Local component:
Coeff. r^n Exp.
3.00000000 -1 5.40423964
16.21271892 1 5.71678458
-11.86640633 0 4.48974455
Non-local component:
Coeff. r^n Exp. Proj.
15.49737620 0 3.43781634 |0><0|');
INSERT INTO "data_tab" VALUES(1,'C','Element Symbol: C
Number of replaced protons: 2
Number of projectors: 1
Pseudopotential data:
Local component:
Coeff. r^n Exp.
4.00000000 -1 8.35973821
33.43895285 1 4.48361888
-19.17537323 0 3.93831258
Non-local component:
Coeff. r^n Exp. Proj.
22.55164191 0 5.02991637 |0><0|');
INSERT INTO "data_tab" VALUES(1,'N','Element Symbol: N
Number of replaced protons: 2
Number of projectors: 1
Pseudopotential data:
Local component:
Coeff. r^n Exp.
5.00000000 -1 9.23501007
46.17505034 1 7.66830008
-30.18893534 0 7.34486070
Non-local component:
Coeff. r^n Exp. Proj.
31.69720409 0 6.99536540 |0><0|');
INSERT INTO "data_tab" VALUES(1,'O','Element Symbol: O
Number of replaced protons: 2
Number of projectors: 1
Pseudopotential data:
Local component:
Coeff. r^n Exp.
6.00000000 -1 9.29793903
55.78763416 1 8.86492204
-38.81978498 0 8.62925665
Non-local component:
Coeff. r^n Exp. Proj.
38.41914135 0 8.71924452 |0><0|');
INSERT INTO "data_tab" VALUES(1,'F','Element Symbol: F
Number of replaced protons: 2
Number of projectors: 1
Pseudopotential data:
Local component:
Coeff. r^n Exp.
7.00000000 -1 11.39210685
79.74474797 1 10.74911370
-49.45159098 0 10.45120693
Non-local component:
Coeff. r^n Exp. Proj.
50.25646328 0 11.30345826 |0><0|');
INSERT INTO "data_tab" VALUES(1,'Ne','Element Symbol: Ne
Number of replaced protons: 2
Number of projectors: 1
Pseudopotential data:
Local component:
Coeff. r^n Exp.
8.00000000 -1 10.74945199
85.99561593 1 10.19801460
-56.79004456 0 10.18694048
Non-local component:
Coeff. r^n Exp. Proj.
55.11144535 0 12.85042963 |0><0|');
INSERT INTO "data_tab" VALUES(1,'Na','Element Symbol: Na
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
1.00000000 -1 5.35838717
5.35838717 1 3.67918975
-2.07764789 0 1.60507673
Non-local component:
Coeff. r^n Exp. Proj.
10.69640234 0 1.32389367 |0><0|
10.11238853 0 1.14052020 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Mg','Element Symbol: Mg
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
2.00000000 -1 4.48537898
8.97075796 1 2.24268949
-7.72153408 0 1.59710474
Non-local component:
Coeff. r^n Exp. Proj.
15.07848048 0 1.57188656 |0><0|
12.37888383 0 1.42757357 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Al','Element Symbol: Al
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
3.00000000 -1 3.07301275
9.21903825 1 9.97055633
-9.65888637 0 2.64134660
Non-local component:
Coeff. r^n Exp. Proj.
17.16680112 0 1.87284747 |0><0|
14.22120694 0 1.79397208 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Si','Element Symbol: Si
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
4.00000000 -1 1.80721061
7.22884246 1 9.99633089
-13.06725590 0 2.50043232
Non-local component:
Coeff. r^n Exp. Proj.
21.20531613 0 2.26686403 |0><0|
15.43693603 0 2.11659661 |1><1|');
INSERT INTO "data_tab" VALUES(1,'P','Element Symbol: P
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
5.00000000 -1 2.02622810
10.13114051 1 9.95970113
-14.94375088 0 2.74841795
Non-local component:
Coeff. r^n Exp. Proj.
23.62479480 0 2.60470698 |0><0|
18.18547203 0 2.54957900 |1><1|');
INSERT INTO "data_tab" VALUES(1,'S','Element Symbol: S
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
6.00000000 -1 2.42178462
14.53070769 1 6.74148698
-17.52965289 0 3.06094751
Non-local component:
Coeff. r^n Exp. Proj.
25.99260928 0 2.94272173 |0><0|
18.93356489 0 2.84566981 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Cl','Element Symbol: Cl
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
7.00000000 -1 2.41079533
16.87556731 1 5.29139158
-18.80917558 0 2.91105513
Non-local component:
Coeff. r^n Exp. Proj.
28.59107316 0 3.34528827 |0><0|
19.37583724 0 3.12408551 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Ar','Element Symbol: Ar
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
8.00000000 -1 3.09403094
24.75224749 1 6.53700323
-20.38446872 0 3.35769859
Non-local component:
Coeff. r^n Exp. Proj.
30.67006675 0 3.68203169 |0><0|
20.84338017 0 3.45735664 |1><1|');
INSERT INTO "data_tab" VALUES(1,'K','Element Symbol: K
Number of replaced protons: 18
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
1.00000000 -1 5.46842198
5.46842198 1 3.45438113
-7.43169289 0 0.86173842
Non-local component:
Coeff. r^n Exp. Proj.
34.84862909 0 0.98715169 |0><0|
33.69024309 0 0.75016011 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Ca','Element Symbol: Ca
Number of replaced protons: 18
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
2.00000000 -1 5.46886760
10.93773521 1 3.45174786
-7.42323874 0 0.93761419
Non-local component:
Coeff. r^n Exp. Proj.
34.65752096 0 1.08463527 |0><0|
34.12044224 0 0.91231229 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Sc','Element Symbol: Sc
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
11.00000000 -1 1.64916555
18.14082106 1 3.06833288
-35.19310566 0 2.41985389
Non-local component:
Coeff. r^n Exp. Proj.
97.62913482 0 7.60319353 |0><0|
33.97033635 0 5.31121835 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Ti','Element Symbol: Ti
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
12.00000000 -1 1.85755224
22.29062683 1 3.30638246
-41.26280518 0 2.70879079
Non-local component:
Coeff. r^n Exp. Proj.
96.94231089 0 8.03040157 |0><0|
38.01641313 0 5.93291405 |1><1|');
INSERT INTO "data_tab" VALUES(1,'V','Element Symbol: V
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
13.00000000 -1 2.16361754
28.12702797 1 4.07901780
-48.27656329 0 3.21436396
Non-local component:
Coeff. r^n Exp. Proj.
96.23226580 0 8.44326050 |0><0|
41.58043539 0 6.53136059 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Cr','Element Symbol: Cr
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
14.00000000 -1 2.94337662
41.20727267 1 3.40750349
-55.51413840 0 3.33587110
Non-local component:
Coeff. r^n Exp. Proj.
103.26274052 0 9.14231779 |0><0|
45.80124572 0 7.21220771 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Mn','Element Symbol: Mn
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
15.00000000 -1 3.29524605
49.42869068 1 3.61599152
-61.66925169 0 3.57405761
Non-local component:
Coeff. r^n Exp. Proj.
112.85037953 0 9.99154195 |0><0|
49.18832867 0 7.80925218 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Fe','Element Symbol: Fe
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
16.00000000 -1 3.72075632
59.53210107 1 3.92321272
-68.75847841 0 3.89595440
Non-local component:
Coeff. r^n Exp. Proj.
112.92561163 0 10.42343546 |0><0|
52.55882759 0 8.41664076 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Co','Element Symbol: Co
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
17.00000000 -1 4.18819322
71.19928469 1 4.42968808
-77.65278252 0 4.39800669
Non-local component:
Coeff. r^n Exp. Proj.
113.90484511 0 10.86075441 |0><0|
56.10698766 0 9.05202771 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Ni','Element Symbol: Ni
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
18.00000000 -1 4.22878924
76.11820629 1 4.46202418
-82.85330412 0 4.44960456
Non-local component:
Coeff. r^n Exp. Proj.
120.84003628 0 11.62700064 |0><0|
58.54148726 0 9.57327085 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Cu','Element Symbol: Cu
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
19.00000000 -1 6.25149628
118.77842940 1 6.72725326
-105.89982403 0 6.61024592
Non-local component:
Coeff. r^n Exp. Proj.
127.35069424 0 12.36568715 |0><0|
71.22984900 0 11.16072939 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Zn','Element Symbol: Zn
Number of replaced protons: 10
Number of projectors: 2
Pseudopotential data:
Local component:
Coeff. r^n Exp.
20.00000000 -1 5.25282726
105.05654526 1 5.85433298
-105.08806248 0 5.80452897
Non-local component:
Coeff. r^n Exp. Proj.
123.87006927 0 12.52174964 |0><0|
72.33499364 0 11.56019052 |1><1|');
INSERT INTO "data_tab" VALUES(1,'Ga','Element Symbol: Ga
Number of replaced protons: 28
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
3.00000000 -1 1.22544253
3.67632759 1 5.71065133
-11.23828733 0 1.33931968
Non-local component:
Coeff. r^n Exp. Proj.
57.88512249 0 2.48772664 |0><0|
43.67871044 0 2.12489462 |1><1|
17.97137628 0 1.27124173 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Ge','Element Symbol: Ge
Number of replaced protons: 28
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
4.00000000 -1 1.88492329
7.53969315 1 9.98137268
-13.13622589 0 1.98008364
Non-local component:
Coeff. r^n Exp. Proj.
61.26369269 0 3.03315885 |0><0|
55.52495744 0 2.76564031 |1><1|
23.49168485 0 1.66026543 |2><2|');
INSERT INTO "data_tab" VALUES(1,'As','Element Symbol: As
Number of replaced protons: 28
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
5.00000000 -1 1.21796059
6.08980295 1 9.96302171
-14.92625816 0 1.86158567
Non-local component:
Coeff. r^n Exp. Proj.
73.75553709 0 3.54546456 |0><0|
68.03580909 0 3.28664249 |1><1|
23.32540737 0 1.95862616 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Se','Element Symbol: Se
Number of replaced protons: 28
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
6.00000000 -1 1.73494732
10.40968393 1 6.91632737
-17.83199463 0 3.10551681
Non-local component:
Coeff. r^n Exp. Proj.
85.94238004 0 4.67503354 |0><0|
78.84838432 0 4.34256579 |1><1|
30.92151589 0 2.61905005 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Br','Element Symbol: Br
Number of replaced protons: 28
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
7.00000000 -1 1.86793881
13.07557164 1 5.30536536
-18.79056037 0 3.32134623
Non-local component:
Coeff. r^n Exp. Proj.
88.58537968 0 5.17694821 |0><0|
79.43718432 0 4.80714881 |1><1|
29.35463757 0 3.03534088 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Kr','Element Symbol: Kr
Number of replaced protons: 28
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
8.00000000 -1 1.72397711
13.79181690 1 6.70510242
-22.77215308 0 2.75463922
Non-local component:
Coeff. r^n Exp. Proj.
92.78570269 0 4.85045356 |0><0|
80.37767796 0 4.52350391 |1><1|
31.36172413 0 3.04556109 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Rb','Element Symbol: Rb
Number of replaced protons: 36
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
1.00000000 -1 4.38377252
4.38377252 1 3.19763054
-14.25890642 0 0.92493011
Non-local component:
Coeff. r^n Exp. Proj.
120.87192624 0 1.14964523 |0><0|
110.02656159 0 0.75675227 |1><1|
22.44443305 0 1.97938422 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Sr','Element Symbol: Sr
Number of replaced protons: 36
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
2.00000000 -1 4.38393876
8.76787753 1 3.19649834
-14.25524313 0 0.98991497
Non-local component:
Coeff. r^n Exp. Proj.
146.87160520 0 1.28683494 |0><0|
140.02622842 0 0.95674483 |1><1|
24.44548775 0 1.94831654 |2><2|');
INSERT INTO "data_tab" VALUES(1,'In','Element Symbol: In
Number of replaced protons: 46
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
3.00000000 -1 0.84508372
2.53525117 1 5.66019931
-7.66579852 0 0.85764327
Non-local component:
Coeff. r^n Exp. Proj.
58.16918845 0 1.87837596 |0><0|
43.63891951 0 1.51547534 |1><1|
17.93363847 0 0.83917399 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Sn','Element Symbol: Sn
Number of replaced protons: 46
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
4.00000000 -1 0.61334103
2.45336413 1 5.67826592
-9.33070594 0 0.86945138
Non-local component:
Coeff. r^n Exp. Proj.
58.04190484 0 2.01380769 |0><0|
43.68447157 0 1.63883815 |1><1|
17.95660523 0 0.92346533 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Sb','Element Symbol: Sb
Number of replaced protons: 46
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
5.00000000 -1 0.74773404
3.73867018 1 5.79307847
-13.88202267 0 1.07909250
Non-local component:
Coeff. r^n Exp. Proj.
57.56076138 0 2.04356327 |0><0|
43.88817098 0 1.70062095 |1><1|
17.82275877 0 1.00414410 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Te','Element Symbol: Te
Number of replaced protons: 46
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
6.00000000 -1 1.28872858
7.73237146 1 5.81046154
-13.65642757 0 1.55436985
Non-local component:
Coeff. r^n Exp. Proj.
57.52907325 0 2.39157624 |0><0|
43.70165092 0 1.96781367 |1><1|
18.64325026 0 1.13849722 |2><2|');
INSERT INTO "data_tab" VALUES(1,'I','Element Symbol: I
Number of replaced protons: 46
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
7.00000000 -1 1.01780923
7.12466460 1 3.60136147
-29.18419372 0 1.68345378
Non-local component:
Coeff. r^n Exp. Proj.
108.68417388 0 2.80278521 |0><0|
99.35380694 0 2.51494051 |1><1|
41.32653157 0 1.56672279 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Xe','Element Symbol: Xe
Number of replaced protons: 46
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
8.00000000 -1 1.22029027
9.76232213 1 5.14625292
-32.16318995 0 2.10563577
Non-local component:
Coeff. r^n Exp. Proj.
110.93633943 0 3.20320603 |0><0|
99.49007680 0 2.86062806 |1><1|
41.81892440 0 1.74523568 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Cs','Element Symbol: Cs
Number of replaced protons: 54
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
1.00000000 -1 4.56591654
4.56591654 1 3.44071226
-10.05901330 0 0.53597632
Non-local component:
Coeff. r^n Exp. Proj.
80.60678398 0 0.75497262 |0><0|
38.41951903 0 0.45697820 |1><1|
34.36421109 0 1.21638889 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Ba','Element Symbol: Ba
Number of replaced protons: 54
Number of projectors: 3
Pseudopotential data:
Local component:
Coeff. r^n Exp.
2.00000000 -1 4.56523866
9.13047731 1 3.44590138
-11.09995006 0 0.73875795
Non-local component:
Coeff. r^n Exp. Proj.
84.60641087 0 0.87685449 |0><0|
42.41929187 0 0.55569714 |1><1|
35.35527943 0 1.68654278 |2><2|');
INSERT INTO "data_tab" VALUES(1,'Tl','Element Symbol: Tl
Number of replaced protons: 78
Number of projectors: 4
Pseudopotential data:
Local component:
Coeff. r^n Exp.
3.00000000 -1 1.34540299
4.03620897 1 1.80622564
-21.68751419 0 0.88272932
Non-local component:
Coeff. r^n Exp. Proj.
91.69120822 0 1.55008480 |0><0|
85.52036047 0 1.34656124 |1><1|
62.00124014 0 0.92193030 |2><2|
45.59483689 0 0.95782546 |3><3|');
INSERT INTO "data_tab" VALUES(1,'Pb','Element Symbol: Pb
Number of replaced protons: 78
Number of projectors: 4
Pseudopotential data:
Local component:
Coeff. r^n Exp.
4.00000000 -1 1.33905502
5.35622008 1 3.57680327
-25.11165802 0 1.08584447
Non-local component:
Coeff. r^n Exp. Proj.
121.69447681 0 1.89957262 |0><0|
114.36466627 0 1.64009233 |1><1|
49.32959048 0 0.93051806 |2><2|
45.59434323 0 1.07638351 |3><3|');
INSERT INTO "data_tab" VALUES(1,'Bi','Element Symbol: Bi
Number of replaced protons: 78
Number of projectors: 4
Pseudopotential data:
Local component:
Coeff. r^n Exp.
5.00000000 -1 1.32337887
6.61689433 1 3.57500194
-25.11872438 0 1.26988106
Non-local component:
Coeff. r^n Exp. Proj.
121.69346942 0 2.17897203 |0><0|
114.36404033 0 1.85248885 |1><1|
59.32816833 0 1.10815262 |2><2|
55.59430794 0 1.13749297 |3><3|');
INSERT INTO "data_tab" VALUES(1,'Po','Element Symbol: Po
Number of replaced protons: 78
Number of projectors: 4
Pseudopotential data:
Local component:
Coeff. r^n Exp.
6.00000000 -1 1.12915012
6.77490072 1 3.57645639
-27.15687722 0 1.43566151
Non-local component:
Coeff. r^n Exp. Proj.
123.68767392 0 2.45166968 |0><0|
116.36047985 0 2.06261668 |1><1|
81.32822758 0 1.28212112 |2><2|
80.59441478 0 1.32883356 |3><3|');
INSERT INTO "data_tab" VALUES(1,'At','Element Symbol: At
Number of replaced protons: 78
Number of projectors: 4
Pseudopotential data:
Local component:
Coeff. r^n Exp.
7.00000000 -1 1.03089969
7.21629781 1 3.69262441
-29.19452905 0 1.45693225
Non-local component:
Coeff. r^n Exp. Proj.
128.68209352 0 2.56609159 |0><0|
119.35613288 0 2.15280178 |1><1|
80.33152157 0 1.37630746 |2><2|
80.58962589 0 1.47835294 |3><3|');
INSERT INTO "data_tab" VALUES(1,'Rn','Element Symbol: Rn
Number of replaced protons: 78
Number of projectors: 4
Pseudopotential data:
Local component:
Coeff. r^n Exp.
8.00000000 -1 1.09138091
8.73104728 1 3.03187023
-29.28994938 0 1.65026568
Non-local component:
Coeff. r^n Exp. Proj.
128.63606138 0 2.92260183 |0><0|
119.39674828 0 2.42144059 |1><1|
80.33096401 0 1.49781359 |2><2|
80.58961959 0 1.47976575 |3><3|');
INSERT INTO "data_tab" VALUES(1,'Rn','Element Symbol: Rn
Number of replaced protons: 78
Number of projectors: 4
Pseudopotential data:
Local component:
Coeff. r^n Exp.
8.00000000 -1 1.09138091
8.73104728 1 3.03187023
-29.28994938 0 1.65026568
Non-local component:
Coeff. r^n Exp. Proj.
128.63606138 0 2.92260183 |0><0|
119.39674828 0 2.42144059 |1><1|
80.33096401 0 1.49781359 |2><2|
80.58961959 0 1.47976575 |3><3|');
CREATE TABLE format_tab(format TEXT);
INSERT INTO "format_tab" VALUES('BFD-Pseudo');
DELETE FROM sqlite_sequence;
INSERT INTO "sqlite_sequence" VALUES('basis_tab',5);
CREATE VIEW output_tab AS
SELECT basis_id,
name,
description,
elt,
data
FROM basis_tab
NATURAL JOIN data_tab;
COMMIT;

717154
data/GAMESS-US.dump Normal file
View File

File diff suppressed because it is too large Load Diff

0
install/Downloads/.gitignore → install/Downloads/.empty
View File

13
install/scripts/install_bats.sh Executable file
View File

@ -0,0 +1,13 @@
#!/bin/bash -x
TARGET=bats
function _install()
{
cp -R ${BUILD} . || exit 1
cd ../bin
ln -s ../install/${TARGET}/libexec/bats . || return 1
cd -
}
source scripts/build.sh

51
install/scripts/install_ocaml.sh
View File

@ -5,16 +5,56 @@ QP_ROOT=$PWD
cd - cd -
# Normal installation # Normal installation
PACKAGES="core cryptokit ocamlfind sexplib" PACKAGES="core cryptokit ocamlfind sexplib ZMQ"
declare -i i # Needed for ZeroMQ
i=$(gcc -dumpversion | cut -d '.' -f 2) export C_INCLUDE_PATH="${QP_ROOT}"/lib:"${C_INCLUDE_PATH}"
if [[ i -lt 6 ]] export LIBRARY_PATH="${QP_ROOT}"/lib:"${LIBRARY_PATH}"
export LD_LIBRARY_PATH="${QP_ROOT}"/lib:"${LD_LIBRARY_PATH}"
# return 0 if program version is equal or greater than check version
check_version () {
if [[ $1 == $2 ]]
then
return 0
fi
local IFS=.
local i ver1=($1) ver2=($2)
# fill empty fields in ver1 with zeros
for ((i=${#ver1[@]}; i<${#ver2[@]}; i++))
do
ver1[i]=0
done
for ((i=0; i<${#ver1[@]}; i++))
do
if [[ -z ${ver2[i]} ]]
then
# fill empty fields in ver2 with zeros
ver2[i]=0
fi
if ((10#${ver1[i]} > 10#${ver2[i]}))
then
return 1
fi
if ((10#${ver1[i]} < 10#${ver2[i]}))
then
return 2
fi
done
return 0
}
i=$(gcc -dumpversion)
check_version 4.6 $i
if [[ $? == 1 ]]
then then
echo "GCC version $(gcc -dumpversion) too old. GCC >= 4.6 required." echo "GCC version $(gcc -dumpversion) too old. GCC >= 4.6 required."
exit 1 exit 1
fi fi
if [[ -d ${HOME}/.opam ]] if [[ -d ${HOME}/.opam ]]
then then
source ${HOME}/.opam/opam-init/init.sh > /dev/null 2> /dev/null || true source ${HOME}/.opam/opam-init/init.sh > /dev/null 2> /dev/null || true
@ -37,6 +77,3 @@ NCPUs=$(cat /proc/cpuinfo | grep -i MHz | wc -l)
${QP_ROOT}/bin/opam install -j ${NCPUs} ${PACKAGES} -y -q || exit 1 ${QP_ROOT}/bin/opam install -j ${NCPUs} ${PACKAGES} -y -q || exit 1
rm -f ../_build/ocaml.log rm -f ../_build/ocaml.log
exit 0

11
install/scripts/install_zeromq.sh
View File

@ -14,12 +14,15 @@ function _install()
cd "${BUILD}" cd "${BUILD}"
./configure --without-libsodium || exit 1 ./configure --without-libsodium || exit 1
make -j 8 || 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 rm -f -- "${QP_ROOT}"/lib/libzmq.a "${QP_ROOT}"/lib/libzmq.so "${QP_ROOT}"/lib/libzmq.so.?
cp .libs/libzmq.a "${QP_ROOT}"/lib # cp .libs/libzmq.a "${QP_ROOT}"/lib
cp .libs/libzmq.so "${QP_ROOT}"/lib/libzmq.so.5 # cp .libs/libzmq.so "${QP_ROOT}"/lib/libzmq.so.5
cp src/.libs/libzmq.a "${QP_ROOT}"/lib
cp src/.libs/libzmq.so "${QP_ROOT}"/lib/libzmq.so.4
cp include/{zmq.h,zmq_utils.h} "${QP_ROOT}"/lib cp include/{zmq.h,zmq_utils.h} "${QP_ROOT}"/lib
cd "${QP_ROOT}"/lib cd "${QP_ROOT}"/lib
ln -s libzmq.so.5 libzmq.so # ln -s libzmq.so.5 libzmq.so
ln -s libzmq.so.4 libzmq.so
cd ${ORIG} cd ${ORIG}
return 0 return 0
} }

9
ocaml/.gitignore vendored
View File

@ -1,6 +1,7 @@
_build _build
ezfio.ml ezfio.ml
.gitignore .gitignore
Git.ml
Input_auto_generated.ml Input_auto_generated.ml
Input_determinants.ml Input_determinants.ml
Input_hartree_fock.ml Input_hartree_fock.ml
@ -16,6 +17,8 @@ qp_edit
qp_edit.ml qp_edit.ml
qp_edit.native qp_edit.native
qp_print qp_print
qp_print_basis
qp_print_basis.native
qp_print.native qp_print.native
qp_run qp_run
qp_run.native qp_run.native
@ -39,9 +42,15 @@ test_excitation
test_excitation.byte test_excitation.byte
test_gto test_gto
test_gto.byte test_gto.byte
test_message
test_message.byte
test_mo_label test_mo_label
test_mo_label.byte test_mo_label.byte
test_molecule test_molecule
test_molecule.byte test_molecule.byte
test_point3d test_point3d
test_point3d.byte test_point3d.byte
test_queuing_system
test_queuing_system.byte
test_task_server
test_task_server.byte

48
ocaml/Address.ml Normal file
View File

@ -0,0 +1,48 @@
open Core.Std
module Tcp : sig
type t
val of_string : string -> t
val to_string : t -> string
end = struct
type t = string
let of_string x =
assert (String.is_prefix ~prefix:"tcp://" x);
x
let to_string x = x
end
module Ipc : sig
type t
val of_string : string -> t
val to_string : t -> string
end = struct
type t = string
let of_string x =
assert (String.is_prefix ~prefix:"ipc://" x);
x
let to_string x = x
end
module Inproc : sig
type t
val of_string : string -> t
val to_string : t -> string
end = struct
type t = string
let of_string x =
assert (String.is_prefix ~prefix:"inproc://" x);
x
let to_string x = x
end
type t =
| Tcp of Tcp.t
| Ipc of Ipc.t
| Inproc of Inproc.t
let to_string = function
| Tcp x -> Tcp.to_string x
| Ipc x -> Ipc.to_string x
| Inproc x -> Inproc.to_string x

20
ocaml/Basis.ml
View File

@ -1,7 +1,7 @@
open Core.Std;; open Core.Std
open Qptypes;; open Qptypes
type t = (Gto.t * Nucl_number.t) list with sexp;; type t = (Gto.t * Nucl_number.t) list with sexp
(** Read all the basis functions of an element *) (** Read all the basis functions of an element *)
let read in_channel at_number = let read in_channel at_number =
@ -12,7 +12,7 @@ let read in_channel at_number =
with with
| Gto.End_Of_Basis -> List.rev result | Gto.End_Of_Basis -> List.rev result
in read [] in read []
;;
(** Find an element in the basis set file *) (** Find an element in the basis set file *)
let find in_channel element = let find in_channel element =
@ -27,13 +27,13 @@ let find in_channel element =
| Element.ElementError _ -> () | Element.ElementError _ -> ()
done ; done ;
!element_read !element_read
;;
(** Read an element from the file *) (** Read an element from the file *)
let read_element in_channel at_number element = let read_element in_channel at_number element =
ignore (find in_channel element) ; ignore (find in_channel element) ;
read in_channel at_number ; read in_channel at_number
;;
let to_string b = let to_string b =
let new_nucleus n = let new_nucleus n =
@ -55,9 +55,9 @@ let to_string b =
in in
do_work [new_nucleus 1] 1 b do_work [new_nucleus 1] 1 b
|> String.concat ~sep:"\n" |> String.concat ~sep:"\n"
;;
include To_md5;;
include To_md5
let to_md5 = to_md5 sexp_of_t let to_md5 = to_md5 sexp_of_t
;;

33
ocaml/Id.ml Normal file
View File

@ -0,0 +1,33 @@
open Core.Std
module Id : sig
type t
val of_int : int -> t
val to_int : t -> int
val of_string : string -> t
val to_string : t -> string
val increment : t -> t
val decrement : t -> t
end
= struct
type t = int
let of_int x =
assert (x>0); x
let to_int x = x
let of_string x =
Int.of_string x
|> of_int
let to_string x =
Int.to_string x
let increment x = x + 1
let decrement x = x - 1
end
module Task = struct
include Id
end
module Client = struct
include Id
end

8
ocaml/Makefile
View File

@ -24,7 +24,7 @@ ALL_EXE=$(patsubst %.ml,%.native,$(wildcard qp_*.ml)) qp_edit.native
default: $(ALL_TESTS) $(ALL_EXE) .gitignore default: $(ALL_TESTS) $(ALL_EXE) .gitignore
.gitignore: $(MLFILES) .gitignore: $(MLFILES)
@for i in .gitignore ezfio.ml Qptypes.ml qptypes_generator.byte _build $(ALL_EXE) $(ALL_TESTS) \ @for i in .gitignore ezfio.ml Qptypes.ml Git.ml qptypes_generator.byte _build $(ALL_EXE) $(ALL_TESTS) \
$(patsubst %.ml,%,$(wildcard test_*.ml)) $(patsubst %.ml,%,$(wildcard qp_*.ml)) \ $(patsubst %.ml,%,$(wildcard test_*.ml)) $(patsubst %.ml,%,$(wildcard qp_*.ml)) \
$(shell grep Input Input_auto_generated.ml | awk '{print $$2 ".ml"}') \ $(shell grep Input Input_auto_generated.ml | awk '{print $$2 ".ml"}') \
qp_edit.ml qp_edit qp_edit.native Input_auto_generated.ml;\ qp_edit.ml qp_edit qp_edit.native Input_auto_generated.ml;\
@ -34,6 +34,7 @@ default: $(ALL_TESTS) $(ALL_EXE) .gitignore
executables: $(QP_ROOT)/data/executables executables: $(QP_ROOT)/data/executables
$(QP_ROOT)/data/executables: remake_executables $(QP_ROOT)/data/executables: remake_executables
$(QP_ROOT)/scripts/module/create_executables_list.sh $(QP_ROOT)/scripts/module/create_executables_list.sh
@ -68,8 +69,9 @@ ezfio.ml: ${QP_ROOT}/install/EZFIO/Ocaml/ezfio.ml
qptypes_generator.byte: qptypes_generator.ml qptypes_generator.byte: qptypes_generator.ml
$(OCAMLBUILD) qptypes_generator.byte -use-ocamlfind $(OCAMLBUILD) qptypes_generator.byte -use-ocamlfind
Qptypes.ml: qptypes_generator.byte Git.ml Qptypes.ml: qptypes_generator.byte
./qptypes_generator.byte > Qptypes.ml ./qptypes_generator.byte > Qptypes.ml
./create_git_sha1.sh
${QP_ROOT}/install/EZFIO/Ocaml/ezfio.ml: ${QP_ROOT}/install/EZFIO/Ocaml/ezfio.ml:
$(NINJA) -C ${QP_ROOT}/install/EZFIO $(NINJA) -C ${QP_ROOT}/install/EZFIO
@ -78,5 +80,5 @@ Input_auto_generated.ml qp_edit.ml:
ei_handler.py ocaml_global ei_handler.py ocaml_global
clean: clean:
rm -rf _build Qptypes.ml Input_auto_generated.ml $(ALL_EXE) $(ALL_TESTS) rm -rf _build Qptypes.ml Git.ml Input_auto_generated.ml $(ALL_EXE) $(ALL_TESTS)

323
ocaml/Message.ml Normal file
View File

@ -0,0 +1,323 @@
open Core.Std
(** New job : Request to create a new multi-tasked job *)
module State : sig
type t
val of_string : string -> t
val to_string : t -> string
end = struct
type t = string
let of_string x = x
let to_string x = x
end
module Newjob_msg : sig
type t =
{ state: State.t;
address_tcp: Address.Tcp.t ;
address_inproc: Address.Inproc.t;
}
val create : address_tcp:string -> address_inproc:string -> state:string -> t
val to_string : t -> string
end = struct
type t =
{ state: State.t;
address_tcp: Address.Tcp.t ;
address_inproc: Address.Inproc.t;
}
let create ~address_tcp ~address_inproc ~state =
{ state = State.of_string state;
address_tcp = Address.Tcp.of_string address_tcp ;
address_inproc = Address.Inproc.of_string address_inproc ;
}
let to_string t =
Printf.sprintf "newjob %s %s %s"
( State.to_string t.state )
( Address.Tcp.to_string t.address_tcp )
( Address.Inproc.to_string t.address_inproc )
end
(** Connect : connect a new client to the task server *)
module Connect_msg : sig
type t = Tcp | Inproc | Ipc
val create : typ:string -> t
val to_string : t -> string
end = struct
type t = Tcp | Inproc | Ipc
let create ~typ =
match typ with
| "tcp" -> Tcp
| "inproc" -> Inproc
| "ipc" -> Ipc
| _ -> assert false
let to_string = function
| Tcp -> "connect tcp"
| Inproc -> "connect inproc"
| Ipc -> "connect ipc"
end
(** ConnectReply : Reply to the connect messsage *)
module ConnectReply_msg : sig
type t =
{ client_id: Id.Client.t ;
state: State.t ;
push_address: Address.t;
}
val create : state:State.t -> client_id:Id.Client.t -> push_address:Address.t -> t
val to_string : t -> string
end = struct
type t =
{ client_id: Id.Client.t ;
state: State.t ;
push_address: Address.t;
}
let create ~state ~client_id ~push_address =
{ client_id ; state ; push_address }
let to_string x =
Printf.sprintf "connect_reply %s %d %s"
(State.to_string x.state)
(Id.Client.to_int x.client_id)
(Address.to_string x.push_address)
end
(** Disconnect : disconnect a client from the task server *)
module Disconnect_msg : sig
type t =
{ client_id: Id.Client.t ;
state: State.t ;
}
val create : state:string -> client_id:string -> t
val to_string : t -> string
end = struct
type t =
{ client_id: Id.Client.t ;
state: State.t ;
}
let create ~state ~client_id =
{ client_id = Id.Client.of_string client_id ; state = State.of_string state }
let to_string x =
Printf.sprintf "disconnect %s %d"
(State.to_string x.state)
(Id.Client.to_int x.client_id)
end
module DisconnectReply_msg : sig
type t =
{ finished: bool ;
state: State.t ;
}
val create : state:State.t -> finished:bool -> t
val to_string : t -> string
end = struct
type t =
{ finished: bool;
state: State.t ;
}
let create ~state ~finished =
{ state ; finished }
let to_string x =
Printf.sprintf "disconnect_reply %s %d"
(State.to_string x.state)
(if x.finished then 1 else 0)
end
(** AddTask : Add a new task to the queue *)
module AddTask_msg : sig
type t =
{ state: State.t;
task: string;
}
val create : state:string -> task:string -> t
val to_string : t -> string
end = struct
type t =
{ state: State.t;
task: string;
}
let create ~state ~task = { state = State.of_string state ; task }
let to_string x =
Printf.sprintf "add_task %s %s" (State.to_string x.state) x.task
end
(** AddTaskReply : Reply to the AddTask message *)
module AddTaskReply_msg : sig
type t
val create : task_id:Id.Task.t -> t
val to_string : t -> string
end = struct
type t = Id.Task.t
let create ~task_id = task_id
let to_string x =
Printf.sprintf "add_task_reply %d" (Id.Task.to_int x)
end
(** GetTask : get a new task to do *)
module GetTask_msg : sig
type t =
{ client_id: Id.Client.t ;
state: State.t ;
}
val create : state:string -> client_id:string -> t
val to_string : t -> string
end = struct
type t =
{ client_id: Id.Client.t ;
state: State.t ;
}
let create ~state ~client_id =
{ client_id = Id.Client.of_string client_id ; state = State.of_string state }
let to_string x =
Printf.sprintf "get_task %s %d"
(State.to_string x.state)
(Id.Client.to_int x.client_id)
end
(** GetTaskReply : Reply to the GetTask message *)
module GetTaskReply_msg : sig
type t
val create : task_id:Id.Task.t -> task:string -> t
val to_string : t -> string
end = struct
type t =
{ task_id: Id.Task.t ;
task : string ;
}
let create ~task_id ~task = { task_id ; task }
let to_string x =
Printf.sprintf "get_task_reply %d %s" (Id.Task.to_int x.task_id) x.task
end
(** TaskDone : Inform the server that a task is finished *)
module TaskDone_msg : sig
type t =
{ client_id: Id.Client.t ;
state: State.t ;
task_id: Id.Task.t;
}
val create : state:string -> client_id:string -> task_id:string -> t
val to_string : t -> string
end = struct
type t =
{ client_id: Id.Client.t ;
state: State.t ;
task_id: Id.Task.t;
}
let create ~state ~client_id ~task_id =
{ client_id = Id.Client.of_string client_id ;
state = State.of_string state ;
task_id = Id.Task.of_string task_id }
let to_string x =
Printf.sprintf "task_done %s %d %d"
(State.to_string x.state)
(Id.Client.to_int x.client_id)
(Id.Task.to_int x.task_id)
end
(** Terminate *)
module Terminate_msg : sig
type t
val create : unit -> t
val to_string : t -> string
end = struct
type t = Terminate
let create () = Terminate
let to_string x = "terminate"
end
(** OK *)
module Ok_msg : sig
type t
val create : unit -> t
val to_string : t -> string
end = struct
type t = Ok
let create () = Ok
let to_string x = "ok"
end
(** Error *)
module Error_msg : sig
type t
val create : string -> t
val to_string : t -> string
end = struct
type t = string
let create x = x
let to_string x =
String.concat ~sep:" " [ "error" ; x ]
end
(** Message *)
type t =
| Newjob of Newjob_msg.t
| Connect of Connect_msg.t
| ConnectReply of ConnectReply_msg.t
| Disconnect of Disconnect_msg.t
| DisconnectReply of DisconnectReply_msg.t
| GetTask of GetTask_msg.t
| GetTaskReply of GetTaskReply_msg.t
| AddTask of AddTask_msg.t
| AddTaskReply of AddTaskReply_msg.t
| TaskDone of TaskDone_msg.t
| Terminate of Terminate_msg.t
| Ok of Ok_msg.t
| Error of Error_msg.t
let of_string s =
let l =
String.split ~on:' ' s
|> List.filter ~f:(fun x -> (String.strip x) <> "")
|> List.map ~f:String.lowercase
in
match l with
| "add_task" :: state :: task ->
AddTask (AddTask_msg.create ~state ~task:(String.concat ~sep:" " task) )
| "get_task" :: state :: client_id :: [] ->
GetTask (GetTask_msg.create ~state ~client_id)
| "task_done" :: state :: client_id :: task_id :: [] ->
TaskDone (TaskDone_msg.create ~state ~client_id ~task_id)
| "disconnect" :: state :: client_id :: [] ->
Disconnect (Disconnect_msg.create ~state ~client_id)
| "connect" :: t :: [] ->
Connect (Connect_msg.create t)
| "new_job" :: state :: push_address_tcp :: push_address_inproc :: [] ->
Newjob (Newjob_msg.create push_address_tcp push_address_inproc state)
| "terminate" :: [] ->
Terminate (Terminate_msg.create () )
| "ok" :: [] ->
Ok (Ok_msg.create ())
| "error" :: rest ->
Error (Error_msg.create (String.concat ~sep:" " rest))
| _ -> failwith "Message not understood"
let to_string = function
| Newjob x -> Newjob_msg.to_string x
| Connect x -> Connect_msg.to_string x
| ConnectReply x -> ConnectReply_msg.to_string x
| Disconnect x -> Disconnect_msg.to_string x
| DisconnectReply x -> DisconnectReply_msg.to_string x
| GetTask x -> GetTask_msg.to_string x
| GetTaskReply x -> GetTaskReply_msg.to_string x
| AddTask x -> AddTask_msg.to_string x
| AddTaskReply x -> AddTaskReply_msg.to_string x
| TaskDone x -> TaskDone_msg.to_string x
| Terminate x -> Terminate_msg.to_string x
| Ok x -> Ok_msg.to_string x
| Error x -> Error_msg.to_string x

1
ocaml/Qpackage.ml
View File

@ -55,7 +55,6 @@ let executables = lazy (
In_channel.input_lines in_channel In_channel.input_lines in_channel
|> List.map ~f:(fun x -> |> List.map ~f:(fun x ->
let e = String.split ~on:' ' x let e = String.split ~on:' ' x
|> List.map ~f:String.strip
|> List.filter ~f:(fun x -> x <> "") |> List.filter ~f:(fun x -> x <> "")
in in
match e with match e with

118
ocaml/Queuing_system.ml Normal file
View File

@ -0,0 +1,118 @@
open Core.Std
type t =
{ queued : Id.Task.t list ;
running : (Id.Task.t, Id.Client.t) Map.Poly.t ;
tasks : (Id.Task.t, string) Map.Poly.t;
clients : Id.Client.t Set.Poly.t;
next_client_id : Id.Client.t;
next_task_id : Id.Task.t;
}
let create () =
{ queued = [] ;
running = Map.Poly.empty ;
tasks = Map.Poly.empty;
clients = Set.Poly.empty;
next_client_id = Id.Client.of_int 1;
next_task_id = Id.Task.of_int 1;
}
let add_task ~task q =
let task_id =
q.next_task_id
in
{ q with
queued = task_id :: q.queued ;
tasks = Map.add q.tasks ~key:task_id ~data:task ;
next_task_id = Id.Task.increment task_id ;
}, task_id
let add_client q =
let client_id =
q.next_client_id
in
{ q with
clients = Set.add q.clients client_id;
next_client_id = Id.Client.increment client_id;
}, client_id
let pop_task ~client_id q =
let { queued ; running ; _ } =
q
in
assert (Set.mem q.clients client_id);
match queued with
| task_id :: new_queue ->
let new_q =
{ q with
queued = new_queue ;
running = Map.add running ~key:task_id ~data:client_id ;
}
in new_q, Some task_id, (Map.find q.tasks task_id)
| [] -> q, None, None
let del_client ~client_id q =
assert (Set.mem q.clients client_id);
{ q with
clients = Set.remove q.clients client_id }
let end_task ~task_id ~client_id q =
let { running ; tasks ; _ } =
q
in
assert (Set.mem q.clients client_id);
let () =
match Map.Poly.find running task_id with
| None -> failwith "Task already finished"
| Some client_id_check -> assert (client_id_check = client_id)
in
{ q with
running = Map.remove running task_id ;
tasks = Map.remove tasks task_id ;
}
let number_of_queued q =
List.length q.queued
let number_of_running q =
Map.length q.running
let to_string { queued ; running ; tasks ; _ } =
let q =
List.map ~f:Id.Task.to_string queued
|> String.concat ~sep:" ; "
and r =
Map.Poly.to_alist running
|> List.map ~f:(fun (t,c) -> "("^(Id.Task.to_string t)^", "
^(Id.Client.to_string c)^")")
|> String.concat ~sep:" ; "
and t =
Map.Poly.to_alist tasks
|> List.map ~f:(fun (t,c) -> "("^(Id.Task.to_string t)^", \""
^c^"\")")
|> String.concat ~sep:" ; "
in
Printf.sprintf "{
queued : { %s }
running : { %s }
tasks : [ %s
]
}" q r t

337
ocaml/TaskServer.ml Normal file
View File

@ -0,0 +1,337 @@
open Core.Std
open Qptypes
(**
The tasks server listens on a REQ socket and accepts the following commands:
* "new_job %s %s %s" state push_address_tcp push_address_inproc -> "OK"
-> "OK"
* "connect %s" ["tcp"|"inproc"]
-> "%d %s %s" id state push_address
* "disconnect %d" id
-> "OK"
* "get_task %d %s" id state
-> "%d %s" task_id task
* "task_done %d task_id %s" id state
-> "%d %s" task_id task
*)
let bind_socket ~socket_type ~socket ~address =
try
ZMQ.Socket.bind socket address
with
| Unix.Unix_error (_, message, f) ->
failwith @@ Printf.sprintf
"\n%s\nUnable to bind the %s socket :\n %s\n%s"
f socket_type address message
| other_exception -> raise other_exception
(** Name of the host on which the server runs *)
let hostname = lazy (
try
Unix.gethostname ()
with
| _ -> "localhost"
)
(** IP address *)
let ip_address = lazy (
match Sys.getenv "QP_NIC" with
| None ->
begin
try
Lazy.force hostname
|> Unix.Inet_addr.of_string_or_getbyname
|> Unix.Inet_addr.to_string
with
| Unix.Unix_error _ ->
failwith "Unable to find IP address from host name."
end
| Some interface ->
begin
try
ok_exn Linux_ext.get_ipv4_address_for_interface interface
with
| Unix.Unix_error _ ->
Lazy.force hostname
|> Unix.Inet_addr.of_string_or_getbyname
|> Unix.Inet_addr.to_string
end
)
let stop ~port =
let zmq_context =
ZMQ.Context.create ()
in
let req_socket =
ZMQ.Socket.create zmq_context ZMQ.Socket.req
and address =
Printf.sprintf "tcp://%s:%d" (Lazy.force ip_address) port
in
ZMQ.Socket.connect req_socket address;
Message.Terminate (Message.Terminate_msg.create ())
|> Message.to_string
|> ZMQ.Socket.send ~block:false req_socket ;
let msg =
ZMQ.Socket.recv req_socket
|> Message.of_string
in
let () =
match msg with
| Message.Ok _ -> ()
| _ -> failwith "Problem in termination"
in
ZMQ.Socket.set_linger_period req_socket 1000;
ZMQ.Socket.close req_socket
(** Run the task server *)
let run ~port =
let zmq_context =
ZMQ.Context.create ()
in
let rep_socket =
ZMQ.Socket.create zmq_context ZMQ.Socket.rep
and address =
Printf.sprintf "tcp://%s:%d" (Lazy.force ip_address) port
in
bind_socket "REP" rep_socket address;
let pollitem =
ZMQ.Poll.mask_of
[| (rep_socket, ZMQ.Poll.In) |]
in
Printf.printf "Task server running : %s\n%!" address;
(** State variables *)
let q = ref
(Queuing_system.create ())
and running =
ref true
and job =
ref None
in
let get_state () =
match !job with
| None -> None
| Some j -> Some j.Message.Newjob_msg.state
in
let get_tcp_address () =
match !job with
| Some j -> Address.Tcp j.Message.Newjob_msg.address_tcp
| None -> assert false
in
let get_inproc_address () =
match !job with
| Some j -> Address.Inproc j.Message.Newjob_msg.address_inproc
| None -> assert false
in
let ok =
Message.Ok (Message.Ok_msg.create ())
in
while ( !running )
do
let state =
get_state ()
and polling =
ZMQ.Poll.poll ~timeout:1000 pollitem
in
let terminate () =
running := false;
Message.to_string ok
|> ZMQ.Socket.send ~block:false rep_socket
and newjob x =
q := Queuing_system.create ();
job := Some x;
Message.to_string ok
|> ZMQ.Socket.send ~block:false rep_socket
and connect state msg =
let push_address =
match msg with
| Message.Connect_msg.Tcp -> get_tcp_address ()
| Message.Connect_msg.Inproc -> get_inproc_address ()
| Message.Connect_msg.Ipc -> assert false
in
let new_q, client_id =
Queuing_system.add_client !q
in
q := new_q;
Message.ConnectReply (Message.ConnectReply_msg.create
~state ~client_id ~push_address)
|> Message.to_string
|> ZMQ.Socket.send ~block:false rep_socket
and disconnect state msg =
let s, c =
msg.Message.Disconnect_msg.state ,
msg.Message.Disconnect_msg.client_id
in
assert (s = state);
let new_q =
Queuing_system.del_client ~client_id:c !q
in
q := new_q;
let finished =
Queuing_system.number_of_queued !q +
Queuing_system.number_of_running !q = 0
in
Message.DisconnectReply (Message.DisconnectReply_msg.create
~state ~finished)
|> Message.to_string
|> ZMQ.Socket.send ~block:false rep_socket
and add_task state msg =
let s, task =
msg.Message.AddTask_msg.state,
msg.Message.AddTask_msg.task
in
assert (s = state);
Message.to_string ok
|> ZMQ.Socket.send ~block:false rep_socket
;
begin
match
String.split ~on:' ' msg.Message.AddTask_msg.task
|> List.filter ~f:(fun x -> x <> "")
with
| "triangle" :: str_l :: [] ->
begin
let l =
Int.of_string str_l
in
for j=1 to l
do
let task =
Printf.sprintf "%d %s" j str_l
in
let new_q, _ =
Queuing_system.add_task ~task !q
in
q := new_q
done
end
| "range" :: str_i :: str_j :: [] ->
begin
let i, j =
Int.of_string str_i,
Int.of_string str_j
in
for k=i to (j+1)
do
let task =
Int.to_string k
in
let new_q, task_id =
Queuing_system.add_task ~task !q
in
q := new_q
done
end
| _ ->
let new_q, task_id =
Queuing_system.add_task ~task !q
in
q := new_q
end
and get_task state msg =
let s, client_id =
msg.Message.GetTask_msg.state,
msg.Message.GetTask_msg.client_id
in
assert (s = state);
let new_q, task_id, task =
Queuing_system.pop_task ~client_id !q
in
q := new_q;
let reply =
match (task, task_id) with
| Some task, Some task_id ->
Message.GetTaskReply (Message.GetTaskReply_msg.create ~task ~task_id)
| _ -> Message.Terminate (Message.Terminate_msg.create ())
in
Message.to_string reply
|> ZMQ.Socket.send ~block:false rep_socket
and task_done state msg =
let s, client_id, task_id =
msg.Message.TaskDone_msg.state,
msg.Message.TaskDone_msg.client_id,
msg.Message.TaskDone_msg.task_id
in
assert (s = state);
let new_q =
Queuing_system.end_task ~task_id ~client_id !q
in
q := new_q;
Message.to_string ok
|> ZMQ.Socket.send ~block:false rep_socket
and error msg =
Message.Error (Message.Error_msg.create msg)
|> Message.to_string
|> ZMQ.Socket.send ~block:false rep_socket
in
if (polling.(0) = Some ZMQ.Poll.In) then
let raw_message =
ZMQ.Socket.recv rep_socket
in
try
let message =
Message.of_string raw_message
in
(*
Printf.printf "%d %d : %s\n%!"
(Queuing_system.number_of_queued !q)
(Queuing_system.number_of_running !q)
(Message.to_string message);
Printf.printf "%s\n%!" (Queuing_system.to_string !q); *)
match (state, message) with
| _ , Message.Terminate _ -> terminate ()
| None , Message.Newjob x -> newjob x
| None , _ -> error "No job is running"
| _ , Message.Newjob _ -> error "A job is already running"
| Some s, Message.Connect x -> connect s x
| Some s, Message.Disconnect x -> disconnect s x
| Some s, Message.AddTask x -> add_task s x
| Some s, Message.GetTask x -> get_task s x
| Some s, Message.TaskDone x -> task_done s x
| _ , _ ->
error ("Invalid message : "^(Message.to_string message))
with
| Failure f -> error (f^" : "^raw_message)
| Assert_failure (f,i,j) -> error (Printf.sprintf "%s:%d:%d : %s" f i j raw_message)
done;
ZMQ.Socket.set_linger_period rep_socket 1000;
ZMQ.Socket.close rep_socket
(*
let () =
Printf.printf "export QP_RUN_ADDRESS=tcp://%s:%d\n%!" (Lazy.force ip_address) (Lazy.force port)
*)

2
ocaml/_tags
View File

@ -1,2 +1,2 @@
true: package(core,sexplib.syntax,cryptokit) true: package(core,sexplib.syntax,cryptokit,ZMQ)
true: thread true: thread

12
ocaml/create_git_sha1.sh Executable file
View File

@ -0,0 +1,12 @@
#!/bin/bash
SHA1=$(git log -1 | head -1 | cut -d ' ' -f 2)
DATE=$(git log -1 | grep Date | cut -d ':' -f 2-)
MESSAGE=$(git log -1 | tail -1)
cat << EOF > Git.ml
open Core.Std
let sha1 = "$SHA1" |> String.strip
let date = "$DATE" |> String.strip
let message = "$MESSAGE" |> String.strip
EOF

6
ocaml/qp_create_ezfio_from_xyz.ml
View File

@ -133,10 +133,6 @@ let run ?o b c d m p xyz_file =
let fetch_channel basis = let fetch_channel basis =
let command = 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 Qpackage.root ^ "/scripts/get_basis.sh \"" ^ temp_filename
^ "." ^ basis ^ "\" \"" ^ basis ^"\"" ^ "." ^ basis ^ "\" \"" ^ basis ^"\""
in in
@ -344,7 +340,7 @@ let run ?o b c d m p xyz_file =
(* Doesn't work... *) (* Doesn't work... *)
if (p) then if (p) then
begin begin
Qpackage.root ^ "/scripts/pseudo/put_pseudo_in_ezfio.py " ^ ezfio_file Qpackage.root ^ "/scripts/pseudo/put_pseudo_in_ezfio.py " ^ ezfio_file ^ " " ^ " BFD-Pseudo"
|> Sys.command_exn |> Sys.command_exn
end; end;

38
ocaml/qp_edit.ml
View File

@ -18,11 +18,11 @@ type keyword =
| Mo_basis | Mo_basis
| Nuclei | Nuclei
| Determinants | Determinants
| Hartree_fock
| Integrals_bielec
| Perturbation | Perturbation
| Properties
| Pseudo | Pseudo
| Integrals_bielec
| Properties
| Hartree_fock
;; ;;
@ -33,11 +33,11 @@ let keyword_to_string = function
| Mo_basis -> "MO basis" | Mo_basis -> "MO basis"
| Nuclei -> "Molecule" | Nuclei -> "Molecule"
| Determinants -> "Determinants" | Determinants -> "Determinants"
| Hartree_fock -> "Hartree_fock"
| Integrals_bielec -> "Integrals_bielec"
| Perturbation -> "Perturbation" | Perturbation -> "Perturbation"
| Properties -> "Properties"
| Pseudo -> "Pseudo" | Pseudo -> "Pseudo"
| Integrals_bielec -> "Integrals_bielec"
| Properties -> "Properties"
| Hartree_fock -> "Hartree_fock"
;; ;;
@ -88,16 +88,16 @@ let get s =
f Determinants_by_hand.(read, to_rst) f Determinants_by_hand.(read, to_rst)
| Determinants -> | Determinants ->
f Determinants.(read, to_rst) f Determinants.(read, to_rst)
| Hartree_fock ->
f Hartree_fock.(read, to_rst)
| Integrals_bielec ->
f Integrals_bielec.(read, to_rst)
| Perturbation -> | Perturbation ->
f Perturbation.(read, to_rst) f Perturbation.(read, to_rst)
| Properties ->
f Properties.(read, to_rst)
| Pseudo -> | Pseudo ->
f Pseudo.(read, to_rst) f Pseudo.(read, to_rst)
| Integrals_bielec ->
f Integrals_bielec.(read, to_rst)
| Properties ->
f Properties.(read, to_rst)
| Hartree_fock ->
f Hartree_fock.(read, to_rst)
end end
with with
| Sys_error msg -> (Printf.eprintf "Info: %s\n%!" msg ; "") | Sys_error msg -> (Printf.eprintf "Info: %s\n%!" msg ; "")
@ -136,11 +136,11 @@ let set str s =
let open Input in let open Input in
match s with match s with
| Determinants -> write Determinants.(of_rst, write) s | 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 | Perturbation -> write Perturbation.(of_rst, write) s
| Properties -> write Properties.(of_rst, write) s
| Pseudo -> write Pseudo.(of_rst, write) s | Pseudo -> write Pseudo.(of_rst, write) s
| Integrals_bielec -> write Integrals_bielec.(of_rst, write) s
| Properties -> write Properties.(of_rst, write) s
| Hartree_fock -> write Hartree_fock.(of_rst, write) s
| Electrons -> write Electrons.(of_rst, write) s | Electrons -> write Electrons.(of_rst, write) s
| Determinants_by_hand -> write Determinants_by_hand.(of_rst, write) s | Determinants_by_hand -> write Determinants_by_hand.(of_rst, write) s
| Nuclei -> write Nuclei.(of_rst, write) s | Nuclei -> write Nuclei.(of_rst, write) s
@ -189,11 +189,11 @@ let run check_only ezfio_filename =
Ao_basis; Ao_basis;
Electrons ; Electrons ;
Determinants ; Determinants ;
Hartree_fock ;
Integrals_bielec ;
Perturbation ; Perturbation ;
Properties ;
Pseudo ; Pseudo ;
Integrals_bielec ;
Properties ;
Hartree_fock ;
Mo_basis; Mo_basis;
Determinants_by_hand ; Determinants_by_hand ;
] ]
@ -212,7 +212,7 @@ let run check_only ezfio_filename =
match check_only with match check_only with
| true -> () | true -> ()
| false -> | false ->
Printf.sprintf "%s %s ; tput sgr0 2> /dev/null" editor temp_filename Printf.sprintf "%s %s" editor temp_filename
|> Sys.command_exn |> Sys.command_exn
; ;

4
ocaml/qp_print.ml
View File

@ -41,7 +41,7 @@ let run_i ~action ezfio_filename =
let action = create_i_action action in let action = create_i_action action in
if (not (Sys.file_exists_exn ezfio_filename)) then if (not (Sys.file_exists_exn ezfio_filename)) then
failwith (ezfio_filename^" does not exists"); failwith (ezfio_filename^" does not exist");
Ezfio.set_file ezfio_filename; Ezfio.set_file ezfio_filename;
@ -133,7 +133,7 @@ let run_i ~action ezfio_filename =
let run_o ~action ezfio_filename = let run_o ~action ezfio_filename =
if (not (Sys.file_exists_exn ezfio_filename)) then if (not (Sys.file_exists_exn ezfio_filename)) then
failwith (ezfio_filename^" does not exists"); failwith (ezfio_filename^" does not exist");
(* Open EZFIO *) (* Open EZFIO *)
Ezfio.set_file ezfio_filename; Ezfio.set_file ezfio_filename;

39
ocaml/qp_print_basis.ml Normal file
View File

@ -0,0 +1,39 @@
open Core.Std
open Qptypes
let basis ()=
let ezfio_filename =
Sys.argv.(1)
in
if (not (Sys.file_exists_exn ezfio_filename)) then
failwith "Error reading EZFIO file";
Ezfio.set_file ezfio_filename;
let basis =
match Input.Ao_basis.read () with
| Some basis -> basis
| _ -> failwith "Error reading basis set"
in
Input.Ao_basis.to_rst basis
|> Rst_string.to_string
|> print_endline
;;
let mo ()=
let ezfio_filename =
Sys.argv.(1)
in
if (not (Sys.file_exists_exn ezfio_filename)) then
failwith "Error reading EZFIO file";
Ezfio.set_file ezfio_filename;
let mo_coef =
match Input.Mo_basis.read () with
| Some mo_coef -> mo_coef
| _ -> failwith "Error reading the mo set"
in
Input.Mo_basis.to_rst mo_coef
|> Rst_string.to_string
|> print_endline
;;
basis ();;
mo ();;

33
ocaml/qp_run.ml
View File

@ -17,13 +17,38 @@ let run exe ezfio_file =
if (not (List.exists ~f:(fun (x,_) -> x = exe) executables)) then if (not (List.exists ~f:(fun (x,_) -> x = exe) executables)) then
failwith ("Executable "^exe^" not found"); failwith ("Executable "^exe^" not found");
Printf.printf "%s\n" (Time.to_string time_start);
Printf.printf "===============\nQuantum Package\n===============\n\n"; Printf.printf "===============\nQuantum Package\n===============\n\n";
Printf.printf "Date : %s\n\n%!" (Time.to_string time_start); Printf.printf "Git Commit: %s\n" Git.message;
Printf.printf "Git Date : %s\n" Git.date;
Printf.printf "Git SHA1 : %s\n" Git.sha1;
Printf.printf "\n\n%!";
(** Check input *)
match (Sys.command ("qp_edit -c "^ezfio_file)) with match (Sys.command ("qp_edit -c "^ezfio_file)) with
| 0 -> () | 0 -> ()
| i -> failwith "Error: Input inconsistent\n"; | i -> failwith "Error: Input inconsistent\n";
; ;
(** Start task server *)
let port_number =
12345
in
let address =
Printf.sprintf "tcp://%s:%d" (Lazy.force TaskServer.ip_address) port_number
in
let task_thread =
let thread =
Thread.create ( fun () ->
TaskServer.run port_number )
in
thread ();
in
Unix.putenv ~key:"QP_RUN_ADDRESS" ~data:address;
(** Run executable *)
let exe = let exe =
match (List.find ~f:(fun (x,_) -> x = exe) executables) with match (List.find ~f:(fun (x,_) -> x = exe) executables) with
| None -> assert false | None -> assert false
@ -34,6 +59,9 @@ let run exe ezfio_file =
| i -> Printf.printf "Program exited with code %d.\n%!" i; | i -> Printf.printf "Program exited with code %d.\n%!" i;
; ;
TaskServer.stop ~port:port_number;
Thread.join task_thread;
let duration = Time.diff (Time.now()) time_start let duration = Time.diff (Time.now()) time_start
|> Core.Span.to_string in |> Core.Span.to_string in
Printf.printf "Wall time : %s\n\n" duration; Printf.printf "Wall time : %s\n\n" duration;
@ -60,6 +88,7 @@ Executes a Quantum Package binary file among these:\n\n"
(fun exe ezfio_file () -> (fun exe ezfio_file () ->
run exe ezfio_file run exe ezfio_file
) )
|> Command.run |> Command.run ~version: Git.sha1 ~build_info: Git.message
;; ;;

89
ocaml/test_message.ml Normal file
View File

@ -0,0 +1,89 @@
open Core.Std
let () =
Message.of_string "new_job tcp://127.0.0.1 inproc://ao_ints:12345 ao_integrals"
|> Message.to_string
|> print_endline
;
Message.of_string "connect tcp"
|> Message.to_string
|> print_endline
;
Message.of_string "connect inproc"
|> Message.to_string
|> print_endline
;
Message.of_string "disconnect 3 mystate"
|> Message.to_string
|> print_endline
;
Message.of_string "get_task 3 mystate"
|> Message.to_string
|> print_endline
;
Message.of_string "task_done 1 mystate 3"
|> Message.to_string
|> print_endline
;
Message.of_string "add_task mystate 1 2 3 4 5 6"
|> Message.to_string
|> print_endline
;
try
Message.of_string "new_job inproc://ao_ints tcp://127.0.0.1:12345 ao_integrals"
|> Message.to_string
|> print_endline
;
failwith "Should have failed"
with
| Assert_failure _ -> print_endline "OK"
;
try
Message.of_string "new_job tcp://ao_ints inproc://ao_ints"
|> Message.to_string
|> print_endline
;
assert false
with
| Failure _ -> print_endline "OK"
;
try
Message.of_string "disconnect -4 mystate"
|> Message.to_string
|> print_endline
;
assert false
with
| Assert_failure _ -> print_endline "OK"
;
try
Message.of_string "disconnect mystate 3"
|> Message.to_string
|> print_endline
;
assert false
with
| Failure _ -> print_endline "OK"
;
try
Message.of_string "connect tcp tcp://127.0.0.1"
|> Message.to_string
|> print_endline
;
assert false
with
| Failure _ -> print_endline "OK"
;

102
ocaml/test_queuing_system.ml Normal file
View File

@ -0,0 +1,102 @@
open Core.Std
let () =
let nclients =
8
in
let q =
Queuing_system.create ()
in
let tasks =
Array.init 20 ~f:(fun i -> Printf.sprintf "Task %d" i)
|> Array.to_list
in
let (q,_) =
List.fold_left tasks ~init:(q, q.Queuing_system.next_task_id)
~f:(fun (q,_) task -> Queuing_system.add_task ~task q)
in
print_endline @@ Queuing_system.to_string q ;
let rec aux q clients = function
| 0 -> q, clients
| i ->
let new_q, client_id =
Queuing_system.add_client q
in
aux new_q (client_id::clients) (i-1)
in
let q, _ =
aux q [] nclients
in
let rec aux q = function
| 0 -> q
| i ->
begin
let c =
Id.Client.of_int i
in
let new_q, task_id, task =
Queuing_system.pop_task ~client_id:c q
in
begin
match task_id, task with
| Some task_id, Some task ->
Printf.printf "Task Running: %d %s\n" (Id.Task.to_int task_id) task
| _ -> Printf.printf "Done!\n"
end;
aux new_q (i-1)
end
in
let rec aux2 q = function
| 0 -> q
| i ->
begin
let task_id =
(Id.Task.of_int i)
in
try
let client_id =
Map.Poly.find_exn q.Queuing_system.running task_id
in
let new_q =
Queuing_system.end_task ~task_id ~client_id q
in
Printf.printf "Task Done : %d\n" (Id.Task.to_int task_id) ;
aux2 new_q (i-1)
with
| _ -> aux2 q 0
end
in
let q =
aux q nclients
in
print_endline @@ Queuing_system.to_string q ;
let q =
aux2 q nclients
in
print_endline @@ Queuing_system.to_string q ;
Printf.printf "Queued : %d\n Running : %d\n"
(Queuing_system.number_of_queued q)
(Queuing_system.number_of_running q)
;
let q =
aux q nclients
in
print_endline @@ Queuing_system.to_string q ;
let q =
aux2 q nclients
in
print_endline @@ Queuing_system.to_string q ;
(*
List.map ~f:Id.Task.to_int tasks
|> List.iter ~f:(fun x -> Printf.printf "%d\n" x)
*)

5
ocaml/test_task_server.ml Normal file
View File

@ -0,0 +1,5 @@
open Core
let () =
TaskServer.run 12345

46
ocaml/test_task_server.py Executable file
View File

@ -0,0 +1,46 @@
#!/usr/bin/python
import zmq
import sys, os
def main():
context = zmq.Context()
socket = context.socket(zmq.REQ)
socket.connect(os.environ["QP_RUN_ADDRESS"])
def send(msg,expected):
print "Send : ", msg
print " -> ", socket.send(msg)
reply = socket.recv()
print "Reply : ", reply
print ""
assert (reply == expected)
send("new_job ao_integrals tcp://130.120.229.139:12345 inproc://ao_integrals",
"ok")
send("new_job ao_integrals tcp://130.120.229.139:12345 inproc://ao_integrals",
"error A job is already running")
send("connect","error Message not understood : connect")
send("connect tcp","connect_reply ao_integrals 1 tcp://130.120.229.139:12345")
send("connect inproc","connect_reply ao_integrals 2 inproc://ao_integrals")
send("disconnect ao_integrals 3","error Queuing_system.ml:65:2 : disconnect ao_integrals 3")
send("disconnect ao_integrals 2","disconnect_reply ao_integrals 1")
send("connect inproc","connect_reply ao_integrals 3 inproc://ao_integrals")
for i in range(10):
send("add_task ao_integrals %d %d"%(i,i+10), "ok")
for i in range(10):
send("get_task ao_integrals 3", "get_task_reply %d %d %d"%(i+1,i,i+10))
send("task_done ao_integrals 3 %d"%(i+1), "ok")
send("get_task ao_integrals 3", "terminate")
send("terminate","ok")
if __name__ == '__main__':
main()

2
plugins/CAS_SD/cas_s_selected.irp.f
View File

@ -30,7 +30,7 @@ program full_ci
endif endif
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_st))) > pt2_max) 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) call H_apply_CAS_S_selected_no_skip(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef PROVIDE psi_coef
PROVIDE psi_det PROVIDE psi_det

1
plugins/CISD/.gitignore vendored
View File

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

1
plugins/CISD_selected/.gitignore vendored
View File

@ -23,6 +23,7 @@ Pseudo
Selectors_full Selectors_full
SingleRefMethod SingleRefMethod
Utils Utils
ZMQ
cisd_selection cisd_selection
ezfio_interface.irp.f ezfio_interface.irp.f
irpf90.make irpf90.make

4
plugins/CISD_selected/README.rst
View File

@ -196,6 +196,10 @@ Documentation
.. by the `update_README.py` script. .. 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 h_apply_cisd
Calls H_apply on the HF determinant and selects all connected single and double 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. excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.

2
plugins/CISD_selected/cisd_selection.irp.f → plugins/CISD_selected/cisd_selected.irp.f
View File

@ -31,9 +31,9 @@ program cisd
print *, 'PT2 = ', pt2(i) print *, 'PT2 = ', pt2(i)
print *, 'E = ', CI_energy(i) print *, 'E = ', CI_energy(i)
print *, 'E_before +PT2 = ', (E_old(i)+pt2(i)) print *, 'E_before +PT2 = ', (E_old(i)+pt2(i))
! print *, 'E+PT2_new= ', (E_old(1)+1.d0*pt2(1)+H_pert_diag(1))/(1.d0 +norm_pert(1))
enddo enddo
E_old = CI_energy E_old = CI_energy
call save_wavefunction
if (abort_all) then if (abort_all) then
exit exit
endif endif

1
plugins/FCIdump/.gitignore vendored
View File

@ -16,6 +16,7 @@ Makefile.depend
Nuclei Nuclei
Pseudo Pseudo
Utils Utils
ZMQ
ezfio_interface.irp.f ezfio_interface.irp.f
fcidump fcidump
irpf90.make irpf90.make

1
plugins/Full_CI/.gitignore vendored
View File

@ -22,6 +22,7 @@ Properties
Pseudo Pseudo
Selectors_full Selectors_full
Utils Utils
ZMQ
ezfio_interface.irp.f ezfio_interface.irp.f
full_ci full_ci
full_ci_no_skip full_ci_no_skip

1
plugins/Hartree_Fock/.gitignore vendored
View File

@ -18,6 +18,7 @@ Nuclei
Pseudo Pseudo
SCF SCF
Utils Utils
ZMQ
ezfio_interface.irp.f ezfio_interface.irp.f
irpf90.make irpf90.make
irpf90_entities irpf90_entities

2
plugins/Hartree_Fock/EZFIO.cfg
View File

@ -2,7 +2,7 @@
type: Threshold type: Threshold
doc: Threshold on the convergence of the Hartree Fock energy doc: Threshold on the convergence of the Hartree Fock energy
interface: ezfio,provider,ocaml interface: ezfio,provider,ocaml
default: 1.e-10 default: 1.e-12
[n_it_scf_max] [n_it_scf_max]
type: Strictly_positive_int type: Strictly_positive_int

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

@ -73,10 +73,6 @@
enddo enddo
endif 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 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)

4
plugins/Hartree_Fock/README.rst
View File

@ -52,7 +52,7 @@ Documentation
Diagonal Fock matrix in the MO basis Diagonal Fock matrix in the MO basis
`diagonal_fock_matrix_mo_sum <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/diagonalize_fock.irp.f#L67>`_ `diagonal_fock_matrix_mo_sum <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/diagonalize_fock.irp.f#L112>`_
diagonal element of the fock matrix calculated as the sum over all the interactions diagonal element of the fock matrix calculated as the sum over all the interactions
with all the electrons in the RHF determinant with all the electrons in the RHF determinant
diagonal_Fock_matrix_mo_sum(i) = sum_{j=1, N_elec} 2 J_ij -K_ij diagonal_Fock_matrix_mo_sum(i) = sum_{j=1, N_elec} 2 J_ij -K_ij
@ -114,7 +114,7 @@ Documentation
.br .br
`fock_mo_to_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L375>`_ `fock_mo_to_ao <http://github.com/LCPQ/quantum_package/tree/master/plugins/Hartree_Fock/Fock_matrix.irp.f#L392>`_
Undocumented Undocumented

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

@ -11,55 +11,55 @@
double precision, allocatable :: work(:), F(:,:), S(:,:) double precision, allocatable :: work(:), F(:,:), S(:,:)
if (mo_tot_num == ao_num) then ! if (mo_tot_num == ao_num) then
! Solve H.C = E.S.C in AO basis set ! ! Solve H.C = E.S.C in AO basis set
!
allocate(F(ao_num_align,ao_num), S(ao_num_align,ao_num) ) ! allocate(F(ao_num_align,ao_num), S(ao_num_align,ao_num) )
do j=1,ao_num ! do j=1,ao_num
do i=1,ao_num ! do i=1,ao_num
S(i,j) = ao_overlap(i,j) ! S(i,j) = ao_overlap(i,j)
F(i,j) = Fock_matrix_ao(i,j) ! F(i,j) = Fock_matrix_ao(i,j)
enddo ! enddo
enddo ! enddo
!
n = ao_num ! n = ao_num
lwork = 1+6*n + 2*n*n ! lwork = 1+6*n + 2*n*n
liwork = 3 + 5*n ! liwork = 3 + 5*n
!
allocate(work(lwork), iwork(liwork) ) ! allocate(work(lwork), iwork(liwork) )
!
lwork = -1 ! lwork = -1
liwork = -1 ! liwork = -1
!
call dsygvd(1,'v','u',ao_num,F,size(F,1),S,size(S,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) ! diagonal_Fock_matrix_mo, work, lwork, iwork, liwork, info)
!
if (info /= 0) then ! if (info /= 0) then
print *, irp_here//' failed : ', info ! print *, irp_here//' failed : ', info
stop 1 ! stop 1
endif ! endif
lwork = int(work(1)) ! lwork = int(work(1))
liwork = iwork(1) ! liwork = iwork(1)
deallocate(work,iwork) ! deallocate(work,iwork)
allocate(work(lwork), iwork(liwork) ) ! allocate(work(lwork), iwork(liwork) )
!
call dsygvd(1,'v','u',ao_num,F,size(F,1),S,size(S,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) ! diagonal_Fock_matrix_mo, work, lwork, iwork, liwork, info)
!
if (info /= 0) then ! if (info /= 0) then
print *, irp_here//' failed : ', info ! print *, irp_here//' failed : ', info
stop 1 ! stop 1
endif ! endif
do j=1,mo_tot_num ! do j=1,mo_tot_num
do i=1,ao_num ! do i=1,ao_num
eigenvectors_Fock_matrix_mo(i,j) = F(i,j) ! eigenvectors_Fock_matrix_mo(i,j) = F(i,j)
enddo ! enddo
enddo ! enddo
!
deallocate(work, iwork, F, S) ! deallocate(work, iwork, F, S)
!
else ! else
!
! Solve H.C = E.C in MO basis set ! Solve H.C = E.C in MO basis set
allocate( F(mo_tot_num_align,mo_tot_num) ) allocate( F(mo_tot_num_align,mo_tot_num) )
@ -69,6 +69,16 @@
enddo enddo
enddo enddo
! Insert level shift here
do i = elec_beta_num+1, elec_alpha_num
F(i,i) += 0.5d0*level_shift
enddo
do i = elec_alpha_num+1, mo_tot_num
F(i,i) += level_shift
enddo
n = mo_tot_num n = mo_tot_num
lwork = 1+6*n + 2*n*n lwork = 1+6*n + 2*n*n
liwork = 3 + 5*n liwork = 3 + 5*n
@ -105,7 +115,8 @@
0.d0, eigenvectors_Fock_matrix_mo, size(eigenvectors_Fock_matrix_mo,1)) 0.d0, eigenvectors_Fock_matrix_mo, size(eigenvectors_Fock_matrix_mo,1))
deallocate(work, iwork, F) deallocate(work, iwork, F)
endif
! endif
END_PROVIDER END_PROVIDER

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

@ -8,8 +8,6 @@ subroutine huckel_guess
double precision :: c double precision :: c
character*(64) :: label character*(64) :: label
mo_coef = ao_ortho_lowdin_coef
TOUCH mo_coef
label = "Guess" label = "Guess"
call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral, & call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral, &
size(mo_mono_elec_integral,1), & size(mo_mono_elec_integral,1), &

2
plugins/MP2/.gitignore vendored
View File

@ -22,8 +22,10 @@ Pseudo
Selectors_full Selectors_full
SingleRefMethod SingleRefMethod
Utils Utils
ZMQ
ezfio_interface.irp.f ezfio_interface.irp.f
irpf90.make irpf90.make
irpf90_entities irpf90_entities
mp2 mp2
mp2_wf
tags tags

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

@ -8,7 +8,7 @@ s.set_perturbation("Moller_plesset")
print s print s
s = H_apply("mp2_selection") s = H_apply("mp2_selection")
s.set_selection_pt2("Moller_plesset") s.set_selection_pt2("Moller_Plesset")
print s print s
END_SHELL END_SHELL

29
plugins/MP2/README.rst
View File

@ -83,6 +83,35 @@ h_apply_mp2_monoexc
Assume N_int is already provided. Assume N_int is already provided.
h_apply_mp2_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_mp2_selection_diexc
Undocumented
h_apply_mp2_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_mp2_selection_diexcp
Undocumented
h_apply_mp2_selection_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>`_ `mp2 <http://github.com/LCPQ/quantum_package/tree/master/plugins/MP2/mp2.irp.f#L1>`_
Undocumented Undocumented
`mp2_wf <http://github.com/LCPQ/quantum_package/tree/master/plugins/MP2/mp2_wf.irp.f#L1>`_
Save the MP2 wave function

3
plugins/MP2/mp2_wf.irp.f
View File

@ -12,8 +12,7 @@ program mp2_wf
allocate (pt2(N_st), norm_pert(N_st), H_pert_diag(N_st)) allocate (pt2(N_st), norm_pert(N_st), H_pert_diag(N_st))
pt2 = 1.d0 pt2 = 1.d0
selection_criterion = 1.e-12 selection_criterion_factor = 0.d0
selection_criterion_min = 1.e-12
TOUCH selection_criterion_min selection_criterion selection_criterion_factor TOUCH selection_criterion_min selection_criterion selection_criterion_factor
call H_apply_mp2_selection(pt2, norm_pert, H_pert_diag, N_st) call H_apply_mp2_selection(pt2, norm_pert, H_pert_diag, N_st)
psi_det = psi_det_sorted psi_det = psi_det_sorted

12
plugins/MRCC_Utils/mrcc_utils.irp.f
View File

@ -23,23 +23,26 @@
call i_h_psi(psi_non_ref(1,1,i), psi_ref_restart, psi_ref_coef_restart, N_int, N_det_ref,& call i_h_psi(psi_non_ref(1,1,i), psi_ref_restart, psi_ref_coef_restart, N_int, N_det_ref,&
size(psi_ref_coef_restart,1), n_states, ihpsi) size(psi_ref_coef_restart,1), n_states, ihpsi)
call i_H_j(psi_non_ref(1,1,i),psi_non_ref(1,1,i),N_int,hii) call i_H_j(psi_non_ref(1,1,i),psi_non_ref(1,1,i),N_int,hii)
! TODO --- Test perturbatif ------
do k=1,N_states do k=1,N_states
lambda_pert(k,i) = 1.d0 / (psi_ref_energy_diagonalized(k)-hii) lambda_pert(k,i) = 1.d0 / (psi_ref_energy_diagonalized(k)-hii)
call i_h_psi(psi_non_ref(1,1,i), psi_ref, psi_ref_coef, N_int, N_det_ref,size(psi_ref_coef,1), n_states, ihpsi_current) call i_h_psi(psi_non_ref(1,1,i), psi_ref, psi_ref_coef, N_int, N_det_ref,size(psi_ref_coef,1), n_states, ihpsi_current)
tmp = psi_non_ref_coef(i,k)/ihpsi_current(k) tmp = psi_non_ref_coef(i,k)/ihpsi_current(k)
i_pert = 1
if((ihpsi(k) * lambda_pert(k,i))/psi_non_ref_coef_restart(i,k) .ge. 0.5d0 &
.and. (ihpsi(k) * lambda_pert(k,i))/psi_non_ref_coef_restart(i,k) > 0.d0 )then ! test on the first order coefficient
i_pert = 0 i_pert = 0
endif ! Perturbation only if 1st order < 0.5 x second order
if((ihpsi(k) * lambda_pert(k,i)) < 0.5d0 * psi_non_ref_coef_restart(i,k) )then
i_pert = 1
else
do j = 1, N_det_ref do j = 1, N_det_ref
call i_H_j(psi_non_ref(1,1,i),psi_ref(1,1,j),N_int,hij) call i_H_j(psi_non_ref(1,1,i),psi_ref(1,1,j),N_int,hij)
! Perturbation diverges when hij*tmp > 0.5
if(dabs(hij * tmp).ge.0.5d0)then if(dabs(hij * tmp).ge.0.5d0)then
i_pert_count +=1 i_pert_count +=1
i_pert = 1 i_pert = 1
exit exit
endif endif
enddo enddo
endif
if( i_pert == 1)then if( i_pert == 1)then
pert_determinants(k,i) = i_pert pert_determinants(k,i) = i_pert
endif endif
@ -50,6 +53,7 @@
lambda_mrcc(k,i) = psi_non_ref_coef(i,k)/ihpsi_current(k) lambda_mrcc(k,i) = psi_non_ref_coef(i,k)/ihpsi_current(k)
endif endif
enddo enddo
! TODO --- Fin test perturbatif ------
enddo enddo
!if(oscillations)then !if(oscillations)then
! print*,'AVERAGING the lambda_mrcc with those of the previous iterations' ! print*,'AVERAGING the lambda_mrcc with those of the previous iterations'

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

@ -2,6 +2,8 @@ BEGIN_SHELL [ /usr/bin/env python ]
import perturbation import perturbation
END_SHELL 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,key_mask,fock_diag_tmp) 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 implicit none
BEGIN_DOC BEGIN_DOC
@ -29,9 +31,18 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c
logical :: fullMatch logical :: fullMatch
logical, external :: is_connected_to logical, external :: is_connected_to
integer(bit_kind), allocatable :: microlist(:,:,:), microlist_zero(:,:,:)
integer, allocatable :: idx_microlist(:), N_microlist(:), ptr_microlist(:), idx_microlist_zero(:)
integer :: mobiles(2), smallerlist
integer(bit_kind), allocatable :: microlist_gen(:,:,:)
integer, allocatable :: idx_microlist_gen(:), N_microlist_gen(:), ptr_microlist_gen(:)
allocate( minilist(Nint,2,N_det_selectors), & allocate( minilist(Nint,2,N_det_selectors), &
minilist_gen(Nint,2,N_det_generators), & minilist_gen(Nint,2,N_det_generators), &
idx_minilist(N_det_selectors) ) idx_minilist(N_det_selectors))
ASSERT (Nint > 0) ASSERT (Nint > 0)
@ -40,14 +51,47 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c
ASSERT (minval(sum_norm_pert) >= 0.d0) ASSERT (minval(sum_norm_pert) >= 0.d0)
ASSERT (N_st > 0) ASSERT (N_st > 0)
call create_minilist_find_previous(key_mask, psi_det_generators, miniList_gen, i_generator-1, N_minilist_gen, fullMatch, Nint) call create_minilist_find_previous(key_mask, psi_det_generators, miniList_gen, i_generator-1, N_minilist_gen, fullMatch, Nint)
if(fullMatch) then if(fullMatch) then
deallocate( minilist, minilist_gen, idx_minilist ) deallocate( minilist, minilist_gen, idx_minilist )
return return
end if end if
call create_minilist(key_mask, psi_selectors, minilist, idx_miniList, N_det_selectors, N_minilist, Nint)
allocate( microlist(Nint,2,N_minilist*4), &
idx_microlist(N_minilist*4), &
ptr_microlist(0:mo_tot_num*2+1), &
N_microlist(0:mo_tot_num*2) )
call create_minilist(key_mask, psi_selectors, miniList, idx_miniList, N_det_selectors, N_minilist, Nint) allocate( microlist_gen(Nint,2,N_minilist_gen*4), &
idx_microlist_gen(N_minilist_gen*4 ), &
ptr_microlist_gen(0:mo_tot_num*2+1), &
N_microlist_gen(0:mo_tot_num*2) )
if(key_mask(1,1) /= 0) then
call create_microlist(minilist, N_minilist, key_mask, microlist, idx_microlist, N_microlist, ptr_microlist, Nint)
call create_microlist(minilist_gen, N_minilist_gen, key_mask, microlist_gen, idx_microlist_gen, N_microlist_gen,ptr_microlist_gen,Nint)
allocate(microlist_zero(Nint,2,N_minilist))
allocate(idx_microlist_zero(N_minilist))
do i=0,mo_tot_num*2
do k=ptr_microlist(i),ptr_microlist(i+1)-1
idx_microlist(k) = idx_minilist(idx_microlist(k))
end do
end do
if(N_microlist(0) > 0) then
microlist_zero(:,:,1:N_microlist(0)) = microlist(:,:,1:N_microlist(0))
idx_microlist_zero(1:N_microlist(0)) = idx_microlist(1:N_microlist(0))
end if
end if
do i=1,buffer_size do i=1,buffer_size
@ -55,14 +99,43 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c
cycle cycle
endif endif
if(key_mask(1,1) /= 0) then
call getMobiles(buffer(:,:,i), key_mask, mobiles, Nint)
if(N_microlist(mobiles(1)) < N_microlist(mobiles(2))) then
smallerlist = mobiles(1)
else
smallerlist = mobiles(2)
end if
if(N_microlist_gen(smallerlist) > 0) then
if(is_connected_to(buffer(1,1,i), microlist_gen(:,:,ptr_microlist_gen(smallerlist):ptr_microlist_gen(smallerlist+1)-1), Nint, N_microlist_gen(smallerlist))) then
cycle
end if
end if
if(N_microlist_gen(0) > 0) then
if(is_connected_to(buffer(1,1,i), microlist_gen(:,:,1:ptr_microlist_gen(1)-1), Nint, N_microlist_gen(0))) then
cycle
end if
end if
if(N_microlist(smallerlist) > 0) then
microlist_zero(:,:,ptr_microlist(1):ptr_microlist(1)+N_microlist(smallerlist)-1) = microlist(:,:,ptr_microlist(smallerlist):ptr_microlist(smallerlist+1)-1)
idx_microlist_zero(ptr_microlist(1):ptr_microlist(1)+N_microlist(smallerlist)-1) = idx_microlist(ptr_microlist(smallerlist):ptr_microlist(smallerlist+1)-1)
! call merdge(microlist(:,:,:,smallerlist), idx_microlist(:,smallerlist), N_microlist(smallerlist), microlist(:,:,:,0), idx_microlist(:,0), N_microlist(0))
end if
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_microlist(smallerlist)+N_microlist(0),n_st,microlist_zero(:,:,:),idx_microlist_zero(:),N_microlist(smallerlist)+N_microlist(0))
else
if(is_connected_to(buffer(1,1,i), miniList_gen, Nint, N_minilist_gen)) then if(is_connected_to(buffer(1,1,i), miniList_gen, Nint, N_minilist_gen)) then
cycle cycle
end if end if
call pt2_$PERT(psi_det_generators(1,1,i_generator),buffer(1,1,i), fock_diag_tmp, & 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) c_pert,e_2_pert,H_pert_diag,Nint,N_minilist,n_st,minilist,idx_minilist,N_minilist)
end if
! 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 do k = 1,N_st
e_2_pert_buffer(k,i) = e_2_pert(k) e_2_pert_buffer(k,i) = e_2_pert(k)
@ -74,11 +147,11 @@ subroutine perturb_buffer_$PERT(i_generator,buffer,buffer_size,e_2_pert_buffer,c
enddo enddo
deallocate( minilist, minilist_gen, idx_minilist ) deallocate( minilist, minilist_gen, idx_minilist )
deallocate( microlist, idx_microlist, N_microlist,ptr_microlist )
deallocate( microlist_gen, idx_microlist_gen,N_microlist_gen,ptr_microlist_gen )
end 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,key_mask,fock_diag_tmp) 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 implicit none
BEGIN_DOC BEGIN_DOC

13
plugins/Perturbation/pt2_equations.irp.f
View File

@ -123,8 +123,8 @@ subroutine pt2_moller_plesset ($arguments)
call get_excitation(ref_bitmask,det_pert,exc,degree,phase,Nint) call get_excitation(ref_bitmask,det_pert,exc,degree,phase,Nint)
if (degree == 2) then if (degree == 2) then
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2) call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
delta_e = Fock_matrix_diag_mo(h1) + Fock_matrix_diag_mo(h2) - & delta_e = (Fock_matrix_diag_mo(h1) - Fock_matrix_diag_mo(p1)) + &
(Fock_matrix_diag_mo(p1) + Fock_matrix_diag_mo(p2)) (Fock_matrix_diag_mo(h2) - Fock_matrix_diag_mo(p2))
delta_e = 1.d0/delta_e delta_e = 1.d0/delta_e
else if (degree == 1) then else if (degree == 1) then
call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2) call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)
@ -134,9 +134,14 @@ subroutine pt2_moller_plesset ($arguments)
delta_e = 0.d0 delta_e = 0.d0
endif endif
call i_H_psi(det_pert,psi_selectors,psi_selectors_coef,Nint,N_det,psi_selectors_size,n_st,i_H_psi_array) if (delta_e /= 0.d0) then
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) h = diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)
do i =1,n_st else
i_H_psi_array(:) = 0.d0
h = 0.d0
endif
do i =1,N_st
H_pert_diag(i) = h H_pert_diag(i) = h
c_pert(i) = i_H_psi_array(i) *delta_e c_pert(i) = i_H_psi_array(i) *delta_e
e_2_pert(i) = c_pert(i) * i_H_psi_array(i) e_2_pert(i) = c_pert(i) * i_H_psi_array(i)

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

@ -20,8 +20,6 @@ subroutine fill_H_apply_buffer_selection(n_selected,det_buffer,e_2_pert_buffer,c
ASSERT (N_int == N_int) ASSERT (N_int == N_int)
ASSERT (N_selected >= 0) ASSERT (N_selected >= 0)
call omp_set_lock(H_apply_buffer_lock(1,iproc)) call omp_set_lock(H_apply_buffer_lock(1,iproc))
smax = selection_criterion
smin = selection_criterion_min
new_size = H_apply_buffer(iproc)%N_det + n_selected new_size = H_apply_buffer(iproc)%N_det + n_selected
if (new_size > h_apply_buffer(iproc)%sze) then if (new_size > h_apply_buffer(iproc)%sze) then
@ -41,8 +39,6 @@ subroutine fill_H_apply_buffer_selection(n_selected,det_buffer,e_2_pert_buffer,c
select_max_out = max(select_max_out,s) select_max_out = max(select_max_out,s)
enddo enddo
if (is_selected) then if (is_selected) then
l = l+1 l = l+1
do j=1,N_int do j=1,N_int
@ -55,8 +51,6 @@ subroutine fill_H_apply_buffer_selection(n_selected,det_buffer,e_2_pert_buffer,c
enddo enddo
ASSERT (sum(popcnt(h_apply_buffer(iproc)%det(:,1,l)) )== elec_alpha_num) ASSERT (sum(popcnt(h_apply_buffer(iproc)%det(:,1,l)) )== elec_alpha_num)
ASSERT (sum(popcnt(h_apply_buffer(iproc)%det(:,2,l))) == elec_beta_num) ASSERT (sum(popcnt(h_apply_buffer(iproc)%det(:,2,l))) == elec_beta_num)
smax = max(s,smax)
smin = min(selection_criterion_min,smin)
endif endif
enddo enddo
H_apply_buffer(iproc)%N_det = l H_apply_buffer(iproc)%N_det = l
@ -65,10 +59,6 @@ subroutine fill_H_apply_buffer_selection(n_selected,det_buffer,e_2_pert_buffer,c
ASSERT (sum(popcnt(h_apply_buffer(iproc)%det(:,2,i))) == elec_beta_num) ASSERT (sum(popcnt(h_apply_buffer(iproc)%det(:,2,i))) == elec_beta_num)
enddo enddo
call omp_unset_lock(H_apply_buffer_lock(1,iproc)) call omp_unset_lock(H_apply_buffer_lock(1,iproc))
!$OMP CRITICAL
selection_criterion = max(selection_criterion,smax)
selection_criterion_min = min(selection_criterion_min,smin)
!$OMP END CRITICAL
end end
BEGIN_PROVIDER [ double precision, selection_criterion ] BEGIN_PROVIDER [ double precision, selection_criterion ]

1
plugins/QmcChem/.gitignore vendored
View File

@ -16,6 +16,7 @@ Makefile.depend
Nuclei Nuclei
Pseudo Pseudo
Utils Utils
ZMQ
ezfio_interface.irp.f ezfio_interface.irp.f
irpf90.make irpf90.make
irpf90_entities irpf90_entities

4
scripts/compilation/qp_create_ninja.py
View File

@ -637,7 +637,7 @@ def ninja_binaries_rule():
# c m d # # c m d #
# ~#~#~ # # ~#~#~ #
l_cmd = ["cd $module_abs/IRPF90_temp", "ninja $out && touch $out"] l_cmd = ["cd $module_abs/IRPF90_temp", "ninja $out && for i in $out ; do [ -x $$i ] && touch $$i ; done"]
# ~#~#~#~#~#~ # # ~#~#~#~#~#~ #
# s t r i n g # # s t r i n g #
@ -911,7 +911,7 @@ if __name__ == "__main__":
if module not in d_binaries: if module not in d_binaries:
l_msg = ["{0} is a root module but does not contain a main file.", l_msg = ["{0} is a root module but does not contain a main file.",
"- Create it in {0}", "- Create it in {0}",
"- Or delete {0} `qp_install_module.py uninstall {0}`", "- Or delete {0} `qp_module.py uninstall {0}`",
"- Or install a module that needs {0} with a main "] "- Or install a module that needs {0} with a main "]
print "\n".join(l_msg).format(module.rel) print "\n".join(l_msg).format(module.rel)

9
scripts/generate_h_apply.py
View File

@ -99,7 +99,8 @@ class H_apply(object):
deallocate(H_jj,iorder) deallocate(H_jj,iorder)
""" """
s["size_max"] = "2048" s["size_max"] = "8192"
s["copy_buffer"] = """call copy_H_apply_buffer_to_wf s["copy_buffer"] = """call copy_H_apply_buffer_to_wf
if (s2_eig) then if (s2_eig) then
call make_s2_eigenfunction call make_s2_eigenfunction
@ -198,7 +199,9 @@ class H_apply(object):
!$ call omp_unset_lock(lck) !$ call omp_unset_lock(lck)
deallocate (e_2_pert_buffer, coef_pert_buffer) deallocate (e_2_pert_buffer, coef_pert_buffer)
""" """
self.data["size_max"] = "2048"
self.data["size_max"] = "8192"
self.data["initialization"] = """ self.data["initialization"] = """
PROVIDE 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
""" """
@ -265,7 +268,7 @@ class H_apply(object):
double precision, intent(inout) :: select_max_out""" double precision, intent(inout) :: select_max_out"""
self.data["params_post"] += ", select_max(min(i_generator,size(select_max,1)))" self.data["params_post"] += ", select_max(min(i_generator,size(select_max,1)))"
self.data["size_max"] = "2048" self.data["size_max"] = "8192"
self.data["copy_buffer"] = """ self.data["copy_buffer"] = """
call copy_H_apply_buffer_to_wf call copy_H_apply_buffer_to_wf
if (s2_eig) then if (s2_eig) then

13
scripts/get_basis.sh
View File

@ -46,15 +46,4 @@ then
exit 1 exit 1
fi fi
pseudo="$1" ${EMSL_API_ROOT}/EMSL_api.py get_basis_data --treat_l --save --path="${tmpfile}" --basis="${basis}"
shift
if [[ -z $pseudo ]]
then
${EMSL_API_ROOT}/EMSL_api.py get_basis_data --treat_l --save --path="${tmpfile}" --basis="${basis}"
else
${EMSL_API_ROOT}/EMSL_api.py get_basis_data --save --path="${tmpfile}" --basis="${basis}" --db_path="${EMSL_API_ROOT}/db/Pseudo.db"
# echo ${EMSL_API_ROOT}/EMSL_api.py get_basis_data --save --path="${tmpfile}" --basis="${basis}" --db_path="${EMSL_API_ROOT}/db/Pseudo.db" 1>&2
# echo $PWD/BASIS
fi

12
scripts/module/qp_install_module.py → scripts/module/qp_module.py
View File

@ -2,13 +2,11 @@
# -*- coding: utf-8 -*- # -*- coding: utf-8 -*-
""" """
Usage: Usage:
qp_install_module.py create -n <name> [<children_modules>...] qp_module.py create -n <name> [<children_modules>...]
qp_install_module.py download -n <name> [<path_folder>...] qp_module.py download -n <name> [<path_folder>...]
qp_install_module.py install <name>... qp_module.py install <name>...
qp_install_module.py list (--installed | --available-local) qp_module.py list (--installed | --available-local)
qp_install_module.py uninstall <name>... qp_module.py uninstall <name>
Options: Options:
list: List all the modules available list: List all the modules available
create: Create a new module create: Create a new module

46
scripts/pseudo/put_pseudo_in_ezfio.py
View File

@ -4,7 +4,7 @@
Create the pseudo potential for a given atom Create the pseudo potential for a given atom
Usage: Usage:
put_pseudo_in_ezfio.py <ezfio_path> put_pseudo_in_ezfio.py <ezfio_path> <pseudo_name> [<db_dump_path>]
Help: Help:
atom is the Abreviation of the atom atom is the Abreviation of the atom
@ -28,7 +28,7 @@ import re
p = re.compile(ur'\|(\d+)><\d+\|') p = re.compile(ur'\|(\d+)><\d+\|')
def get_pseudo_str(l_atom): def get_pseudo_str(db_dump_path,pseudo_name,l_atom):
""" """
Run EMSL_local for getting the str of the pseudo potential Run EMSL_local for getting the str of the pseudo potential
@ -53,7 +53,6 @@ def get_pseudo_str(l_atom):
EMSL_root = "{0}/install/emsl/".format(qpackage_root) EMSL_root = "{0}/install/emsl/".format(qpackage_root)
EMSL_path = "{0}/EMSL_api.py".format(EMSL_root) EMSL_path = "{0}/EMSL_api.py".format(EMSL_root)
db_path = "{0}/db/Pseudo.db".format(EMSL_root)
str_ = "" str_ = ""
@ -63,8 +62,8 @@ def get_pseudo_str(l_atom):
l_cmd_atom = ["--atom", a] l_cmd_atom = ["--atom", a]
l_cmd_head = [EMSL_path, "get_basis_data", l_cmd_head = [EMSL_path, "get_basis_data",
"--db_path", db_path, "--db_dump_path", db_dump_path,
"--basis", "BFD-Pseudo"] "--basis", pseudo_name]
process = Popen(l_cmd_head + l_cmd_atom, stdout=PIPE, stderr=PIPE) process = Popen(l_cmd_head + l_cmd_atom, stdout=PIPE, stderr=PIPE)
@ -180,19 +179,22 @@ def get_zeff_alpha_beta(str_ele):
# #
from elts_num_ele import name_to_elec from elts_num_ele import name_to_elec
from math import ceil, floor
z = name_to_elec[name] z = name_to_elec[name]
z_eff = z - z_remove z_eff = z - z_remove
alpha = (z_remove / 2) alpha = int(ceil(z_remove / 2.))
beta = (z_remove / 2) beta = int(floor(z_remove / 2.))
# Remove more alpha, than beta
# _ # _
# |_) _ _|_ ._ ._ # |_) _ _|_ ._ ._
# | \ (/_ |_ |_| | | | # | \ (/_ |_ |_| | | |
# #
return [z_eff, alpha, beta] return [z_remove, z_eff, alpha, beta]
def add_zero(array, size, type): def add_zero(array, size, type):
@ -224,6 +226,12 @@ def make_it_square(matrix, dim, type=float):
return matrix return matrix
def full_path(path):
path = os.path.expanduser(path)
path = os.path.expandvars(path)
path = os.path.abspath(path)
return path
if __name__ == "__main__": if __name__ == "__main__":
arguments = docopt(__doc__) arguments = docopt(__doc__)
# ___ # ___
@ -235,19 +243,22 @@ if __name__ == "__main__":
# E Z F I O # # E Z F I O #
# ~#~#~#~#~ # # ~#~#~#~#~ #
ezfio_path = arguments["<ezfio_path>"] ezfio_path = full_path(arguments["<ezfio_path>"])
ezfio_path = os.path.expanduser(ezfio_path)
ezfio_path = os.path.expandvars(ezfio_path)
ezfio_path = os.path.abspath(ezfio_path)
ezfio.set_file("{0}".format(ezfio_path)) ezfio.set_file("{0}".format(ezfio_path))
# ~#~#~#~#~#~#~#~#~#~#~ # # ~#~#~#~#~#~#~#~#~#~#~ #
# P s e u d o _ d a t a # # P s e u d o _ d a t a #
# ~#~#~#~#~#~#~#~#~#~#~ # # ~#~#~#~#~#~#~#~#~#~#~ #
if arguments["<db_dump_path>"]:
db_dump_path = full_path(arguments["<db_dump_path>"])
else:
db_dump_path= full_path("{0}/data/BFD-Pseudo.dump".format(qpackage_root))
pseudo_name = arguments["<pseudo_name>"]
l_ele = ezfio.get_nuclei_nucl_label() l_ele = ezfio.get_nuclei_nucl_label()
str_ = get_pseudo_str(l_ele)
str_ = get_pseudo_str(db_dump_path,pseudo_name,l_ele)
# _ # _
# |_) _. ._ _ _ # |_) _. ._ _ _
@ -257,7 +268,7 @@ if __name__ == "__main__":
l_str_ele = [str_ele for str_ele in str_.split("Element Symbol: ") l_str_ele = [str_ele for str_ele in str_.split("Element Symbol: ")
if str_ele] if str_ele]
for i in "l_zeff v_k n_k dz_k v_kl n_kl dz_kl".split(): for i in "l_zeff l_remove v_k n_k dz_k v_kl n_kl dz_kl".split():
exec("{0} = []".format(i)) exec("{0} = []".format(i))
alpha_tot = 0 alpha_tot = 0
@ -296,11 +307,13 @@ if __name__ == "__main__":
# Z _ e f f , a l p h a / b e t a _ e l e c # # Z _ e f f , a l p h a / b e t a _ e l e c #
# ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ # # ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ #
zeff, alpha, beta = get_zeff_alpha_beta(str_ele) zremove, zeff, alpha, beta = get_zeff_alpha_beta(str_ele)
alpha_tot += alpha alpha_tot += alpha
beta_tot += beta beta_tot += beta
l_zeff.append(zeff) l_zeff.append(zeff)
l_remove.append(zremove)
# _ # _
# /\ _| _| _|_ _ _ _ _|_ o _ # /\ _| _| _|_ _ _ _ _|_ o _
# /--\ (_| (_| |_ (_) (/_ /_ | | (_) # /--\ (_| (_| |_ (_) (/_ /_ | | (_)
@ -311,6 +324,7 @@ if __name__ == "__main__":
# ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ # # ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ #
ezfio.nuclei_nucl_charge = l_zeff ezfio.nuclei_nucl_charge = l_zeff
ezfio.pseudo_nucl_charge_remove = l_remove
alpha_tot = ezfio.get_electrons_elec_alpha_num() - alpha_tot alpha_tot = ezfio.get_electrons_elec_alpha_num() - alpha_tot
beta_tot = ezfio.get_electrons_elec_beta_num() - beta_tot beta_tot = ezfio.get_electrons_elec_beta_num() - beta_tot

232
scripts/qp_convert_qmcpack_from_ezfio.py Executable file
View File

@ -0,0 +1,232 @@
#!/usr/bin/python
print "#QP -> QMCPACK"
from ezfio import ezfio
import sys
ezfio_path = sys.argv[1]
ezfio.set_file(ezfio_path)
do_pseudo = ezfio.get_pseudo_do_pseudo()
if do_pseudo:
print "do_pseudo True"
zcore = ezfio.get_pseudo_nucl_charge_remove()
else:
print "do_pseudo False"
n_det =ezfio.get_determinants_n_det()
if n_det == 1:
print "multi_det False"
else:
print "multi_det True"
ao_num = ezfio.get_ao_basis_ao_num()
print "ao_num", ao_num
mo_num = ezfio.get_mo_basis_mo_tot_num()
print "mo_num", mo_num
alpha = ezfio.get_electrons_elec_alpha_num()
beta = ezfio.get_electrons_elec_beta_num()
print "elec_alpha_num", alpha
print "elec_beta_num", beta
print "elec_tot_num", alpha + beta
print "spin_multiplicity", 2*(alpha-beta)+1
l_label = ezfio.get_nuclei_nucl_label()
l_charge = ezfio.get_nuclei_nucl_charge()
l_coord = ezfio.get_nuclei_nucl_coord()
l_coord_str = [" ".join(map(str,i)) for i in l_coord]
print "nucl_num",len(l_label)
print "Atomic coord in Bohr"
for i,t in enumerate(zip(l_label,l_charge,l_coord_str)):
try :
l = (t[0],t[1]+zcore[i],t[1])
except NameError:
l = t
print " ".join(map(str,l))
import subprocess
process = subprocess.Popen(['qp_print_basis', ezfio_path], stdout=subprocess.PIPE)
out, err = process.communicate()
basis_raw, sym_raw, mo_raw = out.split("\n\n\n")
basis_without_header = "\n".join(basis_raw.split("\n")[7:])
for i,l in enumerate(l_label):
basis_without_header=basis_without_header.replace('Atom {0}'.format(i+1),l)
print "BEGIN_BASIS_SET"
print ""
print basis_without_header
print "END_BASIS_SET"
# _
# |\/| / \ _
# | | \_/ _>
#
def same_character(item1):
return item1==item1[0]* len(item1)
def compare_gamess_style(item1, item2):
if len(item1) < len(item2):
return -1
elif len(item1) > len(item2):
return 1
elif same_character(item1) and same_character(item2):
if item1 < item2:
return -1
else:
return 1
elif same_character(item1) and not same_character(item2):
return -1
elif not same_character(item1) and same_character(item2):
return 1
else:
return compare_gamess_style(item1[:-1],item2[:-1])
def expend_and_order_sym(str_):
#Expend
for i,c in enumerate(str_):
try:
n = int(c)
except ValueError:
pass
else:
str_ = str_[:i-1] + str_[i-1]*n + str_[i+1:]
#Order by frequency
return "".join(sorted(str_,key=str_.count,reverse=True))
def get_nb_permutation(str_):
l = len(str_)-1
if l==0:
return 1
else:
return 2*(2*l + 1)
## We will order the symetry
l_sym_without_header = sym_raw.split("\n")[3:-2]
l_l_sym = [i.split() for i in l_sym_without_header]
for l in l_l_sym:
l[2] = expend_and_order_sym(l[2])
l_l_sym_iter = iter(l_l_sym)
for i,l in enumerate(l_l_sym_iter):
n = get_nb_permutation(l[2])
if n !=1:
l_l_sym[i:i+n] = sorted(l_l_sym[i:i+n],key=lambda x : x[2], cmp=compare_gamess_style)
for next_ in range(n-1):
next(l_l_sym_iter)
#Is orderd now
l_block = mo_raw.split("\n\n")[5:-1]
l_block_format=[]
print ""
print "BEGIN_MO"
for block in l_block:
print ""
l_ligne = block.split("\n")
print l_ligne.pop(0)
for l in l_l_sym:
i = int(l[0]) - 1
i_a = int(l[1]) - 1
sym = l[2]
print l_label[i_a],sym,l_ligne[i]
print "END_MO"
if do_pseudo:
print ""
print "BEGIN_PSEUDO"
klocmax = ezfio.get_pseudo_pseudo_klocmax()
kmax = ezfio.get_pseudo_pseudo_kmax()
lmax = ezfio.get_pseudo_pseudo_lmax()
n_k = ezfio.get_pseudo_pseudo_n_k()
v_k = ezfio.get_pseudo_pseudo_v_k()
dz_k = ezfio.get_pseudo_pseudo_dz_k()
n_kl = ezfio.get_pseudo_pseudo_n_kl()
v_kl = ezfio.get_pseudo_pseudo_v_kl()
dz_kl = ezfio.get_pseudo_pseudo_dz_kl()
def list_to_string(l):
return " ".join(map(str,l))
for i,a in enumerate(l_label):
l_str = []
l_dump = []
for k in range(klocmax):
if v_k[k][i]:
l_ = list_to_string([v_k[k][i], n_k[k][i]+2, dz_k[k][i]])
l_dump.append(l_)
l_str.append(l_dump)
for l in range(lmax+1):
l_dump = []
for k in range(kmax):
if v_kl[l][k][i]:
l_ = list_to_string([v_kl[l][k][i], n_kl[l][k][i]+2, dz_kl[l][k][i]])
l_dump.append(l_)
if l_dump:
l_str.append(l_dump)
str_ = "PARAMETERS FOR {0} ON ATOM {1} WITH ZCORE {2} AND LMAX {3} ARE"
print str_.format(a,i+1,zcore[i],len(l_str))
for i, l in enumerate(l_str):
str_ = "FOR L= {0} COEFF N ZETA"
print str_.format(len(l_str)-i-1)
for ii, ll in enumerate(l):
print " ",ii+1, ll
str_ = "THE ECP RUN REMOVES {0} CORE ELECTRONS, AND THE SAME NUMBER OF PROTONS."
print str_.format(sum(zcore))
print "END_PSEUDO"
print ""
print "BEGIN_DET"
print ""
print "mo_num", mo_num
print "det_num", n_det
print ""
psi_det = ezfio.get_determinants_psi_det()
psi_coef = ezfio.get_determinants_psi_coef()[0]
for c, (l_det_bit_alpha, l_det_bit_beta) in zip(psi_coef,psi_det):
print c
for det in l_det_bit_alpha:
bin_det_raw = "{0:b}".format(det)[::-1]
bin_det = bin_det_raw+"0"*(mo_num-len(bin_det_raw))
print bin_det
for det in l_det_bit_beta:
bin_det_raw = "{0:b}".format(det)[::-1]
bin_det = bin_det_raw+"0"*(mo_num-len(bin_det_raw))
print bin_det
print ""
print "END_DET"

2
src/AO_Basis/README.rst
View File

@ -117,7 +117,7 @@ Documentation
:math:`\int \chi_i(r) \chi_j(r) dr)` :math:`\int \chi_i(r) \chi_j(r) dr)`
`ao_overlap_abs <http://github.com/LCPQ/quantum_package/tree/master/src/AO_Basis/ao_overlap.irp.f#L65>`_ `ao_overlap_abs <http://github.com/LCPQ/quantum_package/tree/master/src/AO_Basis/ao_overlap.irp.f#L66>`_
Overlap between absolute value of atomic basis functions: Overlap between absolute value of atomic basis functions:
:math:`\int |\chi_i(r)| |\chi_j(r)| dr)` :math:`\int |\chi_i(r)| |\chi_j(r)| dr)`

16
src/AO_Basis/aos.irp.f
View File

@ -24,13 +24,14 @@ BEGIN_PROVIDER [ double precision, ao_coef_normalized, (ao_num_align,ao_prim_num
BEGIN_DOC BEGIN_DOC
! Coefficients including the AO normalization ! Coefficients including the AO normalization
END_DOC END_DOC
double precision :: norm, norm2,overlap_x,overlap_y,overlap_z,C_A(3) double precision :: norm, norm2,overlap_x,overlap_y,overlap_z,C_A(3), c
integer :: l, powA(3), nz integer :: l, powA(3), nz
integer :: i,j integer :: i,j,k
nz=100 nz=100
C_A(1) = 0.d0 C_A(1) = 0.d0
C_A(2) = 0.d0 C_A(2) = 0.d0
C_A(3) = 0.d0 C_A(3) = 0.d0
ao_coef_normalized = 0.d0
do i=1,ao_num do i=1,ao_num
powA(1) = ao_power(i,1) powA(1) = ao_power(i,1)
powA(2) = ao_power(i,2) powA(2) = ao_power(i,2)
@ -39,6 +40,17 @@ BEGIN_PROVIDER [ double precision, ao_coef_normalized, (ao_num_align,ao_prim_num
call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j),powA,powA,overlap_x,overlap_y,overlap_z,norm,nz) call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j),powA,powA,overlap_x,overlap_y,overlap_z,norm,nz)
ao_coef_normalized(i,j) = ao_coef(i,j)/sqrt(norm) ao_coef_normalized(i,j) = ao_coef(i,j)/sqrt(norm)
enddo enddo
! Normalization of the contracted basis functions
norm = 0.d0
do j=1,ao_prim_num(i)
do k=1,ao_prim_num(i)
call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,k),powA,powA,overlap_x,overlap_y,overlap_z,c,nz)
norm = norm+c*ao_coef_normalized(i,j)*ao_coef_normalized(i,k)
enddo
enddo
do j=1,ao_prim_num(i)
ao_coef_normalized(i,j) = ao_coef_normalized(i,j)/sqrt(norm)
enddo
enddo enddo
END_PROVIDER END_PROVIDER

14
src/Bitmask/README.rst
View File

@ -80,7 +80,7 @@ Documentation
Bitmask to include all possible single excitations from Hartree-Fock Bitmask to include all possible single excitations from Hartree-Fock
`core_bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L350>`_ `core_bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L344>`_
Reunion of the inactive, active and virtual bitmasks Reunion of the inactive, active and virtual bitmasks
@ -142,7 +142,7 @@ Documentation
Hartree Fock bit mask Hartree Fock bit mask
`i_bitmask_gen <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L364>`_ `i_bitmask_gen <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L358>`_
Current bitmask for the generators Current bitmask for the generators
@ -150,7 +150,7 @@ Documentation
Bitmasks for the inactive orbitals that are excited in post CAS method Bitmasks for the inactive orbitals that are excited in post CAS method
`inact_virt_bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L338>`_ `inact_virt_bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L332>`_
Reunion of the inactive and virtual bitmasks Reunion of the inactive and virtual bitmasks
@ -158,7 +158,7 @@ Documentation
Undocumented Undocumented
`list_inact <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L304>`_ `list_inact <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L298>`_
Undocumented Undocumented
@ -167,7 +167,7 @@ Documentation
occupations "list(N_int*bit_kind_size,2) occupations "list(N_int*bit_kind_size,2)
`list_virt <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L305>`_ `list_virt <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L299>`_
Undocumented Undocumented
@ -219,11 +219,11 @@ Documentation
Reference bit mask, used in Slater rules, chosen as Hartree-Fock bitmask Reference bit mask, used in Slater rules, chosen as Hartree-Fock bitmask
`reunion_of_bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L325>`_ `reunion_of_bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L319>`_
Reunion of the inactive, active and virtual bitmasks Reunion of the inactive, active and virtual bitmasks
`unpaired_alpha_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L373>`_ `unpaired_alpha_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks.irp.f#L367>`_
Bitmask reprenting the unpaired alpha electrons in the HF_bitmask Bitmask reprenting the unpaired alpha electrons in the HF_bitmask

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

@ -141,6 +141,19 @@ BEGIN_PROVIDER [ integer(bit_kind), generators_bitmask_restart, (N_int,2,6,N_gen
enddo enddo
endif endif
integer :: i
do k=1,N_generators_bitmask
do ispin=1,2
do i=1,N_int
generators_bitmask_restart(i,ispin,s_hole ,k) = iand(full_ijkl_bitmask(i,d_hole1),generators_bitmask_restart(i,ispin,s_hole,k) )
generators_bitmask_restart(i,ispin,s_part ,k) = iand(full_ijkl_bitmask(i,d_part1),generators_bitmask_restart(i,ispin,s_part,k) )
generators_bitmask_restart(i,ispin,d_hole1,k) = iand(full_ijkl_bitmask(i,d_hole1),generators_bitmask_restart(i,ispin,d_hole1,k) )
generators_bitmask_restart(i,ispin,d_part1,k) = iand(full_ijkl_bitmask(i,d_part1),generators_bitmask_restart(i,ispin,d_part1,k) )
generators_bitmask_restart(i,ispin,d_hole2,k) = iand(full_ijkl_bitmask(i,d_hole2),generators_bitmask_restart(i,ispin,d_hole2,k) )
generators_bitmask_restart(i,ispin,d_part2,k) = iand(full_ijkl_bitmask(i,d_part2),generators_bitmask_restart(i,ispin,d_part2,k) )
enddo
enddo
enddo
END_PROVIDER END_PROVIDER
@ -172,7 +185,7 @@ BEGIN_PROVIDER [ integer(bit_kind), generators_bitmask, (N_int,2,6,N_generators_
if (exists) then if (exists) then
call ezfio_get_bitmasks_generators(generators_bitmask) call ezfio_get_bitmasks_generators(generators_bitmask)
else else
integer :: k, ispin integer :: k, ispin, i
do k=1,N_generators_bitmask do k=1,N_generators_bitmask
do ispin=1,2 do ispin=1,2
generators_bitmask(:,ispin,s_hole ,k) = full_ijkl_bitmask(:,d_hole1) generators_bitmask(:,ispin,s_hole ,k) = full_ijkl_bitmask(:,d_hole1)
@ -185,6 +198,18 @@ BEGIN_PROVIDER [ integer(bit_kind), generators_bitmask, (N_int,2,6,N_generators_
enddo enddo
endif endif
do k=1,N_generators_bitmask
do ispin=1,2
do i=1,N_int
generators_bitmask(i,ispin,s_hole ,k) = iand(full_ijkl_bitmask(i,d_hole1),generators_bitmask(i,ispin,s_hole,k) )
generators_bitmask(i,ispin,s_part ,k) = iand(full_ijkl_bitmask(i,d_part1),generators_bitmask(i,ispin,s_part,k) )
generators_bitmask(i,ispin,d_hole1,k) = iand(full_ijkl_bitmask(i,d_hole1),generators_bitmask(i,ispin,d_hole1,k) )
generators_bitmask(i,ispin,d_part1,k) = iand(full_ijkl_bitmask(i,d_part1),generators_bitmask(i,ispin,d_part1,k) )
generators_bitmask(i,ispin,d_hole2,k) = iand(full_ijkl_bitmask(i,d_hole2),generators_bitmask(i,ispin,d_hole2,k) )
generators_bitmask(i,ispin,d_part2,k) = iand(full_ijkl_bitmask(i,d_part2),generators_bitmask(i,ispin,d_part2,k) )
enddo
enddo
enddo
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ integer, N_cas_bitmask ] BEGIN_PROVIDER [ integer, N_cas_bitmask ]
@ -223,7 +248,7 @@ BEGIN_PROVIDER [ integer(bit_kind), cas_bitmask, (N_int,2,N_cas_bitmask) ]
! Bitmasks for CAS reference determinants. (N_int, alpha/beta, CAS reference) ! Bitmasks for CAS reference determinants. (N_int, alpha/beta, CAS reference)
END_DOC END_DOC
logical :: exists logical :: exists
integer :: i,i_part,i_gen,j integer :: i,i_part,i_gen,j,k
PROVIDE ezfio_filename PROVIDE ezfio_filename
call ezfio_has_bitmasks_cas(exists) call ezfio_has_bitmasks_cas(exists)
@ -240,14 +265,21 @@ BEGIN_PROVIDER [ integer(bit_kind), cas_bitmask, (N_int,2,N_cas_bitmask) ]
else else
i_part = 2 i_part = 2
i_gen = 1 i_gen = 1
do j = 1, N_cas_bitmask do j=1, N_cas_bitmask
do i = 1, N_int do i=1, N_int
cas_bitmask(i,1,j) = generators_bitmask_restart(i,1,i_part,i_gen) cas_bitmask(i,1,j) = generators_bitmask_restart(i,1,i_part,i_gen)
cas_bitmask(i,2,j) = generators_bitmask_restart(i,2,i_part,i_gen) cas_bitmask(i,2,j) = generators_bitmask_restart(i,2,i_part,i_gen)
enddo enddo
enddo enddo
endif endif
endif endif
do i=1,N_cas_bitmask
do j = 1, N_cas_bitmask
do k=1,N_int
cas_bitmask(k,j,i) = iand(cas_bitmask(k,j,i),full_ijkl_bitmask(k,j))
enddo
enddo
enddo
END_PROVIDER END_PROVIDER

2
src/Determinants/EZFIO.cfg
View File

@ -44,7 +44,7 @@ default: False
type: Threshold type: Threshold
doc: Thresholds of Davidson's algorithm doc: Thresholds of Davidson's algorithm
interface: ezfio,provider,ocaml interface: ezfio,provider,ocaml
default: 1.e-8 default: 1.e-12
[threshold_generators] [threshold_generators]
type: Threshold type: Threshold

18
src/Determinants/H_apply.template.f
View File

@ -97,24 +97,27 @@ end subroutine
subroutine $subroutine_diexcP(key_in, fs1, fh1, particl_1, fs2, fh2, particl_2, fock_diag_tmp, i_generator, iproc_in $parameters ) subroutine $subroutine_diexcP(key_in, fs1, fh1, particl_1, fs2, fh2, particl_2, fock_diag_tmp, i_generator, iproc_in $parameters )
implicit none
integer(bit_kind), intent(in) :: key_in(N_int, 2), particl_1(N_int, 2), particl_2(N_int, 2) integer(bit_kind), intent(in) :: key_in(N_int, 2), particl_1(N_int, 2), particl_2(N_int, 2)
double precision, intent(in) :: fock_diag_tmp(2,mo_tot_num+1) double precision, intent(in) :: fock_diag_tmp(2,mo_tot_num+1)
integer(bit_kind) :: p1_mask(N_int, 2), p2_mask(N_int, 2), key_mask(N_int, 2) integer(bit_kind) :: p1_mask(N_int, 2), p2_mask(N_int, 2), key_mask(N_int, 2)
integer,intent(in) :: fh1,fh2,fs1,fs2,i_generator,iproc_in integer,intent(in) :: fs1,fs2,i_generator,iproc_in, fh1,fh2
integer(bit_kind) :: miniList(N_int, 2, N_det) integer(bit_kind) :: miniList(N_int, 2, N_det)
integer :: n_minilist, n_alpha, n_beta, deg(2), i, ni integer :: n_minilist, n_alpha, n_beta, deg(2), i, ni
$declarations $declarations
integer(bit_kind), parameter :: one = 1_bit_kind
p1_mask(:,:) = 0_bit_kind p1_mask(:,:) = 0_bit_kind
p2_mask(:,:) = 0_bit_kind p2_mask(:,:) = 0_bit_kind
p1_mask(ishft(fh1-1,-bit_kind_shift) + 1, fs1) = ishft(1_8,iand(fh1-1,bit_kind_size-1))
p2_mask(ishft(fh2-1,-bit_kind_shift) + 1, fs2) = ishft(1_8,iand(fh2-1,bit_kind_size-1)) p1_mask(ishft(fh1-1,-bit_kind_shift) + 1, fs1) = ishft(one,iand(fh1-1,bit_kind_size-1))
p2_mask(ishft(fh2-1,-bit_kind_shift) + 1, fs2) = ishft(one,iand(fh2-1,bit_kind_size-1))
key_mask(:,:) = key_in(:,:) key_mask(:,:) = key_in(:,:)
key_mask(ishft(fh1-1,-bit_kind_shift) + 1, fs1) -= ishft(1_8,iand(fh1-1,bit_kind_size-1)) key_mask(ishft(fh1-1,-bit_kind_shift) + 1, fs1) -= ishft(one,iand(fh1-1,bit_kind_size-1))
key_mask(ishft(fh2-1,-bit_kind_shift) + 1, fs2) -= ishft(1_8,iand(fh2-1,bit_kind_size-1)) key_mask(ishft(fh2-1,-bit_kind_shift) + 1, fs2) -= ishft(one,iand(fh2-1,bit_kind_size-1))
call $subroutine_diexcOrg(key_in, key_mask, p1_mask, particl_1, p2_mask, particl_2, fock_diag_tmp, i_generator, iproc_in $parameters ) call $subroutine_diexcOrg(key_in, key_mask, p1_mask, particl_1, p2_mask, particl_2, fock_diag_tmp, i_generator, iproc_in $parameters )
end subroutine end subroutine
@ -261,6 +264,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl
! Build array of the non-zero integrals of second excitation ! Build array of the non-zero integrals of second excitation
$filter_integrals $filter_integrals
if (ispin == 1) then if (ispin == 1) then
integer :: jjj integer :: jjj
@ -269,7 +273,7 @@ subroutine $subroutine_diexcOrg(key_in,key_mask,hole_1,particl_1,hole_2, particl
i_b = occ_hole_tmp(kk,other_spin) i_b = occ_hole_tmp(kk,other_spin)
ASSERT (i_b > 0) ASSERT (i_b > 0)
ASSERT (i_b <= mo_tot_num) ASSERT (i_b <= mo_tot_num)
do jjj=1,N_elec_in_key_part_2(other_spin) ! particule do jjj=1,N_elec_in_key_part_2(other_spin) ! particle
j_b = occ_particle_tmp(jjj,other_spin) j_b = occ_particle_tmp(jjj,other_spin)
ASSERT (j_b > 0) ASSERT (j_b > 0)
ASSERT (j_b <= mo_tot_num) ASSERT (j_b <= mo_tot_num)

2
src/Determinants/README.rst
View File

@ -316,7 +316,7 @@ Documentation
idx(0) is the number of determinants that interact with key1 idx(0) is the number of determinants that interact with key1
`filter_connected_i_h_psi0_sc2 <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/filter_connected.irp.f#L205>`_ `filter_connected_i_h_psi0_sc2 <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/filter_connected.irp.f#L197>`_
standard filter_connected_i_H_psi but returns in addition standard filter_connected_i_H_psi but returns in addition
.br .br
the array of the index of the non connected determinants to key1 the array of the index of the non connected determinants to key1

14
src/Determinants/density_matrix.irp.f
View File

@ -10,7 +10,7 @@
double precision :: ck, cl, ckl double precision :: ck, cl, ckl
double precision :: phase double precision :: phase
integer :: h1,h2,p1,p2,s1,s2, degree integer :: h1,h2,p1,p2,s1,s2, degree
integer :: exc(0:2,2,2),n_occ_alpha integer :: exc(0:2,2,2),n_occ(2)
double precision, allocatable :: tmp_a(:,:), tmp_b(:,:) double precision, allocatable :: tmp_a(:,:), tmp_b(:,:)
if(only_single_double_dm)then if(only_single_double_dm)then
@ -22,7 +22,7 @@
one_body_dm_mo_beta = 0.d0 one_body_dm_mo_beta = 0.d0
!$OMP PARALLEL DEFAULT(NONE) & !$OMP PARALLEL DEFAULT(NONE) &
!$OMP PRIVATE(j,k,l,m,occ,ck, cl, ckl,phase,h1,h2,p1,p2,s1,s2, degree,exc, & !$OMP PRIVATE(j,k,l,m,occ,ck, cl, ckl,phase,h1,h2,p1,p2,s1,s2, degree,exc, &
!$OMP tmp_a, tmp_b, n_occ_alpha)& !$OMP tmp_a, tmp_b, n_occ)&
!$OMP SHARED(psi_det,psi_coef,N_int,N_states,state_average_weight,elec_alpha_num,& !$OMP SHARED(psi_det,psi_coef,N_int,N_states,state_average_weight,elec_alpha_num,&
!$OMP elec_beta_num,one_body_dm_mo_alpha,one_body_dm_mo_beta,N_det,mo_tot_num_align,& !$OMP elec_beta_num,one_body_dm_mo_alpha,one_body_dm_mo_beta,N_det,mo_tot_num_align,&
!$OMP mo_tot_num) !$OMP mo_tot_num)
@ -31,8 +31,7 @@
tmp_b = 0.d0 tmp_b = 0.d0
!$OMP DO SCHEDULE(dynamic) !$OMP DO SCHEDULE(dynamic)
do k=1,N_det do k=1,N_det
call bitstring_to_list(psi_det(1,1,k), occ(1,1), n_occ_alpha, N_int) call bitstring_to_list_ab(psi_det(1,1,k), occ, n_occ, N_int)
call bitstring_to_list(psi_det(1,2,k), occ(1,2), n_occ_alpha, N_int)
do m=1,N_states do m=1,N_states
ck = psi_coef(k,m)*psi_coef(k,m) * state_average_weight(m) ck = psi_coef(k,m)*psi_coef(k,m) * state_average_weight(m)
do l=1,elec_alpha_num do l=1,elec_alpha_num
@ -182,13 +181,10 @@ subroutine set_natural_mos
END_DOC END_DOC
character*(64) :: label character*(64) :: label
double precision, allocatable :: tmp(:,:) double precision, allocatable :: tmp(:,:)
allocate(tmp(size(one_body_dm_mo,1),size(one_body_dm_mo,2)))
! Negation to have the occupied MOs first after the diagonalization
tmp = one_body_dm_mo
label = "Natural" label = "Natural"
call mo_as_eigvectors_of_mo_matrix(tmp,size(tmp,1),size(tmp,2),label,-1) ! call mo_as_eigvectors_of_mo_matrix(one_body_dm_mo,size(one_body_dm_mo,1),mo_tot_num,label,-1)
deallocate(tmp) call mo_as_svd_vectors_of_mo_matrix(one_body_dm_mo,size(one_body_dm_mo,1),mo_tot_num,mo_tot_num,label)
end end
subroutine save_natural_mos subroutine save_natural_mos

124
src/Determinants/filter_connected.irp.f
View File

@ -98,6 +98,130 @@ subroutine filter_connected(key1,key2,Nint,sze,idx)
end end
subroutine getMobiles(key,key_mask, mobiles,Nint)
use bitmasks
integer(bit_kind),intent(in) :: key(Nint,2), key_mask(Nint,2)
integer,intent(out) :: mobiles(2)
integer,intent(in) :: Nint
integer(bit_kind) :: mobileMask(Nint,2)
integer :: list(Nint*bit_kind_size), nel
do j=1,Nint
mobileMask(j,1) = xor(key(j,1), key_mask(j,1))
mobileMask(j,2) = xor(key(j,2), key_mask(j,2))
end do
call bitstring_to_list(mobileMask(:,1), list(:), nel, Nint)
if(nel == 2) then
mobiles(1) = list(1)
mobiles(2) = list(2)
else if(nel == 1) then
mobiles(1) = list(1)
call bitstring_to_list(mobileMask(:,2), list(:), nel, Nint)
mobiles(2) = list(1) + mo_tot_num
else
call bitstring_to_list(mobileMask(:,2), list(:), nel, Nint)
mobiles(1) = list(1) + mo_tot_num
mobiles(2) = list(2) + mo_tot_num
end if
end subroutine
subroutine create_microlist(minilist, N_minilist, key_mask, microlist, idx_microlist, N_microlist, ptr_microlist, Nint)
use bitmasks
integer, intent(in) :: Nint, N_minilist
integer(bit_kind), intent(in) :: minilist(Nint,2,N_minilist), key_mask(Nint,2)
integer, intent(out) :: N_microlist(0:mo_tot_num*2), ptr_microlist(0:mo_tot_num*2+1), idx_microlist(N_minilist*4)
integer(bit_kind), intent(out) :: microlist(Nint,2,N_minilist*4)
integer :: i,j,k,nt,n_element(2)
integer :: list(Nint*bit_kind_size,2), cur_microlist(0:mo_tot_num*2+1)
integer(bit_kind) :: key_mask_neg(Nint,2), mobileMask(Nint,2)
do i=1,Nint
key_mask_neg(i,1) = not(key_mask(i,1))
key_mask_neg(i,2) = not(key_mask(i,2))
end do
N_microlist(:) = 0
do i=1, N_minilist
do j=1,Nint
mobileMask(j,1) = iand(key_mask_neg(j,1), minilist(j,1,i))
mobileMask(j,2) = iand(key_mask_neg(j,2), minilist(j,2,i))
end do
call bitstring_to_list(mobileMask(:,1), list(:,1), n_element(1), Nint)
call bitstring_to_list(mobileMask(:,2), list(:,2), n_element(2), Nint)
if(n_element(1) + n_element(2) /= 4) then
N_microlist(0) = N_microlist(0) + 1
else
do j=1,n_element(1)
nt = list(j,1)
N_microlist(nt) = N_microlist(nt) + 1
end do
do j=1,n_element(2)
nt = list(j,2) + mo_tot_num
N_microlist(nt) = N_microlist(nt) + 1
end do
end if
end do
ptr_microlist(0) = 1
do i=1,mo_tot_num*2+1
ptr_microlist(i) = ptr_microlist(i-1) + N_microlist(i-1)
end do
cur_microlist(:) = ptr_microlist(:)
do i=1, N_minilist
do j=1,Nint
mobileMask(j,1) = iand(key_mask_neg(j,1), minilist(j,1,i))
mobileMask(j,2) = iand(key_mask_neg(j,2), minilist(j,2,i))
end do
call bitstring_to_list(mobileMask(:,1), list(:,1), n_element(1), Nint)
call bitstring_to_list(mobileMask(:,2), list(:,2), n_element(2), Nint)
if(n_element(1) + n_element(2) /= 4) then
idx_microlist(cur_microlist(0)) = i
microlist(:,:,cur_microlist(0)) = minilist(:,:,i)
cur_microlist(0) = cur_microlist(0) + 1
else
do j=1,n_element(1)
nt = list(j,1)
idx_microlist(cur_microlist(nt)) = i
microlist(:,:,cur_microlist(nt)) = minilist(:,:,i)
cur_microlist(nt) = cur_microlist(nt) + 1
end do
do j=1,n_element(2)
nt = list(j,2) + mo_tot_num
idx_microlist(cur_microlist(nt)) = i
microlist(:,:,cur_microlist(nt)) = minilist(:,:,i)
cur_microlist(nt) = cur_microlist(nt) + 1
end do
end if
end do
end subroutine
subroutine merdge(mic, idx_mic, N_mic, mic0, idx_mic0, N_mic0, Nint)
use bitmasks
integer(bit_kind) :: mic(Nint,2,N_mic), mic0(Nint,2,*)
integer :: idx_mic(N_mic), idx_mic0(N_mic0), N_mic, N_mic0
mic0(:,:,N_mic0+1:N_mic0+N_mic) = mic(:,:,:)
idx_mic0(N_mic0+1:N_mic0+N_mic) = idx_mic(:)
end subroutine
subroutine filter_connected_i_H_psi0(key1,key2,Nint,sze,idx) subroutine filter_connected_i_H_psi0(key1,key2,Nint,sze,idx)
use bitmasks use bitmasks
BEGIN_DOC BEGIN_DOC

36
src/Ezfio_files/README.rst
View File

@ -203,6 +203,14 @@ output_bitmask
Output file for Bitmask Output file for Bitmask
output_cisd
Output file for CISD
output_cisd_selected
Output file for CISD_selected
`output_cpu_time_0 <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files/output.irp.f#L2>`_ `output_cpu_time_0 <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files/output.irp.f#L2>`_
Initial CPU and wall times when printing in the output files Initial CPU and wall times when printing in the output files
@ -219,6 +227,10 @@ output_ezfio_files
Output file for Ezfio_files Output file for Ezfio_files
output_fcidump
Output file for FCIdump
output_full_ci output_full_ci
Output file for Full_CI Output file for Full_CI
@ -247,12 +259,8 @@ output_moguess
Output file for MOGuess Output file for MOGuess
output_mrcc_cassd output_mp2
Output file for MRCC_CASSD Output file for MP2
output_mrcc_utils
Output file for MRCC_Utils
output_nuclei output_nuclei
@ -271,18 +279,14 @@ output_pseudo
Output file for Pseudo Output file for Pseudo
output_psiref_cas
Output file for Psiref_CAS
output_psiref_utils
Output file for Psiref_Utils
output_selectors_full output_selectors_full
Output file for Selectors_full Output file for Selectors_full
output_singlerefmethod
Output file for SingleRefMethod
output_utils output_utils
Output file for Utils Output file for Utils
@ -291,6 +295,10 @@ output_utils
Initial CPU and wall times when printing in the output files Initial CPU and wall times when printing in the output files
output_zmq
Output file for ZMQ
`write_bool <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files/output.irp.f#L88>`_ `write_bool <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files/output.irp.f#L88>`_
Write an logical value in output Write an logical value in output

1
src/Integrals_Bielec/.gitignore vendored
View File

@ -13,6 +13,7 @@ Makefile.depend
Nuclei Nuclei
Pseudo Pseudo
Utils Utils
ZMQ
ezfio_interface.irp.f ezfio_interface.irp.f
irpf90.make irpf90.make
irpf90_entities irpf90_entities

77
src/Integrals_Bielec/README.rst
View File

@ -31,6 +31,7 @@ Needed Modules
* `Pseudo <http://github.com/LCPQ/quantum_package/tree/master/src/Pseudo>`_ * `Pseudo <http://github.com/LCPQ/quantum_package/tree/master/src/Pseudo>`_
* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_ * `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
* `ZMQ <http://github.com/LCPQ/quantum_package/tree/master/src/ZMQ>`_
Documentation Documentation
============= =============
@ -47,7 +48,7 @@ Documentation
i(r1) j(r1) 1/r12 k(r2) l(r2) i(r1) j(r1) 1/r12 k(r2) l(r2)
`ao_bielec_integral_schwartz <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L501>`_ `ao_bielec_integral_schwartz <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L442>`_
Needed to compute Schwartz inequalities Needed to compute Schwartz inequalities
@ -61,6 +62,14 @@ Documentation
i(r1) j(r2) 1/r12 k(r1) l(r2) i(r1) j(r2) 1/r12 k(r1) l(r2)
`ao_bielec_integrals_in_map_collector <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bielec_integrals_in_map_slave.irp.f#L57>`_
Collects results from the AO integral calculation
`ao_bielec_integrals_in_map_slave <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bielec_integrals_in_map_slave.irp.f#L1>`_
Computes a buffer of integrals
`ao_integrals_map <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L6>`_ `ao_integrals_map <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L6>`_
AO integrals AO integrals
@ -89,7 +98,7 @@ Documentation
Frees the memory of the AO map Frees the memory of the AO map
`clear_mo_map <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L422>`_ `clear_mo_map <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L508>`_
Frees the memory of the MO map Frees the memory of the MO map
@ -97,6 +106,10 @@ Documentation
Compute AO 1/r12 integrals for all i and fixed j,k,l Compute AO 1/r12 integrals for all i and fixed j,k,l
`compute_ao_integrals_jl <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L1202>`_
Parallel client for AO integrals
`disk_access_ao_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ezfio_interface.irp.f#L28>`_ `disk_access_ao_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ezfio_interface.irp.f#L28>`_
Read/Write AO integrals from/to disk [ Write | Read | None ] Read/Write AO integrals from/to disk [ Write | Read | None ]
@ -109,15 +122,15 @@ Documentation
Compute integrals on the fly Compute integrals on the fly
`dump_ao_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f_template_567#L3>`_ `dump_ao_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f_template_562#L3>`_
Save to disk the $ao integrals Save to disk the $ao integrals
`dump_mo_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f_template_567#L137>`_ `dump_mo_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f_template_562#L137>`_
Save to disk the $ao integrals Save to disk the $ao integrals
`eri <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L664>`_ `eri <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L605>`_
ATOMIC PRIMTIVE bielectronic integral between the 4 primitives :: ATOMIC PRIMTIVE bielectronic integral between the 4 primitives ::
primitive_1 = x1**(a_x) y1**(a_y) z1**(a_z) exp(-alpha * r1**2) primitive_1 = x1**(a_x) y1**(a_y) z1**(a_z) exp(-alpha * r1**2)
primitive_2 = x1**(b_x) y1**(b_y) z1**(b_z) exp(- beta * r1**2) primitive_2 = x1**(b_x) y1**(b_y) z1**(b_z) exp(- beta * r1**2)
@ -139,7 +152,7 @@ Documentation
t_w(i,2,k) = t(i) t_w(i,2,k) = t(i)
`general_primitive_integral <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L526>`_ `general_primitive_integral <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L467>`_
Computes the integral <pq|rs> where p,q,r,s are Gaussian primitives Computes the integral <pq|rs> where p,q,r,s are Gaussian primitives
@ -161,11 +174,11 @@ Documentation
Returns the number of elements in the AO map Returns the number of elements in the AO map
`get_mo_bielec_integral <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L281>`_ `get_mo_bielec_integral <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L280>`_
Returns one integral <ij|kl> in the MO basis Returns one integral <ij|kl> in the MO basis
`get_mo_bielec_integral_schwartz <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L299>`_ `get_mo_bielec_integral_schwartz <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L298>`_
Returns one integral <ij|kl> in the MO basis Returns one integral <ij|kl> in the MO basis
@ -174,47 +187,47 @@ Documentation
i for j,k,l fixed. i for j,k,l fixed.
`get_mo_bielec_integrals_existing_ik <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L364>`_ `get_mo_bielec_integrals_ij <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L364>`_
Returns multiple integrals <ij|kl> in the MO basis, all Returns multiple integrals <ij|kl> in the MO basis, all
i(1)j(1) 1/r12 k(2)l(2) i(1)j(2) 1/r12 k(1)l(2)
i for j,k,l fixed. i, j for k,l fixed.
`get_mo_map_size <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L414>`_ `get_mo_map_size <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L418>`_
Return the number of elements in the MO map Return the number of elements in the MO map
`give_polynom_mult_center_x <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L874>`_ `give_polynom_mult_center_x <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L819>`_
subroutine that returns the explicit polynom in term of the "t" subroutine that returns the explicit polynom in term of the "t"
variable of the following polynomw : variable of the following polynomw :
I_x1(a_x, d_x,p,q) * I_x1(a_y, d_y,p,q) * I_x1(a_z, d_z,p,q) I_x1(a_x, d_x,p,q) * I_x1(a_y, d_y,p,q) * I_x1(a_z, d_z,p,q)
`i_x1_new <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L795>`_ `i_x1_new <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L738>`_
recursive function involved in the bielectronic integral recursive function involved in the bielectronic integral
`i_x1_pol_mult <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L937>`_ `i_x1_pol_mult <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L882>`_
recursive function involved in the bielectronic integral recursive function involved in the bielectronic integral
`i_x1_pol_mult_a1 <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L1057>`_ `i_x1_pol_mult_a1 <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L1002>`_
recursive function involved in the bielectronic integral recursive function involved in the bielectronic integral
`i_x1_pol_mult_a2 <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L1111>`_ `i_x1_pol_mult_a2 <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L1056>`_
recursive function involved in the bielectronic integral recursive function involved in the bielectronic integral
`i_x1_pol_mult_recurs <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L971>`_ `i_x1_pol_mult_recurs <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L916>`_
recursive function involved in the bielectronic integral recursive function involved in the bielectronic integral
`i_x2_new <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L830>`_ `i_x2_new <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L774>`_
recursive function involved in the bielectronic integral recursive function involved in the bielectronic integral
`i_x2_pol_mult <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L1173>`_ `i_x2_pol_mult <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L1118>`_
recursive function involved in the bielectronic integral recursive function involved in the bielectronic integral
@ -222,21 +235,21 @@ Documentation
Create new entry into AO map Create new entry into AO map
`insert_into_mo_integrals_map <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L265>`_ `insert_into_mo_integrals_map <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f#L264>`_
Create new entry into MO map, or accumulate in an existing entry Create new entry into MO map, or accumulate in an existing entry
`integrale_new <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L721>`_ `integrale_new <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L662>`_
calculate the integral of the polynom :: calculate the integral of the polynom ::
I_x1(a_x+b_x, c_x+d_x,p,q) * I_x1(a_y+b_y, c_y+d_y,p,q) * I_x1(a_z+b_z, c_z+d_z,p,q) I_x1(a_x+b_x, c_x+d_x,p,q) * I_x1(a_y+b_y, c_y+d_y,p,q) * I_x1(a_z+b_z, c_z+d_z,p,q)
between ( 0 ; 1) between ( 0 ; 1)
`load_ao_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f_template_567#L89>`_ `load_ao_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f_template_562#L89>`_
Read from disk the $ao integrals Read from disk the $ao integrals
`load_mo_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f_template_567#L223>`_ `load_mo_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/map_integrals.irp.f_template_562#L223>`_
Read from disk the $ao integrals Read from disk the $ao integrals
@ -244,43 +257,43 @@ Documentation
Returns one integral <ij|kl> in the MO basis Returns one integral <ij|kl> in the MO basis
`mo_bielec_integral_jj <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L465>`_ `mo_bielec_integral_jj <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L464>`_
mo_bielec_integral_jj(i,j) = J_ij mo_bielec_integral_jj(i,j) = J_ij
mo_bielec_integral_jj_exchange(i,j) = K_ij mo_bielec_integral_jj_exchange(i,j) = K_ij
mo_bielec_integral_jj_anti(i,j) = J_ij - K_ij mo_bielec_integral_jj_anti(i,j) = J_ij - K_ij
`mo_bielec_integral_jj_anti <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L467>`_ `mo_bielec_integral_jj_anti <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L466>`_
mo_bielec_integral_jj(i,j) = J_ij mo_bielec_integral_jj(i,j) = J_ij
mo_bielec_integral_jj_exchange(i,j) = K_ij mo_bielec_integral_jj_exchange(i,j) = K_ij
mo_bielec_integral_jj_anti(i,j) = J_ij - K_ij mo_bielec_integral_jj_anti(i,j) = J_ij - K_ij
`mo_bielec_integral_jj_anti_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L327>`_ `mo_bielec_integral_jj_anti_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L326>`_
mo_bielec_integral_jj_from_ao(i,j) = J_ij mo_bielec_integral_jj_from_ao(i,j) = J_ij
mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij
mo_bielec_integral_jj_anti_from_ao(i,j) = J_ij - K_ij mo_bielec_integral_jj_anti_from_ao(i,j) = J_ij - K_ij
`mo_bielec_integral_jj_exchange <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L466>`_ `mo_bielec_integral_jj_exchange <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L465>`_
mo_bielec_integral_jj(i,j) = J_ij mo_bielec_integral_jj(i,j) = J_ij
mo_bielec_integral_jj_exchange(i,j) = K_ij mo_bielec_integral_jj_exchange(i,j) = K_ij
mo_bielec_integral_jj_anti(i,j) = J_ij - K_ij mo_bielec_integral_jj_anti(i,j) = J_ij - K_ij
`mo_bielec_integral_jj_exchange_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L326>`_ `mo_bielec_integral_jj_exchange_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L325>`_
mo_bielec_integral_jj_from_ao(i,j) = J_ij mo_bielec_integral_jj_from_ao(i,j) = J_ij
mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij
mo_bielec_integral_jj_anti_from_ao(i,j) = J_ij - K_ij mo_bielec_integral_jj_anti_from_ao(i,j) = J_ij - K_ij
`mo_bielec_integral_jj_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L325>`_ `mo_bielec_integral_jj_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L324>`_
mo_bielec_integral_jj_from_ao(i,j) = J_ij mo_bielec_integral_jj_from_ao(i,j) = J_ij
mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij
mo_bielec_integral_jj_anti_from_ao(i,j) = J_ij - K_ij mo_bielec_integral_jj_anti_from_ao(i,j) = J_ij - K_ij
`mo_bielec_integral_schwartz <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L492>`_ `mo_bielec_integral_schwartz <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/mo_bi_integrals.irp.f#L491>`_
Needed to compute Schwartz inequalities Needed to compute Schwartz inequalities
@ -304,7 +317,7 @@ Documentation
Aligned n_pt_max_integrals Aligned n_pt_max_integrals
`n_pt_sup <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L860>`_ `n_pt_sup <http://github.com/LCPQ/quantum_package/tree/master/src/Integrals_Bielec/ao_bi_integrals.irp.f#L805>`_
Returns the upper boundary of the degree of the polynomial involved in the Returns the upper boundary of the degree of the polynomial involved in the
bielctronic integral : bielctronic integral :
Ix(a_x,b_x,c_x,d_x) * Iy(a_y,b_y,c_y,d_y) * Iz(a_z,b_z,c_z,d_z) Ix(a_x,b_x,c_x,d_x) * Iy(a_y,b_y,c_y,d_y) * Iz(a_z,b_z,c_z,d_z)

116
src/Integrals_Bielec/ao_bi_integrals.irp.f
View File

@ -40,24 +40,22 @@ double precision function ao_bielec_integral(i,j,k,l)
L_center(p) = nucl_coord(num_l,p) L_center(p) = nucl_coord(num_l,p)
enddo enddo
double precision :: coef1, coef2, coef3, coef4
double precision :: p_inv,q_inv
double precision :: general_primitive_integral
do p = 1, ao_prim_num(i) do p = 1, ao_prim_num(i)
double precision :: coef1
coef1 = ao_coef_normalized_ordered_transp(p,i) coef1 = ao_coef_normalized_ordered_transp(p,i)
do q = 1, ao_prim_num(j) do q = 1, ao_prim_num(j)
double precision :: coef2
coef2 = coef1*ao_coef_normalized_ordered_transp(q,j) coef2 = coef1*ao_coef_normalized_ordered_transp(q,j)
double precision :: p_inv,q_inv
call give_explicit_poly_and_gaussian(P_new,P_center,pp,fact_p,iorder_p,& call give_explicit_poly_and_gaussian(P_new,P_center,pp,fact_p,iorder_p,&
ao_expo_ordered_transp(p,i),ao_expo_ordered_transp(q,j), & ao_expo_ordered_transp(p,i),ao_expo_ordered_transp(q,j), &
I_power,J_power,I_center,J_center,dim1) I_power,J_power,I_center,J_center,dim1)
p_inv = 1.d0/pp p_inv = 1.d0/pp
do r = 1, ao_prim_num(k) do r = 1, ao_prim_num(k)
double precision :: coef3
coef3 = coef2*ao_coef_normalized_ordered_transp(r,k) coef3 = coef2*ao_coef_normalized_ordered_transp(r,k)
do s = 1, ao_prim_num(l) do s = 1, ao_prim_num(l)
double precision :: coef4
coef4 = coef3*ao_coef_normalized_ordered_transp(s,l) coef4 = coef3*ao_coef_normalized_ordered_transp(s,l)
double precision :: general_primitive_integral
call give_explicit_poly_and_gaussian(Q_new,Q_center,qq,fact_q,iorder_q,& call give_explicit_poly_and_gaussian(Q_new,Q_center,qq,fact_q,iorder_q,&
ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l), & ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l), &
K_power,L_power,K_center,L_center,dim1) K_power,L_power,K_center,L_center,dim1)
@ -65,7 +63,7 @@ double precision function ao_bielec_integral(i,j,k,l)
integral = general_primitive_integral(dim1, & integral = general_primitive_integral(dim1, &
P_new,P_center,fact_p,pp,p_inv,iorder_p, & P_new,P_center,fact_p,pp,p_inv,iorder_p, &
Q_new,Q_center,fact_q,qq,q_inv,iorder_q) Q_new,Q_center,fact_q,qq,q_inv,iorder_q)
ao_bielec_integral += coef4 * integral ao_bielec_integral = ao_bielec_integral + coef4 * integral
enddo ! s enddo ! s
enddo ! r enddo ! r
enddo ! q enddo ! q
@ -94,7 +92,7 @@ double precision function ao_bielec_integral(i,j,k,l)
I_power(1),J_power(1),K_power(1),L_power(1), & I_power(1),J_power(1),K_power(1),L_power(1), &
I_power(2),J_power(2),K_power(2),L_power(2), & I_power(2),J_power(2),K_power(2),L_power(2), &
I_power(3),J_power(3),K_power(3),L_power(3)) I_power(3),J_power(3),K_power(3),L_power(3))
ao_bielec_integral += coef4 * integral ao_bielec_integral = ao_bielec_integral + coef4 * integral
enddo ! s enddo ! s
enddo ! r enddo ! r
enddo ! q enddo ! q
@ -129,8 +127,8 @@ double precision function ao_bielec_integral_schwartz_accel(i,j,k,l)
num_k = ao_nucl(k) num_k = ao_nucl(k)
num_l = ao_nucl(l) num_l = ao_nucl(l)
ao_bielec_integral_schwartz_accel = 0.d0 ao_bielec_integral_schwartz_accel = 0.d0
double precision :: thresh double precision :: thr
thresh = ao_integrals_threshold*ao_integrals_threshold thr = ao_integrals_threshold*ao_integrals_threshold
allocate(schwartz_kl(0:ao_prim_num(l),0:ao_prim_num(k))) allocate(schwartz_kl(0:ao_prim_num(l),0:ao_prim_num(k)))
@ -181,18 +179,18 @@ double precision function ao_bielec_integral_schwartz_accel(i,j,k,l)
P_new,P_center,fact_p,pp,p_inv,iorder_p, & P_new,P_center,fact_p,pp,p_inv,iorder_p, &
P_new,P_center,fact_p,pp,p_inv,iorder_p) * & P_new,P_center,fact_p,pp,p_inv,iorder_p) * &
coef2*coef2 coef2*coef2
if (schwartz_kl(0,0)*schwartz_ij < thresh) then if (schwartz_kl(0,0)*schwartz_ij < thr) then
cycle cycle
endif endif
do r = 1, ao_prim_num(k) do r = 1, ao_prim_num(k)
if (schwartz_kl(0,r)*schwartz_ij < thresh) then if (schwartz_kl(0,r)*schwartz_ij < thr) then
cycle cycle
endif endif
double precision :: coef3 double precision :: coef3
coef3 = coef2*ao_coef_normalized_ordered_transp(r,k) coef3 = coef2*ao_coef_normalized_ordered_transp(r,k)
do s = 1, ao_prim_num(l) do s = 1, ao_prim_num(l)
double precision :: coef4 double precision :: coef4
if (schwartz_kl(s,r)*schwartz_ij < thresh) then if (schwartz_kl(s,r)*schwartz_ij < thr) then
cycle cycle
endif endif
coef4 = coef3*ao_coef_normalized_ordered_transp(s,l) coef4 = coef3*ao_coef_normalized_ordered_transp(s,l)
@ -246,16 +244,16 @@ double precision function ao_bielec_integral_schwartz_accel(i,j,k,l)
I_power(1),J_power(1),I_power(1),J_power(1), & I_power(1),J_power(1),I_power(1),J_power(1), &
I_power(2),J_power(2),I_power(2),J_power(2), & I_power(2),J_power(2),I_power(2),J_power(2), &
I_power(3),J_power(3),I_power(3),J_power(3))*coef2*coef2 I_power(3),J_power(3),I_power(3),J_power(3))*coef2*coef2
if (schwartz_kl(0,0)*schwartz_ij < thresh) then if (schwartz_kl(0,0)*schwartz_ij < thr) then
cycle cycle
endif endif
do r = 1, ao_prim_num(k) do r = 1, ao_prim_num(k)
if (schwartz_kl(0,r)*schwartz_ij < thresh) then if (schwartz_kl(0,r)*schwartz_ij < thr) then
cycle cycle
endif endif
coef3 = coef2*ao_coef_normalized_ordered_transp(r,k) coef3 = coef2*ao_coef_normalized_ordered_transp(r,k)
do s = 1, ao_prim_num(l) do s = 1, ao_prim_num(l)
if (schwartz_kl(s,r)*schwartz_ij < thresh) then if (schwartz_kl(s,r)*schwartz_ij < thr) then
cycle cycle
endif endif
coef4 = coef3*ao_coef_normalized_ordered_transp(s,l) coef4 = coef3*ao_coef_normalized_ordered_transp(s,l)
@ -295,11 +293,10 @@ subroutine compute_ao_bielec_integrals(j,k,l,sze,buffer_value)
! Compute AO 1/r12 integrals for all i and fixed j,k,l ! Compute AO 1/r12 integrals for all i and fixed j,k,l
END_DOC END_DOC
include 'Utils/constants.include.F'
integer, intent(in) :: j,k,l,sze integer, intent(in) :: j,k,l,sze
real(integral_kind), intent(out) :: buffer_value(sze) real(integral_kind), intent(out) :: buffer_value(sze)
double precision :: ao_bielec_integral double precision :: ao_bielec_integral
double precision :: thresh
thresh = ao_integrals_threshold
integer :: i integer :: i
@ -339,8 +336,7 @@ BEGIN_PROVIDER [ logical, ao_bielec_integrals_in_map ]
integer :: i,j,k,l integer :: i,j,k,l
double precision :: ao_bielec_integral,cpu_1,cpu_2, wall_1, wall_2 double precision :: ao_bielec_integral,cpu_1,cpu_2, wall_1, wall_2
double precision :: integral, wall_0 double precision :: integral, wall_0
double precision :: thresh include 'Utils/constants.include.F'
thresh = ao_integrals_threshold
! For integrals file ! For integrals file
integer(key_kind),allocatable :: buffer_i(:) integer(key_kind),allocatable :: buffer_i(:)
@ -348,7 +344,7 @@ BEGIN_PROVIDER [ logical, ao_bielec_integrals_in_map ]
real(integral_kind),allocatable :: buffer_value(:) real(integral_kind),allocatable :: buffer_value(:)
integer :: n_integrals, rc integer :: n_integrals, rc
integer :: jl_pairs(2,ao_num*(ao_num+1)/2), kk, m, j1, i1, lmax integer :: kk, m, j1, i1, lmax
integral = ao_bielec_integral(1,1,1,1) integral = ao_bielec_integral(1,1,1,1)
@ -368,55 +364,21 @@ BEGIN_PROVIDER [ logical, ao_bielec_integrals_in_map ]
call wall_time(wall_1) call wall_time(wall_1)
call cpu_time(cpu_1) call cpu_time(cpu_1)
integer(ZMQ_PTR) :: zmq_socket_rep_inproc, zmq_socket_push_inproc integer(ZMQ_PTR) :: zmq_to_qp_run_socket
zmq_socket_rep_inproc = f77_zmq_socket(zmq_context, ZMQ_REP) call new_parallel_job(zmq_to_qp_run_socket,'ao_integrals')
rc = f77_zmq_bind(zmq_socket_rep_inproc, 'inproc://req_rep')
if (rc /= 0) then
stop 'Unable to connect zmq_socket_rep_inproc'
endif
integer(ZMQ_PTR) :: thread(0:nproc)
external :: ao_bielec_integrals_in_map_slave, ao_bielec_integrals_in_map_collector character*(32) :: task
rc = pthread_create( thread(0), ao_bielec_integrals_in_map_collector )
! Create client threads do l=1,ao_num
do i=1,nproc write(task,*) 'triangle', l
rc = pthread_create( thread(i), ao_bielec_integrals_in_map_slave ) call add_task_to_taskserver(zmq_to_qp_run_socket,task)
enddo enddo
character*(64) :: message_string external :: ao_bielec_integrals_in_map_slave_inproc, ao_bielec_integrals_in_map_collector
call new_parallel_threads(ao_bielec_integrals_in_map_slave_inproc, ao_bielec_integrals_in_map_collector)
do l = ao_num, 1, -1 call end_parallel_job(zmq_to_qp_run_socket,'ao_integrals')
rc = f77_zmq_recv( zmq_socket_rep_inproc, message_string, 64, 0)
print *, l
! TODO : error handling
ASSERT (rc >= 0)
ASSERT (message == 'get_ao_integrals')
rc = f77_zmq_send( zmq_socket_rep_inproc, l, 4, 0)
enddo
do i=1,nproc
rc = f77_zmq_recv( zmq_socket_rep_inproc, message_string, 64, 0)
! TODO : error handling
ASSERT (rc >= 0)
ASSERT (message == 'get_ao_integrals')
rc = f77_zmq_send( zmq_socket_rep_inproc, 0, 4, 0)
enddo
! TODO terminer thread(0)
rc = f77_zmq_unbind(zmq_socket_rep_inproc, 'inproc://req_rep')
do i=1,nproc
rc = pthread_join( thread(i) )
enddo
zmq_socket_push_inproc = f77_zmq_socket(zmq_context, ZMQ_PUSH)
rc = f77_zmq_connect(zmq_socket_push_inproc, 'inproc://push_pull')
if (rc /= 0) then
stop 'Unable to connect zmq_socket_push_inproc'
endif
rc = f77_zmq_send( zmq_socket_push_inproc, -1, 4, ZMQ_SNDMORE)
rc = f77_zmq_send( zmq_socket_push_inproc, 0_key_kind, key_kind, ZMQ_SNDMORE)
rc = f77_zmq_send( zmq_socket_push_inproc, 0_integral_kind, integral_kind, 0)
rc = pthread_join( thread(0) )
print*, 'Sorting the map' print*, 'Sorting the map'
call map_sort(ao_integrals_map) call map_sort(ao_integrals_map)
@ -513,11 +475,11 @@ double precision function general_primitive_integral(dim, &
enddo enddo
n_Ix = 0 n_Ix = 0
do ix = 0, iorder_p(1) do ix = 0, iorder_p(1)
if (abs(P_new(ix,1)) < 1.d-8) cycle if (abs(P_new(ix,1)) < thresh) cycle
a = P_new(ix,1) a = P_new(ix,1)
do jx = 0, iorder_q(1) do jx = 0, iorder_q(1)
d = a*Q_new(jx,1) d = a*Q_new(jx,1)
if (abs(d) < 1.d-8) cycle if (abs(d) < thresh) cycle
!DEC$ FORCEINLINE !DEC$ FORCEINLINE
call give_polynom_mult_center_x(P_center(1),Q_center(1),ix,jx,p,q,iorder,pq_inv,pq_inv_2,p10_1,p01_1,p10_2,p01_2,dx,nx) call give_polynom_mult_center_x(P_center(1),Q_center(1),ix,jx,p,q,iorder,pq_inv,pq_inv_2,p10_1,p01_1,p10_2,p01_2,dx,nx)
!DEC$ FORCEINLINE !DEC$ FORCEINLINE
@ -534,11 +496,11 @@ double precision function general_primitive_integral(dim, &
enddo enddo
n_Iy = 0 n_Iy = 0
do iy = 0, iorder_p(2) do iy = 0, iorder_p(2)
if (abs(P_new(iy,2)) > 1.d-8) then if (abs(P_new(iy,2)) > thresh) then
b = P_new(iy,2) b = P_new(iy,2)
do jy = 0, iorder_q(2) do jy = 0, iorder_q(2)
e = b*Q_new(jy,2) e = b*Q_new(jy,2)
if (abs(e) < 1.d-8) cycle if (abs(e) < thresh) cycle
!DEC$ FORCEINLINE !DEC$ FORCEINLINE
call give_polynom_mult_center_x(P_center(2),Q_center(2),iy,jy,p,q,iorder,pq_inv,pq_inv_2,p10_1,p01_1,p10_2,p01_2,dy,ny) call give_polynom_mult_center_x(P_center(2),Q_center(2),iy,jy,p,q,iorder,pq_inv,pq_inv_2,p10_1,p01_1,p10_2,p01_2,dy,ny)
!DEC$ FORCEINLINE !DEC$ FORCEINLINE
@ -556,11 +518,11 @@ double precision function general_primitive_integral(dim, &
enddo enddo
n_Iz = 0 n_Iz = 0
do iz = 0, iorder_p(3) do iz = 0, iorder_p(3)
if (abs(P_new(iz,3)) > 1.d-8) then if (abs(P_new(iz,3)) > thresh) then
c = P_new(iz,3) c = P_new(iz,3)
do jz = 0, iorder_q(3) do jz = 0, iorder_q(3)
f = c*Q_new(jz,3) f = c*Q_new(jz,3)
if (abs(f) < 1.d-8) cycle if (abs(f) < thresh) cycle
!DEC$ FORCEINLINE !DEC$ FORCEINLINE
call give_polynom_mult_center_x(P_center(3),Q_center(3),iz,jz,p,q,iorder,pq_inv,pq_inv_2,p10_1,p01_1,p10_2,p01_2,dz,nz) call give_polynom_mult_center_x(P_center(3),Q_center(3),iz,jz,p,q,iorder,pq_inv,pq_inv_2,p10_1,p01_1,p10_2,p01_2,dz,nz)
!DEC$ FORCEINLINE !DEC$ FORCEINLINE
@ -1214,10 +1176,10 @@ subroutine compute_ao_integrals_jl(j,l,n_integrals,buffer_i,buffer_value)
integer :: i,k integer :: i,k
double precision :: ao_bielec_integral,cpu_1,cpu_2, wall_1, wall_2 double precision :: ao_bielec_integral,cpu_1,cpu_2, wall_1, wall_2
double precision :: integral, wall_0 double precision :: integral, wall_0
double precision :: thresh double precision :: thr
integer :: kk, m, j1, i1 integer :: kk, m, j1, i1
thresh = ao_integrals_threshold thr = ao_integrals_threshold
n_integrals = 0 n_integrals = 0
@ -1232,15 +1194,15 @@ subroutine compute_ao_integrals_jl(j,l,n_integrals,buffer_i,buffer_value)
if (i1 > j1) then if (i1 > j1) then
exit exit
endif endif
if (ao_overlap_abs(i,k)*ao_overlap_abs(j,l) < thresh) then if (ao_overlap_abs(i,k)*ao_overlap_abs(j,l) < thr) then
cycle cycle
endif endif
if (ao_bielec_integral_schwartz(i,k)*ao_bielec_integral_schwartz(j,l) < thresh ) then if (ao_bielec_integral_schwartz(i,k)*ao_bielec_integral_schwartz(j,l) < thr ) then
cycle cycle
endif endif
!DIR$ FORCEINLINE !DIR$ FORCEINLINE
integral = ao_bielec_integral(i,k,j,l) integral = ao_bielec_integral(i,k,j,l)
if (abs(integral) < thresh) then if (abs(integral) < thr) then
cycle cycle
endif endif
n_integrals += 1 n_integrals += 1

119
src/Integrals_Bielec/ao_bielec_integrals_in_map_slave.irp.f
View File

@ -1,4 +1,20 @@
subroutine ao_bielec_integrals_in_map_slave subroutine ao_bielec_integrals_in_map_slave_tcp
implicit none
BEGIN_DOC
! Computes a buffer of integrals
END_DOC
call ao_bielec_integrals_in_map_slave(0)
end
subroutine ao_bielec_integrals_in_map_slave_inproc
implicit none
BEGIN_DOC
! Computes a buffer of integrals
END_DOC
call ao_bielec_integrals_in_map_slave(1)
end
subroutine ao_bielec_integrals_in_map_slave(thread)
use map_module use map_module
use f77_zmq use f77_zmq
implicit none implicit none
@ -6,51 +22,65 @@ subroutine ao_bielec_integrals_in_map_slave
! Computes a buffer of integrals ! Computes a buffer of integrals
END_DOC END_DOC
integer, intent(in) :: thread
integer :: j,l,n_integrals integer :: j,l,n_integrals
integer :: rc integer :: rc
character*(8), external :: zmq_port
integer(ZMQ_PTR) :: zmq_socket_req_inproc, zmq_socket_push_inproc
real(integral_kind), allocatable :: buffer_value(:) real(integral_kind), allocatable :: buffer_value(:)
integer(key_kind), allocatable :: buffer_i(:) integer(key_kind), allocatable :: buffer_i(:)
integer :: worker_id, task_id
character*(512) :: task
integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket
integer(ZMQ_PTR) :: zmq_to_qp_run_socket
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
integer(ZMQ_PTR) :: zmq_socket_push
! zmq_socket_push = f77_zmq_socket(zmq_context, ZMQ_PUSH)
zmq_socket_push = f77_zmq_socket(zmq_context, ZMQ_REQ )
if (thread == 1) then
rc = f77_zmq_connect(zmq_socket_push, trim(zmq_socket_pull_inproc_address))
else
rc = f77_zmq_connect(zmq_socket_push, trim(zmq_socket_push_tcp_address))
endif
if (rc /= 0) then
stop 'Unable to connect zmq_socket_push_tcp'
endif
allocate ( buffer_i(ao_num*ao_num), buffer_value(ao_num*ao_num) ) allocate ( buffer_i(ao_num*ao_num), buffer_value(ao_num*ao_num) )
! Sockets call connect_to_taskserver(zmq_to_qp_run_socket,worker_id,thread)
zmq_socket_req_inproc = f77_zmq_socket(zmq_context, ZMQ_REQ)
rc = f77_zmq_connect(zmq_socket_req_inproc, 'inproc://req_rep') do
if (rc /= 0) then call get_task_from_taskserver(zmq_to_qp_run_socket,worker_id, task_id, task)
stop 'Unable to connect zmq_socket_req_inproc' if (task_id == 0) then
endif exit
zmq_socket_push_inproc = f77_zmq_socket(zmq_context, ZMQ_PUSH)
rc = f77_zmq_connect(zmq_socket_push_inproc, 'inproc://push_pull')
if (rc /= 0) then
stop 'Unable to connect zmq_socket_push_inproc'
endif
rc = f77_zmq_send( zmq_socket_req_inproc, 'get_ao_integrals', 16, 0)
rc = f77_zmq_recv( zmq_socket_req_inproc, l, 4, 0)
do while (l > 0)
rc = f77_zmq_send( zmq_socket_req_inproc, 'get_ao_integrals', 16, 0)
do j = 1, l
if (ao_overlap_abs(j,l) < ao_integrals_threshold) then
cycle
endif endif
read(task,*) j, l
call compute_ao_integrals_jl(j,l,n_integrals,buffer_i,buffer_value) call compute_ao_integrals_jl(j,l,n_integrals,buffer_i,buffer_value)
rc = f77_zmq_send( zmq_socket_push_inproc, n_integrals, 4, ZMQ_SNDMORE) rc = f77_zmq_send( zmq_socket_push, n_integrals, 4, ZMQ_SNDMORE)
rc = f77_zmq_send( zmq_socket_push_inproc, buffer_i, key_kind*n_integrals, ZMQ_SNDMORE) rc = f77_zmq_send( zmq_socket_push, buffer_i, key_kind*n_integrals, ZMQ_SNDMORE)
rc = f77_zmq_send( zmq_socket_push_inproc, buffer_value, integral_kind*n_integrals, 0) rc = f77_zmq_send( zmq_socket_push, buffer_value, integral_kind*n_integrals, 0)
enddo call task_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id)
rc = f77_zmq_recv( zmq_socket_req_inproc, l, 4, 0) character*(2) :: ok
rc = f77_zmq_recv( zmq_socket_push, ok, 2, 0)
enddo enddo
deallocate( buffer_i, buffer_value ) deallocate( buffer_i, buffer_value )
rc = f77_zmq_disconnect(zmq_socket_req_inproc, 'inproc://req_rep') integer :: finished
call disconnect_from_taskserver(zmq_to_qp_run_socket,worker_id,finished)
if (finished /= 0) then
rc = f77_zmq_send( zmq_socket_push, -1, 4, 0)
rc = f77_zmq_recv( zmq_socket_push, ok, 2, ZMQ_NOBLOCK)
endif
rc = f77_zmq_disconnect(zmq_socket_push,trim(zmq_socket_push_tcp_address))
rc = f77_zmq_close(zmq_socket_push)
end end
@ -64,36 +94,27 @@ subroutine ao_bielec_integrals_in_map_collector
integer :: j,l,n_integrals integer :: j,l,n_integrals
integer :: rc integer :: rc
character*(8), external :: zmq_port
integer(ZMQ_PTR) :: zmq_socket_pull_inproc
real(integral_kind), allocatable :: buffer_value(:) real(integral_kind), allocatable :: buffer_value(:)
integer(key_kind), allocatable :: buffer_i(:) integer(key_kind), allocatable :: buffer_i(:)
allocate ( buffer_i(ao_num*ao_num), buffer_value(ao_num*ao_num) ) allocate ( buffer_i(ao_num*ao_num), buffer_value(ao_num*ao_num) )
zmq_socket_pull_inproc = f77_zmq_socket(zmq_context, ZMQ_PULL)
rc = f77_zmq_bind(zmq_socket_pull_inproc, 'inproc://push_pull')
if (rc /= 0) then
stop 'Unable to connect zmq_socket_pull_inproc'
endif
n_integrals = 0 n_integrals = 0
do while (n_integrals >= 0) do while (n_integrals >= 0)
rc = f77_zmq_recv( zmq_socket_pull_inproc, n_integrals, 4, 0) rc = f77_zmq_recv( zmq_socket_pull, n_integrals, 4, 0)
if (n_integrals > -1) then if (n_integrals >= 0) then
rc = f77_zmq_recv( zmq_socket_pull_inproc, buffer_i, key_kind*n_integrals, 0) rc = f77_zmq_recv( zmq_socket_pull, buffer_i, key_kind*n_integrals, 0)
rc = f77_zmq_recv( zmq_socket_pull_inproc, buffer_value, integral_kind*n_integrals, 0) rc = f77_zmq_recv( zmq_socket_pull, buffer_value, integral_kind*n_integrals, 0)
rc = f77_zmq_send( zmq_socket_pull, 'ok', 2, 0)
call insert_into_ao_integrals_map(n_integrals,buffer_i,buffer_value) call insert_into_ao_integrals_map(n_integrals,buffer_i,buffer_value)
else else
rc = f77_zmq_recv( zmq_socket_pull_inproc, buffer_i, key_kind, 0) rc = f77_zmq_send( zmq_socket_pull, 'ok', 2, 0)
rc = f77_zmq_recv( zmq_socket_pull_inproc, buffer_value, integral_kind, 0)
endif endif
enddo enddo
rc = f77_zmq_unbind(zmq_socket_pull_inproc, 'inproc://push_pull')
deallocate( buffer_i, buffer_value ) deallocate( buffer_i, buffer_value )
end end

2
src/Integrals_Bielec/mo_bi_integrals.irp.f
View File

@ -70,7 +70,7 @@ subroutine add_integrals_to_map(mask_ijkl)
integer :: i1,j1,k1,l1, ii1, kmax, thread_num integer :: i1,j1,k1,l1, ii1, kmax, thread_num
integer :: i2,i3,i4 integer :: i2,i3,i4
double precision,parameter :: thr_coef = 0.d0 double precision,parameter :: thr_coef = 1.d-10
PROVIDE ao_bielec_integrals_in_map PROVIDE ao_bielec_integrals_in_map

20
src/Integrals_Bielec/qp_ao_ints.irp.f Normal file
View File

@ -0,0 +1,20 @@
program qp_ao_ints
implicit none
BEGIN_DOC
! Increments a running calculation to compute AO integrals
END_DOC
! Set the state of the ZMQ
zmq_state = 'ao_integrals'
! Provide everything needed
double precision :: integral, ao_bielec_integral
integral = ao_bielec_integral(1,1,1,1)
!$OMP PARALLEL DEFAULT(PRIVATE)
call ao_bielec_integrals_in_map_slave_tcp
!$OMP END PARALLEL
print *, 'Done'
end

34
src/Integrals_Monoelec/pot_ao_ints.irp.f
View File

@ -26,7 +26,7 @@
n_pt_in = n_pt_max_integrals n_pt_in = n_pt_max_integrals
!$OMP DO SCHEDULE (guided) !$OMP DO SCHEDULE (dynamic)
do j = 1, ao_num do j = 1, ao_num
num_A = ao_nucl(j) num_A = ao_nucl(j)
@ -81,23 +81,17 @@
integer :: power_A(3),power_B(3) integer :: power_A(3),power_B(3)
integer :: i,j,k,l,n_pt_in,m integer :: i,j,k,l,n_pt_in,m
double precision ::overlap_x,overlap_y,overlap_z,overlap,dx,NAI_pol_mult double precision ::overlap_x,overlap_y,overlap_z,overlap,dx,NAI_pol_mult
! Important for OpenMP
ao_nucl_elec_integral_per_atom = 0.d0 ao_nucl_elec_integral_per_atom = 0.d0
do k = 1, nucl_num
C_center(1) = nucl_coord(k,1)
C_center(2) = nucl_coord(k,2)
C_center(3) = nucl_coord(k,3)
!$OMP PARALLEL & !$OMP PARALLEL &
!$OMP DEFAULT (NONE) & !$OMP DEFAULT (NONE) &
!$OMP PRIVATE (i,j,l,m,alpha,beta,A_center,B_center,power_A,power_B, & !$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,power_A,power_B, &
!$OMP num_A,num_B,c,n_pt_in) & !$OMP num_A,num_B,c,n_pt_in,C_center) &
!$OMP SHARED (k,ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp, & !$OMP SHARED (ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp, &
!$OMP n_pt_max_integrals,ao_nucl_elec_integral_per_atom,nucl_num,C_center) !$OMP n_pt_max_integrals,ao_nucl_elec_integral_per_atom,nucl_num)
n_pt_in = n_pt_max_integrals n_pt_in = n_pt_max_integrals
!$OMP DO SCHEDULE (guided) !$OMP DO SCHEDULE (dynamic)
double precision :: c double precision :: c
do j = 1, ao_num do j = 1, ao_num
@ -108,6 +102,10 @@
A_center(1) = nucl_coord(num_A,1) A_center(1) = nucl_coord(num_A,1)
A_center(2) = nucl_coord(num_A,2) A_center(2) = nucl_coord(num_A,2)
A_center(3) = nucl_coord(num_A,3) A_center(3) = nucl_coord(num_A,3)
do k = 1, nucl_num
C_center(1) = nucl_coord(k,1)
C_center(2) = nucl_coord(k,2)
C_center(3) = nucl_coord(k,3)
do i = 1, ao_num do i = 1, ao_num
power_B(1)= ao_power(i,1) power_B(1)= ao_power(i,1)
power_B(2)= ao_power(i,2) power_B(2)= ao_power(i,2)
@ -128,9 +126,9 @@
ao_nucl_elec_integral_per_atom(i,j,k) = -c ao_nucl_elec_integral_per_atom(i,j,k) = -c
enddo enddo
enddo enddo
enddo
!$OMP END DO !$OMP END DO
!$OMP END PARALLEL !$OMP END PARALLEL
enddo
END_PROVIDER END_PROVIDER
@ -173,7 +171,7 @@ include 'Utils/constants.include.F'
enddo enddo
const_factor = dist*rho const_factor = dist*rho
const = p * dist_integral const = p * dist_integral
if(const_factor.ge.80.d0)then if(const_factor > 80.d0)then
NAI_pol_mult = 0.d0 NAI_pol_mult = 0.d0
return return
endif endif
@ -378,7 +376,7 @@ recursive subroutine I_x1_pol_mult_mono_elec(a,c,R1x,R1xp,R2x,d,nd,n_pt_in)
enddo enddo
call I_x2_pol_mult_mono_elec(c-1,R1x,R1xp,R2x,X,nx,n_pt_in) call I_x2_pol_mult_mono_elec(c-1,R1x,R1xp,R2x,X,nx,n_pt_in)
do ix=0,nx do ix=0,nx
X(ix) *= c X(ix) *= dble(c)
enddo enddo
call multiply_poly(X,nx,R2x,2,d,nd) call multiply_poly(X,nx,R2x,2,d,nd)
ny=0 ny=0
@ -393,7 +391,7 @@ recursive subroutine I_x1_pol_mult_mono_elec(a,c,R1x,R1xp,R2x,d,nd,n_pt_in)
call I_x1_pol_mult_mono_elec(a-2,c,R1x,R1xp,R2x,X,nx,n_pt_in) call I_x1_pol_mult_mono_elec(a-2,c,R1x,R1xp,R2x,X,nx,n_pt_in)
! print*,'nx a-2,c= ',nx ! print*,'nx a-2,c= ',nx
do ix=0,nx do ix=0,nx
X(ix) *= a-1 X(ix) *= dble(a-1)
enddo enddo
call multiply_poly(X,nx,R2x,2,d,nd) call multiply_poly(X,nx,R2x,2,d,nd)
! print*,'nd out = ',nd ! print*,'nd out = ',nd
@ -405,7 +403,7 @@ recursive subroutine I_x1_pol_mult_mono_elec(a,c,R1x,R1xp,R2x,d,nd,n_pt_in)
call I_x1_pol_mult_mono_elec(a-1,c-1,R1x,R1xp,R2x,X,nx,n_pt_in) call I_x1_pol_mult_mono_elec(a-1,c-1,R1x,R1xp,R2x,X,nx,n_pt_in)
! print*,'nx a-1,c-1 = ',nx ! print*,'nx a-1,c-1 = ',nx
do ix=0,nx do ix=0,nx
X(ix) *= c X(ix) *= dble(c)
enddo enddo
call multiply_poly(X,nx,R2x,2,d,nd) call multiply_poly(X,nx,R2x,2,d,nd)
ny=0 ny=0
@ -446,7 +444,7 @@ recursive subroutine I_x2_pol_mult_mono_elec(c,R1x,R1xp,R2x,d,nd,dim)
call I_x1_pol_mult_mono_elec(0,c-2,R1x,R1xp,R2x,X,nx,dim) call I_x1_pol_mult_mono_elec(0,c-2,R1x,R1xp,R2x,X,nx,dim)
! print*,'nx 0,c-2 = ',nx ! print*,'nx 0,c-2 = ',nx
do ix=0,nx do ix=0,nx
X(ix) *= c-1 X(ix) *= dble(c-1)
enddo enddo
call multiply_poly(X,nx,R2x,2,d,nd) call multiply_poly(X,nx,R2x,2,d,nd)
! print*,'nd = ',nd ! print*,'nd = ',nd

90
src/Integrals_Monoelec/pot_ao_pseudo_ints.irp.f
View File

@ -15,6 +15,7 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
BEGIN_DOC BEGIN_DOC
! Local pseudo-potential ! Local pseudo-potential
END_DOC END_DOC
include 'Utils/constants.include.F'
double precision :: alpha, beta, gama, delta double precision :: alpha, beta, gama, delta
integer :: num_A,num_B integer :: num_A,num_B
double precision :: A_center(3),B_center(3),C_center(3) double precision :: A_center(3),B_center(3),C_center(3)
@ -24,16 +25,10 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
double precision :: cpu_1, cpu_2, wall_1, wall_2, wall_0 double precision :: cpu_1, cpu_2, wall_1, wall_2, wall_0
integer :: thread_num integer :: thread_num
!$ integer :: omp_get_thread_num integer :: omp_get_thread_num
ao_pseudo_integral_local = 0.d0 ao_pseudo_integral_local = 0.d0
!! Dump array
integer, allocatable :: n_k_dump(:)
double precision, allocatable :: v_k_dump(:), dz_k_dump(:)
allocate(n_k_dump(1:pseudo_klocmax), v_k_dump(1:pseudo_klocmax), dz_k_dump(1:pseudo_klocmax))
print*, 'Providing the nuclear electron pseudo integrals (local)' print*, 'Providing the nuclear electron pseudo integrals (local)'
call wall_time(wall_1) call wall_time(wall_1)
@ -44,11 +39,10 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
!$OMP DEFAULT (NONE) & !$OMP DEFAULT (NONE) &
!$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,C_center,power_A,power_B,& !$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,C_center,power_A,power_B,&
!$OMP num_A,num_B,Z,c,n_pt_in, & !$OMP num_A,num_B,Z,c,n_pt_in, &
!$OMP v_k_dump,n_k_dump, dz_k_dump, &
!$OMP wall_0,wall_2,thread_num) & !$OMP wall_0,wall_2,thread_num) &
!$OMP SHARED (ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp,& !$OMP SHARED (ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp,&
!$OMP ao_pseudo_integral_local,nucl_num,nucl_charge, & !$OMP ao_pseudo_integral_local,nucl_num,nucl_charge, &
!$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_v_k,pseudo_n_k, pseudo_dz_k,& !$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_v_k_transp,pseudo_n_k_transp, pseudo_dz_k_transp,&
!$OMP wall_1) !$OMP wall_1)
!$ thread_num = omp_get_thread_num() !$ thread_num = omp_get_thread_num()
@ -75,7 +69,7 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
c = 0.d0 c = 0.d0
if (dabs(ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i))& if (dabs(ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i))&
< 1.d-10) then < thresh) then
cycle cycle
endif endif
do k = 1, nucl_num do k = 1, nucl_num
@ -84,11 +78,10 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
C_center(1:3) = nucl_coord(k,1:3) C_center(1:3) = nucl_coord(k,1:3)
v_k_dump = pseudo_v_k(k,1:pseudo_klocmax) c = c + Vloc(pseudo_klocmax, &
n_k_dump = pseudo_n_k(k,1:pseudo_klocmax) pseudo_v_k_transp (1,k), &
dz_k_dump = pseudo_dz_k(k,1:pseudo_klocmax) pseudo_n_k_transp (1,k), &
pseudo_dz_k_transp(1,k), &
c = c + Vloc(pseudo_klocmax, v_k_dump,n_k_dump, dz_k_dump,&
A_center,power_A,alpha,B_center,power_B,beta,C_center) A_center,power_A,alpha,B_center,power_B,beta,C_center)
enddo enddo
@ -112,8 +105,6 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
!$OMP END PARALLEL !$OMP END PARALLEL
deallocate(n_k_dump,v_k_dump, dz_k_dump)
END_PROVIDER END_PROVIDER
@ -122,25 +113,20 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
BEGIN_DOC BEGIN_DOC
! Local pseudo-potential ! Local pseudo-potential
END_DOC END_DOC
include 'Utils/constants.include.F'
double precision :: alpha, beta, gama, delta double precision :: alpha, beta, gama, delta
integer :: num_A,num_B integer :: num_A,num_B
double precision :: A_center(3),B_center(3),C_center(3) double precision :: A_center(3),B_center(3),C_center(3)
integer :: power_A(3),power_B(3) integer :: power_A(3),power_B(3)
integer :: i,j,k,l,n_pt_in,m integer :: i,j,k,l,n_pt_in,m
double precision :: Vloc, Vpseudo double precision :: Vloc, Vpseudo
!$ integer :: omp_get_thread_num integer :: omp_get_thread_num
double precision :: cpu_1, cpu_2, wall_1, wall_2, wall_0 double precision :: cpu_1, cpu_2, wall_1, wall_2, wall_0
integer :: thread_num integer :: thread_num
ao_pseudo_integral_non_local = 0.d0 ao_pseudo_integral_non_local = 0.d0
!! Dump array
integer, allocatable :: n_kl_dump(:,:)
double precision, allocatable :: v_kl_dump(:,:), dz_kl_dump(:,:)
allocate(n_kl_dump(pseudo_kmax,0:pseudo_lmax), v_kl_dump(pseudo_kmax,0:pseudo_lmax), dz_kl_dump(pseudo_kmax,0:pseudo_lmax))
print*, 'Providing the nuclear electron pseudo integrals (non-local)' print*, 'Providing the nuclear electron pseudo integrals (non-local)'
call wall_time(wall_1) call wall_time(wall_1)
@ -150,14 +136,14 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
!$OMP DEFAULT (NONE) & !$OMP DEFAULT (NONE) &
!$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,C_center,power_A,power_B,& !$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,C_center,power_A,power_B,&
!$OMP num_A,num_B,Z,c,n_pt_in, & !$OMP num_A,num_B,Z,c,n_pt_in, &
!$OMP n_kl_dump, v_kl_dump, dz_kl_dump, &
!$OMP wall_0,wall_2,thread_num) & !$OMP wall_0,wall_2,thread_num) &
!$OMP SHARED (ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp,& !$OMP SHARED (ao_num,ao_prim_num,ao_expo_ordered_transp,ao_power,ao_nucl,nucl_coord,ao_coef_normalized_ordered_transp,&
!$OMP ao_pseudo_integral_non_local,nucl_num,nucl_charge,& !$OMP ao_pseudo_integral_non_local,nucl_num,nucl_charge,&
!$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_n_kl, pseudo_v_kl, pseudo_dz_kl,& !$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_n_kl_transp, pseudo_v_kl_transp, pseudo_dz_kl_transp,&
!$OMP wall_1) !$OMP wall_1)
!$ thread_num = omp_get_thread_num() !$ thread_num = omp_get_thread_num()
!$OMP DO SCHEDULE (guided) !$OMP DO SCHEDULE (guided)
do j = 1, ao_num do j = 1, ao_num
@ -181,7 +167,7 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
c = 0.d0 c = 0.d0
if (dabs(ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i))& if (dabs(ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i))&
< 1.d-10) then < thresh) then
cycle cycle
endif endif
@ -191,12 +177,11 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
C_center(1:3) = nucl_coord(k,1:3) C_center(1:3) = nucl_coord(k,1:3)
n_kl_dump = pseudo_n_kl(k,1:pseudo_kmax,0:pseudo_lmax) c = c + Vpseudo(pseudo_lmax,pseudo_kmax, &
v_kl_dump = pseudo_v_kl(k,1:pseudo_kmax,0:pseudo_lmax) pseudo_v_kl_transp(1,0,k), &
dz_kl_dump = pseudo_dz_kl(k,1:pseudo_kmax,0:pseudo_lmax) pseudo_n_kl_transp(1,0,k), &
pseudo_dz_kl_transp(1,0,k), &
c = c + Vpseudo(pseudo_lmax,pseudo_kmax,v_kl_dump,n_kl_dump,dz_kl_dump,A_center,power_A,alpha,B_center,power_B,beta,C_center) A_center,power_A,alpha,B_center,power_B,beta,C_center)
enddo enddo
ao_pseudo_integral_non_local(i,j) = ao_pseudo_integral_non_local(i,j) +& ao_pseudo_integral_non_local(i,j) = ao_pseudo_integral_non_local(i,j) +&
ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i)*c ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i)*c
@ -215,13 +200,48 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral_local, (ao_num_align,ao_nu
enddo enddo
!$OMP END DO !$OMP END DO
!$OMP END PARALLEL !$OMP END PARALLEL
deallocate(n_kl_dump,v_kl_dump, dz_kl_dump)
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [ double precision, pseudo_v_k_transp, (pseudo_klocmax,nucl_num) ]
&BEGIN_PROVIDER [ integer , pseudo_n_k_transp, (pseudo_klocmax,nucl_num) ]
&BEGIN_PROVIDER [ double precision, pseudo_dz_k_transp, (pseudo_klocmax,nucl_num)]
implicit none
BEGIN_DOC
! Transposed arrays for pseudopotentials
END_DOC
integer :: i,j
do j=1,nucl_num
do i=1,pseudo_klocmax
pseudo_v_k_transp (i,j) = pseudo_v_k (j,i)
pseudo_n_k_transp (i,j) = pseudo_n_k (j,i)
pseudo_dz_k_transp(i,j) = pseudo_dz_k(j,i)
enddo
enddo
END_PROVIDER
BEGIN_PROVIDER [ double precision, pseudo_v_kl_transp, (pseudo_kmax,0:pseudo_lmax,nucl_num) ]
&BEGIN_PROVIDER [ integer , pseudo_n_kl_transp, (pseudo_kmax,0:pseudo_lmax,nucl_num) ]
&BEGIN_PROVIDER [ double precision, pseudo_dz_kl_transp, (pseudo_kmax,0:pseudo_lmax,nucl_num)]
implicit none
BEGIN_DOC
! Transposed arrays for pseudopotentials
END_DOC
integer :: i,j,l
do j=1,nucl_num
do l=0,pseudo_lmax
do i=1,pseudo_kmax
pseudo_v_kl_transp (i,l,j) = pseudo_v_kl (j,i,l)
pseudo_n_kl_transp (i,l,j) = pseudo_n_kl (j,i,l)
pseudo_dz_kl_transp(i,l,j) = pseudo_dz_kl(j,i,l)
enddo
enddo
enddo
END_PROVIDER

12
src/Integrals_Monoelec/spread_dipole_ao.irp.f
View File

@ -55,11 +55,11 @@
beta = ao_expo_ordered_transp(l,i) beta = ao_expo_ordered_transp(l,i)
call overlap_gaussian_xyz(A_center,B_center,alpha,beta,power_A,power_B,overlap_x,overlap_y,overlap_z,overlap,dim1) call overlap_gaussian_xyz(A_center,B_center,alpha,beta,power_A,power_B,overlap_x,overlap_y,overlap_z,overlap,dim1)
call overlap_bourrin_spread(A_center(1),B_center(1),alpha,beta,power_A(1),power_B(1),tmp,lower_exp_val,dx,dim1) call overlap_bourrin_spread(A_center(1),B_center(1),alpha,beta,power_A(1),power_B(1),tmp,lower_exp_val,dx,dim1)
accu_x += c*(tmp*overlap_y*overlap_z) accu_x += c*tmp*overlap_y*overlap_z
call overlap_bourrin_spread(A_center(2),B_center(2),alpha,beta,power_A(2),power_B(2),tmp,lower_exp_val,dx,dim1) call overlap_bourrin_spread(A_center(2),B_center(2),alpha,beta,power_A(2),power_B(2),tmp,lower_exp_val,dx,dim1)
accu_y += c*(tmp*overlap_x*overlap_z) accu_y += c*tmp*overlap_x*overlap_z
call overlap_bourrin_spread(A_center(3),B_center(3),alpha,beta,power_A(3),power_B(3),tmp,lower_exp_val,dx,dim1) call overlap_bourrin_spread(A_center(3),B_center(3),alpha,beta,power_A(3),power_B(3),tmp,lower_exp_val,dx,dim1)
accu_z += c*(tmp*overlap_y*overlap_x) accu_z += c*tmp*overlap_y*overlap_x
enddo enddo
enddo enddo
ao_spread_x(i,j) = accu_x ao_spread_x(i,j) = accu_x
@ -130,11 +130,11 @@
call overlap_gaussian_xyz(A_center,B_center,alpha,beta,power_A,power_B,overlap_x,overlap_y,overlap_z,overlap,dim1) call overlap_gaussian_xyz(A_center,B_center,alpha,beta,power_A,power_B,overlap_x,overlap_y,overlap_z,overlap,dim1)
call overlap_bourrin_dipole(A_center(1),B_center(1),alpha,beta,power_A(1),power_B(1),tmp,lower_exp_val,dx,dim1) call overlap_bourrin_dipole(A_center(1),B_center(1),alpha,beta,power_A(1),power_B(1),tmp,lower_exp_val,dx,dim1)
accu_x = accu_x + c*(tmp*overlap_y*overlap_z) accu_x = accu_x + c*tmp*overlap_y*overlap_z
call overlap_bourrin_dipole(A_center(2),B_center(2),alpha,beta,power_A(2),power_B(2),tmp,lower_exp_val,dx,dim1) call overlap_bourrin_dipole(A_center(2),B_center(2),alpha,beta,power_A(2),power_B(2),tmp,lower_exp_val,dx,dim1)
accu_y = accu_y + c*(tmp*overlap_x*overlap_z) accu_y = accu_y + c*tmp*overlap_x*overlap_z
call overlap_bourrin_dipole(A_center(3),B_center(3),alpha,beta,power_A(3),power_B(3),tmp,lower_exp_val,dx,dim1) call overlap_bourrin_dipole(A_center(3),B_center(3),alpha,beta,power_A(3),power_B(3),tmp,lower_exp_val,dx,dim1)
accu_z = accu_z + c*(tmp*overlap_y*overlap_x) accu_z = accu_z + c*tmp*overlap_y*overlap_x
enddo enddo
enddo enddo
ao_dipole_x(i,j) = accu_x ao_dipole_x(i,j) = accu_x

0
src/Integrals_Monoelec/var_pt2_ratio.irp.f
View File

2
src/MOGuess/.gitignore vendored
View File

@ -15,6 +15,8 @@ Nuclei
Pseudo Pseudo
Utils Utils
ezfio_interface.irp.f ezfio_interface.irp.f
guess_overlap
irpf90.make irpf90.make
irpf90_entities irpf90_entities
tags tags
truncate_mos

2
src/MOGuess/H_CORE_guess.irp.f
View File

@ -5,8 +5,6 @@ program H_CORE_guess
END_DOC END_DOC
implicit none implicit none
character*(64) :: label character*(64) :: label
mo_coef = ao_ortho_lowdin_coef
TOUCH mo_coef
label = "Guess" label = "Guess"
call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral, & call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral, &
size(mo_mono_elec_integral,1), & size(mo_mono_elec_integral,1), &

8
src/MOGuess/README.rst
View File

@ -43,6 +43,10 @@ Documentation
supposed to be the Identity supposed to be the Identity
`guess_mimi <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess/guess_overlap.irp.f#L1>`_
Produce `H_core` MO orbital
`h_core_guess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess/H_CORE_guess.irp.f#L1>`_ `h_core_guess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess/H_CORE_guess.irp.f#L1>`_
Produce `H_core` MO orbital Produce `H_core` MO orbital
output: mo_basis.mo_tot_num mo_basis.mo_label mo_basis.ao_md5 mo_basis.mo_coef mo_basis.mo_occ output: mo_basis.mo_tot_num mo_basis.mo_label mo_basis.ao_md5 mo_basis.mo_coef mo_basis.mo_occ
@ -51,3 +55,7 @@ Documentation
`hcore_guess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess/h_core_guess_routine.irp.f#L1>`_ `hcore_guess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess/h_core_guess_routine.irp.f#L1>`_
Produce `H_core` MO orbital Produce `H_core` MO orbital
`prog_truncate_mo <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess/truncate_mos.irp.f#L1>`_
Truncate MO set

2
src/MOGuess/guess_overlap.irp.f
View File

@ -5,8 +5,6 @@ program guess_mimi
implicit none implicit none
character*(64) :: label character*(64) :: label
mo_coef = ao_ortho_lowdin_coef
TOUCH mo_coef
label = "Guess" label = "Guess"
call mo_as_eigvectors_of_mo_matrix(ao_overlap, & call mo_as_eigvectors_of_mo_matrix(ao_overlap, &
size(ao_overlap,1), & size(ao_overlap,1), &

2
src/MOGuess/h_core_guess_routine.irp.f
View File

@ -4,8 +4,6 @@ subroutine hcore_guess
END_DOC END_DOC
implicit none implicit none
character*(64) :: label character*(64) :: label
mo_coef = ao_ortho_lowdin_coef
TOUCH mo_coef
label = "Guess" label = "Guess"
call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral, & call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral, &
size(mo_mono_elec_integral,1), & size(mo_mono_elec_integral,1), &

12
src/MOGuess/mo_ortho_lowdin.irp.f
View File

@ -1,4 +1,3 @@
BEGIN_PROVIDER [double precision, ao_ortho_lowdin_coef, (ao_num_align,ao_num)] BEGIN_PROVIDER [double precision, ao_ortho_lowdin_coef, (ao_num_align,ao_num)]
implicit none implicit none
BEGIN_DOC BEGIN_DOC
@ -8,13 +7,9 @@ BEGIN_PROVIDER [double precision, ao_ortho_lowdin_coef, (ao_num_align,ao_num)]
END_DOC END_DOC
integer :: i,j,k,l integer :: i,j,k,l
double precision :: tmp_matrix(ao_num_align,ao_num),accu double precision :: tmp_matrix(ao_num_align,ao_num),accu
tmp_matrix(:,:) = 0.d0
do j=1, ao_num do j=1, ao_num
do i=1, ao_num tmp_matrix(j,j) = 1.d0
tmp_matrix(i,j) = 0.d0
enddo
enddo
do i=1, ao_num
tmp_matrix(i,i) = 1.d0
enddo enddo
call ortho_lowdin(ao_overlap,ao_num_align,ao_num,tmp_matrix,ao_num_align,ao_num) call ortho_lowdin(ao_overlap,ao_num_align,ao_num,tmp_matrix,ao_num_align,ao_num)
do i=1, ao_num do i=1, ao_num
@ -23,6 +18,7 @@ BEGIN_PROVIDER [double precision, ao_ortho_lowdin_coef, (ao_num_align,ao_num)]
enddo enddo
enddo enddo
END_PROVIDER END_PROVIDER
BEGIN_PROVIDER [double precision, ao_ortho_lowdin_overlap, (ao_num_align,ao_num)] BEGIN_PROVIDER [double precision, ao_ortho_lowdin_overlap, (ao_num_align,ao_num)]
implicit none implicit none
BEGIN_DOC BEGIN_DOC
@ -40,7 +36,7 @@ BEGIN_PROVIDER [double precision, ao_ortho_lowdin_overlap, (ao_num_align,ao_num)
do j=1, ao_num do j=1, ao_num
c = 0.d0 c = 0.d0
do l=1, ao_num do l=1, ao_num
c = ao_ortho_lowdin_coef(j,l) * ao_overlap(k,l) c += ao_ortho_lowdin_coef(j,l) * ao_overlap(k,l)
enddo enddo
do i=1, ao_num do i=1, ao_num
ao_ortho_lowdin_overlap(i,j) += ao_ortho_lowdin_coef(i,k) * c ao_ortho_lowdin_overlap(i,j) += ao_ortho_lowdin_coef(i,k) * c

25
src/MOGuess/pot_mo_ortho_canonical_ints.irp.f Normal file
View File

@ -0,0 +1,25 @@
BEGIN_PROVIDER [double precision, ao_ortho_canonical_nucl_elec_integral, (mo_tot_num_align,mo_tot_num)]
implicit none
integer :: i1,j1,i,j
double precision :: c_i1,c_j1
ao_ortho_canonical_nucl_elec_integral = 0.d0
!$OMP PARALLEL DO DEFAULT(none) &
!$OMP PRIVATE(i,j,i1,j1,c_j1,c_i1) &
!$OMP SHARED(mo_tot_num,ao_num,ao_ortho_canonical_coef, &
!$OMP ao_ortho_canonical_nucl_elec_integral, ao_nucl_elec_integral)
do i = 1, mo_tot_num
do j = 1, mo_tot_num
do i1 = 1,ao_num
c_i1 = ao_ortho_canonical_coef(i1,i)
do j1 = 1,ao_num
c_j1 = c_i1*ao_ortho_canonical_coef(j1,j)
ao_ortho_canonical_nucl_elec_integral(j,i) = ao_ortho_canonical_nucl_elec_integral(j,i) + &
c_j1 * ao_nucl_elec_integral(j1,i1)
enddo
enddo
enddo
enddo
!$OMP END PARALLEL DO
END_PROVIDER

10
src/MO_Basis/README.rst
View File

@ -76,11 +76,11 @@ Documentation
by convention, the '-' MO is in the greater index (max(j,k)) by convention, the '-' MO is in the greater index (max(j,k))
`mo_as_eigvectors_of_mo_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis/utils.irp.f#L24>`_ `mo_as_eigvectors_of_mo_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis/utils.irp.f#L47>`_
Undocumented Undocumented
`mo_as_eigvectors_of_mo_matrix_sort_by_observable <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis/utils.irp.f#L62>`_ `mo_as_eigvectors_of_mo_matrix_sort_by_observable <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis/utils.irp.f#L103>`_
Undocumented Undocumented
@ -116,7 +116,7 @@ Documentation
Undocumented Undocumented
`mo_sort_by_observable <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis/utils.irp.f#L144>`_ `mo_sort_by_observable <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis/utils.irp.f#L185>`_
Undocumented Undocumented
@ -143,3 +143,7 @@ Documentation
`save_mos <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis/utils.irp.f#L1>`_ `save_mos <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis/utils.irp.f#L1>`_
Undocumented Undocumented
`save_mos_truncated <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis/utils.irp.f#L24>`_
Undocumented

41
src/MO_Basis/ao_ortho_canonical.irp.f Normal file
View File

@ -0,0 +1,41 @@
BEGIN_PROVIDER [ double precision, ao_ortho_canonical_coef, (ao_num_align,ao_num)]
&BEGIN_PROVIDER [ integer, ao_ortho_canonical_num ]
implicit none
BEGIN_DOC
! matrix of the coefficients of the mos generated by the
! orthonormalization by the S^{-1/2} canonical transformation of the aos
! ao_ortho_canonical_coef(i,j) = coefficient of the ith ao on the jth ao_ortho_canonical orbital
END_DOC
integer :: i
ao_ortho_canonical_coef(:,:) = 0.d0
do i=1,ao_num
ao_ortho_canonical_coef(i,i) = 1.d0
enddo
call ortho_canonical(ao_overlap,ao_num_align,ao_num,ao_ortho_canonical_coef,ao_num_align,ao_ortho_canonical_num)
END_PROVIDER
BEGIN_PROVIDER [double precision, ao_ortho_canonical_overlap, (ao_ortho_canonical_num,ao_ortho_canonical_num)]
implicit none
BEGIN_DOC
! overlap matrix of the ao_ortho_canonical
! supposed to be the Identity
END_DOC
integer :: i,j,k,l
double precision :: c
do j=1, ao_ortho_canonical_num
do i=1, ao_ortho_canonical_num
ao_ortho_canonical_overlap(i,j) = 0.d0
enddo
enddo
do j=1, ao_ortho_canonical_num
do k=1, ao_num
c = 0.d0
do l=1, ao_num
c += ao_ortho_canonical_coef(l,j) * ao_overlap(l,k)
enddo
do i=1, ao_ortho_canonical_num
ao_ortho_canonical_overlap(i,j) += ao_ortho_canonical_coef(k,i) * c
enddo
enddo
enddo
END_PROVIDER

2
src/MO_Basis/mo_overlap.irp.f
View File

@ -4,7 +4,7 @@ BEGIN_PROVIDER [ double precision, mo_overlap,(mo_tot_num_align,mo_tot_num)]
integer :: i,j,n,l integer :: i,j,n,l
double precision :: f double precision :: f
integer :: lmax integer :: lmax
lmax = iand(ao_num,-4) lmax = (ao_num/4) * 4
!$OMP PARALLEL DO SCHEDULE(STATIC) DEFAULT(NONE) & !$OMP PARALLEL DO SCHEDULE(STATIC) DEFAULT(NONE) &
!$OMP PRIVATE(i,j,n,l) & !$OMP PRIVATE(i,j,n,l) &
!$OMP SHARED(mo_overlap,mo_coef,ao_overlap, & !$OMP SHARED(mo_overlap,mo_coef,ao_overlap, &

12
src/MO_Basis/mos.irp.f
View File

@ -9,8 +9,9 @@ BEGIN_PROVIDER [ integer, mo_tot_num ]
if (exists) then if (exists) then
call ezfio_get_mo_basis_mo_tot_num(mo_tot_num) call ezfio_get_mo_basis_mo_tot_num(mo_tot_num)
else else
mo_tot_num = ao_num mo_tot_num = ao_ortho_canonical_num
endif endif
call write_int(6,mo_tot_num,'mo_tot_num')
ASSERT (mo_tot_num > 0) ASSERT (mo_tot_num > 0)
END_PROVIDER END_PROVIDER
@ -56,7 +57,14 @@ END_PROVIDER
deallocate(buffer) deallocate(buffer)
else else
! Orthonormalized AO basis ! Orthonormalized AO basis
mo_coef = 0. do i=1,mo_tot_num
do j=1,ao_num
mo_coef(j,i) = ao_ortho_canonical_coef(j,i)
enddo
do j=ao_num+1,ao_num_align
mo_coef(j,i) = 0.d0
enddo
enddo
endif endif
END_PROVIDER END_PROVIDER

47
src/MO_Basis/utils.irp.f
View File

@ -100,6 +100,53 @@ subroutine mo_as_eigvectors_of_mo_matrix(matrix,n,m,label,sign)
mo_label = label mo_label = label
end end
subroutine mo_as_svd_vectors_of_mo_matrix(matrix,lda,m,n,label)
implicit none
integer,intent(in) :: lda,m,n
character*(64), intent(in) :: label
double precision, intent(in) :: matrix(lda,n)
integer :: i,j
double precision, allocatable :: mo_coef_new(:,:), U(:,:),D(:), A(:,:), Vt(:,:), work(:)
!DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: mo_coef_new, U, Vt, A
call write_time(output_mo_basis)
if (m /= mo_tot_num) then
print *, irp_here, ': Error : m/= mo_tot_num'
stop 1
endif
allocate(A(lda,n),U(lda,n),mo_coef_new(ao_num_align,m),D(m),Vt(lda,n))
do j=1,n
do i=1,m
A(i,j) = matrix(i,j)
enddo
enddo
mo_coef_new = mo_coef
call svd(A,lda,U,lda,D,Vt,lda,m,n)
write (output_mo_basis,'(A)'), 'MOs are now **'//trim(label)//'**'
write (output_mo_basis,'(A)'), ''
write (output_mo_basis,'(A)'), 'Eigenvalues'
write (output_mo_basis,'(A)'), '-----------'
write (output_mo_basis,'(A)'), ''
write (output_mo_basis,'(A)'), '======== ================'
do i=1,m
write (output_mo_basis,'(I8,X,F16.10)'), i,D(i)
enddo
write (output_mo_basis,'(A)'), '======== ================'
write (output_mo_basis,'(A)'), ''
call dgemm('N','N',ao_num,m,m,1.d0,mo_coef_new,size(mo_coef_new,1),U,size(U,1),0.d0,mo_coef,size(mo_coef,1))
deallocate(A,mo_coef_new,U,Vt,D)
call write_time(output_mo_basis)
mo_label = label
end
subroutine mo_as_eigvectors_of_mo_matrix_sort_by_observable(matrix,observable,n,m,label) subroutine mo_as_eigvectors_of_mo_matrix_sort_by_observable(matrix,observable,n,m,label)
implicit none implicit none
integer,intent(in) :: n,m integer,intent(in) :: n,m

7
src/Pseudo/EZFIO.cfg
View File

@ -1,3 +1,10 @@
[nucl_charge_remove]
doc: Nuclear charges removed
type:double precision
size: (nuclei.nucl_num)
interface: ezfio, provider
[pseudo_klocmax] [pseudo_klocmax]
doc: test doc: test
type:integer type:integer

90
src/Utils/LinearAlgebra.irp.f
View File

@ -42,11 +42,90 @@ subroutine svd(A,LDA,U,LDU,D,Vt,LDVt,m,n)
end end
subroutine ortho_canonical(overlap,LDA,N,C,LDC,m)
implicit none
BEGIN_DOC
! Compute C_new=C_old.U.s^-1/2 canonical orthogonalization.
!
! overlap : overlap matrix
!
! LDA : leftmost dimension of overlap array
!
! N : Overlap matrix is NxN (array is (LDA,N) )
!
! C : Coefficients of the vectors to orthogonalize. On exit,
! orthogonal vectors
!
! LDC : leftmost dimension of C
!
! m : Coefficients matrix is MxN, ( array is (LDC,N) )
!
END_DOC
integer, intent(in) :: lda, ldc, n
integer, intent(out) :: m
double precision, intent(in) :: overlap(lda,n)
double precision, intent(inout) :: C(ldc,n)
double precision, allocatable :: U(:,:)
double precision, allocatable :: Vt(:,:)
double precision, allocatable :: D(:)
double precision, allocatable :: S_half(:,:)
!DEC$ ATTRIBUTES ALIGN : 64 :: U, Vt, D
integer :: info, i, j
allocate (U(ldc,n), Vt(lda,n), D(n), S_half(lda,n))
call svd(overlap,lda,U,ldc,D,Vt,lda,n,n)
m=n
do i=1,n
if ( D(i) >= 1.d-6 ) then
D(i) = 1.d0/dsqrt(D(i))
else
m = i-1
print *, 'Removed Linear dependencies below:', 1.d0/D(m)
exit
endif
enddo
do i=m+1,n
print *, D(i)
D(i) = 0.d0
enddo
do i=1,m
if ( D(i) >= 1.d5 ) then
print *, 'Warning: Basis set may have linear dependence problems'
endif
enddo
!$OMP PARALLEL DEFAULT(NONE) &
!$OMP SHARED(S_half,U,D,Vt,n,C,m) &
!$OMP PRIVATE(i,j)
!$OMP DO
do j=1,n
do i=1,n
S_half(i,j) = U(i,j)*D(j)
enddo
do i=1,n
U(i,j) = C(i,j)
enddo
enddo
!$OMP END DO
!$OMP END PARALLEL
call dgemm('N','N',n,m,n,1.d0,U,size(U,1),S_half,size(S_half,1),0.d0,C,size(C,1))
deallocate (U, Vt, D, S_half)
end
subroutine ortho_lowdin(overlap,LDA,N,C,LDC,m) subroutine ortho_lowdin(overlap,LDA,N,C,LDC,m)
implicit none implicit none
BEGIN_DOC BEGIN_DOC
! Compute C_new=C_old.S^-1/2 canonical orthogonalization. ! Compute C_new=C_old.S^-1/2 orthogonalization.
! !
! overlap : overlap matrix ! overlap : overlap matrix
! !
@ -70,9 +149,8 @@ subroutine ortho_lowdin(overlap,LDA,N,C,LDC,m)
double precision :: Vt(lda,n) double precision :: Vt(lda,n)
double precision :: D(n) double precision :: D(n)
double precision :: S_half(lda,n) double precision :: S_half(lda,n)
double precision,allocatable :: work(:) !DEC$ ATTRIBUTES ALIGN : 64 :: U, Vt, D
!DEC$ ATTRIBUTES ALIGN : 64 :: U, Vt, D, work integer :: info, i, j, k
integer :: info, lwork, i, j, k
call svd(overlap,lda,U,ldc,D,Vt,lda,m,n) call svd(overlap,lda,U,ldc,D,Vt,lda,m,n)
@ -82,7 +160,7 @@ subroutine ortho_lowdin(overlap,LDA,N,C,LDC,m)
!$OMP DO !$OMP DO
do i=1,n do i=1,n
if ( D(i) < 1.d-12 ) then if ( D(i) < 1.d-6 ) then
D(i) = 0.d0 D(i) = 0.d0
else else
D(i) = 1.d0/dsqrt(D(i)) D(i) = 1.d0/dsqrt(D(i))
@ -94,6 +172,7 @@ subroutine ortho_lowdin(overlap,LDA,N,C,LDC,m)
!$OMP END DO !$OMP END DO
do k=1,n do k=1,n
if (D(k) /= 0.d0) then
!$OMP DO !$OMP DO
do j=1,n do j=1,n
do i=1,n do i=1,n
@ -101,6 +180,7 @@ subroutine ortho_lowdin(overlap,LDA,N,C,LDC,m)
enddo enddo
enddo enddo
!$OMP END DO NOWAIT !$OMP END DO NOWAIT
endif
enddo enddo
!$OMP BARRIER !$OMP BARRIER

27
src/Utils/README.rst
View File

@ -36,7 +36,7 @@ Documentation
Compute 1st dimension such that it is aligned for vectorization. Compute 1st dimension such that it is aligned for vectorization.
`apply_rotation <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L168>`_ `apply_rotation <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L196>`_
Apply the rotation found by find_rotation Apply the rotation found by find_rotation
@ -122,7 +122,7 @@ Documentation
1/n! 1/n!
`find_rotation <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L149>`_ `find_rotation <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L177>`_
Find A.C = B Find A.C = B
@ -148,7 +148,7 @@ Documentation
Undocumented Undocumented
`get_pseudo_inverse <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L95>`_ `get_pseudo_inverse <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L123>`_
Find C = A^-1 Find C = A^-1
@ -420,7 +420,7 @@ Documentation
contains the new order of the elements. contains the new order of the elements.
`lapack_diag <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L247>`_ `lapack_diag <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L275>`_
Diagonalize matrix H Diagonalize matrix H
.br .br
H is untouched between input and ouptut H is untouched between input and ouptut
@ -431,7 +431,7 @@ Documentation
.br .br
`lapack_diag_s2 <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L310>`_ `lapack_diag_s2 <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L338>`_
Diagonalize matrix H Diagonalize matrix H
.br .br
H is untouched between input and ouptut H is untouched between input and ouptut
@ -442,7 +442,7 @@ Documentation
.br .br
`lapack_diagd <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L180>`_ `lapack_diagd <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L208>`_
Diagonalize matrix H Diagonalize matrix H
.br .br
H is untouched between input and ouptut H is untouched between input and ouptut
@ -453,7 +453,7 @@ Documentation
.br .br
`lapack_partial_diag <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L376>`_ `lapack_partial_diag <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L404>`_
Diagonalize matrix H Diagonalize matrix H
.br .br
H is untouched between input and ouptut H is untouched between input and ouptut
@ -482,7 +482,7 @@ Documentation
Number of current OpenMP threads Number of current OpenMP threads
`ortho_lowdin <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L1>`_ `ortho_lowdin <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L46>`_
Compute C_new=C_old.S^-1/2 canonical orthogonalization. Compute C_new=C_old.S^-1/2 canonical orthogonalization.
.br .br
overlap : overlap matrix overlap : overlap matrix
@ -597,7 +597,7 @@ Documentation
to be in integer*8 format to be in integer*8 format
`set_zero_extra_diag <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L433>`_ `set_zero_extra_diag <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L461>`_
Undocumented Undocumented
@ -615,6 +615,15 @@ Documentation
Stop the progress bar Stop the progress bar
`svd <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/LinearAlgebra.irp.f#L1>`_
Compute A = U.D.Vt
.br
LDx : leftmost dimension of x
.br
Dimsneion of A is m x n
.br
`trap_signals <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/abort.irp.f#L19>`_ `trap_signals <http://github.com/LCPQ/quantum_package/tree/master/src/Utils/abort.irp.f#L19>`_
What to do when a signal is caught. Here, trap Ctrl-C and call the control_C subroutine. What to do when a signal is caught. Here, trap Ctrl-C and call the control_C subroutine.

3
src/Utils/constants.include.F
View File

@ -1,4 +1,4 @@
integer, parameter :: max_dim = 255 integer, parameter :: max_dim = 511
integer, parameter :: SIMD_vector = 32 integer, parameter :: SIMD_vector = 32
double precision, parameter :: pi = dacos(-1.d0) double precision, parameter :: pi = dacos(-1.d0)
@ -8,3 +8,4 @@ double precision, parameter :: dfour_pi = 4.d0*dacos(-1.d0)
double precision, parameter :: dtwo_pi = 2.d0*dacos(-1.d0) double precision, parameter :: dtwo_pi = 2.d0*dacos(-1.d0)
double precision, parameter :: inv_sq_pi = 1.d0/dsqrt(dacos(-1.d0)) double precision, parameter :: inv_sq_pi = 1.d0/dsqrt(dacos(-1.d0))
double precision, parameter :: inv_sq_pi_2 = 0.5d0/dsqrt(dacos(-1.d0)) double precision, parameter :: inv_sq_pi_2 = 0.5d0/dsqrt(dacos(-1.d0))
double precision, parameter :: thresh = 1.d-15

8
src/Utils/integration.irp.f
View File

@ -78,7 +78,7 @@ subroutine give_explicit_poly_and_gaussian(P_new,P_center,p,fact_k,iorder,alpha,
!DEC$ FORCEINLINE !DEC$ FORCEINLINE
call gaussian_product(alpha,A_center,beta,B_center,fact_k,p,P_center) call gaussian_product(alpha,A_center,beta,B_center,fact_k,p,P_center)
if (fact_k < 1.d-8) then if (fact_k < thresh) then
fact_k = 0.d0 fact_k = 0.d0
return return
endif endif
@ -210,7 +210,7 @@ subroutine gaussian_product(a,xa,b,xb,k,p,xp)
xab(3) = xa(3)-xb(3) xab(3) = xa(3)-xb(3)
ab = ab*p_inv ab = ab*p_inv
k = ab*(xab(1)*xab(1)+xab(2)*xab(2)+xab(3)*xab(3)) k = ab*(xab(1)*xab(1)+xab(2)*xab(2)+xab(3)*xab(3))
if (k > 20.d0) then if (k > 40.d0) then
k=0.d0 k=0.d0
return return
endif endif
@ -249,7 +249,7 @@ subroutine gaussian_product_x(a,xa,b,xb,k,p,xp)
xab = xa-xb xab = xa-xb
ab = ab*p_inv ab = ab*p_inv
k = ab*xab*xab k = ab*xab*xab
if (k > 20.d0) then if (k > 40.d0) then
k=0.d0 k=0.d0
return return
endif endif
@ -580,7 +580,7 @@ double precision function rint_large_n(n,rho)
enddo enddo
t2=0.d0 t2=0.d0
do l=0,k do l=0,k
t2=t2+(-1.d0)**l/fact(l+1)/fact(k-l) t2=t2+(-1.d0)**l/(fact(l+1)*fact(k-l))
enddo enddo
alpha_k=t2*fact(k+1)*fact(k)*(-1.d0)**k alpha_k=t2*fact(k+1)*fact(k)*(-1.d0)**k
alpha_k= alpha_k/t1 alpha_k= alpha_k/t1

6
src/Utils/one_e_integration.irp.f
View File

@ -150,19 +150,19 @@ subroutine overlap_gaussian_xyz(A_center,B_center,alpha,beta,power_A,&
integer :: i integer :: i
do i = 1,iorder_p(1) do i = 1,iorder_p(1)
overlap_x += P_new(i,1) * F_integral_tab(i) overlap_x = overlap_x + P_new(i,1) * F_integral_tab(i)
enddo enddo
call gaussian_product_x(alpha,A_center(1),beta,B_center(1),fact_p,p,P_center(1)) call gaussian_product_x(alpha,A_center(1),beta,B_center(1),fact_p,p,P_center(1))
overlap_x *= fact_p overlap_x *= fact_p
do i = 1,iorder_p(2) do i = 1,iorder_p(2)
overlap_y += P_new(i,2) * F_integral_tab(i) overlap_y = overlap_y + P_new(i,2) * F_integral_tab(i)
enddo enddo
call gaussian_product_x(alpha,A_center(2),beta,B_center(2),fact_p,p,P_center(2)) call gaussian_product_x(alpha,A_center(2),beta,B_center(2),fact_p,p,P_center(2))
overlap_y *= fact_p overlap_y *= fact_p
do i = 1,iorder_p(3) do i = 1,iorder_p(3)
overlap_z += P_new(i,3) * F_integral_tab(i) overlap_z = overlap_z + P_new(i,3) * F_integral_tab(i)
enddo enddo
call gaussian_product_x(alpha,A_center(3),beta,B_center(3),fact_p,p,P_center(3)) call gaussian_product_x(alpha,A_center(3),beta,B_center(3),fact_p,p,P_center(3))
overlap_z *= fact_p overlap_z *= fact_p

2
src/ZMQ/NEEDED_CHILDREN_MODULES
View File

@ -1 +1 @@
Utils

35
src/ZMQ/README.rst
View File

@ -5,11 +5,38 @@ ZMQ
Socket address : defined as an environment variable : QP_RUN_ADDRESS Socket address : defined as an environment variable : QP_RUN_ADDRESS
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
Documentation Documentation
============= =============
.. Do not edit this section It was auto-generated .. Do not edit this section It was auto-generated
.. by the `update_README.py` script. .. by the `update_README.py` script.
`qp_run_address <http://github.com/LCPQ/quantum_package/tree/master/src/ZMQ/zmq.irp.f#L13>`_
Address of the qp_run socket
Example : tcp://130.120.229.139:12345
`zmq_context <http://github.com/LCPQ/quantum_package/tree/master/src/ZMQ/zmq.irp.f#L4>`_
Context for the ZeroMQ library
`zmq_port <http://github.com/LCPQ/quantum_package/tree/master/src/ZMQ/zmq.irp.f#L38>`_
Undocumented
`zmq_port_start <http://github.com/LCPQ/quantum_package/tree/master/src/ZMQ/zmq.irp.f#L14>`_
Address of the qp_run socket
Example : tcp://130.120.229.139:12345
`zmq_socket_pull <http://github.com/LCPQ/quantum_package/tree/master/src/ZMQ/zmq.irp.f#L87>`_
Socket which pulls the results (2)
`zmq_socket_push <http://github.com/LCPQ/quantum_package/tree/master/src/ZMQ/zmq.irp.f#L70>`_
Socket on which to push the results (1)
`zmq_to_qp_run_socket <http://github.com/LCPQ/quantum_package/tree/master/src/ZMQ/zmq.irp.f#L47>`_
Socket on which the qp_run process replies

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