-
Notifications
You must be signed in to change notification settings - Fork 156
/
1e_mints-helper.ipynb
211 lines (211 loc) · 10 KB
/
1e_mints-helper.ipynb
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# MintsHelper: Generating 1- and 2-electron Integrals with <span style='font-variant: small-caps'> Psi4 </span>\n",
"\n",
"In all of quantum chemistry, one process which is common to nearly every method is the evaluation of one-\n",
"and two-electron integrals. Fortunately, we can leverage infrastructure in <span style='font-variant: small-caps'> \n",
"Psi4 </span> to perform this task for us. This tutorial will discuss the [``psi4.core.MintsHelper``](http://psicode.org/psi4manual/master/psi4api.html#psi4.core.MintsHelper \"Go to API\") class, which is an\n",
"interface for the powerful Psi4 ``libmints`` library which wraps the `libint` library, where these integrals are actually computed. \n",
"\n",
"## MintsHelper Overview\n",
"In order to compute 1- and 2-electron integrals, we first need a molecule and basis set with which to work. So, \n",
"before diving into `MintsHelper`, we need to build these objects. In the cell below, we have imported\n",
"<span style='font-variant: small-caps'> Psi4 </span> and NumPy, defined a water molecule, and set the basis to\n",
"cc-pVDZ. We've also set the memory available to <span style='font-variant: small-caps'> Psi4</span>, as well as\n",
"defined a variable `numpy_memory` which we will discuss later."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"# ==> Setup <==\n",
"# Import statements\n",
"import psi4\n",
"import numpy as np\n",
"\n",
"# Memory & Output file\n",
"psi4.set_memory(int(2e9))\n",
"numpy_memory = 2\n",
"psi4.core.set_output_file('output.dat', False)\n",
"\n",
"# Molecule definition\n",
"h2o = psi4.geometry(\"\"\"\n",
"O\n",
"H 1 0.96\n",
"H 1 0.96 2 104.5\n",
"\"\"\")\n",
"\n",
"# Basis Set\n",
"psi4.set_options({'basis': 'cc-pvdz'})"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now, we are ready to create an instance of the `MintsHelper` class. To do this, we need to pass a `BasisSet`\n",
"object to the `MintsHelper` initializer. Fortunately, from the previous tutorial on the `Wavefunction` class, we know\n",
"that we can obtain such an object from an existing wavefunction. So, let's build a new wavefunction for our molecule,\n",
"get the basis set object, and build an instance of `MintsHelper`:"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"# ==> Build MintsHelper Instance <==\n",
"# Build new wavefunction\n",
"wfn = psi4.core.Wavefunction.build(h2o, psi4.core.get_global_option('basis'))\n",
"\n",
"# Initialize MintsHelper with wavefunction's basis set\n",
"mints = psi4.core.MintsHelper(wfn.basisset())"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Below are summarized several commonly computed quantities and how to obtain them using a `MintsHelper` class method:\n",
"\n",
"| Quantity | Function | Description |\n",
"|----------|----------|-------------|\n",
"| AO Overlap integrals | [mints.ao_overlap()](http://psicode.org/psi4manual/master/psi4api.html#psi4.core.MintsHelper.ao_overlap \"Go to Documentation\") | Returns AO overlap matrix as a `psi4.core.Matrix` object |\n",
"| AO Kinetic Energy | [mints.ao_kinetic()](http://psicode.org/psi4manual/master/psi4api.html#psi4.core.MintsHelper.ao_kinetic \"Go to Documentation\") | Returns AO kinetic energy matrix as a `psi4.core.Matrix` object |\n",
"| AO Potential Energy | [mints.ao_potential()](http://psicode.org/psi4manual/master/psi4api.html#psi4.core.MintsHelper.ao_potential \"Go to Documentation\") | Returns AO potential energy matrix as a `psi4.core.Matrix` object |\n",
"| AO Electron Repulsion Integrals | [mints.ao_eri()](http://psicode.org/psi4manual/master/psi4api.html#psi4.core.MintsHelper.ao_eri \"Go to Documentation\") | Returns AO electron repulsion integrals as a `psi4.core.Matrix` object \n",
"\n",
"As discussed previously, any of these `psi4.core.Matrix` objects can be accessed as NumPy arrays, which is the preferred \n",
"method in Psi4NumPy. For a Psi4 matrix `A`, we can access a NumPy view using `np.asarray(A)` or `A.np`, or we can make a\n",
"copy of the matrix using `np.array(A)`. This works as one would expect, converting square matrices into rank-2 NumPy \n",
"arrays, for the overlap (S), kinetic energy (T), and potential energy (V) matrices. In Psi4, the electron repulsion integrals \n",
"(ERIs) are handled somewhat differently; `mints.ao_eri()` returns the rank-4 ERI tensor packed into a 2D matrix. If the \n",
"four indices of the ERI are p, q, r, s, then this element of the Psi4 Matrix can be accessed by first computing composite \n",
"indices `pq = p * nbf + q` and `rs = r * nbf + s`, and then accessing element `A.get(pq,rs)`. However, for convenience, \n",
"the NumPy view is a rank-4 tensor, and a particular ERI is more simply accessed like this:\n",
"~~~python\n",
"I = np.asarray(mints.ao_eri())\n",
"val = I[p][q][r][s]\n",
"~~~\n",
"\n",
"In addition to these methods, another which is worth mentioning is the [`MintsHelper.mo_eri()`](http://psicode.org\n",
"/psi4manual/master/psi4api.html#psi4.core.MintsHelper.mo_eri \"Go to Documentation\") function, which can transform \n",
"the four-index, two-electron repulsion integrals from the atomic orbital (AO) to the molecular orbital (MO) basis,\n",
"which will be important in MP2 theory. \n",
"\n",
"## Memory Considerations\n",
"\n",
"Before moving forward to computing any 1- or 2-electron integrals, we must first discuss the memory requirements of\n",
"these objects. Whenever these quantities are computed, they are stored directly in memory (a.k.a. RAM,\n",
"*not* on the hard drive) which, for a typical laptop or personal computer, usually tops out at around 16 GB of \n",
"space. The storage space required by the two-index AO overlap integrals and four-index ERIs scales as ${\\cal O}(N^2)$ \n",
"and ${\\cal O}(N^4)$, respectively, where $N$ is the number of AO basis functions. This means that for a\n",
"system with 500 AO basis functions, while the AO overlap integrals will only require 1 MB of memory to store,\n",
"the ERIs will require a staggering **500 GB** of memory!! This can be reduced to **62.5 GB** of memory if integral permutational symmetry is used. \n",
"However, this complicates the bookkeeping, and is not used in the `mints` functions discussed above. For this reason, as well as the steep computational \n",
"scaling of many of the methods demonstrated here, we limit ourselves to small systems ($\\sim50$ basis functions)\n",
"which should not require such egregious amounts of memory. Additionally, we will employ a \"memory check\" to catch\n",
"any case which could potentially try to use more memory than is available:\n",
"~~~python\n",
"# Memory check for ERI tensor\n",
"I_size = (nbf**4) * 8.e-9\n",
"print('Size of the ERI tensor will be %4.2f GB.' % (I_size))\n",
"memory_footprint = I_size * 1.5\n",
"if I_size > numpy_memory:\n",
" psi4.core.clean()\n",
" raise Exception(\"Estimated memory utilization (%4.2f GB) exceeds allotted memory \\\n",
" limit of %4.2f GB.\" % (memory_footprint, numpy_memory))\n",
"~~~\n",
"In this example, we have somewhat arbitrarily assumed that whatever other matrices we may need, in total their memory\n",
"requirement will not exceed 50% of the size of the ERIs (hence, the total memory footprint of `I_size * 1.5`)\n",
"Using the `numpy_memory` variable, we are able to control whether the ERIs will be computed, based on the amount of\n",
"memory required to store them. \n",
"\n",
"<font color=\"red\">**NOTE: DO NOT EXCEED YOUR SYSTEM'S MEMORY. THIS MAY RESULT IN YOUR PROGRAM AND/OR COMPUTER CRASHING!**</font>\n",
"\n",
"## Examples: AO Overlap, AO ERIs, Core Hamiltonian\n",
"The cell below demonstrates obtaining the AO overlap integrals, conducting the\n",
"above memory check, and computing the ERIs and core Hamiltonian matrix for our water molecule."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# ==> Integrals galore! <==\n",
"# AO Overlap\n",
"S = np.asarray(mints.ao_overlap())\n",
"\n",
"# Number of basis functions\n",
"nbf = S.shape[0]\n",
"\n",
"# Memory check\n",
"I_size = (nbf ** 4) * 8.e-9\n",
"print('Size of the ERI tensor will be %4.2f GB.' % (I_size))\n",
"memory_footprint = I_size * 1.5\n",
"if I_size > numpy_memory:\n",
" psi4.core.clean()\n",
" raise Exception(\"Estimated memory utilization (%4.2f GB) exceeds allotted memory \\\n",
" limit of %4.2f GB.\" % (memory_footprint, numpy_memory))\n",
"\n",
"# Compute AO-basis ERIs\n",
"I = mints.ao_eri()\n",
"\n",
"# Compute AO Core Hamiltonian\n",
"T = np.asarray(mints.ao_kinetic())\n",
"V = np.asarray(mints.ao_potential())\n",
"H = T + V"
]
}
],
"metadata": {
"anaconda-cloud": {},
"kernelspec": {
"display_name": "Python 2",
"language": "python",
"name": "python2"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.12"
},
"latex_envs": {
"LaTeX_envs_menu_present": true,
"bibliofile": "biblio.bib",
"cite_by": "apalike",
"current_citInitial": 1,
"eqLabelWithNumbers": true,
"eqNumInitial": 1,
"labels_anchors": false,
"latex_user_defs": false,
"report_style_numbering": false,
"user_envs_cfg": false
}
},
"nbformat": 4,
"nbformat_minor": 0
}