{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "b0184a3a-d155-4ea7-82bb-d88b29ddee16",
   "metadata": {},
   "source": [
    "## The Original EigenDecompression.EigenDecompose\n",
    "\n",
    "See\n",
    "\n",
    "* https://mobile.twitter.com/realize_ss/status/1615160291108745216\n",
    "* https://qiita.com/lelele/items/8408410a94f5c6b8f76e"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "9c72fb52-e0c2-4b6e-ada5-b0fba6982ae5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Eigen{ComplexF64, ComplexF64, Matrix{ComplexF64}, Vector{ComplexF64}}\n",
       "values:\n",
       "100-element Vector{ComplexF64}:\n",
       " -2.6713385149783693 - 0.5238483074100407im\n",
       " -2.6713385149783693 + 0.5238483074100407im\n",
       "  -2.380694749857707 - 1.5819003151310398im\n",
       "  -2.380694749857707 + 1.5819003151310398im\n",
       " -2.3445049148115737 - 0.8348768899082996im\n",
       " -2.3445049148115737 + 0.8348768899082996im\n",
       " -2.1226846073905064 - 0.22379425662578242im\n",
       " -2.1226846073905064 + 0.22379425662578242im\n",
       "  -2.049294063421365 - 1.6320269189857641im\n",
       "  -2.049294063421365 + 1.6320269189857641im\n",
       " -2.0178436003856652 + 0.0im\n",
       " -1.9229326770807857 - 1.7886924334849146im\n",
       " -1.9229326770807857 + 1.7886924334849146im\n",
       "                     ⋮\n",
       "  1.9543269536275725 + 0.5040016799489412im\n",
       "  1.9954076975615094 - 1.0426916311345158im\n",
       "  1.9954076975615094 + 1.0426916311345158im\n",
       "  2.0560135421866663 - 0.48888998497523045im\n",
       "  2.0560135421866663 + 0.48888998497523045im\n",
       "   2.416364338786246 + 0.0im\n",
       "   2.486887950907518 - 1.224341365689113im\n",
       "   2.486887950907518 + 1.224341365689113im\n",
       "  2.6754719890150778 - 0.39289719212777263im\n",
       "  2.6754719890150778 + 0.39289719212777263im\n",
       "  3.0130359977740344 + 0.0im\n",
       "    50.0991251975706 + 0.0im\n",
       "vectors:\n",
       "100×100 Matrix{ComplexF64}:\n",
       "    0.0928227+0.0727977im   …  0.00480172+0.0im  -0.0897557+0.0im\n",
       "   -0.0220687+0.0949872im      -0.0683635+0.0im   -0.107816+0.0im\n",
       "    -0.062205-0.0221392im        0.150972+0.0im   -0.104096+0.0im\n",
       " -0.000929597-0.0558679im         0.21009+0.0im   -0.105664+0.0im\n",
       "   -0.0858376+0.0687551im       0.0157017+0.0im  -0.0975072+0.0im\n",
       "    0.0620079-0.0456313im   …    0.120079+0.0im   -0.104935+0.0im\n",
       "    0.0910325-0.0859799im       0.0523112+0.0im  -0.0967437+0.0im\n",
       "   -0.0549273+0.0307295im       0.0888397+0.0im   -0.106674+0.0im\n",
       "     0.142685+0.0105723im       0.0982881+0.0im  -0.0831259+0.0im\n",
       "    0.0213531-0.0756764im      -0.0218936+0.0im  -0.0993371+0.0im\n",
       "   -0.0338946-0.108939im    …   0.0366303+0.0im  -0.0905221+0.0im\n",
       "    0.0784134+0.105462im        0.0501874+0.0im  -0.0885164+0.0im\n",
       "   -0.0326154+0.101677im        0.0967074+0.0im   -0.103145+0.0im\n",
       "             ⋮              ⋱                    \n",
       "    0.0293922-0.0230606im      -0.0531536+0.0im   -0.102552+0.0im\n",
       "    0.0929652+0.0773336im        0.073754+0.0im  -0.0988448+0.0im\n",
       "    0.0585327-0.0135158im   …   0.0340458+0.0im   -0.100852+0.0im\n",
       "    0.0205471-0.0148881im       -0.205876+0.0im   -0.107353+0.0im\n",
       "     0.174649-0.0368603im      -0.0668649+0.0im   -0.108233+0.0im\n",
       "   -0.0567686-0.0222343im       0.0213794+0.0im   -0.100761+0.0im\n",
       "      0.10441-0.0925521im       0.0107639+0.0im  -0.0950483+0.0im\n",
       "    0.0780944+0.102672im    …   -0.166945+0.0im   -0.102682+0.0im\n",
       "   -0.0634175-0.00904586im      0.0309103+0.0im   -0.102103+0.0im\n",
       "     -0.11412-0.0948267im       0.0745764+0.0im  -0.0954268+0.0im\n",
       "    0.0400251+0.0790165im       0.0617443+0.0im   -0.095069+0.0im\n",
       "     0.068089-0.00059643im     -0.0551766+0.0im  -0.0948394+0.0im"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "using LinearAlgebra\n",
    "M = rand(100,100)#対角化したい行列\n",
    "E, P = eigen(M)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "f69a793a-abb9-4c47-ab07-912b4bd711de",
   "metadata": {},
   "outputs": [
    {
     "ename": "LoadError",
     "evalue": "MethodError: no method matching exp(::Eigen{ComplexF64, ComplexF64, Matrix{ComplexF64}, Vector{ComplexF64}})\n\n\u001b[0mClosest candidates are:\n\u001b[0m  exp(\u001b[91m::Union{Float16, Float32, Float64}\u001b[39m)\n\u001b[0m\u001b[90m   @\u001b[39m \u001b[90mBase\u001b[39m \u001b[90mspecial\\\u001b[39m\u001b[90m\u001b[4mexp.jl:325\u001b[24m\u001b[39m\n\u001b[0m  exp(\u001b[91m::Adjoint{T, <:AbstractMatrix} where T\u001b[39m)\n\u001b[0m\u001b[90m   @\u001b[39m \u001b[35mLinearAlgebra\u001b[39m \u001b[90mD:\\Julia-1.9.0-beta2\\share\\julia\\stdlib\\v1.9\\LinearAlgebra\\src\\\u001b[39m\u001b[90m\u001b[4mdense.jl:595\u001b[24m\u001b[39m\n\u001b[0m  exp(\u001b[91m::Transpose{T, <:AbstractMatrix} where T\u001b[39m)\n\u001b[0m\u001b[90m   @\u001b[39m \u001b[35mLinearAlgebra\u001b[39m \u001b[90mD:\\Julia-1.9.0-beta2\\share\\julia\\stdlib\\v1.9\\LinearAlgebra\\src\\\u001b[39m\u001b[90m\u001b[4mdense.jl:596\u001b[24m\u001b[39m\n\u001b[0m  ...\n",
     "output_type": "error",
     "traceback": [
      "MethodError: no method matching exp(::Eigen{ComplexF64, ComplexF64, Matrix{ComplexF64}, Vector{ComplexF64}})\n\n\u001b[0mClosest candidates are:\n\u001b[0m  exp(\u001b[91m::Union{Float16, Float32, Float64}\u001b[39m)\n\u001b[0m\u001b[90m   @\u001b[39m \u001b[90mBase\u001b[39m \u001b[90mspecial\\\u001b[39m\u001b[90m\u001b[4mexp.jl:325\u001b[24m\u001b[39m\n\u001b[0m  exp(\u001b[91m::Adjoint{T, <:AbstractMatrix} where T\u001b[39m)\n\u001b[0m\u001b[90m   @\u001b[39m \u001b[35mLinearAlgebra\u001b[39m \u001b[90mD:\\Julia-1.9.0-beta2\\share\\julia\\stdlib\\v1.9\\LinearAlgebra\\src\\\u001b[39m\u001b[90m\u001b[4mdense.jl:595\u001b[24m\u001b[39m\n\u001b[0m  exp(\u001b[91m::Transpose{T, <:AbstractMatrix} where T\u001b[39m)\n\u001b[0m\u001b[90m   @\u001b[39m \u001b[35mLinearAlgebra\u001b[39m \u001b[90mD:\\Julia-1.9.0-beta2\\share\\julia\\stdlib\\v1.9\\LinearAlgebra\\src\\\u001b[39m\u001b[90m\u001b[4mdense.jl:596\u001b[24m\u001b[39m\n\u001b[0m  ...\n",
      "",
      "Stacktrace:",
      " [1] top-level scope",
      "   @ In[2]:1"
     ]
    }
   ],
   "source": [
    "exp(eigen(M))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "0ceab85b-f2a5-43b9-a45a-19bebb655d60",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Main.EigenDecompression"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "module EigenDecompression\n",
    "\n",
    "export EigenDecompose, eigDecomp\n",
    "using LinearAlgebra\n",
    "import Base.*, Base./\n",
    "\n",
    "#対角化された行列型\n",
    "struct EigenDecompose{T<:Number} <: AbstractMatrix{T}\n",
    "    P::AbstractMatrix{T}\n",
    "    D::Diagonal{T}\n",
    "    invP::AbstractMatrix{T}\n",
    "end\n",
    "\n",
    "#普通のMatrixを対角化する\n",
    "function eigDecomp(mat::AbstractMatrix)\n",
    "    E, P = eigen(mat)\n",
    "    EigenDecompose(P, Diagonal(E), inv(P))\n",
    "end\n",
    "\n",
    "#EigenDecompose型に対する関数\n",
    "Base.exp(eig::EigenDecompose) = EigenDecompose(eig.P, exp(eig.D), eig.invP)\n",
    "*(eig::EigenDecompose, vec::AbstractVector) = eig.P * eig.D * eig.invP * vec\n",
    "*(eig::EigenDecompose, sc::Number) = EigenDecompose(eig.P, eig.D*sc, eig.invP)\n",
    "/(eig::EigenDecompose, sc::Number) = EigenDecompose(eig.P, eig.D/sc, eig.invP)\n",
    "\n",
    "#普通のMatrixに戻す\n",
    "Base.Array(eig::EigenDecompose) = eig.P * eig.D * eig.invP\n",
    "\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "a6f2726f-780f-4cbc-a980-e13b78f159a6",
   "metadata": {},
   "outputs": [],
   "source": [
    "using .EigenDecompression\n",
    "\n",
    "M = rand(100, 100)\n",
    "eM = eigDecomp(M)\n",
    "for i in 1:100\n",
    "    v = rand(100)\n",
    "    rnd = rand()\n",
    "    @assert exp(M*rnd)*v ≈ exp(eM*rnd)*v\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "b3566c65-ea4a-48b2-9736-00e7904ebe3b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "bench2 (generic function with 1 method)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "using BenchmarkTools\n",
    "\n",
    "M = rand(100, 100);\n",
    "\n",
    "#普通な方\n",
    "function bench1(M)\n",
    "    for i in 1:100\n",
    "        v = rand(100)\n",
    "        exp(M*rand())*v\n",
    "    end\n",
    "end\n",
    "\n",
    "#今回実装した方\n",
    "function bench2(M)\n",
    "    eM = eigDecomp(M)\n",
    "    for i in 1:100\n",
    "        v = rand(100)\n",
    "        exp(eM*rand())*v\n",
    "    end\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "7ce7e7d9-dc0b-4384-8296-a0c0c884deae",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "BenchmarkTools.Trial: 58 samples with 1 evaluation.\n",
       " Range \u001b[90m(\u001b[39m\u001b[36m\u001b[1mmin\u001b[22m\u001b[39m … \u001b[35mmax\u001b[39m\u001b[90m):  \u001b[39m\u001b[36m\u001b[1m78.395 ms\u001b[22m\u001b[39m … \u001b[35m100.653 ms\u001b[39m  \u001b[90m┊\u001b[39m GC \u001b[90m(\u001b[39mmin … max\u001b[90m): \u001b[39m2.43% … 3.65%\n",
       " Time  \u001b[90m(\u001b[39m\u001b[34m\u001b[1mmedian\u001b[22m\u001b[39m\u001b[90m):     \u001b[39m\u001b[34m\u001b[1m86.035 ms               \u001b[22m\u001b[39m\u001b[90m┊\u001b[39m GC \u001b[90m(\u001b[39mmedian\u001b[90m):    \u001b[39m2.82%\n",
       " Time  \u001b[90m(\u001b[39m\u001b[32m\u001b[1mmean\u001b[22m\u001b[39m ± \u001b[32mσ\u001b[39m\u001b[90m):   \u001b[39m\u001b[32m\u001b[1m87.298 ms\u001b[22m\u001b[39m ± \u001b[32m  4.854 ms\u001b[39m  \u001b[90m┊\u001b[39m GC \u001b[90m(\u001b[39mmean ± σ\u001b[90m):  \u001b[39m3.62% ± 1.60%\n",
       "\n",
       "  \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m▂\u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m█\u001b[39m \u001b[39m \u001b[39m \u001b[34m \u001b[39m\u001b[39m \u001b[39m \u001b[39m \u001b[32m \u001b[39m\u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \n",
       "  \u001b[39m▄\u001b[39m▄\u001b[39m▁\u001b[39m▁\u001b[39m▄\u001b[39m▄\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m█\u001b[39m▁\u001b[39m▁\u001b[39m█\u001b[39m▄\u001b[39m▆\u001b[39m▁\u001b[39m█\u001b[39m█\u001b[39m█\u001b[39m▄\u001b[34m▆\u001b[39m\u001b[39m▄\u001b[39m▄\u001b[39m▁\u001b[32m█\u001b[39m\u001b[39m█\u001b[39m▁\u001b[39m▄\u001b[39m▄\u001b[39m▆\u001b[39m▄\u001b[39m▁\u001b[39m▄\u001b[39m▁\u001b[39m▁\u001b[39m▆\u001b[39m▁\u001b[39m▄\u001b[39m▄\u001b[39m▆\u001b[39m▁\u001b[39m▁\u001b[39m▄\u001b[39m▆\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▄\u001b[39m▆\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▄\u001b[39m \u001b[39m▁\n",
       "  78.4 ms\u001b[90m         Histogram: frequency by time\u001b[39m         98.7 ms \u001b[0m\u001b[1m<\u001b[22m\n",
       "\n",
       " Memory estimate\u001b[90m: \u001b[39m\u001b[33m53.80 MiB\u001b[39m, allocs estimate\u001b[90m: \u001b[39m\u001b[33m1900\u001b[39m."
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "@benchmark bench1(M)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "8441c537-41ff-425b-a83b-bc7c455a1ff4",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "BenchmarkTools.Trial: 690 samples with 1 evaluation.\n",
       " Range \u001b[90m(\u001b[39m\u001b[36m\u001b[1mmin\u001b[22m\u001b[39m … \u001b[35mmax\u001b[39m\u001b[90m):  \u001b[39m\u001b[36m\u001b[1m6.103 ms\u001b[22m\u001b[39m … \u001b[35m15.411 ms\u001b[39m  \u001b[90m┊\u001b[39m GC \u001b[90m(\u001b[39mmin … max\u001b[90m): \u001b[39m0.00% … 0.00%\n",
       " Time  \u001b[90m(\u001b[39m\u001b[34m\u001b[1mmedian\u001b[22m\u001b[39m\u001b[90m):     \u001b[39m\u001b[34m\u001b[1m6.824 ms              \u001b[22m\u001b[39m\u001b[90m┊\u001b[39m GC \u001b[90m(\u001b[39mmedian\u001b[90m):    \u001b[39m0.00%\n",
       " Time  \u001b[90m(\u001b[39m\u001b[32m\u001b[1mmean\u001b[22m\u001b[39m ± \u001b[32mσ\u001b[39m\u001b[90m):   \u001b[39m\u001b[32m\u001b[1m7.234 ms\u001b[22m\u001b[39m ± \u001b[32m 1.191 ms\u001b[39m  \u001b[90m┊\u001b[39m GC \u001b[90m(\u001b[39mmean ± σ\u001b[90m):  \u001b[39m0.48% ± 2.20%\n",
       "\n",
       "  \u001b[39m \u001b[39m \u001b[39m \u001b[39m█\u001b[39m▅\u001b[39m▂\u001b[39m▁\u001b[34m \u001b[39m\u001b[39m \u001b[39m \u001b[39m \u001b[32m \u001b[39m\u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \u001b[39m \n",
       "  \u001b[39m▃\u001b[39m▅\u001b[39m▇\u001b[39m█\u001b[39m█\u001b[39m█\u001b[39m█\u001b[34m▇\u001b[39m\u001b[39m█\u001b[39m▇\u001b[39m▆\u001b[32m▇\u001b[39m\u001b[39m▄\u001b[39m▃\u001b[39m▃\u001b[39m▃\u001b[39m▃\u001b[39m▃\u001b[39m▃\u001b[39m▃\u001b[39m▃\u001b[39m▃\u001b[39m▃\u001b[39m▃\u001b[39m▂\u001b[39m▃\u001b[39m▂\u001b[39m▃\u001b[39m▂\u001b[39m▃\u001b[39m▂\u001b[39m▂\u001b[39m▂\u001b[39m▂\u001b[39m▂\u001b[39m▃\u001b[39m▂\u001b[39m▂\u001b[39m▁\u001b[39m▂\u001b[39m▂\u001b[39m▂\u001b[39m▂\u001b[39m▂\u001b[39m▂\u001b[39m▁\u001b[39m▂\u001b[39m▂\u001b[39m▁\u001b[39m▂\u001b[39m▂\u001b[39m▂\u001b[39m▂\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▁\u001b[39m▂\u001b[39m \u001b[39m▃\n",
       "  6.1 ms\u001b[90m         Histogram: frequency by time\u001b[39m        12.3 ms \u001b[0m\u001b[1m<\u001b[22m\n",
       "\n",
       " Memory estimate\u001b[90m: \u001b[39m\u001b[33m1.82 MiB\u001b[39m, allocs estimate\u001b[90m: \u001b[39m\u001b[33m1728\u001b[39m."
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "@benchmark bench2(M)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a9ecc6b7-4601-40b9-937f-05854d56f58b",
   "metadata": {},
   "source": [
    "## EigenDecomposedMatrices.EigenDecomposed"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "57691738-e725-481e-ba43-9d7fea5b8c60",
   "metadata": {},
   "outputs": [],
   "source": [
    "using LinearAlgebra\n",
    "using BenchmarkTools"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "a3672e5b-5030-4952-b50a-da95ef2eb11e",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Main.EigenDecomposedMatrices"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "module EigenDecomposedMatrices\n",
    "\n",
    "export EigenDecomposed\n",
    "\n",
    "using LinearAlgebra\n",
    "using Memoization\n",
    "\n",
    "struct EigenDecomposed{\n",
    "        T,\n",
    "        TE<:AbstractVector{T},\n",
    "        TP<:AbstractMatrix{T},\n",
    "        TinvP<:AbstractMatrix{T}\n",
    "    } <: AbstractMatrix{T}\n",
    "    E::TE\n",
    "    P::TP\n",
    "    invP::TinvP\n",
    "end\n",
    "\n",
    "function EigenDecomposed(A::AbstractMatrix)\n",
    "    E, P = eigen(A)\n",
    "    invP = ishermitian(A) ? P' : inv(P)\n",
    "    EigenDecomposed(E, P, invP)\n",
    "end\n",
    "\n",
    "LinearAlgebra.eigvals(ed::EigenDecomposed) = ed.E\n",
    "LinearAlgebra.eigvecs(ed::EigenDecomposed) = ed.P\n",
    "inveigvecs(ed::EigenDecomposed) = ed.invP\n",
    "@memoize Base.parent(ed::EigenDecomposed) = eigvecs(ed) * Diagonal(eigvals(ed)) * inveigvecs(ed)\n",
    "Base.convert(::Type{Array}, ed::EigenDecomposed) = convert(Array, parent(ed))\n",
    "for op in (:eltype, :size)\n",
    "    @eval Base.$op(ed::EigenDecomposed) = $op(eigvecs(ed))\n",
    "end\n",
    "Base.getindex(ed::EigenDecomposed, I...) = getindex(parent(ed), I...)\n",
    "\n",
    "Base.:*(c::Number, ed::EigenDecomposed) = EigenDecomposed(c*eigvals(ed), eigvecs(ed), inveigvecs(ed))\n",
    "Base.:*(ed::EigenDecomposed, c::Number) = EigenDecomposed(eigvals(ed)*c, eigvecs(ed), inveigvecs(ed))\n",
    "Base.:\\(c::Number, ed::EigenDecomposed) = EigenDecomposed(c\\eigvals(ed), eigvecs(ed), inveigvecs(ed))\n",
    "Base.:/(ed::EigenDecomposed, c::Number) = EigenDecomposed(eigvals(ed)/c, eigvecs(ed), inveigvecs(ed))\n",
    "for T in (AbstractVector, AbstractMatrix)\n",
    "    @eval function Base.:*(ed::EigenDecomposed, v::$T)\n",
    "        E, P, invP = eigvals(ed), eigvecs(ed), inveigvecs(ed)\n",
    "        P * (Diagonal(E) * (invP * v))\n",
    "    end\n",
    "end\n",
    "\n",
    "function exp_old(ed::EigenDecomposed)\n",
    "    E, P, invP = eigvals(ed), eigvecs(ed), inveigvecs(ed)\n",
    "    expE = exp.(E)\n",
    "    expA = P * Diagonal(expE) * invP \n",
    "    EigenDecomposed(expE, P, invP)\n",
    "end\n",
    "\n",
    "LinearAlgebra.lmul!(c::Number, ed::EigenDecomposed) = lmul!(c, eigvals(ed))\n",
    "LinearAlgebra.rmul!(ed::EigenDecomposed, c::Number) = rmul!(eigvals(ed), c)\n",
    "LinearAlgebra.ldiv!(c::Number, ed::EigenDecomposed) = ldiv!(c, eigvals(ed))\n",
    "LinearAlgebra.rdiv!(ed::EigenDecomposed, c::Number) = rdiv!(eigvals(ed), c)\n",
    "\n",
    "for op in (:exp, :log, :sin, :cos)\n",
    "    opE = Symbol(op, \"E\")\n",
    "    op_eigendecomposed = Symbol(op, \"_eigendecomposed\")\n",
    "    op_eigendecomposed! = Symbol(op_eigendecomposed, \"!\")\n",
    "    op_eigendecomposed!_doc =\n",
    "        \"\"\"\n",
    "        $op_eigendecomposed!(Y, ed::EigenDecomposed, $opE=similar(ed.E), tmpY=similar(Y))\n",
    "\n",
    "        returns the `$op` of `ed` and stores the result in `Y`, overwriting the existing value of `Y`. \n",
    "        It does not overwrite `ed` and uses `$opE` and `tmpY` as workspaces.\n",
    "        \"\"\"\n",
    "    @eval begin\n",
    "        @doc $op_eigendecomposed!_doc\n",
    "        function $op_eigendecomposed!(Y, ed::EigenDecomposed, $opE=similar(ed.E), tmpY=similar(Y))\n",
    "            E, P, invP = eigvals(ed), eigvecs(ed), inveigvecs(ed)\n",
    "            @. $opE = $op.(E)\n",
    "            mul!(tmpY, P, Diagonal($opE))\n",
    "            mul!(Y, tmpY, invP)\n",
    "        end\n",
    "        $op_eigendecomposed(ed::EigenDecomposed) = $op_eigendecomposed!(similar(eigvecs(ed)), ed)\n",
    "        Base.$op(ed::EigenDecomposed) = $op_eigendecomposed(ed)\n",
    "    end\n",
    "end\n",
    "\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "f4402c98-850b-44e0-9407-6088e6421f19",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "exp\\_eigendecomposed!(Y, ed::EigenDecomposed, expE=similar(ed.E), tmpY=similar(Y))\n",
       "\n",
       "returns the \\texttt{exp} of \\texttt{ed} and stores the result in \\texttt{Y}, overwriting the existing value of \\texttt{Y}.  It does not overwrite \\texttt{ed} and uses \\texttt{expE} and \\texttt{tmpY} as workspaces.\n",
       "\n"
      ],
      "text/markdown": [
       "exp_eigendecomposed!(Y, ed::EigenDecomposed, expE=similar(ed.E), tmpY=similar(Y))\n",
       "\n",
       "returns the `exp` of `ed` and stores the result in `Y`, overwriting the existing value of `Y`.  It does not overwrite `ed` and uses `expE` and `tmpY` as workspaces.\n"
      ],
      "text/plain": [
       "  exp_eigendecomposed!(Y, ed::EigenDecomposed, expE=similar(ed.E),\n",
       "  tmpY=similar(Y))\n",
       "\n",
       "  returns the \u001b[36mexp\u001b[39m of \u001b[36med\u001b[39m and stores the result in \u001b[36mY\u001b[39m, overwriting the existing\n",
       "  value of \u001b[36mY\u001b[39m. It does not overwrite \u001b[36med\u001b[39m and uses \u001b[36mexpE\u001b[39m and \u001b[36mtmpY\u001b[39m as workspaces."
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "?EigenDecomposedMatrices.exp_eigendecomposed!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "beb73095-1445-4318-b010-1f4c495edec6",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "log\\_eigendecomposed!(Y, ed::EigenDecomposed, logE=similar(ed.E), tmpY=similar(Y))\n",
       "\n",
       "returns the \\texttt{log} of \\texttt{ed} and stores the result in \\texttt{Y}, overwriting the existing value of \\texttt{Y}.  It does not overwrite \\texttt{ed} and uses \\texttt{logE} and \\texttt{tmpY} as workspaces.\n",
       "\n"
      ],
      "text/markdown": [
       "log_eigendecomposed!(Y, ed::EigenDecomposed, logE=similar(ed.E), tmpY=similar(Y))\n",
       "\n",
       "returns the `log` of `ed` and stores the result in `Y`, overwriting the existing value of `Y`.  It does not overwrite `ed` and uses `logE` and `tmpY` as workspaces.\n"
      ],
      "text/plain": [
       "  log_eigendecomposed!(Y, ed::EigenDecomposed, logE=similar(ed.E),\n",
       "  tmpY=similar(Y))\n",
       "\n",
       "  returns the \u001b[36mlog\u001b[39m of \u001b[36med\u001b[39m and stores the result in \u001b[36mY\u001b[39m, overwriting the existing\n",
       "  value of \u001b[36mY\u001b[39m. It does not overwrite \u001b[36med\u001b[39m and uses \u001b[36mlogE\u001b[39m and \u001b[36mtmpY\u001b[39m as workspaces."
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "?EigenDecomposedMatrices.log_eigendecomposed!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "312bddb6-0a2d-4caf-8b60-d421be7143f3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "# 2 methods for type constructor:<ul><li> Main.EigenDecomposedMatrices.EigenDecomposed(E::<b>TE</b>, P::<b>TP</b>, invP::<b>TinvP</b>)<i> where {T, TE<:AbstractVector{T}, TP<:AbstractMatrix{T}, TinvP<:AbstractMatrix{T}}</i> in Main.EigenDecomposedMatrices at In[9]:14</li> <li> Main.EigenDecomposedMatrices.EigenDecomposed(A::<b>AbstractMatrix</b>) in Main.EigenDecomposedMatrices at In[9]:19</li> </ul>"
      ],
      "text/plain": [
       "# 2 methods for type constructor:\n",
       " [1] Main.EigenDecomposedMatrices.EigenDecomposed(\u001b[90mE\u001b[39m::\u001b[1mTE\u001b[22m, \u001b[90mP\u001b[39m::\u001b[1mTP\u001b[22m, \u001b[90minvP\u001b[39m::\u001b[1mTinvP\u001b[22m) where {T, TE<:AbstractVector{T}, TP<:AbstractMatrix{T}, TinvP<:AbstractMatrix{T}}\n",
       "\u001b[90m     @\u001b[39m \u001b[90m\u001b[4mIn[9]:14\u001b[24m\u001b[39m\n",
       " [2] Main.EigenDecomposedMatrices.EigenDecomposed(\u001b[90mA\u001b[39m::\u001b[1mAbstractMatrix\u001b[22m)\n",
       "\u001b[90m     @\u001b[39m \u001b[90m\u001b[4mIn[9]:19\u001b[24m\u001b[39m"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "methods(EigenDecomposedMatrices.EigenDecomposed)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "880f55dc-77ea-46e3-8bdf-d4bb540a005a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "# 1 method for type constructor:<ul><li> (var\"#ctor-self#\"::<b>Type{Main.EigenDecomposedMatrices.EigenDecomposed{T, TE, TP, TinvP}} where {T, TE<:AbstractVector{T}, TP<:AbstractMatrix{T}, TinvP<:AbstractMatrix{T}}</b>)(E, P, invP) in Main.EigenDecomposedMatrices at In[9]:14</li> </ul>"
      ],
      "text/plain": [
       "# 1 method for type constructor:\n",
       " [1] (var\"#ctor-self#\"::Type{Main.EigenDecomposedMatrices.EigenDecomposed{T, TE, TP, TinvP}} where {T, TE<:AbstractVector{T}, TP<:AbstractMatrix{T}, TinvP<:AbstractMatrix{T}})(\u001b[90mE\u001b[39m, \u001b[90mP\u001b[39m, \u001b[90minvP\u001b[39m)\n",
       "\u001b[90m     @\u001b[39m \u001b[90m\u001b[4mIn[9]:14\u001b[24m\u001b[39m"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "methods(EigenDecomposedMatrices.EigenDecomposed{Float64, Vector{Float64}, Matrix{Float64}, Matrix{Float64}})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "3b092734-4748-40e2-9081-9e4309d93e8d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "21-element Vector{Method}:<ul><li> *(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, c::<b>Number</b>) in Main.EigenDecomposedMatrices at In[9]:36<li> *(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, v::<b>AbstractVector</b>) in Main.EigenDecomposedMatrices at In[9]:40<li> *(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, v::<b>AbstractMatrix</b>) in Main.EigenDecomposedMatrices at In[9]:40<li> *(c::<b>Number</b>, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:35<li> /(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, c::<b>Number</b>) in Main.EigenDecomposedMatrices at In[9]:38<li> \\(c::<b>Number</b>, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:37<li> convert(::<b>Type{Array}</b>, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:29<li> cos(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:78<li> eltype(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:31<li> exp(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:78<li> getindex(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, I...) in Main.EigenDecomposedMatrices at In[9]:33<li> log(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:78<li> parent(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at <a href=\"file://D:/.julia/packages/Memoization/ut5GT/src/Memoization.jl\" target=\"_blank\">D:\\.julia\\packages\\Memoization\\ut5GT\\src\\Memoization.jl:162</a><li> sin(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:78<li> size(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:31<li> eigvals(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:25<li> eigvecs(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:26<li> ldiv!(c::<b>Number</b>, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:55<li> lmul!(c::<b>Number</b>, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:53<li> rdiv!(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, c::<b>Number</b>) in Main.EigenDecomposedMatrices at In[9]:56<li> rmul!(ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, c::<b>Number</b>) in Main.EigenDecomposedMatrices at In[9]:54</ul>"
      ],
      "text/plain": [
       "[1] *(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mc\u001b[39m::\u001b[1mNumber\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:36\u001b[24m\u001b[39m\n",
       "[2] *(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mv\u001b[39m::\u001b[1mAbstractVector\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:40\u001b[24m\u001b[39m\n",
       "[3] *(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mv\u001b[39m::\u001b[1mAbstractMatrix\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:40\u001b[24m\u001b[39m\n",
       "[4] *(\u001b[90mc\u001b[39m::\u001b[1mNumber\u001b[22m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:35\u001b[24m\u001b[39m\n",
       "[5] /(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mc\u001b[39m::\u001b[1mNumber\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:38\u001b[24m\u001b[39m\n",
       "[6] \\(\u001b[90mc\u001b[39m::\u001b[1mNumber\u001b[22m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:37\u001b[24m\u001b[39m\n",
       "[7] convert(::\u001b[1mType\u001b[22m\u001b[0m{Array}, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:29\u001b[24m\u001b[39m\n",
       "[8] cos(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:78\u001b[24m\u001b[39m\n",
       "[9] eltype(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:31\u001b[24m\u001b[39m\n",
       "[10] exp(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:78\u001b[24m\u001b[39m\n",
       "[11] getindex(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mI...\u001b[39m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:33\u001b[24m\u001b[39m\n",
       "[12] log(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:78\u001b[24m\u001b[39m\n",
       "[13] parent(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90mD:\\.julia\\packages\\Memoization\\ut5GT\\src\\\u001b[39m\u001b[90m\u001b[4mMemoization.jl:162\u001b[24m\u001b[39m\n",
       "[14] sin(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:78\u001b[24m\u001b[39m\n",
       "[15] size(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:31\u001b[24m\u001b[39m\n",
       "[16] eigvals(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:25\u001b[24m\u001b[39m\n",
       "[17] eigvecs(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:26\u001b[24m\u001b[39m\n",
       "[18] ldiv!(\u001b[90mc\u001b[39m::\u001b[1mNumber\u001b[22m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:55\u001b[24m\u001b[39m\n",
       "[19] lmul!(\u001b[90mc\u001b[39m::\u001b[1mNumber\u001b[22m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:53\u001b[24m\u001b[39m\n",
       "[20] rdiv!(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mc\u001b[39m::\u001b[1mNumber\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:56\u001b[24m\u001b[39m\n",
       "[21] rmul!(\u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mc\u001b[39m::\u001b[1mNumber\u001b[22m)\u001b[90m @\u001b[39m \u001b[90mMain.EigenDecomposedMatrices\u001b[39m \u001b[90m\u001b[4mIn[9]:54\u001b[24m\u001b[39m"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "methodswith(EigenDecomposedMatrices.EigenDecomposed)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "d00504b2-cf67-4d6a-abdf-1f71f61d9c30",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "# 3 methods for generic function <b>exp_eigendecomposed!</b> from \u001b[36mMain.EigenDecomposedMatrices\u001b[39m:<ul><li> exp_eigendecomposed!(Y, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:71</li> <li> exp_eigendecomposed!(Y, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, expE) in Main.EigenDecomposedMatrices at In[9]:71</li> <li> exp_eigendecomposed!(Y, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, expE, tmpY) in Main.EigenDecomposedMatrices at In[9]:71</li> </ul>"
      ],
      "text/plain": [
       "# 3 methods for generic function \"exp_eigendecomposed!\" from \u001b[36mMain.EigenDecomposedMatrices\u001b[39m:\n",
       " [1] exp_eigendecomposed!(\u001b[90mY\u001b[39m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\n",
       "\u001b[90m     @\u001b[39m \u001b[90m\u001b[4mIn[9]:71\u001b[24m\u001b[39m\n",
       " [2] exp_eigendecomposed!(\u001b[90mY\u001b[39m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mexpE\u001b[39m)\n",
       "\u001b[90m     @\u001b[39m \u001b[90m\u001b[4mIn[9]:71\u001b[24m\u001b[39m\n",
       " [3] exp_eigendecomposed!(\u001b[90mY\u001b[39m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mexpE\u001b[39m, \u001b[90mtmpY\u001b[39m)\n",
       "\u001b[90m     @\u001b[39m \u001b[90m\u001b[4mIn[9]:71\u001b[24m\u001b[39m"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "methods(EigenDecomposedMatrices.exp_eigendecomposed!)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "e2e20647-3c4e-4fd2-8c18-4a2f8e439fd0",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "# 3 methods for generic function <b>log_eigendecomposed!</b> from \u001b[36mMain.EigenDecomposedMatrices\u001b[39m:<ul><li> log_eigendecomposed!(Y, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>) in Main.EigenDecomposedMatrices at In[9]:71</li> <li> log_eigendecomposed!(Y, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, logE) in Main.EigenDecomposedMatrices at In[9]:71</li> <li> log_eigendecomposed!(Y, ed::<b>Main.EigenDecomposedMatrices.EigenDecomposed</b>, logE, tmpY) in Main.EigenDecomposedMatrices at In[9]:71</li> </ul>"
      ],
      "text/plain": [
       "# 3 methods for generic function \"log_eigendecomposed!\" from \u001b[36mMain.EigenDecomposedMatrices\u001b[39m:\n",
       " [1] log_eigendecomposed!(\u001b[90mY\u001b[39m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m)\n",
       "\u001b[90m     @\u001b[39m \u001b[90m\u001b[4mIn[9]:71\u001b[24m\u001b[39m\n",
       " [2] log_eigendecomposed!(\u001b[90mY\u001b[39m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mlogE\u001b[39m)\n",
       "\u001b[90m     @\u001b[39m \u001b[90m\u001b[4mIn[9]:71\u001b[24m\u001b[39m\n",
       " [3] log_eigendecomposed!(\u001b[90mY\u001b[39m, \u001b[90med\u001b[39m::\u001b[1mMain.EigenDecomposedMatrices.EigenDecomposed\u001b[22m, \u001b[90mlogE\u001b[39m, \u001b[90mtmpY\u001b[39m)\n",
       "\u001b[90m     @\u001b[39m \u001b[90m\u001b[4mIn[9]:71\u001b[24m\u001b[39m"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "methods(EigenDecomposedMatrices.log_eigendecomposed!)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "6f34704c-58d3-46ac-9333-19417c300f84",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3×3 Main.EigenDecomposedMatrices.EigenDecomposed{Float64, Vector{Float64}, Matrix{Float64}, Adjoint{Float64, Matrix{Float64}}}:\n",
       "  2.0          -1.0  -3.33067e-16\n",
       " -1.0           2.0  -1.0\n",
       " -3.33067e-16  -1.0   2.0"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = [\n",
    "    2 -1 0\n",
    "    -1 2 -1\n",
    "    0 -1 2\n",
    "]\n",
    "\n",
    "edA = EigenDecomposedMatrices.EigenDecomposed(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "04d9e6cc-253b-4340-a448-6050a71bc349",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3×3 Matrix{Float64}:\n",
       "  0.51986   -0.623225  -0.173287\n",
       " -0.623225   0.346574  -0.623225\n",
       " -0.173287  -0.623225   0.51986"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "log(edA)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "070a1c06-ed0a-4d68-9716-168bcb1c46d6",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3×3 Matrix{Float64}:\n",
       "  0.51986   -0.623225  -0.173287\n",
       " -0.623225   0.346574  -0.623225\n",
       " -0.173287  -0.623225   0.51986"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "log(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "a64fdfe1-f568-459a-b9c2-5d35aa060134",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "log(edA) ≈ log(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "a2d8c089-f82a-43c6-b324-370b6268ddd5",
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 2^8\n",
    "M = 5I + randn(n, n)\n",
    "v = randn(n)\n",
    "c = randn()\n",
    "\n",
    "edM = EigenDecomposedMatrices.EigenDecomposed(M)\n",
    "\n",
    "Y = similar(eigvecs(edM))\n",
    "expE = similar(eigvals(edM))\n",
    "tmpY = similar(Y)\n",
    "\n",
    "y = similar(eigvals(edM))\n",
    "alpha = randn()\n",
    "beta = randn();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "5fe62250-1d58-4adb-b9a3-352e56f2e075",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "256×256 Main.EigenDecomposedMatrices.EigenDecomposed{ComplexF64, Vector{ComplexF64}, Matrix{ComplexF64}, Matrix{ComplexF64}}:\n",
       "    3.56334-1.65165e-13im  …  -0.0799542+7.35922e-14im\n",
       "   -1.29471-3.33504e-14im        2.79704-1.73821e-15im\n",
       "   0.188711-1.38904e-14im      -0.146065-7.04386e-16im\n",
       "  -0.139691-1.2809e-13im        -0.11464+1.72234e-13im\n",
       "   0.242261+2.90866e-14im       -1.44222+1.48799e-13im\n",
       "   0.476316+3.58411e-14im  …    0.497164-2.26801e-13im\n",
       "   0.550097-3.61273e-14im      -0.134351+1.41119e-13im\n",
       "   0.717261+1.67825e-13im       0.895324-9.16857e-14im\n",
       "   -1.04636-8.49634e-14im      -0.190769-2.40075e-13im\n",
       "     1.6971-2.75795e-14im       0.770272-6.1418e-14im\n",
       "    1.26714+5.92264e-14im  …   -0.567488-1.78777e-13im\n",
       " -0.0854576-1.88943e-14im      -0.822014+4.62595e-14im\n",
       "   0.826429+6.51037e-14im       -1.53237-1.375e-13im\n",
       "           ⋮               ⋱            ⋮\n",
       "  -0.997608-2.67526e-14im      0.0412252+1.99661e-13im\n",
       " -0.0798052-1.49411e-14im  …    0.144672-1.3304e-13im\n",
       "   -1.17006+1.4458e-14im       -0.431936+1.83672e-13im\n",
       " -0.0696381+8.73041e-14im      -0.825939-1.10264e-14im\n",
       "   0.553973+1.47842e-14im        1.41854-6.53152e-14im\n",
       "    1.05755+4.38409e-14im      -0.994764-9.63611e-14im\n",
       "   0.143704+9.43705e-15im  …   -0.252374-1.1153e-13im\n",
       "   0.089135-3.98563e-14im       -1.87436-1.5556e-14im\n",
       "   -1.30653-7.54455e-14im        1.53995-1.05616e-13im\n",
       "   0.642478-5.84129e-14im        1.36199-7.11765e-14im\n",
       "    0.24546-2.13334e-14im       0.401485-2.4925e-14im\n",
       "    -1.3641+7.4302e-14im   …     5.51035+9.54081e-14im"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "edM"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "e613f081-b46d-406a-97dd-9fbff8168a51",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Main.EigenDecomposedMatrices.EigenDecomposed{ComplexF64, Vector{ComplexF64}, Matrix{ComplexF64}, Matrix{ComplexF64}}\n",
      "  E: Array{ComplexF64}((256,)) ComplexF64[-10.529011730788033 + 0.0im, -10.280415425445305 - 5.0635332402051345im, -10.280415425445305 + 5.0635332402051345im, -9.365531653942693 - 7.310921121421777im, -9.365531653942693 + 7.310921121421777im, -9.066423960276218 - 0.8998268797512045im, -9.066423960276218 + 0.8998268797512045im, -8.732971181835605 - 1.4609873383902336im, -8.732971181835605 + 1.4609873383902336im, -8.46124963643927 - 8.33273467746558im  …  18.632067916841073 - 4.578055223585328im, 18.632067916841073 + 4.578055223585328im, 18.93064239638064 - 7.073434597909069im, 18.93064239638064 + 7.073434597909069im, 19.20170503001391 - 3.658032128271094im, 19.20170503001391 + 3.658032128271094im, 20.930683959675996 - 0.7030285754956695im, 20.930683959675996 + 0.7030285754956695im, 21.184754820465372 - 3.0347453281790697im, 21.184754820465372 + 3.0347453281790697im]\n",
      "  P: Array{ComplexF64}((256, 256)) ComplexF64[0.016739760214174452 + 0.0im 0.0479671595323605 - 0.01662350156735747im … 0.08231100940410303 + 0.03720618909538042im 0.08231100940410303 - 0.03720618909538042im; -0.05678105463615567 + 0.0im -0.02753855144836508 - 0.03255712650078475im … 0.01670703557346129 + 0.024130937819675107im 0.01670703557346129 - 0.024130937819675107im; … ; 0.023221059785412463 + 0.0im -0.017275722649348676 - 0.009015532372298877im … 0.049235288759355164 + 0.03697492625403338im 0.049235288759355164 - 0.03697492625403338im; -0.0578255784977335 + 0.0im -0.002679638667629752 - 0.007116991095201977im … 0.02694655853750013 + 0.0004816998404853853im 0.02694655853750013 - 0.0004816998404853853im]\n",
      "  invP: Array{ComplexF64}((256, 256)) ComplexF64[0.3043258839874642 + 4.5105492045176234e-15im 0.11647828660858314 + 8.454429333244175e-15im … 0.15219311634828636 + 1.4654943925052066e-14im -0.17563710572420543 + 2.1649348980190553e-15im; 0.08912104015473599 + 0.06654733900857696im -0.08069004242424725 + 0.18453572669503732im … -0.2156137951599397 + 0.002174603423839006im 0.02083977906498702 - 0.0641592478208322im; … ; 0.0385888730719014 - 0.09966379283416266im 0.11941897834644964 + 0.033343849341449926im … -0.10247759199925907 + 0.07210367195681458im 0.048176002794752926 - 0.18025959146290277im; 0.03858887307190356 + 0.09966379283416041im 0.1194189783464544 - 0.03334384934145064im … -0.1024775919992489 - 0.07210367195681436im 0.048176002794748166 + 0.18025959146290615im]\n"
     ]
    }
   ],
   "source": [
    "dump(edM)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "283dbc81-a65a-4cfa-8b7b-2b4e3444a526",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "M ≈ parent(edM) == Matrix(edM)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "id": "38622469-86c8-4542-97a3-e97138c061fd",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "M ≈ edM"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "id": "4de76f66-fa04-4293-8d77-0b00c3ecb182",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "c*M ≈ c*edM ≈ edM*c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "id": "73f06380-1bc3-45a6-ab28-3c30fd64c4c6",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "c\\M ≈ c\\edM ≈ edM/c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "id": "98def48b-63dc-4dd3-88b0-0a8332c409b5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(\n",
    "    exp(M)\n",
    "    ≈ exp(edM)\n",
    "    ≈ EigenDecomposedMatrices.exp_old(edM)\n",
    "    ≈ EigenDecomposedMatrices.exp_eigendecomposed!(Y, edM)\n",
    "    ≈ EigenDecomposedMatrices.exp_eigendecomposed!(Y, edM, expE, tmpY)\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "id": "23c99064-3d0e-4b04-9810-c176a273316b",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "typeof(y) = Vector{ComplexF64}\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "@show typeof(y)\n",
    "(\n",
    "    exp(M) * v\n",
    "    ≈ exp(edM) * v\n",
    "    ≈ EigenDecomposedMatrices.exp_old(edM) * v\n",
    "    ≈ mul!(y, EigenDecomposedMatrices.exp_eigendecomposed!(Y, edM), v)\n",
    "    ≈ mul!(y, EigenDecomposedMatrices.exp_eigendecomposed!(Y, edM, expE, tmpY), v)\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "id": "bd540d02-5314-4565-8ed1-5eaa48307f7a",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  40.167 ms (26 allocations: 3.53 MiB)\n"
     ]
    }
   ],
   "source": [
    "@btime edM = EigenDecomposedMatrices.EigenDecomposed(M);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "id": "ffb45dc7-2f11-44e1-8acd-e4114f7449a8",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  10.090 ms (16 allocations: 3.01 MiB)\n",
      "  2.712 ms (7 allocations: 2.01 MiB)\n",
      "  2.696 ms (9 allocations: 2.02 MiB)\n",
      "  2.540 ms (3 allocations: 1.00 MiB)\n",
      "  2.447 ms (0 allocations: 0 bytes)\n"
     ]
    }
   ],
   "source": [
    "@btime exp($M) * $v\n",
    "@btime exp($edM) * $v\n",
    "@btime EigenDecomposedMatrices.exp_old($edM) * $v\n",
    "@btime mul!($y, EigenDecomposedMatrices.exp_eigendecomposed!($Y, $edM), $v)\n",
    "@btime mul!($y, EigenDecomposedMatrices.exp_eigendecomposed!($Y, $edM, $expE, $tmpY), $v);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "id": "5d55c783-7d12-4e09-8121-ec398884a7e4",
   "metadata": {},
   "outputs": [],
   "source": [
    "n2 = 2^8\n",
    "M2 = Symmetric(5I + randn(n2, n2))\n",
    "v2 = randn(n)\n",
    "c2 = randn()\n",
    "\n",
    "edM2 = EigenDecomposedMatrices.EigenDecomposed(M2)\n",
    "\n",
    "Y2 = similar(eigvecs(edM2))\n",
    "expE2 = similar(eigvals(edM2))\n",
    "tmpY2 = similar(Y2)\n",
    "\n",
    "y2 = similar(eigvals(edM2))\n",
    "alpha2 = randn()\n",
    "beta2 = randn();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "id": "4ca2d50f-c879-410c-8e8b-0b815fa54a66",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "256×256 Main.EigenDecomposedMatrices.EigenDecomposed{Float64, Vector{Float64}, Matrix{Float64}, Adjoint{Float64, Matrix{Float64}}}:\n",
       "  6.03077    -0.43245     0.599883    …  -0.872836      0.980907\n",
       " -0.43245     7.07582    -0.965762       -2.46577       1.27889\n",
       "  0.599883   -0.965762    6.28325        -0.897865      0.902268\n",
       "  0.554185   -1.61176    -1.88607         1.52239      -0.800269\n",
       " -0.0394612  -0.472386   -2.16329        -0.666488      1.60583\n",
       " -0.684983    0.673263    1.44873     …  -1.11132      -1.06482\n",
       "  0.544451    1.8252      0.567267        0.568297      1.60087\n",
       " -0.60012    -1.6835     -0.332066        0.151455     -1.76098\n",
       " -1.32028    -1.71558    -0.503773       -0.276013      1.46591\n",
       " -0.619846    0.546905   -1.79987         0.501719      0.114942\n",
       " -1.96161     0.575018   -0.648195    …   1.16492      -1.17966\n",
       " -0.303545   -0.670053    0.967682       -0.000749713   0.373673\n",
       "  1.33041    -0.0295437  -0.717992        0.155946      1.85579\n",
       "  ⋮                                   ⋱                 ⋮\n",
       " -0.590561    0.333331    0.250405        0.90949      -1.53465\n",
       " -0.748456    1.29363     0.1356      …  -1.60434       0.683966\n",
       "  1.53178     0.157997    0.00791121     -1.8972        0.444282\n",
       "  0.269089   -0.167004    1.49271        -0.630714     -2.68816\n",
       "  2.1885     -1.43587    -1.18202        -1.17287       1.17593\n",
       " -0.787258   -0.512251   -0.294453        2.0117        0.0495256\n",
       " -0.550615   -0.564627    1.8215      …   0.254222     -0.881607\n",
       "  1.39278     0.444108    0.341001       -0.411732      0.491344\n",
       "  0.960274   -0.737275    0.302913        1.16593      -0.539727\n",
       " -1.638      -0.849312    0.108656       -0.763322      1.67947\n",
       " -0.872836   -2.46577    -0.897865        6.45844       0.754344\n",
       "  0.980907    1.27889     0.902268    …   0.754344      5.11205"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "edM2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "id": "a4d5952e-473f-4f9a-b051-819bbc67a055",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Main.EigenDecomposedMatrices.EigenDecomposed{Float64, Vector{Float64}, Matrix{Float64}, Adjoint{Float64, Matrix{Float64}}}\n",
      "  E: Array{Float64}((256,)) [-26.281246714455474, -26.079893599626686, -25.32277485066867, -25.0801104625073, -23.954875700929986, -23.74285195174835, -23.168019607432065, -22.60925439400774, -22.454572564359246, -21.705343122758094  …  31.569737139516377, 32.20706026446964, 32.53684276189552, 33.124400141791504, 33.974663809226215, 34.44625078387939, 34.499597457606136, 34.84570025054353, 36.07654166938045, 36.899884804719484]\n",
      "  P: Array{Float64}((256, 256)) [0.1449279306214216 0.005919536496859418 … 0.00196042283590344 0.005370199356810865; 0.0977684655406359 -0.04973701370055792 … -0.0878992312227967 -0.05778779773302181; … ; 0.04308743597446929 0.04324139883354045 … -0.019127423566981328 -0.04138098814577151; -0.13325192873325928 -0.028583547120122572 … 0.12277657998369695 -0.056353074728779456]\n",
      "  invP: Adjoint{Float64, Matrix{Float64}}\n",
      "    parent: Array{Float64}((256, 256)) [0.1449279306214216 0.005919536496859418 … 0.00196042283590344 0.005370199356810865; 0.0977684655406359 -0.04973701370055792 … -0.0878992312227967 -0.05778779773302181; … ; 0.04308743597446929 0.04324139883354045 … -0.019127423566981328 -0.04138098814577151; -0.13325192873325928 -0.028583547120122572 … 0.12277657998369695 -0.056353074728779456]\n"
     ]
    }
   ],
   "source": [
    "dump(edM2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "id": "48b70d58-74b6-4a9e-a95b-9973cc775118",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "M2 ≈ parent(edM2) == Matrix(edM2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "id": "0dde64d3-71f4-40de-a2a8-c5bd79e24d64",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "M2 ≈ edM2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "id": "f0a53088-b448-4880-b082-265b9f435b06",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "c2*M2 ≈ c2*edM2 ≈ edM2*c2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "id": "b3cf000b-7d36-4354-ad55-324218517899",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 38,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "c2\\M2 ≈ c2\\edM2 ≈ edM2/c2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "id": "8a9f4b9d-964e-4ddf-950e-7ed6eae067e6",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(\n",
    "    exp(M2)\n",
    "    ≈ exp(edM2)\n",
    "    ≈ EigenDecomposedMatrices.exp_old(edM2)\n",
    "    ≈ EigenDecomposedMatrices.exp_eigendecomposed!(Y2, edM2)\n",
    "    ≈ EigenDecomposedMatrices.exp_eigendecomposed!(Y2, edM2, expE2, tmpY2)\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "id": "2ba1d808-754a-4075-8c5a-122fa6360e62",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "typeof(y2) = Vector{Float64}\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "@show typeof(y2)\n",
    "(\n",
    "    exp(M2) * v2\n",
    "    ≈ exp(edM2) * v2\n",
    "    ≈ EigenDecomposedMatrices.exp_old(edM2) * v2\n",
    "    ≈ mul!(y2, EigenDecomposedMatrices.exp_eigendecomposed!(Y2, edM2), v2)\n",
    "    ≈ mul!(y2, EigenDecomposedMatrices.exp_eigendecomposed!(Y2, edM2, expE2, tmpY2), v2)\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "id": "a1925b60-427c-45ed-9053-2d79a9d412af",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  6.472 ms (14 allocations: 1.59 MiB)\n"
     ]
    }
   ],
   "source": [
    "@btime edM2 = EigenDecomposedMatrices.EigenDecomposed(M2);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "id": "2292854e-259d-4615-8551-8c2c6060b7b5",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  7.427 ms (19 allocations: 2.60 MiB)\n",
      "  726.000 μs (6 allocations: 1.00 MiB)\n",
      "  719.700 μs (8 allocations: 1.01 MiB)\n",
      "  690.300 μs (3 allocations: 514.17 KiB)\n",
      "  662.900 μs (0 allocations: 0 bytes)\n"
     ]
    }
   ],
   "source": [
    "@btime exp($M2) * $v2\n",
    "@btime exp($edM2) * $v2\n",
    "@btime EigenDecomposedMatrices.exp_old($edM2) * $v2\n",
    "@btime mul!($y2, EigenDecomposedMatrices.exp_eigendecomposed!($Y2, $edM2), $v2)\n",
    "@btime mul!($y2, EigenDecomposedMatrices.exp_eigendecomposed!($Y2, $edM2, $expE2, $tmpY2), $v2);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ba4c485f-0f08-4cb5-be8f-990f48c6c6f7",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "jupytext": {
   "encoding": "# -*- coding: utf-8 -*-",
   "formats": "ipynb,jl:hydrogen"
  },
  "kernelspec": {
   "display_name": "Julia 1.9.0-beta3",
   "language": "julia",
   "name": "julia-1.9"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.9.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}