Skip to content

Commit

Permalink
#130: Adding tests and fixing problems
Browse files Browse the repository at this point in the history
  • Loading branch information
@giovanni-br
giovanni-br committed Jul 17, 2024
1 parent bd006fa commit 38f6323
Show file tree
Hide file tree
Showing 4 changed files with 137 additions and 125 deletions.
2 changes: 1 addition & 1 deletion examples/Ising.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -453,7 +453,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.5"
"version": "3.10.11"
}
},
"nbformat": 4,
Expand Down
171 changes: 105 additions & 66 deletions examples/TJUV.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@
},
{
"cell_type": "code",
"execution_count": 1,
"execution_count": 2,
"metadata": {},
"outputs": [
{
Expand Down Expand Up @@ -78,7 +78,95 @@
},
{
"cell_type": "code",
"execution_count": 2,
"execution_count": 3,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Numerical energies: [-1. -1. -1. -1.]\n",
"Analytical energies: [-2.0000000e+00 -1.2246468e-16 3.6739404e-16 2.0000000e+00]\n"
]
}
],
"source": [
"import numpy as np\n",
"\n",
"def test_tjuv_energy_spectrum():\n",
" # Define parameters for a simple 1D chain\n",
" N = 4\n",
" t = 1.0\n",
" alpha = -t\n",
" beta = 0.0\n",
" u_onsite = np.zeros(N)\n",
" gamma = None\n",
" charges = None\n",
" sym = None\n",
" J_eq = 0.0\n",
" J_ax = 0.0\n",
"\n",
" # Connectivity matrix for a 1D chain with periodic boundary conditions\n",
" connectivity = np.zeros((N, N))\n",
" for i in range(N):\n",
" connectivity[i, (i + 1) % N] = 1\n",
" connectivity[(i + 1) % N, i] = 1\n",
"\n",
" # Initialize the HamTJUV object\n",
" tjuv_hamiltonian = HamTJUV(connectivity=connectivity,\n",
" alpha=alpha,\n",
" beta=beta,\n",
" u_onsite=u_onsite,\n",
" gamma=gamma,\n",
" charges=charges,\n",
" sym=sym,\n",
" J_eq=J_eq,\n",
" J_ax=J_ax)\n",
"\n",
" # Generate the one-body integral (Hamiltonian matrix)\n",
" tjuv_one_body = tjuv_hamiltonian.generate_one_body_integral(basis='spatial basis', dense=True)\n",
"\n",
" # Calculate the eigenvalues (energy spectrum)\n",
" energies, _ = np.linalg.eigh(tjuv_one_body)\n",
"\n",
" # Analytical energy levels for comparison\n",
" k_vals = np.arange(N)\n",
" analytical_energies = -2 * t * np.cos(2 * np.pi * k_vals / N)\n",
"\n",
" # Sort the energies for comparison\n",
" energies_sorted = np.sort(energies)\n",
" analytical_energies_sorted = np.sort(analytical_energies)\n",
"\n",
" # Print the energy spectrum and analytical values\n",
" print(\"Numerical energies:\", energies_sorted)\n",
" print(\"Analytical energies:\", analytical_energies_sorted)\n",
"\n",
"# Outside of the function, call the test function\n",
"test_tjuv_energy_spectrum()\n"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"ename": "ModuleNotFoundError",
"evalue": "No module named 'pyscf'",
"output_type": "error",
"traceback": [
"\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[1;31mModuleNotFoundError\u001b[0m Traceback (most recent call last)",
"Cell \u001b[1;32mIn[5], line 3\u001b[0m\n\u001b[0;32m 1\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mnumpy\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01mnp\u001b[39;00m\n\u001b[0;32m 2\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mmatplotlib\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mpyplot\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01mplt\u001b[39;00m\n\u001b[1;32m----> 3\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mpyscf\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m fci\n\u001b[0;32m 5\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mtest_tjuv_energy_spectrum\u001b[39m():\n\u001b[0;32m 6\u001b[0m \u001b[38;5;66;03m# Define parameters\u001b[39;00m\n\u001b[0;32m 7\u001b[0m N_values \u001b[38;5;241m=\u001b[39m np\u001b[38;5;241m.\u001b[39mlinspace(\u001b[38;5;241m4\u001b[39m, \u001b[38;5;241m10\u001b[39m, \u001b[38;5;241m4\u001b[39m)\u001b[38;5;241m.\u001b[39mastype(\u001b[38;5;28mint\u001b[39m) \u001b[38;5;66;03m# Varying number of orbitals\u001b[39;00m\n",
"\u001b[1;31mModuleNotFoundError\u001b[0m: No module named 'pyscf'"
]
}
],
"source": []
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
Expand Down Expand Up @@ -129,74 +217,25 @@
"print(\"One body integrals in spin basis: \\n\", h1_spin)\n",
"print(\"Shape of two body integral in spinorbital basis: \", h2_spin.shape)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Energy Spectrum: [-1. -1. -1. -1.]\n",
"Analytical Energy Spectrum: [-2.0000000e+00 -1.2246468e-16 2.0000000e+00 3.6739404e-16]\n"
]
}
],
"source": [
"import numpy as np\n",
"import scipy.linalg\n",
"\n",
"\n",
"# Define parameters for a simple 1D chain\n",
"N = 4\n",
"t = 1.0\n",
"alpha = -t\n",
"beta = 0.0\n",
"u_onsite = np.zeros(N)\n",
"gamma = None\n",
"charges = None\n",
"sym = None\n",
"J_eq = 0.0\n",
"J_ax = 0.0 \n",
"\n",
"# Connectivity matrix for a 1D chain with periodic boundary conditions\n",
"connectivity = np.zeros((N, N))\n",
"for i in range(N):\n",
" connectivity[i, (i+1)%N] = 1\n",
" connectivity[(i+1)%N, i] = 1\n",
"\n",
"# Initialize the HamTJUV object\n",
"tjuv_hamiltonian = HamTJUV(connectivity=connectivity,\n",
" alpha=alpha,\n",
" beta=beta,\n",
" u_onsite=u_onsite,\n",
" gamma=gamma,\n",
" charges=charges,\n",
" sym=sym,\n",
" J_eq=J_eq,\n",
" J_ax=J_ax)\n",
"\n",
"# Generate the one-body integral (Hamiltonian matrix)\n",
"tjuv_one_body = tjuv_hamiltonian.generate_one_body_integral(basis='spatial basis', dense=True)\n",
"\n",
"# Calculate the eigenvalues (energy spectrum)\n",
"energies, _ = np.linalg.eigh(tjuv_one_body)\n",
"\n",
"# Print the energy spectrum\n",
"print(\"Energy Spectrum:\", energies)\n",
"\n",
"# Analytical energy levels for comparison\n",
"k_vals = np.arange(N)\n",
"analytical_energies = -2 * t * np.cos(2 * np.pi * k_vals / N)\n",
"print(\"Analytical Energy Spectrum:\", analytical_energies)\n"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"name": "python"
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.10.11"
}
},
"nbformat": 4,
Expand Down
13 changes: 5 additions & 8 deletions moha/hamiltonians.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -559,9 +559,9 @@ def generate_one_body_integral(self,
return self.one_body.todense() if dense else self.one_body

def generate_two_body_integral(self,
sym=1,
basis='spinorbital basis',
dense=False,
basis='spinorbital basis'):
sym=1):
r"""Generate two body integral in spatial or spinorbital basis.
Parameters
Expand Down Expand Up @@ -732,13 +732,10 @@ def __init__(self,
Parameters
----------
connectivity_ppp: Union[list, np.ndarray]
connectivity: Union[list, np.ndarray]
List of tuples specifying sites and bonds
np.ndarray of shape (n_sites, n_sites)
between sites for the PPP model.
connectivity_heisenberg: np.ndarray
Numpy array that specifies the connectivity beetwen sites
Heisenberg model.
between sites for the TJUV model.
alpha: float
Specifies the site energy if all sites are equivalent.
Default value is the 2p-pi orbital of Carbon.
Expand Down Expand Up @@ -871,6 +868,6 @@ def generate_two_body_integral(self, basis: str, dense: bool, sym=1):
two_body_ppp = self.ocupation_part.generate_two_body_integral(
basis, dense, sym)
two_body_heisenberg = self.spin_part.generate_two_body_integral(
sym, dense, basis)
basis, dense, sym)
self.two_body = two_body_ppp + two_body_heisenberg
return self.two_body
76 changes: 26 additions & 50 deletions moha/test/test_tjuv.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -50,9 +50,7 @@ def test_tjuv_consistency_zero_body():
connectivity=connectivity,
alpha=alpha,
beta=beta,
gamma=np.zeros(
(2,
2)),
gamma=None,
charges=None,
sym=None,
u_onsite=u_onsite)
Expand Down Expand Up @@ -138,7 +136,8 @@ def test_tjuv_consistency_two_body():
u_onsite = np.array([1, 1, 1, 1, 1, 1])
gamma = None
charges = 1
sym = 8
sym = 8 # Use an integer value for symmetry

J_eq = 1
J_ax = 1

Expand All @@ -153,52 +152,29 @@ def test_tjuv_consistency_two_body():
J_eq=J_eq,
J_ax=J_ax)

# Generate the two-body integral (not implemented in this example)


def test_tjuv_energy_spectrum():
# Define parameters for a simple 1D chain
N = 4
t = 1.0
alpha = -t
beta = 0.0
u_onsite = np.zeros(N)
gamma = None
charges = None
sym = None
J_eq = None
J_ax = None

# Connectivity matrix for a 1D chain with periodic boundary conditions
connectivity = np.zeros((N, N))
for i in range(N):
connectivity[i, (i + 1) % N] = 1
connectivity[(i + 1) % N, i] = 1
# Generate the two-body integral
tjuv_two_body = tjuv_hamiltonian.generate_two_body_integral(
basis='spatial basis', dense=True, sym=sym)

# Initialize the HamTJUV object
tjuv_hamiltonian = HamTJUV(connectivity=connectivity,
alpha=alpha,
beta=beta,
u_onsite=u_onsite,
gamma=gamma,
charges=charges,
sym=sym,
J_eq=J_eq,
J_ax=J_ax)

# Generate the one-body integral (Hamiltonian matrix)
tjuv_one_body = tjuv_hamiltonian.generate_one_body_integral()

# Calculate the eigenvalues (energy spectrum)
energies, _ = np.linalg.eigh(tjuv_one_body)

# Analytical energy levels for comparison
k_vals = np.arange(N)
analytical_energies = -2 * t * np.cos(2 * np.pi * k_vals / N)
heisenberg = HamHeisenberg(
J_eq=J_eq,
J_ax=J_ax,
mu=0,
connectivity=connectivity)
heisenberg_two_body = heisenberg.generate_two_body_integral(
basis='spatial basis', dense=True, sym=sym)

# Sort the energies for comparison
energies_sorted = np.sort(energies)
analytical_energies_sorted = np.sort(analytical_energies)
hpp = HamPPP(
connectivity=connectivity,
alpha=alpha,
beta=beta,
gamma=None,
charges=None,
sym=sym,
u_onsite=u_onsite)
hpp_two_body = hpp.generate_two_body_integral(
basis='spatial basis', dense=True, sym=sym)

# Assert that the numerical and analytical energies are close
assert_allclose(energies_sorted, analytical_energies_sorted, atol=1e-10)
# Assert that the TJUV two-body integral is close to the sum of Heisenberg
# and PPP two-body integrals
assert_allclose(tjuv_two_body, heisenberg_two_body + hpp_two_body)

0 comments on commit 38f6323

Please sign in to comment.