{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "import numpy as np\n",
    "from numpy.linalg import norm, eig\n",
    "from EigenFactorization import find_eigenvalues, find_eigenvectors\n",
    "\n",
    "%precision 3\n",
    "np.set_printoptions(precision=3, suppress=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def verify(values, vectors):\n",
    "    for i in range(len(values)):\n",
    "        value = values[i]\n",
    "        vector = vectors[:, i]\n",
    "        \n",
    "        print()\n",
    "        print(f'norm({vector}) = {norm(vector):.3f}')\n",
    "        print(f'Av = {A@vector}')\n",
    "        print(f'𝜆v = {value*vector}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def test_factorization(A, x_index=math.nan, tolerance=0.0001):\n",
    "    n, RQ, values = find_eigenvalues(A, tolerance)\n",
    "    \n",
    "    print(f'Upper-triangular matrix after {n} iterations:')\n",
    "    print(RQ)\n",
    "    print()\n",
    "    \n",
    "    print('Eigenvalues:')\n",
    "    print(values)\n",
    "    print()\n",
    "    print('Eigenvectors:')\n",
    "    vectors = find_eigenvectors(A, values, x_index)\n",
    "    print(vectors)\n",
    "    \n",
    "    verify(values, vectors)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def test_eig(A):\n",
    "    values, vectors = eig(A)\n",
    "    \n",
    "    print('Eigenvalues:')\n",
    "    print(values)\n",
    "    print()\n",
    "    print('Eigenvectors:')\n",
    "    print(vectors)    \n",
    "    \n",
    "    verify(values, vectors)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Symmetric matrix\n",
    "\n",
    "A = np.array([[-5,  2],\n",
    "              [ 2, -2]])\n",
    "\n",
    "test_factorization(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Symmetric matrix\n",
    "\n",
    "A = np.array([[1, 3, 0, 0],\n",
    "              [3, 2, 1, 0],\n",
    "              [0, 1, 3, 4],\n",
    "              [0, 0, 4, 1]])\n",
    "\n",
    "test_factorization(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Symmetric matrix\n",
    "\n",
    "A = np.array([[ 89,  75,  22, 102],\n",
    "              [ 75, 116,  27, 120],\n",
    "              [ 22,  27,  33,  62],\n",
    "              [102, 120,  62, 200]])\n",
    "\n",
    "test_factorization(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.array([[5, 1], \n",
    "              [3, 3]])\n",
    "\n",
    "test_factorization(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.array([[1, 2,  3],\n",
    "              [3, 2,  1],\n",
    "              [1, 0, -1]])\n",
    "\n",
    "test_factorization(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.array([[1, 0, 0,   0],\n",
    "              [0, 1, 5, -10],\n",
    "              [1, 0, 2,   0],\n",
    "              [1, 0, 0,  3]])\n",
    "\n",
    "test_factorization(A, tolerance=1.0e-5)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# must have x[1] = 1\n",
    "# need high tolerance = 1.0e-12\n",
    "\n",
    "A = np.array([[1, -3, 3],\n",
    "              [3, -5, 3],\n",
    "              [6, -6, 4]])\n",
    "\n",
    "test_factorization(A, tolerance=1.0e-12)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# wrong eigenvectors unless x[1] = 1\n",
    "\n",
    "A = np.array([[ 92, -32, -15],\n",
    "              [-64,  34,  39],\n",
    "              [176, -68, -99]])\n",
    "\n",
    "test_factorization(A, x_index=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# sigular unless x[1] = 1\n",
    "\n",
    "A = np.array([[1, -2, 1],\n",
    "              [0,  0, 0],\n",
    "              [0,  1, 1]])\n",
    "\n",
    "test_factorization(A, x_index=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# singular unless x[1] = 1\n",
    "\n",
    "A = np.array([[6, 3, -9],\n",
    "              [0, 9, 15],\n",
    "              [0, 0, 15]])\n",
    "\n",
    "test_factorization(A, x_index=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# singular always\n",
    "\n",
    "A = np.array([[2,  0],\n",
    "              [0, -1]])\n",
    "\n",
    "test_factorization(A, x_index=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# singular always\n",
    "\n",
    "A = np.array([[1, 0, 0], \n",
    "              [1, 1, 0], \n",
    "              [0, 0, 1]])\n",
    "\n",
    "test_factorization(A, x_index=3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_eig(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.7.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}