|
37 | 37 | "from qiskit_addon_sqd.fermion import rotate_integrals\n",
|
38 | 38 | "\n",
|
39 | 39 | "# Specify molecule properties\n",
|
40 |
| - "open_shell = False\n", |
41 | 40 | "spin_sq = 0\n",
|
42 | 41 | "\n",
|
43 | 42 | "# Build N2 molecule\n",
|
44 | 43 | "mol = pyscf.gto.Mole()\n",
|
45 | 44 | "mol.build(\n",
|
46 | 45 | " atom=[[\"N\", (0, 0, 0)], [\"N\", (1.0, 0, 0)]],\n",
|
47 |
| - " basis=\"6-31g\",\n", |
| 46 | + " # basis=\"6-31g\",\n", |
| 47 | + " basis=\"sto-6g\",\n", |
48 | 48 | " symmetry=\"Dooh\",\n",
|
49 | 49 | ")\n",
|
50 | 50 | "\n",
|
|
104 | 104 | "metadata": {},
|
105 | 105 | "outputs": [],
|
106 | 106 | "source": [
|
107 |
| - "from qiskit_addon_sqd.counts import counts_to_arrays, generate_counts_uniform\n", |
| 107 | + "from qiskit_addon_sqd.counts import generate_bit_array_uniform\n", |
108 | 108 | "\n",
|
109 | 109 | "# Generate random samples\n",
|
110 |
| - "counts_dict = generate_counts_uniform(10_000, num_orbitals * 2, rand_seed=rng)\n", |
111 |
| - "\n", |
112 |
| - "# Convert counts into bitstring and probability arrays\n", |
113 |
| - "bitstring_matrix_full, probs_arr_full = counts_to_arrays(counts_dict)" |
| 110 | + "bit_array = generate_bit_array_uniform(10_000, num_orbitals * 2, rand_seed=rng)" |
114 | 111 | ]
|
115 | 112 | },
|
116 | 113 | {
|
|
131 | 128 | "name": "stdout",
|
132 | 129 | "output_type": "stream",
|
133 | 130 | "text": [
|
134 |
| - "Starting configuration recovery iteration 0\n", |
135 |
| - "Starting configuration recovery iteration 1\n", |
136 |
| - "Starting configuration recovery iteration 2\n", |
137 |
| - "Starting configuration recovery iteration 3\n", |
138 |
| - "Starting configuration recovery iteration 4\n" |
| 131 | + "Iteration 1\n", |
| 132 | + "\tSubsample 0\n", |
| 133 | + "\t\tEnergy: -108.51514815848884\n", |
| 134 | + "\t\tSubspace dimension: 900\n", |
| 135 | + "\tSubsample 1\n", |
| 136 | + "\t\tEnergy: -107.74386023793795\n", |
| 137 | + "\t\tSubspace dimension: 784\n", |
| 138 | + "\tSubsample 2\n", |
| 139 | + "\t\tEnergy: -107.7411531716499\n", |
| 140 | + "\t\tSubspace dimension: 841\n", |
| 141 | + "Iteration 2\n", |
| 142 | + "\tSubsample 0\n", |
| 143 | + "\t\tEnergy: -108.53322549223714\n", |
| 144 | + "\t\tSubspace dimension: 625\n", |
| 145 | + "\tSubsample 1\n", |
| 146 | + "\t\tEnergy: -108.57716824917271\n", |
| 147 | + "\t\tSubspace dimension: 676\n", |
| 148 | + "\tSubsample 2\n", |
| 149 | + "\t\tEnergy: -108.508931350816\n", |
| 150 | + "\t\tSubspace dimension: 576\n", |
| 151 | + "Iteration 3\n", |
| 152 | + "\tSubsample 0\n", |
| 153 | + "\t\tEnergy: -108.55902215621444\n", |
| 154 | + "\t\tSubspace dimension: 625\n", |
| 155 | + "\tSubsample 1\n", |
| 156 | + "\t\tEnergy: -108.5776045917047\n", |
| 157 | + "\t\tSubspace dimension: 529\n", |
| 158 | + "\tSubsample 2\n", |
| 159 | + "\t\tEnergy: -108.56103466207057\n", |
| 160 | + "\t\tSubspace dimension: 676\n", |
| 161 | + "Iteration 4\n", |
| 162 | + "\tSubsample 0\n", |
| 163 | + "\t\tEnergy: -108.56619527937809\n", |
| 164 | + "\t\tSubspace dimension: 529\n", |
| 165 | + "\tSubsample 1\n", |
| 166 | + "\t\tEnergy: -108.56463009757981\n", |
| 167 | + "\t\tSubspace dimension: 625\n", |
| 168 | + "\tSubsample 2\n", |
| 169 | + "\t\tEnergy: -108.56630478170408\n", |
| 170 | + "\t\tSubspace dimension: 784\n", |
| 171 | + "Iteration 5\n", |
| 172 | + "\tSubsample 0\n", |
| 173 | + "\t\tEnergy: -108.56883735605814\n", |
| 174 | + "\t\tSubspace dimension: 529\n", |
| 175 | + "\tSubsample 1\n", |
| 176 | + "\t\tEnergy: -108.58150746468066\n", |
| 177 | + "\t\tSubspace dimension: 625\n", |
| 178 | + "\tSubsample 2\n", |
| 179 | + "\t\tEnergy: -108.55836297172482\n", |
| 180 | + "\t\tSubspace dimension: 529\n" |
139 | 181 | ]
|
140 | 182 | }
|
141 | 183 | ],
|
142 | 184 | "source": [
|
143 |
| - "from qiskit_addon_sqd.configuration_recovery import recover_configurations\n", |
144 |
| - "from qiskit_addon_sqd.fermion import solve_fermion\n", |
145 |
| - "from qiskit_addon_sqd.subsampling import postselect_and_subsample\n", |
| 185 | + "from functools import partial\n", |
| 186 | + "\n", |
| 187 | + "from qiskit_addon_sqd.fermion import SCIResult, diagonalize_fermionic_hamiltonian, solve_sci_batch\n", |
146 | 188 | "\n",
|
147 |
| - "# SQSD options\n", |
148 |
| - "iterations = 5\n", |
| 189 | + "# SQD options\n", |
| 190 | + "energy_tol = 1e-3\n", |
| 191 | + "occupancies_tol = 1e-3\n", |
| 192 | + "max_iterations = 5\n", |
149 | 193 | "\n",
|
150 | 194 | "# Eigenstate solver options\n",
|
151 |
| - "n_batches = 3\n", |
152 |
| - "samples_per_batch = 100\n", |
153 |
| - "max_davidson_cycles = 200\n", |
| 195 | + "num_batches = 3\n", |
| 196 | + "samples_per_batch = 20\n", |
| 197 | + "symmetrize_spin = True\n", |
| 198 | + "carryover_threshold = 1e-1\n", |
| 199 | + "max_cycle = 200\n", |
154 | 200 | "\n",
|
155 |
| - "# Self-consistent configuration recovery loop\n", |
156 |
| - "e_hist = np.zeros((iterations, n_batches)) # energy history\n", |
157 |
| - "s_hist = np.zeros((iterations, n_batches)) # spin history\n", |
158 |
| - "occupancy_hist = []\n", |
159 |
| - "avg_occupancy = None\n", |
160 |
| - "for i in range(iterations):\n", |
161 |
| - " print(f\"Starting configuration recovery iteration {i}\")\n", |
162 |
| - " # On the first iteration, we have no orbital occupancy information from the\n", |
163 |
| - " # solver, so we just post-select from the full bitstring set based on hamming weight.\n", |
164 |
| - " if avg_occupancy is None:\n", |
165 |
| - " bs_mat_tmp = bitstring_matrix_full\n", |
166 |
| - " probs_arr_tmp = probs_arr_full\n", |
| 201 | + "# Pass options to the built-in eigensolver. If you just want to use the defaults,\n", |
| 202 | + "# you can omit this step, in which case you would not specify the sci_solver argument\n", |
| 203 | + "# in the call to diagonalize_fermionic_hamiltonian below.\n", |
| 204 | + "sci_solver = partial(solve_sci_batch, spin_sq=0.0, max_cycle=max_cycle)\n", |
167 | 205 | "\n",
|
168 |
| - " # In following iterations, we use both the occupancy info and the target hamming\n", |
169 |
| - " # weight to refine bitstrings.\n", |
170 |
| - " else:\n", |
171 |
| - " bs_mat_tmp, probs_arr_tmp = recover_configurations(\n", |
172 |
| - " bitstring_matrix_full,\n", |
173 |
| - " probs_arr_full,\n", |
174 |
| - " avg_occupancy,\n", |
175 |
| - " num_elec_a,\n", |
176 |
| - " num_elec_b,\n", |
177 |
| - " rand_seed=rng,\n", |
178 |
| - " )\n", |
| 206 | + "# List to capture intermediate results\n", |
| 207 | + "result_history = []\n", |
179 | 208 | "\n",
|
180 |
| - " # Throw out samples with incorrect hamming weight and create batches of subsamples.\n", |
181 |
| - " batches = postselect_and_subsample(\n", |
182 |
| - " bs_mat_tmp,\n", |
183 |
| - " probs_arr_tmp,\n", |
184 |
| - " hamming_right=num_elec_a,\n", |
185 |
| - " hamming_left=num_elec_b,\n", |
186 |
| - " samples_per_batch=samples_per_batch,\n", |
187 |
| - " num_batches=n_batches,\n", |
188 |
| - " rand_seed=rng,\n", |
189 |
| - " )\n", |
190 | 209 | "\n",
|
191 |
| - " # Run eigenstate solvers in a loop. This loop should be parallelized for larger problems.\n", |
192 |
| - " int_e = np.zeros(n_batches)\n", |
193 |
| - " int_s = np.zeros(n_batches)\n", |
194 |
| - " int_occs = []\n", |
195 |
| - " cs = []\n", |
196 |
| - " for j in range(n_batches):\n", |
197 |
| - " energy_sci, coeffs_sci, avg_occs, spin = solve_fermion(\n", |
198 |
| - " batches[j],\n", |
199 |
| - " hcore_rot,\n", |
200 |
| - " eri_rot,\n", |
201 |
| - " open_shell=open_shell,\n", |
202 |
| - " spin_sq=spin_sq,\n", |
203 |
| - " max_cycle=max_davidson_cycles,\n", |
204 |
| - " )\n", |
205 |
| - " energy_sci += nuclear_repulsion_energy\n", |
206 |
| - " int_e[j] = energy_sci\n", |
207 |
| - " int_s[j] = spin\n", |
208 |
| - " int_occs.append(avg_occs)\n", |
209 |
| - " cs.append(coeffs_sci)\n", |
| 210 | + "def callback(results: list[SCIResult]):\n", |
| 211 | + " result_history.append(results)\n", |
| 212 | + " iteration = len(result_history)\n", |
| 213 | + " print(f\"Iteration {iteration}\")\n", |
| 214 | + " for i, result in enumerate(results):\n", |
| 215 | + " print(f\"\\tSubsample {i}\")\n", |
| 216 | + " print(f\"\\t\\tEnergy: {result.energy + nuclear_repulsion_energy}\")\n", |
| 217 | + " print(f\"\\t\\tSubspace dimension: {np.prod(result.sci_state.amplitudes.shape)}\")\n", |
210 | 218 | "\n",
|
211 |
| - " # Combine batch results\n", |
212 |
| - " avg_occupancy = tuple(np.mean(int_occs, axis=0))\n", |
213 | 219 | "\n",
|
214 |
| - " # Track optimization history\n", |
215 |
| - " e_hist[i, :] = int_e\n", |
216 |
| - " s_hist[i, :] = int_s\n", |
217 |
| - " occupancy_hist.append(avg_occupancy)" |
| 220 | + "result = diagonalize_fermionic_hamiltonian(\n", |
| 221 | + " hcore,\n", |
| 222 | + " eri,\n", |
| 223 | + " bit_array,\n", |
| 224 | + " samples_per_batch=samples_per_batch,\n", |
| 225 | + " norb=num_orbitals,\n", |
| 226 | + " nelec=(num_elec_a, num_elec_b),\n", |
| 227 | + " num_batches=num_batches,\n", |
| 228 | + " energy_tol=energy_tol,\n", |
| 229 | + " occupancies_tol=occupancies_tol,\n", |
| 230 | + " max_iterations=max_iterations,\n", |
| 231 | + " sci_solver=sci_solver,\n", |
| 232 | + " symmetrize_spin=symmetrize_spin,\n", |
| 233 | + " carryover_threshold=carryover_threshold,\n", |
| 234 | + " callback=callback,\n", |
| 235 | + " seed=rng,\n", |
| 236 | + ")" |
218 | 237 | ]
|
219 | 238 | },
|
220 | 239 | {
|
|
224 | 243 | "source": [
|
225 | 244 | "### Refine the subspace\n",
|
226 | 245 | "\n",
|
227 |
| - "To refine the subspace, we will take the CI strings of the batch with the lowest energy\n", |
228 |
| - "from the last configuration recovery step. Other strategies may be used, like taking the union \n", |
229 |
| - "of the CI strings of the batches in the last configuration recovery iteration." |
| 246 | + "To refine the subspace, we will take the CI strings of the SCI state returned by SQD. Other strategies may be used, like taking the union of the CI strings of the batches in the last configuration recovery iteration. You can save information from the recovery iterations by passing a callback function to `run_sqd` (see [Improving energy estimation of a chemistry Hamiltonian with SQD](../tutorials/01_chemistry_hamiltonian.ipynb))." |
230 | 247 | ]
|
231 | 248 | },
|
232 | 249 | {
|
|
239 | 256 | "name": "stdout",
|
240 | 257 | "output_type": "stream",
|
241 | 258 | "text": [
|
242 |
| - "Subspace dimension: 32400\n", |
243 |
| - "Energy of that batch from SQD: -107.55130314653121\n" |
| 259 | + "Subspace dimension: 625\n", |
| 260 | + "Energy of that batch from SQD: -108.58150746468066\n" |
244 | 261 | ]
|
245 | 262 | }
|
246 | 263 | ],
|
247 | 264 | "source": [
|
248 |
| - "from qiskit_addon_sqd.fermion import bitstring_matrix_to_ci_strs\n", |
249 |
| - "\n", |
250 |
| - "best_batch = batches[np.argmin(e_hist[-1])]\n", |
251 |
| - "ci_strs_up, ci_strs_dn = bitstring_matrix_to_ci_strs(best_batch, open_shell=open_shell)\n", |
252 |
| - "print(f\"Subspace dimension: {len(ci_strs_up) * len(ci_strs_dn)}\")\n", |
253 |
| - "print(f\"Energy of that batch from SQD: {e_hist[-1, np.argmin(e_hist[-1])]}\")\n", |
254 |
| - "\n", |
255 |
| - "# Union strategy\n", |
256 |
| - "\n", |
257 |
| - "# batches_union = np.concatenate((batches[0], batches[1]), axis = 0)\n", |
258 |
| - "# for i in range(n_batches-2):\n", |
259 |
| - "# batches_union = np.concatenate((batches_union, batches[ i+ 2]))\n", |
260 |
| - "# ci_strs_up, ci_strs_dn = bitstring_matrix_to_ci_strs(\n", |
261 |
| - "# batches_union, open_shell=open_shell\n", |
262 |
| - "# )\n", |
263 |
| - "# print (f\"Subspace dimension: {len(ci_strs_up) * len(ci_strs_dn)}\")" |
| 265 | + "print(f\"Subspace dimension: {np.prod(result.sci_state.amplitudes.shape)}\")\n", |
| 266 | + "print(f\"Energy of that batch from SQD: {result.energy + nuclear_repulsion_energy}\")" |
264 | 267 | ]
|
265 | 268 | },
|
266 | 269 | {
|
|
316 | 319 | "learning_rate = 0.05\n",
|
317 | 320 | "\n",
|
318 | 321 | "e_improved, k_flat, orbital_occupancies = optimize_orbitals(\n",
|
319 |
| - " best_batch,\n", |
| 322 | + " (result.sci_state.ci_strs_a, result.sci_state.ci_strs_b),\n", |
320 | 323 | " hcore_rot,\n",
|
321 | 324 | " eri_rot,\n",
|
322 | 325 | " k_flat,\n",
|
323 |
| - " open_shell=open_shell,\n", |
| 326 | + " open_shell=False,\n", |
324 | 327 | " spin_sq=spin_sq,\n",
|
325 | 328 | " num_iters=num_iters,\n",
|
326 | 329 | " num_steps_grad=num_steps_grad,\n",
|
327 | 330 | " learning_rate=learning_rate,\n",
|
328 |
| - " max_cycle=max_davidson_cycles,\n", |
| 331 | + " max_cycle=max_cycle,\n", |
329 | 332 | ")"
|
330 | 333 | ]
|
331 | 334 | },
|
|
347 | 350 | "name": "stdout",
|
348 | 351 | "output_type": "stream",
|
349 | 352 | "text": [
|
350 |
| - "Exact energy: -109.04667177808027\n", |
351 |
| - "SQD energy: -107.55130314653121\n", |
352 |
| - "Energy after OO: -107.39122543051276\n" |
| 353 | + "Exact energy: -108.59598735098602\n", |
| 354 | + "SQD energy: -108.58150746468066\n", |
| 355 | + "Energy after OO: -108.49008547897714\n" |
353 | 356 | ]
|
354 | 357 | }
|
355 | 358 | ],
|
356 | 359 | "source": [
|
357 | 360 | "print(f\"Exact energy: {exact_energy}\")\n",
|
358 |
| - "print(f\"SQD energy: {np.min(e_hist[-1])}\")\n", |
| 361 | + "print(f\"SQD energy: {result.energy + nuclear_repulsion_energy}\")\n", |
359 | 362 | "print(f\"Energy after OO: {e_improved + nuclear_repulsion_energy}\")"
|
360 | 363 | ]
|
361 | 364 | }
|
|
0 commit comments