Final: Fix hopping and QR tasks
This commit is contained in:
@ -606,7 +606,7 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 7,
|
||||
"execution_count": 17,
|
||||
"metadata": {
|
||||
"deletable": false,
|
||||
"nbgrader": {
|
||||
@ -625,7 +625,16 @@
|
||||
"source": [
|
||||
"def hopping(i, j, n):\n",
|
||||
" \"\"\"\n",
|
||||
" Calculates hopping matrix elements t_ij for sigma = 0.25 in a 1D\n",
|
||||
" chain of n atoms at distance a = 1 from eachother.\n",
|
||||
" \n",
|
||||
" Args:\n",
|
||||
" i: origin site index\n",
|
||||
" j: destination site index\n",
|
||||
" n: number of atoms in the chain\n",
|
||||
"\n",
|
||||
" Returns:\n",
|
||||
" Hopping parameter t_ij.\n",
|
||||
" \"\"\"\n",
|
||||
" \n",
|
||||
" sigma = .25\n",
|
||||
@ -633,10 +642,11 @@
|
||||
" \n",
|
||||
" # This 'infinity' is large enough, as the Gaussians decay quite quickly\n",
|
||||
" # away from the atomic positions, which we already saw in the overlap\n",
|
||||
" # above.\n",
|
||||
" # TODO: Remember https://www.wolframalpha.com/input?i=integrate+1%2F%28.25*sqrt%282pi%29%29*exp%28-.5*%28x%2F.25%29%5E2%29+for+x%3D-1..1.\n",
|
||||
" infty = 1e2\n",
|
||||
" h = 1e-6\n",
|
||||
" # above. In fact, 99.7% of all probability mass is under the integral\n",
|
||||
" # for x radius of 3*sigma from the centers x_i.\n",
|
||||
" # TODO: Make integration bounds relative to n and a, and maybe sigma.\n",
|
||||
" infty = 5e1\n",
|
||||
" h = 1e-5\n",
|
||||
" x = np.arange(-infty, infty, h)\n",
|
||||
" \n",
|
||||
" def V(x):\n",
|
||||
@ -653,21 +663,12 @@
|
||||
" #V = lambda x: np.sum( -1/( np.abs(np.subtract.outer(x, positions)) + 0.001 ), axis=1 )\n",
|
||||
" \n",
|
||||
" integrand = lambda x: atomic_basis(x, positions[i], sigma)*V(x)*atomic_basis(x, positions[j], sigma)\n",
|
||||
" return integrate(integrand, x)\n",
|
||||
"\n",
|
||||
"# Yann had\n",
|
||||
"# 1) a plot of the peaks as before\n",
|
||||
"# 2) a weird line under it with peaks downward at these x\n",
|
||||
"# and slowly rising energy values outside the chain\n",
|
||||
"# 3)\n",
|
||||
"# For i = 1 ...\n",
|
||||
"# t_{i,i-1} = -0.1203\n",
|
||||
"# t_{i,i+1} = ..."
|
||||
" return integrate(integrand, x)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 8,
|
||||
"execution_count": 27,
|
||||
"metadata": {
|
||||
"deletable": false,
|
||||
"nbgrader": {
|
||||
@ -687,15 +688,88 @@
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"-3.08808170797426e-06 -0.13849169656144708\n"
|
||||
"For i = 0 ...\n",
|
||||
"\tt_{i,i+1} = -0.13849171527765194\n",
|
||||
"\tt_{i,i+2} = -3.0880816778755327e-06\n",
|
||||
"\tt_{i,i+3} = -1.8833562169250255e-15\n",
|
||||
"\n",
|
||||
"For i = 1 ...\n",
|
||||
"\tt_{i,i-1} = -0.13849171527765194\n",
|
||||
"\tt_{i,i+1} = -0.14871538196685558\n",
|
||||
"\tt_{i,i+2} = -3.1307050912510583e-06\n",
|
||||
"\tt_{i,i+3} = -1.9456306931892115e-15\n",
|
||||
"\n",
|
||||
"For i = 2 ...\n",
|
||||
"\tt_{i,i-1} = -0.14871538196685555\n",
|
||||
"\tt_{i,i+1} = -0.1536327589565012\n",
|
||||
"\tt_{i,i-2} = -3.0880816778755327e-06\n",
|
||||
"\tt_{i,i+2} = -3.152251398227753e-06\n",
|
||||
"\tt_{i,i+3} = -1.976347912754562e-15\n",
|
||||
"\n",
|
||||
"For i = 3 ...\n",
|
||||
"\tt_{i,i-1} = -0.1536327589565012\n",
|
||||
"\tt_{i,i+1} = -0.15605828154196247\n",
|
||||
"\tt_{i,i-2} = -3.1307050912510583e-06\n",
|
||||
"\tt_{i,i+2} = -3.1616580016478533e-06\n",
|
||||
"\tt_{i,i-3} = -1.8833562169250255e-15\n",
|
||||
"\tt_{i,i+3} = -1.98575211950378e-15\n",
|
||||
"\n",
|
||||
"For i = 4 ...\n",
|
||||
"\tt_{i,i-1} = -0.15605828154196244\n",
|
||||
"\tt_{i,i+1} = -0.15680086554505168\n",
|
||||
"\tt_{i,i-2} = -3.152251398227753e-06\n",
|
||||
"\tt_{i,i+2} = -3.161664327360904e-06\n",
|
||||
"\tt_{i,i-3} = -1.9456306931892115e-15\n",
|
||||
"\tt_{i,i+3} = -1.976348151505985e-15\n",
|
||||
"\n",
|
||||
"For i = 5 ...\n",
|
||||
"\tt_{i,i-1} = -0.15680086554505168\n",
|
||||
"\tt_{i,i+1} = -0.156058300394482\n",
|
||||
"\tt_{i,i-2} = -3.1616580016478533e-06\n",
|
||||
"\tt_{i,i+2} = -3.152251395817713e-06\n",
|
||||
"\tt_{i,i-3} = -1.9763479127545625e-15\n",
|
||||
"\tt_{i,i+3} = -1.945630454073864e-15\n",
|
||||
"\n",
|
||||
"For i = 6 ...\n",
|
||||
"\tt_{i,i-1} = -0.15605830039448199\n",
|
||||
"\tt_{i,i+1} = -0.153632740075246\n",
|
||||
"\tt_{i,i-2} = -3.161664327360904e-06\n",
|
||||
"\tt_{i,i+2} = -3.1306987607180875e-06\n",
|
||||
"\tt_{i,i-3} = -1.9857521195037804e-15\n",
|
||||
"\tt_{i,i+3} = -1.8833562165610968e-15\n",
|
||||
"\n",
|
||||
"For i = 7 ...\n",
|
||||
"\tt_{i,i-1} = -0.153632740075246\n",
|
||||
"\tt_{i,i+1} = -0.14871538193811829\n",
|
||||
"\tt_{i,i-2} = -3.152251395817713e-06\n",
|
||||
"\tt_{i,i+2} = -3.088087998768686e-06\n",
|
||||
"\tt_{i,i-3} = -1.976348151505985e-15\n",
|
||||
"\n",
|
||||
"For i = 8 ...\n",
|
||||
"\tt_{i,i-1} = -0.14871538193811829\n",
|
||||
"\tt_{i,i+1} = -0.1384917341014353\n",
|
||||
"\tt_{i,i-2} = -3.1306987607180875e-06\n",
|
||||
"\tt_{i,i-3} = -1.945630454073864e-15\n",
|
||||
"\n",
|
||||
"For i = 9 ...\n",
|
||||
"\tt_{i,i-1} = -0.1384917341014353\n",
|
||||
"\tt_{i,i-2} = -3.088087998768686e-06\n",
|
||||
"\tt_{i,i-3} = -1.8833562165610968e-15\n",
|
||||
"\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"print(hopping(0, 2, 10), hopping(0, 1, 10))\n",
|
||||
"n = 10\n",
|
||||
"\n",
|
||||
"# YOUR CODE HERE\n",
|
||||
"#raise NotImplementedError()"
|
||||
"for i in range(n):\n",
|
||||
" print(\"For i =\", i, \"...\")\n",
|
||||
" for r in range(1, 4):\n",
|
||||
" if i - r >= 0:\n",
|
||||
" print(\"\\tt_{{i,i-{}}} = {}\".format(r, hopping(i, i - r, n)))\n",
|
||||
" if i + r < n:\n",
|
||||
" print(\"\\tt_{{i,i+{}}} = {}\".format(r, hopping(i, i + r, n)))\n",
|
||||
" print()"
|
||||
]
|
||||
},
|
||||
{
|
||||
@ -723,7 +797,7 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"execution_count": 54,
|
||||
"metadata": {
|
||||
"deletable": false,
|
||||
"nbgrader": {
|
||||
@ -740,10 +814,30 @@
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# YOUR CODE HERE\n",
|
||||
"raise NotImplementedError()\n",
|
||||
"\n",
|
||||
"# Yann said we just have to be able to solve tri-diagonal matrix diagonalization, but more general is ok."
|
||||
"def QREig(T, eps=1e-6, k_max=10000):\n",
|
||||
" \"\"\"\n",
|
||||
" Follows the method of the QR decomposition based diagonalization routine\n",
|
||||
" for tridiagonal matrices. The matrix T is diagonalized, resulting in\n",
|
||||
" all diagonal elements being an eigenvalue.\n",
|
||||
" \n",
|
||||
" Args:\n",
|
||||
" T: a tridiagonaliz matrix.\n",
|
||||
" eps: the desired accuracy.\n",
|
||||
" k_max: maximum number of iterations after which to cut off\n",
|
||||
" \n",
|
||||
" Returns:\n",
|
||||
" A one dimensional array with the eigenvalues of the matrix T.\n",
|
||||
" \"\"\"\n",
|
||||
" \n",
|
||||
" e = eps + 1\n",
|
||||
" k = 0\n",
|
||||
" while e > eps and k < k_max:\n",
|
||||
" k += 1\n",
|
||||
" Q, R = np.linalg.qr(T)\n",
|
||||
" T = np.matmul(R,Q)\n",
|
||||
" e = np.sum(np.abs(np.diag(T, k=1)))\n",
|
||||
" print(k)\n",
|
||||
" return np.diag(T)"
|
||||
]
|
||||
},
|
||||
{
|
||||
@ -769,7 +863,7 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"execution_count": 75,
|
||||
"metadata": {
|
||||
"deletable": false,
|
||||
"nbgrader": {
|
||||
@ -784,10 +878,46 @@
|
||||
"task": false
|
||||
}
|
||||
},
|
||||
"outputs": [],
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"10000\n",
|
||||
"1000\n",
|
||||
"[-3.35311057 0.14162732 3. 4.21148325]\n",
|
||||
"10000\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"# YOUR CODE HERE\n",
|
||||
"raise NotImplementedError()"
|
||||
"def test_QREig():\n",
|
||||
" # Test case one\n",
|
||||
" T = np.array([\n",
|
||||
" [1,4,0,0],\n",
|
||||
" [3,4,1,0],\n",
|
||||
" [0,2,3,4],\n",
|
||||
" [0,0,1,3]\n",
|
||||
" ])\n",
|
||||
" # Eigenvalues are roots of λ^4 - 11*λ^3 + 25*λ^2 + 31*λ - 46.\n",
|
||||
" eigenvalues_of_T = np.array([-1.45350244, 1., 4.65531023, 6.79819221])\n",
|
||||
" assert np.allclose(np.sort(QREig(T)), eigenvalues_of_T)\n",
|
||||
" #print(QREig(T))\n",
|
||||
" \n",
|
||||
" # Test case two\n",
|
||||
" T = np.array([\n",
|
||||
" [1,4,0,0],\n",
|
||||
" [3,0,1,0],\n",
|
||||
" [0,2,0,4],\n",
|
||||
" [0,0,0,3]\n",
|
||||
" ])\n",
|
||||
" #eigenvalues_of_T = np.array([-3.35311, .141627, 3. , 4.21148])\n",
|
||||
" eigenvalues_of_T = np.sort(eig(T)[0])\n",
|
||||
" #print(QREig(T+np.eye(len(T)), k_max=1000000)-1)\n",
|
||||
" print(np.sort(QREig(T, k_max=1000)))\n",
|
||||
" assert np.allclose(np.sort(QREig(T)), eigenvalues_of_T)\n",
|
||||
"\n",
|
||||
"test_QREig()"
|
||||
]
|
||||
},
|
||||
{
|
||||
@ -815,7 +945,7 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"execution_count": 77,
|
||||
"metadata": {
|
||||
"deletable": false,
|
||||
"nbgrader": {
|
||||
@ -833,8 +963,14 @@
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"def TBHamiltonian(n):\n",
|
||||
" # YOUR CODE HERE\n",
|
||||
" raise NotImplementedError()"
|
||||
" \"\"\"\n",
|
||||
" \n",
|
||||
" \"\"\"\n",
|
||||
" \n",
|
||||
" # TODO: Comment on the weird 20% differences in hopping parameters.\n",
|
||||
" \n",
|
||||
" t = hopping(i, j, n)\n",
|
||||
" H = np.eye(n) + np.eye(n, -1) + np.eye(n, 1)"
|
||||
]
|
||||
},
|
||||
{
|
||||
@ -856,6 +992,9 @@
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"for n in [10, 20, 40, 80, 100]:\n",
|
||||
" pass\n",
|
||||
"\n",
|
||||
"# Do your plotting here ...\n",
|
||||
"\n",
|
||||
"# YOUR CODE HERE\n",
|
||||
|
||||
Reference in New Issue
Block a user