prb10_p.ipynb (8477B)
1 { 2 "cells": [ 3 { 4 "cell_type": "code", 5 "execution_count": 1, 6 "id": "e8c4b067", 7 "metadata": {}, 8 "outputs": [], 9 "source": [ 10 "import numpy as np" 11 ] 12 }, 13 { 14 "cell_type": "code", 15 "execution_count": 59, 16 "id": "35139a3d", 17 "metadata": {}, 18 "outputs": [], 19 "source": [ 20 "def lanczos_alg(A: np.ndarray, b: np.array, N: int=None) -> np.ndarray:\n", 21 " if N == None:\n", 22 " N = A.shape[0]\n", 23 " v_1 = b/np.linalg.norm(b)\n", 24 " V = np.zeros([A.shape[0], N])\n", 25 " \n", 26 " V[:,0] = v_1\n", 27 " w̃_j = A @ V[:,0]\n", 28 " α_j = np.conjugate(w̃_j) @ V[:,0]\n", 29 " w_j = w̃_j - α_j*V[:,0]\n", 30 " \n", 31 " ritz = [np.conjugate(V[:,0]) @ A @ V[:,0]]\n", 32 " for j in range(1, N):\n", 33 " if (β_j:=np.linalg.norm(w_j)) != 0:\n", 34 " V[:,j] = w_j/β_j\n", 35 " else:\n", 36 " break\n", 37 " w̃_j = A @ V[:,j]\n", 38 " α_j = np.conjugate(w̃_j) @ V[:,j]\n", 39 " w_j = w̃_j - α_j*V[:,j] - β_j * V[:,j-1]\n", 40 " \n", 41 " ritz.append(np.conjugate(V[:,j]) @ A @ V[:,j])\n", 42 " return V, ritz" 43 ] 44 }, 45 { 46 "cell_type": "code", 47 "execution_count": 73, 48 "id": "aa9ac15b", 49 "metadata": {}, 50 "outputs": [ 51 { 52 "name": "stdout", 53 "output_type": "stream", 54 "text": [ 55 "[[ 0.57900658 -0.32435077 -0.56269734 0.49286881 -0.42996386 -0.82840343\n", 56 " 0.29946173 0.19799369 0.44790496 0.62677608 0.30550101 0.55964456\n", 57 " -0.44827518 -0.65545559 0.49439402 0.35355848 0.44829467 0.65322912\n", 58 " 0.4691473 0.38546294 -0.45234048 -0.65328039 -0.29575937 0.3119645\n", 59 " 0.62324105 0.65328105 -0.41041207 -0.5310114 -0.35045076 -0.65327485\n", 60 " 0.44828037 0.57780436 0.19626784 0.65314259 -0.44838273 0.18210056\n", 61 " -0.62814265 -0.60853106 0.44940642 0.75811541 0.39629793 -0.21773652\n", 62 " -0.46989337 -0.65860905 0.63005841 -0.30019911 -0.28132012 0.48236852\n", 63 " -0.62832911 -0.5187234 -0.32163515 0.29769028 0.73425597 -0.47967195\n", 64 " -0.37704813 -0.43765746 -0.66045568 -0.40777599 -0.45326236 0.44812517\n", 65 " 0.65380112 -0.30669785 -0.6051936 -0.3358828 -0.65334514 -0.02324727\n", 66 " 0.45896652 0.60195637 0.65303739 -0.45292671 0.23926742 0.55576214\n", 67 " -0.65478006 0.4482482 0.31696558 0.52987942 0.64640066 -0.44829888\n", 68 " 0.38549514 0.44541757 -0.67232128 0.4483082 0.46400403 0.38569858\n", 69 " -0.66016044 -0.4473534 0.54359334 -0.61275386 -0.39239591 -0.41840701\n", 70 " -0.54758134 0.15479221 -0.82132777 -0.0401838 0.32731399 -0.62451906\n", 71 " 0.64802289 0.2879372 -0.45982671 0.5518971 ]\n", 72 " [ 0.58805175 -0.05423688 -0.06798499 -0.80413403 -0.5455362 0.01059911\n", 73 " -0.83007415 0.11512958 0.54666835 0.26787692 0.09236427 -0.78794957\n", 74 " -0.54682618 -0.26728689 -0.78877515 -0.08586063 0.54683546 0.27063279\n", 75 " -0.284458 -0.74509647 -0.53914128 -0.27059833 0.05982443 0.73244018\n", 76 " -0.62187536 0.27059779 -0.57079841 0.30321086 0.71346104 -0.27059759\n", 77 " 0.54686975 -0.77908366 0.14408408 0.27049265 -0.54687059 -0.57656055\n", 78 " 0.52268989 -0.3078303 0.54736985 -0.34213788 0.43095611 -0.81346418\n", 79 " 0.18839992 -0.24489227 0.33655805 0.62300199 0.66228768 0.08814893\n", 80 " -0.60383867 0.38361569 0.6931433 0.47604162 0.38162238 0.30471331\n", 81 " 0.73136448 -0.5423859 -0.27941799 0.17456422 0.77353787 0.54552404\n", 82 " 0.2712379 0.00992259 0.67047551 -0.69075046 -0.27060662 0.04018302\n", 83 " 0.60116859 -0.7507227 0.27091894 -0.54786238 0.06259372 -0.78982004\n", 84 " -0.2685389 0.5468414 -0.0277812 -0.78831559 0.28062614 -0.54683543\n", 85 " -0.13969699 -0.79016108 -0.23895029 0.54683931 -0.27550642 -0.7469929\n", 86 " -0.25977012 -0.54645887 -0.44674577 -0.35911728 0.56980911 -0.58887229\n", 87 " 0.46580066 -0.34670594 -0.33969343 -0.73988723 -0.62453819 -0.33669113\n", 88 " 0.27858351 -0.64728847 -0.53769625 -0.45259002]\n", 89 " [ 0.30509728 0.93431284 -0.07857044 0.16673885 -0.55200807 0.19652033\n", 90 " 0.25878826 -0.76791633 0.5469731 -0.20466397 -0.78185586 0.21825329\n", 91 " -0.54684481 0.26920035 0.36451884 -0.70399848 0.54683972 -0.27053318\n", 92 " -0.58773768 0.52572562 -0.55225742 0.2705957 0.75410176 0.40915197\n", 93 " 0.43669568 -0.27059837 -0.58485955 0.50059999 -0.57802549 0.27060699\n", 94 " 0.54682743 0.14591815 -0.77885009 -0.27039935 -0.54680038 0.71500562\n", 95 " 0.36081765 0.24502887 0.54638668 -0.49172198 -0.10982104 -0.0510733\n", 96 " -0.86228786 0.29659175 0.32148198 0.62758918 -0.64406094 -0.42645651\n", 97 " -0.4363783 0.50865706 -0.60743461 0.59215986 -0.14697125 0.57280492\n", 98 " -0.5473977 -0.55114273 0.26170673 0.65832555 -0.44016397 0.54770734\n", 99 " -0.26996275 0.73353648 -0.39697567 -0.48077268 0.27054053 0.92054875\n", 100 " 0.24150568 0.14500457 -0.27062647 -0.53278517 -0.77219061 0.21652209\n", 101 " 0.27014777 0.5469465 -0.72882184 0.30754943 -0.27401056 -0.54683426\n", 102 " -0.6772041 0.42098046 0.25523824 0.54683146 -0.59400045 0.52829137\n", 103 " 0.26196833 -0.54721502 -0.44984181 0.35763158 -0.64622441 -0.50213147\n", 104 " 0.42716543 -0.4564788 -0.40350926 0.66808463 -0.62051456 -0.33380815\n", 105 " -0.29503483 0.64535901 -0.53733834 -0.46604763]\n", 106 " [ 0.47524959 -0.13752981 0.82010757 0.28748383 -0.46131163 0.52441889\n", 107 " 0.39284631 0.59820414 0.4487263 -0.70250008 0.53557389 0.13527362\n", 108 " -0.44832047 0.6530443 0.02311554 0.60992662 0.44829587 -0.65334631\n", 109 " 0.59459976 -0.14093311 -0.44690773 0.65328343 -0.58333349 0.44587465\n", 110 " 0.18476611 -0.65328189 -0.40458672 -0.61277215 0.18451055 0.6532846\n", 111 " 0.44828334 0.19462443 0.57802357 -0.65354623 -0.44821293 -0.3509762\n", 112 " -0.44949161 0.68912358 0.44708132 0.25770561 -0.80322109 -0.53689696\n", 113 " -0.01285863 0.64675747 0.62161436 -0.3576043 0.25965269 -0.76005606\n", 114 " -0.22395394 -0.57011668 0.21708613 0.57802385 -0.5418782 -0.59072763\n", 115 " 0.15263178 -0.45881933 0.64593611 -0.60815586 0.04947804 0.44900337\n", 116 " -0.65275892 -0.606434 -0.16313681 -0.42297067 0.65323809 -0.38785937\n", 117 " 0.60796469 0.23029901 -0.65338075 -0.45917771 0.58528181 0.14293683\n", 118 " 0.65281706 0.44820483 0.60628356 0.05656569 -0.65447183 -0.44829835\n", 119 " 0.6109606 0.00491072 0.65248778 0.44828772 0.59662378 -0.11893881\n", 120 " 0.65427844 -0.44923642 0.55006076 0.60635562 0.32208235 -0.47542505\n", 121 " -0.54837407 -0.80464991 0.21727711 -0.06793428 0.34318986 -0.62063447\n", 122 " -0.64452464 -0.28569487 -0.45903121 0.52284938]]\n", 123 "\n", 124 "\n", 125 "[[ 10.36 2.58 -0. ... -0.08 -10.67 -0.05]\n", 126 " [ 2.58 1.94 1.05 ... 0.31 -3.03 -0.68]\n", 127 " [ -0. 1.05 2.91 ... -0.3 -0.39 0.11]\n", 128 " ...\n", 129 " [ -0.08 0.31 -0.3 ... 0.59 0. 0. ]\n", 130 " [-10.67 -3.03 -0.39 ... 0. 11.09 0.21]\n", 131 " [ -0.05 -0.68 0.11 ... 0. 0.21 0.91]]\n", 132 "\n", 133 "\n", 134 "[ 3.41 0.59 11.09 0.91]\n", 135 "[11.1 3.41 0.9 0.59]\n" 136 ] 137 } 138 ], 139 "source": [ 140 "A = np.array([[4, 3, 2, 1],\n", 141 " [3, 4, 3 ,2], \n", 142 " [2, 3, 4, 3], \n", 143 " [1, 2, 3, 4]])\n", 144 "b = np.random.rand(4)\n", 145 "V, ritz = lanczos_alg(A, b, N = 100)\n", 146 "\n", 147 "print(V)\n", 148 "print('\\n')\n", 149 "print(np.round(np.conjugate(V).T @ A @ V, 2))\n", 150 "print('\\n')\n", 151 "print(np.round(ritz, 2)[-4:])\n", 152 "print(np.round(np.linalg.eig(A)[0], 2))" 153 ] 154 } 155 ], 156 "metadata": { 157 "kernelspec": { 158 "display_name": "Python 3", 159 "language": "python", 160 "name": "python3" 161 }, 162 "language_info": { 163 "codemirror_mode": { 164 "name": "ipython", 165 "version": 3 166 }, 167 "file_extension": ".py", 168 "mimetype": "text/x-python", 169 "name": "python", 170 "nbconvert_exporter": "python", 171 "pygments_lexer": "ipython3", 172 "version": "3.10.4" 173 } 174 }, 175 "nbformat": 4, 176 "nbformat_minor": 5 177 }