<jupyter_start><jupyter_text>Running `H.stark_map()` is essentially equivilent to:<jupyter_code>field_au = 100.0 * field * e * a0 / En_h  # atomic units

stark_map = np.array(
    [linalg.eigvalsh(H.matrix(Fz=f).toarray()) for f in tqdm(field_au)]
)

# plot
for level in stark_map.T:
    plt.plot(field, level, c="k", alpha=0.1)

plt.xlabel("electric field (V / cm)")
plt.ylabel("energy (a. u.)")
plt.show()
%timeit linalg.eigvalsh(H.matrix(Fz=2.5e-10).toarray())<jupyter_output>304 ms ± 56.8 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)<jupyter_text>Sometimes a large basis is required to get accurate energy levels, even if only a fraction of the eigenvalues are of interest. Fortunately, it is possible to compute a subset of the eignevalues, which is often a little faster.  Moreover, the total Hamiltonian matrix can remain sparse, potentially saving a lot of memory.<jupyter_code>subset = list(H.basis.argwhere(lambda x: x.n == 36))
k = len(subset)
sigma = np.median(H.e0[subset])

sparse_map = np.array(
    [
        sp.linalg.eigsh(H.matrix(Fz=f), k=k, sigma=sigma, return_eigenvectors=False)
        for f in tqdm(field_au)
    ]
)

# plot
for level in sparse_map.T:
    plt.plot(field, level, c="k", alpha=0.5)

plt.xlabel("electric field (V / cm)")
plt.ylabel("energy (a. u.)")
plt.show()<jupyter_output>100%|██████████| 20/20 [00:03<00:00,  6.61it/s]<jupyter_text>(The above calculation can probably be sped up by moving some of the shift-inversion steps outside of the loop.) Compare sparse and dense calculations:<jupyter_code># dense
for i in subset:
    plt.plot(field, stark_map[:, i], c="k")

# sparse
for level in sparse_map.T:
    plt.plot(-field, level, c="r")

plt.ylim(np.min(sparse_map), np.max(sparse_map))
plt.xlabel("electric field (V / cm)")
plt.ylabel("energy (a. u.)")
plt.show()<jupyter_output><empty_output><jupyter_text>If the subset contains a partial $n$ manifold, exactly which eigenvalues are computed can vary abruptly from one field to the next.  This can get confusing when attempting to trace states through crossings.<jupyter_code>k = 12
sparse_map_12 = np.array(
    [
        sp.linalg.eigsh(H.matrix(Fz=f), k=k, sigma=sigma, return_eigenvectors=False)
        for f in tqdm(field_au)
    ]
)


# dense
for i in subset:
    plt.plot(field, stark_map[:, i], c="k", alpha=0.2)

# sparse
for level in sparse_map_12.T:
    plt.plot(field, level, c="r", ls="", marker="o", alpha=0.5)

plt.ylim(np.min(sparse_map_12), np.max(sparse_map_12))
plt.xlabel("electric field (V / cm)")
plt.ylabel("energy (a. u.)")
plt.show()
%timeit sp.linalg.eigsh(H.matrix(Fz=2.5e-10), k=k, sigma=sigma, return_eigenvectors=False)<jupyter_output>21.4 ms ± 1.33 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
{'fraction_non_alphanumeric': 0.14764855997083484, 'fraction_numerical': 0.025154939846882975, 'mean_word_length': 4.444444444444445, 'pattern_counts': {'":': 0, '<': 15, '<?xml version=': 0, '>': 14, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '55949', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 1123, 'n_tokens_neox': 1042, 'n_words': 290, 'path': 'notebooks/Sparse.ipynb', 'repo_name': 'ad3ller/pypositronium', 'size': '2743'}