diff --git a/exp-321.ipynb b/exp-321.ipynb new file mode 100644 index 0000000000000000000000000000000000000000..0ba74f06c7412fc45a79f622636d2b6c457d12bd --- /dev/null +++ b/exp-321.ipynb @@ -0,0 +1,108 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### IPython notebook for Example 3.2.1 from the lecture" + ] + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "import numpy as np\n", + "from IPython.display import display, Math\n", + "\n", + "\n", + "def kappa(z):\n", + " return np.power(z, -1.5)*np.exp(-1/(4*z)) / (2*np.sqrt(np.pi))\n", + "\n", + "\n", + "n = 20\n", + "h = 1/n\n", + "\n", + "display(Math(r'\\text{Create }A,b_1,b_2\\text{ from Example 1.7.1}'))\n", + "A = np.zeros((n, n))\n", + "for i in range(n):\n", + " for j in range(i+1):\n", + " A[i, j] = h * kappa((i+1)*h - (j+0.5)*h)\n", + "\n", + "b_1 = A[:, 0]\n", + "b_2 = A[:, 1]" + ], + "outputs": [], + "execution_count": null + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "display(Math(r'U^TAV=\\Sigma'))\n", + "U, sigma, VT = np.linalg.svd(A)" + ], + "outputs": [], + "execution_count": null + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "print('sigma_1 = {:.4f}'.format(sigma[0]))\n", + "print('sigma_19 = {:.4f}'.format(sigma[18]))\n", + "print('sigma_20 = {:.2e}'.format(sigma[19]))" + ], + "outputs": [], + "execution_count": null + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "display(Math(r'\\Sigma_{\\rm reg}^\\dagger:= {\\rm diag}(\\sigma_1^{-1}, \\ldots, \\sigma_{19}^{-1},0)'))\n", + "sigma_dagger = 1 / sigma\n", + "sigma_dagger[-1] = 0\n", + "display(Math(r'A_{\\rm reg}^\\dagger:=V \\Sigma_{\\rm reg}^\\dagger U^T'))\n", + "A_dagger_reg = np.matmul(np.matmul(VT.T, np.diag(sigma_dagger)), U.T)\n", + "display(Math(r'\\tilde z_{\\rm reg}= A_{\\rm reg}^\\dagger(b_1+b_2)'))\n", + "z_reg_tilde = np.matmul(A_dagger_reg, b_1+b_2)\n", + "z = np.zeros(n)" + ], + "outputs": [], + "execution_count": null + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "z[0:2] = 1\n", + "print('||z_reg_tilde-z||_2 = {:.1e}'.format(np.linalg.norm(z_reg_tilde-z)))\n", + "display(Math(r'\\text{In Example 1.7.1, we had }\\|\\tilde z\\|_2=9.5\\cdot10^{24}\\text{ for }\\tilde z=A\\setminus (b_1+b_2)'))\n" + ], + "outputs": [], + "execution_count": null + } + ], + "metadata": { + "anaconda-cloud": {}, + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.6.1" + } + }, + "nbformat": 4, + "nbformat_minor": 1 +} \ No newline at end of file