notes

uni notes
git clone git://popovic.xyz/notes.git
Log | Files | Refs

prb9_p.ipynb (4665B)


      1 {
      2  "cells": [
      3   {
      4    "cell_type": "code",
      5    "execution_count": 1,
      6    "id": "0ae2f006",
      7    "metadata": {},
      8    "outputs": [],
      9    "source": [
     10     "import numpy as np"
     11    ]
     12   },
     13   {
     14    "cell_type": "code",
     15    "execution_count": 184,
     16    "id": "6946cd07",
     17    "metadata": {},
     18    "outputs": [],
     19    "source": [
     20     "def fourier_frame(n: int, d: int) -> (np.ndarray, np.ndarray):\n",
     21     "    ω = lambda i, j : np.exp(-2*np.pi *1j/n * (i-1) * (j-1))\n",
     22     "    H = np.zeros([n, d], dtype='complex')\n",
     23     "    for i in range(n):\n",
     24     "        for j in range(d):\n",
     25     "            H[i, j] = 1/np.sqrt(n) * ω(i, j)\n",
     26     "    S = np.conjugate(H.T) @ H \n",
     27     "    return H, S \n",
     28     "\n",
     29     "def reconstruct(v: np.array, S: np.ndarray, H: np.ndarray) -> np.array:\n",
     30     "    s = 0\n",
     31     "    for j in range(H.shape[0]):\n",
     32     "        s += v @ np.conjugate(H.T)[:,j] * invs(S) @ H[j]\n",
     33     "    return np.round(s, 2)\n",
     34     "\n",
     35     "def qr_eig(A: np.ndarray, eps: float=1e-5) -> np.ndarray:\n",
     36     "    A_k = np.copy(A)\n",
     37     "    while True:\n",
     38     "        Q, R = np.linalg.qr(A_k)\n",
     39     "        A_1 = np.copy(A_k)\n",
     40     "        A_k = R@Q\n",
     41     "        if np.linalg.norm(A - A_1)/np.linalg.norm(A_1) < eps:\n",
     42     "            break\n",
     43     "    return np.diagonal(A_k)\n",
     44     "\n",
     45     "# in our particular case them matrix S is unitary\n",
     46     "def invs(A: np.ndarray) -> np.ndarray:\n",
     47     "    return np.linalg.inv(A)"
     48    ]
     49   },
     50   {
     51    "cell_type": "code",
     52    "execution_count": 185,
     53    "id": "75be68d5",
     54    "metadata": {},
     55    "outputs": [
     56     {
     57      "name": "stdout",
     58      "output_type": "stream",
     59      "text": [
     60       "(1-0j) (1-0j) (1+0j) (1+0j) (1+0j) (1-0j)\n"
     61      ]
     62     }
     63    ],
     64    "source": [
     65     "d = 6; n = 10\n",
     66     "H, S =  fourier_frame(n, d)\n",
     67     "print(*np.round(qr_eig(S), 2))"
     68    ]
     69   },
     70   {
     71    "cell_type": "code",
     72    "execution_count": 186,
     73    "id": "fbe1004e",
     74    "metadata": {},
     75    "outputs": [
     76     {
     77      "name": "stdout",
     78      "output_type": "stream",
     79      "text": [
     80       "(1-0j) (2-0j) (3+0j) (4+0j) (5+0j) (6-0j)\n"
     81      ]
     82     }
     83    ],
     84    "source": [
     85     "v = np.arange(1, d+1) \n",
     86     "print(*reconstruct(v, S, H))"
     87    ]
     88   },
     89   {
     90    "cell_type": "code",
     91    "execution_count": 187,
     92    "id": "24f82555",
     93    "metadata": {},
     94    "outputs": [
     95     {
     96      "name": "stdout",
     97      "output_type": "stream",
     98      "text": [
     99       "(1.44-1.29j) (1.56-1.48j) (2.29+0.38j) (4+1.71j) (5.71+0.68j) (6.44-1.29j)\n"
    100      ]
    101     }
    102    ],
    103    "source": [
    104     "H̃ = H[:-1]\n",
    105     "S̃ = np.conjugate(H̃.T) @ H̃ \n",
    106     "print(*reconstruct(v, S̃, H̃))"
    107    ]
    108   },
    109   {
    110    "cell_type": "code",
    111    "execution_count": null,
    112    "id": "c5bdeb99",
    113    "metadata": {},
    114    "outputs": [],
    115    "source": [
    116     "H, S = fourier_frame(5, 5)"
    117    ]
    118   },
    119   {
    120    "cell_type": "code",
    121    "execution_count": 14,
    122    "id": "ba0f89b5",
    123    "metadata": {},
    124    "outputs": [
    125     {
    126      "data": {
    127       "text/plain": [
    128        "array([[ 1.+0.j,  0.-0.j,  0.+0.j, -0.+0.j, -0.-0.j],\n",
    129        "       [ 0.+0.j,  1.+0.j,  0.-0.j,  0.+0.j, -0.+0.j],\n",
    130        "       [ 0.-0.j,  0.+0.j,  1.+0.j,  0.-0.j,  0.+0.j],\n",
    131        "       [-0.-0.j,  0.-0.j,  0.+0.j,  1.+0.j,  0.-0.j],\n",
    132        "       [-0.+0.j, -0.-0.j,  0.-0.j,  0.+0.j,  1.+0.j]])"
    133       ]
    134      },
    135      "execution_count": 14,
    136      "metadata": {},
    137      "output_type": "execute_result"
    138     }
    139    ],
    140    "source": [
    141     "np.round(np.conjugate(H.T) @ H, 2)"
    142    ]
    143   },
    144   {
    145    "cell_type": "code",
    146    "execution_count": 19,
    147    "id": "85380980",
    148    "metadata": {},
    149    "outputs": [
    150     {
    151      "data": {
    152       "text/plain": [
    153        "array([[ 0.-0.j,  0.+0.j,  1.+0.j,  0.-0.j,  0.+0.j],\n",
    154        "       [ 0.+0.j,  1.+0.j,  0.-0.j,  0.+0.j, -0.+0.j],\n",
    155        "       [ 1.+0.j,  0.-0.j,  0.+0.j, -0.+0.j, -0.-0.j],\n",
    156        "       [ 0.-0.j,  0.+0.j, -0.+0.j, -0.-0.j,  1.+0.j],\n",
    157        "       [ 0.+0.j, -0.+0.j, -0.-0.j,  1.+0.j,  0.-0.j]])"
    158       ]
    159      },
    160      "execution_count": 19,
    161      "metadata": {},
    162      "output_type": "execute_result"
    163     }
    164    ],
    165    "source": [
    166     "np.round(H@H, 2)"
    167    ]
    168   },
    169   {
    170    "cell_type": "code",
    171    "execution_count": null,
    172    "id": "4e2bedda",
    173    "metadata": {},
    174    "outputs": [],
    175    "source": []
    176   }
    177  ],
    178  "metadata": {
    179   "kernelspec": {
    180    "display_name": "Python 3 (ipykernel)",
    181    "language": "python",
    182    "name": "python3"
    183   },
    184   "language_info": {
    185    "codemirror_mode": {
    186     "name": "ipython",
    187     "version": 3
    188    },
    189    "file_extension": ".py",
    190    "mimetype": "text/x-python",
    191    "name": "python",
    192    "nbconvert_exporter": "python",
    193    "pygments_lexer": "ipython3",
    194    "version": "3.10.4"
    195   }
    196  },
    197  "nbformat": 4,
    198  "nbformat_minor": 5
    199 }