diff --git a/Lilly1966/.ipynb_checkpoints/Untitled3-checkpoint.ipynb b/Lilly1966/.ipynb_checkpoints/Untitled3-checkpoint.ipynb
new file mode 100644
index 0000000..363fcab
--- /dev/null
+++ b/Lilly1966/.ipynb_checkpoints/Untitled3-checkpoint.ipynb
@@ -0,0 +1,6 @@
+{
+ "cells": [],
+ "metadata": {},
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/Lilly1966/OSE_solve.ipynb b/Lilly1966/OSE_solve.ipynb
new file mode 100644
index 0000000..653c929
--- /dev/null
+++ b/Lilly1966/OSE_solve.ipynb
@@ -0,0 +1,735 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "d7f053a4",
+   "metadata": {},
+   "source": [
+    "# Orr-Sommerfeld 方程式固有値解析 (Chebyshev 多項式展開)\n",
+    "[Orszag (1971)](https://doi.org/10.1017/S0022112071002842) による Orr-Sommerfeld 方程式の Chebyshev 多項式展開による固有値解析\n",
+    "\n",
+    "ただし, 切断波数$N$が大きい場合, 4 階微分項の評価が $O(N^7)$となり, 固有値計算時の丸め誤差が無視できなくなる可能性がある ([Dongarra et al. 1996](https://doi.org/10.1016/S0168-9274(96)00049-9)). そのため, 中間変数$\\chi $を導入して 2 階微分のみで計算する. \n",
+    "\n",
+    "\\begin{equation}\n",
+    "\\dfrac{d^2\\upsilon}{dy^2} -\\alpha ^2\\upsilon \\equiv \\chi , \\label{eq:Orszag1971-AA5-1}\n",
+    "\\end{equation}\n",
+    "となる変数$\\chi $を導入すると, Orr-Sommerfeld 方程式［Orszag (1971) の (2) 式］は\n",
+    "\\begin{equation}\n",
+    "\\dfrac{d^2\\chi }{dy^2} -\\alpha ^2\\chi -i\\alpha R\\left[\\overline{u} \\chi -\\dfrac{d^2\\overline{u}}{dy^2} \\upsilon \\right] +i\\alpha R\\lambda \\chi =0. \\label{eq:Orszag1971-AA5-2}\n",
+    "\\end{equation}\n",
+    "$\\upsilon ,\\; \\chi ,\\; \\overline{u} $についてそれぞれ\n",
+    "\\begin{equation}\n",
+    "\\upsilon =\\sum^{\\infty}_{n=0}{a_{\\upsilon ,n}T_n}, \\quad \\chi =\\sum^{\\infty}_{n=0}{a_{\\chi ,n}T_n} ,\\quad \\overline{u} =\\sum^{\\infty}_{n=0}{b_nT_n} \\label{eq:Orszag1971-AA5-3}\n",
+    "\\end{equation}\n",
+    "という展開を行うと, 支配方程式は <!-- (\\ref{eq:Orszag1971-AA5-1}), (\\ref{eq:Orszag1971-AA5-2}) 式は\n",
+    "\\begin{align}\n",
+    "&\\; \\sum^{\\infty}_{n=0}{\\left\\{a^{(2)}_{\\upsilon ,n}-\\alpha ^2a_{\\upsilon ,n}-a_{\\chi ,n} \\right\\} T_n} =0, \\label{eq:Orszag1971-AA5-4} \\\\\n",
+    "&\\; \\sum^{\\infty}_{n=0}{\\left\\{a^{(2)}_{\\chi ,n}-\\alpha ^2a_{\\chi ,n}-\\dfrac{i\\alpha R}{2} \\left[e_{u,\\chi ,n}-f_{u,\\upsilon ,n} \\right] +i\\alpha R\\lambda a_{\\chi ,n} \\right\\} T_n} =0, \\label{eq:Orszag1971-AA5-5}\n",
+    "\\end{align}\n",
+    "ここで, 非線形項については Orszag (1971) の (A7) 式を用いて\n",
+    "\\begin{equation}\n",
+    "\\begin{split}\n",
+    "e_{u,\\chi ,n}\\equiv &\\; \\dfrac{1}{c_n} \\sum^{\\infty}_{m=-\\infty}{\\overline{a}_{\\chi ,m}\\overline{b}_{n-m}} , \\\\\n",
+    "f_{u,\\upsilon ,n}\\equiv &\\; \\dfrac{1}{c_n} \\sum^{\\infty}_{m=-\\infty}{\\overline{a}_{\\upsilon ,m}\\sum^{\\infty}_{\\substack{p=|n-m|+2\\\\ p\\equiv |n-m|\\; (\\mathrm{mod} \\; 2)}}{p(p^2-(n-m)^2)b_p}} .\n",
+    "\\end{split}\n",
+    " \\label{eq:Orszag1971-AA5-6}\n",
+    "\\end{equation}\n",
+    "(\\ref{eq:Orszag1971-A5}) 式から (\\ref{eq:Orszag1971-AA5-4}), (\\ref{eq:Orszag1971-AA5-5}) 式は\n",
+    "\\begin{align}\n",
+    "&\\; \\sum^{\\infty}_{n=0}{\\left\\{c^{-1}_n\\sum^{\\infty}_{\\substack{p=n+2\\\\ p\\equiv n\\; (\\mathrm{mod} \\; 2)}}{p(p^2-n^2)a_{\\upsilon ,p}} -\\alpha ^2a_{\\upsilon ,n}-a_{\\chi ,n} \\right\\} T_n} =0, \\label{eq:Orszag1971-AA5-7} \\\\\n",
+    "&\\; \\sum^{\\infty}_{n=0}{\\left\\{c^{-1}_n\\sum^{\\infty}_{\\substack{p=n+2\\\\ p\\equiv n\\; (\\mathrm{mod} \\; 2)}}{p(p^2-n^2)a_{\\chi ,p}} -\\alpha ^2a_{\\chi ,n}-\\dfrac{i\\alpha R}{2} \\left[e_{u,\\chi ,n}-f_{u,\\upsilon ,n} \\right] +i\\alpha R\\lambda a_{\\chi ,n} \\right\\} T_n} =0. \\label{eq:Orszag1971-AA5-8}\n",
+    "\\end{align}\n",
+    "よって, $T_n$の係数についての関係式は以下の 2 つになる: -->\n",
+    "$n$次チェビシェフ多項式$T_n$の$0\\leq n\\leq N-2$について, 以下の 2 式が得られる:\n",
+    "\\begin{align}\n",
+    "&\\; \\sum^{\\infty}_{\\substack{p=n+2\\\\ p\\equiv n\\; (\\mathrm{mod} \\; 2)}}{p(p^2-n^2)a_{\\upsilon ,p}} -c_n\\alpha ^2a_{\\upsilon ,n}-c_na_{\\chi ,n} =0, \\label{eq:Orszag1971-AA5-7} \\\\\n",
+    "&\\; \\sum^{\\infty}_{\\substack{p=n+2\\\\ p\\equiv n\\; (\\mathrm{mod} \\; 2)}}{p(p^2-n^2)a_{\\chi ,p}} -c_n\\alpha ^2a_{\\chi ,n} -\\dfrac{i\\alpha R}{2} \\left[\\sum^{\\infty}_{m=-\\infty}{\\overline{a}_{\\chi ,m}\\overline{b}_{n-m}} -\\sum^{\\infty}_{m=-\\infty}{\\overline{a}_{\\upsilon ,m}\\sum^{\\infty}_{\\substack{p=|n-m|+2\\\\ p\\equiv |n-m|\\; (\\mathrm{mod} \\; 2)}}{p(p^2-(n-m)^2)b_p}} \\right] =-i\\alpha R\\lambda c_na_{\\chi ,n} . \\label{eq:Orszag1971-AA5-8}\n",
+    "\\end{align}\n",
+    "境界条件［Orszag (1971) の (15)-(16) 式］から得られる係数の関係式:\n",
+    "$$\n",
+    "\\sum^{N_{\\mathrm{e}}}_{\\substack{p=0\\\\ p=0\\; (\\mathrm{mod}\\; 2)}}{a_{\\upsilon ,p}} =0,\\quad \n",
+    "\\sum^{N_{\\mathrm{e}}}_{\\substack{p=0\\\\ p=0\\; (\\mathrm{mod}\\; 2)}}{p^2a_{\\upsilon ,p}} =0,\\quad \n",
+    "\\sum^{N_{\\mathrm{o}}}_{\\substack{p=1\\\\ p=1\\; (\\mathrm{mod}\\; 2)}}{a_{\\upsilon ,p}} =0,\\quad \n",
+    "\\sum^{N_{\\mathrm{o}}}_{\\substack{p=1\\\\ p=1\\; (\\mathrm{mod}\\; 2)}}{p^2a_{\\upsilon ,p}} =0.\n",
+    "$$\n",
+    "上の 2 つの支配方程式および境界条件から, 求める固有値方程式は\n",
+    "$$A\\textbf{x} =\\lambda B\\textbf{x} ,$$\n",
+    "ただし, 区分行列$A_{1,2}$が支配方程式 1 において, 変数 2 の係数を表す$N-1$次正方行列とすると, \n",
+    "$$A\\equiv \\left(\\begin{array}{cc}\n",
+    "A_{1,\\upsilon} & 0 & 0 & A_{1,\\chi} & 0 & 0 \\\\\n",
+    "\\mathrm{BC1}   &   &   &      0     & 0 & 0 \\\\\n",
+    "\\mathrm{BC2}   &   &   &      0     & 0 & 0 \\\\\n",
+    "A_{2,\\upsilon} & 0 & 0 & A_{2,\\chi} & 0 & 0 \\\\\n",
+    "\\mathrm{BC3}   &   &   &      0     & 0 & 0 \\\\\n",
+    "\\mathrm{BC4}   &   &   &      0     & 0 & 0 \\\\\n",
+    "\\end{array} \\right) ,\\quad B\\equiv \\left(\\begin{array}{cc}\n",
+    "\\textbf{0} & 0 & 0 & \\textbf{0} & 0 & 0 \\\\\n",
+    "0          &   &   &      0     & 0 & 0 \\\\\n",
+    "0          &   &   &      0     & 0 & 0 \\\\\n",
+    "\\textbf{0} & 0 & 0 & B_{2,\\chi} & 0 & 0 \\\\\n",
+    "0          &   &   &      0     & 0 & 0 \\\\\n",
+    "0          &   &   &      0     & 0 & 0 \\\\\n",
+    "\\end{array} \\right) ,\\quad \\textbf{x} \\equiv \\left[a_{\\upsilon ,0},a_{\\upsilon ,1},\\cdots ,a_{\\upsilon ,N},a_{\\chi ,0},a_{\\chi ,1},\\cdots ,a_{\\chi ,N} \\right] ^T.$$"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 33,
+   "id": "7e7d6392",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "Nemax = 200, Nomax = 199\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): 0.02883820761178906 - 7578.242545455161im\n",
+      "Schur(Eig): 0.009405756684565517 - 1314.4376797941554im\n",
+      "Schur(Eig): 0.014244895720953492 - 537.8574465367766im\n",
+      "Schur(Eig): 0.016450157282747912 - 288.51348477771575im\n",
+      "Schur(Eig): 0.02103346120333471 - 180.39666028187585im\n",
+      "Schur(Eig): 0.025848848456281236 - 123.19195139204454im\n",
+      "Schur(Eig): 0.03205497200377114 - 89.63935804312915im\n",
+      "Schur(Eig): 0.038929374231420434 - 68.13201026223997im\n",
+      "Schur(Eig): 0.04694294410795094 - 53.61348178339329im\n",
+      "Schur(Eig): 0.05575204045304491 - 43.30296211449413im\n",
+      "Schur(Eig): 0.06559442227475659 - 35.75067193808814im\n",
+      "Schur(Eig): 0.07627540357947211 - 30.033554949848465im\n",
+      "Schur(Eig): 0.08792850011508294 - 25.61607108907338im\n",
+      "Schur(Eig): 0.10043070899048237 - 22.122849306497507im\n",
+      "Schur(Eig): 0.11386077327612609 - 19.320140971145047im\n",
+      "Schur(Eig): 0.1281349105002331 - 17.03250325944365im\n",
+      "Schur(Eig): 0.1432998698583759 - 15.145068407245297im\n",
+      "Schur(Eig): 0.15929548384870762 - 13.567075557511187im\n",
+      "Schur(Eig): 0.17614847010943072 - 12.236839362139095im\n",
+      "Schur(Eig): 0.19381421003017424 - 11.103631184207575im\n",
+      "Schur(Eig): 0.2123065069883877 - 10.13196095432081im\n",
+      "Schur(Eig): 0.23159223133775889 - 9.291696984526405im\n",
+      "Schur(Eig): 0.25167732522671404 - 8.56126472534193im\n",
+      "Schur(Eig): 0.27253885904938296 - 7.92186980103505im\n",
+      "Schur(Eig): 0.2941796917483501 - 7.359796586150784im\n",
+      "Schur(Eig): 0.3165886521483531 - 6.862846950645903im\n",
+      "Schur(Eig): 0.3397718683313401 - 6.421976647284318im\n",
+      "Schur(Eig): 0.36373573381017626 - 6.029010383316242im\n",
+      "Schur(Eig): 0.38850049383070856 - 5.6778068207963015im\n",
+      "Schur(Eig): 0.41410382288308695 - 5.3627488870422075im\n",
+      "Schur(Eig): 0.4406013358733947 - 5.079574494956819im\n",
+      "Schur(Eig): 0.46809249340979986 - 4.824359741766496im\n",
+      "Schur(Eig): 0.4967146295049181 - 4.594112202357917im\n",
+      "Schur(Eig): 0.5267049464373432 - 4.386087682420243im\n",
+      "Schur(Eig): 0.5584293662940205 - 4.1982417641691345im\n",
+      "Schur(Eig): 0.5926782370272254 - 4.028910498644979im\n",
+      "Schur(Eig): 0.6318873363430142 - 3.878163779238189im\n",
+      "Schur(Eig): 0.6774124487451393 - 3.763500053650113im\n",
+      "Schur(Eig): 0.6724188114836468 - 3.669503887191289im\n",
+      "Schur(Eig): 0.6666346567856843 - 3.5475978950393094im\n",
+      "Schur(Eig): 0.6667369112552747 - 3.429129755491025im\n",
+      "Schur(Eig): 0.6668393733825629 - 3.3134147800032405im\n",
+      "Schur(Eig): 0.666856311072032 - 3.199687877896379im\n",
+      "Schur(Eig): 0.6668698179540508 - 3.0879154146174312im\n",
+      "Schur(Eig): 0.6668850993839182 - 2.9781031043927317im\n",
+      "Schur(Eig): 0.6669018359263257 - 2.870250403683462im\n",
+      "Schur(Eig): 0.6669201648031444 - 2.7643562438989986im\n",
+      "Schur(Eig): 0.6669402701188141 - 2.6604194605606213im\n",
+      "Schur(Eig): 0.6669623585767975 - 2.5584387838874467im\n",
+      "Schur(Eig): 0.6669866659509853 - 2.4584128226652346im\n",
+      "Schur(Eig): 0.6670134605322281 - 2.3603400490888857im\n",
+      "Schur(Eig): 0.667043048541936 - 2.2642187830663874im\n",
+      "Schur(Eig): 0.6670757804235424 - 2.1700471755920914im\n",
+      "Schur(Eig): 0.6671120580366842 - 2.0778231866826116im\n",
+      "Schur(Eig): 0.6671523431785931 - 1.9875445633351365im\n",
+      "Schur(Eig): 0.6671971677312326 - 1.8992088120201795im\n",
+      "Schur(Eig): 0.6672471464171399 - 1.8128131702983779im\n",
+      "Schur(Eig): 0.6673029901943451 - 1.7283545703860157im\n",
+      "Schur(Eig): 0.667365524784893 - 1.6458296030286281im\n",
+      "Schur(Eig): 0.6674357106791179 - 1.5652344670759777im\n",
+      "Schur(Eig): 0.6675146678396625 - 1.486564926666811im\n",
+      "Schur(Eig): 0.23752648881857785 + 0.0037396706222145603im\n",
+      "Schur(Eig): 0.6676037066353783 - 1.409816245299302im\n",
+      "Schur(Eig): 0.6677043627182395 - 1.3349831229686595im\n",
+      "Schur(Eig): 0.19005924936950183 - 0.1828219254116369im\n",
+      "Schur(Eig): 0.667818440284551 - 1.2620596189975923im\n",
+      "Schur(Eig): 0.667948064841178 - 1.1910390667367665im\n",
+      "Schur(Eig): 0.6680957458771399 - 1.1219139721740363im\n",
+      "Schur(Eig): 0.9646425098876398 - 0.0351865839248333im\n",
+      "Schur(Eig): 0.34910682012256655 - 0.12450197754871109im\n",
+      "Schur(Eig): 0.6682644469478218 - 1.0546759064547755im\n",
+      "Schur(Eig): 0.6684576869017671 - 0.989315378367853im\n",
+      "Schur(Eig): 0.9363517811244114 - 0.06325156893621248im\n",
+      "Schur(Eig): 0.36849847834757754 - 0.23882483171782218im\n",
+      "Schur(Eig): 0.6686796236220521 - 0.9258216904581784im\n",
+      "Schur(Eig): 0.9080563344935584 - 0.09131286178046318im\n",
+      "Schur(Eig): 0.668935211379344 - 0.8641827952264582im\n",
+      "Schur(Eig): 0.4749011869127081 - 0.20873122005003544im\n",
+      "Schur(Eig): 0.8797556958156186 - 0.11937073101026466im\n",
+      "Schur(Eig): 0.6692302988188383 - 0.8043850926397907im\n",
+      "Schur(Eig): 0.8514493818764631 - 0.14742560074344907im\n",
+      "Schur(Eig): 0.5129162046302708 - 0.28662504241366654im\n",
+      "Schur(Eig): 0.8231369612514765 - 0.17547807357377554im\n",
+      "Schur(Eig): 0.7948183878980362 - 0.20352914421365847im\n",
+      "Schur(Eig): 0.5872129303355416 - 0.26716171254068716im\n",
+      "Schur(Eig): 0.7664940772247112 - 0.23158507171718073im\n",
+      "Schur(Eig): 0.6695718470141386 - 0.7464132968048778im\n",
+      "Schur(Eig): 0.6699680740685054 - 0.6902501370427058im\n",
+      "Schur(Eig): 0.7381150003547364 - 0.25969189546075594im\n",
+      "Schur(Eig): 0.6704286339966932 - 0.6358762450837643im\n",
+      "Schur(Eig): 0.7088747432491715 - 0.28765535594742253im\n",
+      "Schur(Eig): 0.6361048559307513 - 0.3251971610110621im\n",
+      "Schur(Eig): 0.6828633658739772 - 0.3076116006474047im\n",
+      "Schur(Eig): 0.6709649733030304 - 0.5832698557777029im\n",
+      "Schur(Eig): 0.6715904009498063 - 0.5324065147651961im\n",
+      "Schur(Eig): 0.6776436663762307 - 0.3437328375681226im\n",
+      "Schur(Eig): 0.6723226907333578 - 0.48325970758858783im\n",
+      "Schur(Eig): 0.6745109420659848 - 0.3898261303001523im\n",
+      "Schur(Eig): 0.6732077446999831 - 0.43579516755895803im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n"
+     ]
+    },
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): NaN + NaN*im\n",
+      "Schur(Eig): 0.028851580470429898 - 7427.412242041225im\n",
+      "Schur(Eig): 0.009533588970052825 - 1288.2787233690856im\n",
+      "Schur(Eig): 0.014396958205195097 - 527.1566441913569im\n",
+      "Schur(Eig): 0.016628580243660488 - 282.7769998739129im\n",
+      "Schur(Eig): 0.02125292187361028 - 176.81237598404454im\n",
+      "Schur(Eig): 0.026119154949576815 - 120.74641332223527im\n",
+      "Schur(Eig): 0.032385535531883815 - 87.86174245066204im\n",
+      "Schur(Eig): 0.039329601038889495 - 66.7825557976176im\n",
+      "Schur(Eig): 0.04742195916392423 - 52.553069738534845im\n",
+      "Schur(Eig): 0.05631896904887011 - 42.44783633451269im\n",
+      "Schur(Eig): 0.06625806836511919 - 35.04592966338191im\n",
+      "Schur(Eig): 0.0770445496186914 - 29.442672481245136im\n",
+      "Schur(Eig): 0.08881157642882291 - 25.113181360813634im\n",
+      "Schur(Eig): 0.1014361042011465 - 21.689560561534627im\n",
+      "Schur(Eig): 0.11499649821234285 - 18.942711695841773im\n",
+      "Schur(Eig): 0.1294089198181688 - 16.700686600350803im\n",
+      "Schur(Eig): 0.14471971920455637 - 14.850901929562848im\n",
+      "Schur(Eig): 0.16086866736631758 - 13.304405276029732im\n",
+      "Schur(Eig): 0.17788206988671296 - 12.000737477058372im\n",
+      "Schur(Eig): 0.195715259165424 - 10.890181166870757im\n",
+      "Schur(Eig): 0.21438163512998343 - 9.937951376838607im\n",
+      "Schur(Eig): 0.23384807257382156 - 9.114517601005865im\n",
+      "Schur(Eig): 0.2541201752826265 - 8.39873379878657im\n",
+      "Schur(Eig): 0.2751751671241231 - 7.772181247029363im\n",
+      "Schur(Eig): 0.29701576942380037 - 7.221416953352625im\n",
+      "Schur(Eig): 0.3196313288656588 - 6.73448758325907im\n",
+      "Schur(Eig): 0.3430283427264351 - 6.302529352368175im\n",
+      "Schur(Eig): 0.3672145572259243 - 5.917532503421068im\n",
+      "Schur(Eig): 0.39221176530007096 - 5.57347953513074im\n",
+      "Schur(Eig): 0.4180608555520649 - 5.264869231416005im\n",
+      "Schur(Eig): 0.4448216679242136 - 4.9875274895461im\n",
+      "Schur(Eig): 0.4726011259137812 - 4.737614362022132im\n",
+      "Schur(Eig): 0.501547449209317 - 4.5122027838334775im\n",
+      "Schur(Eig): 0.531918178839164 - 4.308614522888535im\n",
+      "Schur(Eig): 0.5641228380998303 - 4.124868378717572im\n",
+      "Schur(Eig): 0.5991030757694588 - 3.9594161375064765im\n",
+      "Schur(Eig): 0.640132839257879 - 3.8130758046841033im\n",
+      "Schur(Eig): 0.6818692445843058 - 3.70965326550958im\n",
+      "Schur(Eig): 0.6699810259372195 - 3.609524425068888im\n",
+      "Schur(Eig): 0.666593083128247 - 3.4881139930364884im\n",
+      "Schur(Eig): 0.6668188394688124 - 3.3709146755875516im\n",
+      "Schur(Eig): 0.666893879320586 - 3.2562057758916922im\n",
+      "Schur(Eig): 0.6669098098466136 - 3.143458355138187im\n",
+      "Schur(Eig): 0.6669255853951719 - 3.032666881924751im\n",
+      "Schur(Eig): 0.6669431586289576 - 2.9238354325428593im\n",
+      "Schur(Eig): 0.66696231652267 - 2.81696315510865im\n",
+      "Schur(Eig): 0.6669832230117946 - 2.71204894108598im\n",
+      "Schur(Eig): 0.6670060760038876 - 2.6090915843046405im\n",
+      "Schur(Eig): 0.6670310976186595 - 2.508089766507401im\n",
+      "Schur(Eig): 0.6670585401646524 - 2.4090420414631817im\n",
+      "Schur(Eig): 0.6670886906235786 - 2.31194682097146im\n",
+      "Schur(Eig): 0.6671218762052807 - 2.216802357546517im\n",
+      "Schur(Eig): 0.6671584705378463 - 2.1236067257038025im\n",
+      "Schur(Eig): 0.6671989015460608 - 2.032357799685315im\n",
+      "Schur(Eig): 0.6672436603264309 - 1.9430532304452381im\n",
+      "Schur(Eig): 0.6672933118803402 - 1.8556904157206815im\n",
+      "Schur(Eig): 0.6673485079369096 - 1.7702664692964718im\n",
+      "Schur(Eig): 0.667410001879772 - 1.686778184264008im\n",
+      "Schur(Eig): 0.6674786675328377 - 1.6052219899119025im\n",
+      "Schur(Eig): 0.6675555199095552 - 1.5255939061026371im\n",
+      "Schur(Eig): 0.6676417426648599 - 1.4478894860899223im\n",
+      "Schur(Eig): 0.6677387180326144 - 1.3721037552372932im\n",
+      "Schur(Eig): 0.27720434380851733 - 0.05089872725658317im\n",
+      "Schur(Eig): 0.6678480672909648 - 1.2982311373576947im\n",
+      "Schur(Eig): 0.21272578235378325 - 0.19936069475434753im\n",
+      "Schur(Eig): 0.6679716919722942 - 1.2262653745565242im\n",
+      "Schur(Eig): 0.6681118358356443 - 1.1561994303389072im\n",
+      "Schur(Eig): 0.6682711424754932 - 1.0880253863300817im\n",
+      "Schur(Eig): 0.9646309154505646 - 0.03516727763104578im\n",
+      "Schur(Eig): 0.6684527468496686 - 1.0217343142581585im\n",
+      "Schur(Eig): 0.4163510155763265 - 0.1382265253000814im\n",
+      "Schur(Eig): 0.6686603562446971 - 0.9573161471892272im\n",
+      "Schur(Eig): 0.9363165358812776 - 0.06320149583980876im\n",
+      "Schur(Eig): 0.3839876109074957 - 0.26510649960903276im\n",
+      "Schur(Eig): 0.9079830546296593 - 0.09122273543365653im\n",
+      "Schur(Eig): 0.6688983789719299 - 0.8947595190105614im\n",
+      "Schur(Eig): 0.8796272922085633 - 0.11923285262208567im\n",
+      "Schur(Eig): 0.6691720559744425 - 0.834051590850208im\n",
+      "Schur(Eig): 0.5320452088016702 - 0.2064652191117368im\n",
+      "Schur(Eig): 0.8512458401120699 - 0.1472339290642039im\n",
+      "Schur(Eig): 0.6694875891020631 - 0.7751778696700894im\n",
+      "Schur(Eig): 0.8228350407533278 - 0.1752286787491913im\n",
+      "Schur(Eig): 0.5222281555545437 - 0.3142970585337267im\n",
+      "Schur(Eig): 0.7943883846246975 - 0.20322066451492837im\n",
+      "Schur(Eig): 0.6698523995623807 - 0.7181220051500524im\n",
+      "Schur(Eig): 0.7658768124412472 - 0.2311859878020367im\n",
+      "Schur(Eig): 0.6367193825465925 - 0.25988571005387817im\n",
+      "Schur(Eig): 0.7374157502565554 - 0.2587170852291249im\n",
+      "Schur(Eig): 0.6702751958863984 - 0.6628655095746344im\n",
+      "Schur(Eig): 0.7123159247889363 - 0.2855146556769195im\n",
+      "Schur(Eig): 0.6707663761310181 - 0.6093876651212601im\n",
+      "Schur(Eig): 0.6713381853349581 - 0.5576650860648068im\n",
+      "Schur(Eig): 0.6941568719197801 - 0.31949859450281015im\n",
+      "Schur(Eig): 0.6424176412961435 - 0.35254526646010687im\n",
+      "Schur(Eig): 0.6720075259566941 - 0.5076720092943916im\n",
+      "Schur(Eig): 0.6784240677887209 - 0.36421754924582495im\n",
+      "Schur(Eig): 0.6728223857272042 - 0.4593698713921976im\n",
+      "Schur(Eig): 0.6740837771202557 - 0.41247342802481496im\n",
+      "Schur(Eig): NaN + NaN*im\n"
+     ]
+    }
+   ],
+   "source": [
+    "using OffsetArrays  # 配列要素番号を任意開始にするモジュール (行列計算時は 1 始まりの通常 matrix に copy する)\n",
+    "using LinearAlgebra\n",
+    "\n",
+    "#-- Drawing parameters    \n",
+    "N = 200  # wavenumber in Chebyshev\n",
+    "Ncal = N + 1\n",
+    "Ntot = 2 * (N + 1)\n",
+    "\n",
+    "# initialize working arrays\n",
+    "# A: non-lambda coefficients, B: lambda coefficients\n",
+    "A_J = reshape(zeros(Complex{Float64},Ntot,Ntot),Ntot,Ntot)  # Matrix for calculation\n",
+    "B_J = reshape(zeros(Complex{Float64},Ntot,Ntot),Ntot,Ntot)  # Matrix for calculation\n",
+    "### --- From here, supposing matrix element as row, col as in mathematical notation\n",
+    "Acal = reshape(zeros(Complex{Float64},Ntot,Ntot),Ntot,Ntot)  # Matrix for calculation\n",
+    "Bcal = reshape(zeros(Complex{Float64},Ntot,Ntot),Ntot,Ntot)  # Matrix for calculation\n",
+    "A_vv = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for v in v Chebyshev\n",
+    "A_vchi = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for chi in v Chebyshev\n",
+    "A_chiv = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for chi in chi Chebyshev\n",
+    "A_chichi = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for chi in chi Chebyshev\n",
+    "B_chichi = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for chi in chi Chebyshev\n",
+    "### --- From here, supposing matrix element as col, row as in julia notation (only using these arrays for matrix calc)\n",
+    "b = OffsetArray(zeros(Float64,Ncal),0:N)  # A vector composed of Chebyshev coefficients for the basic state flow\n",
+    "c = OffsetArray(zeros(Float64,Ncal),0:N)  # factorized vector (c_0=2, c_n=1, (n>0))\n",
+    "\n",
+    "# setting parameters\n",
+    "α = 1.0\n",
+    "Re = 10000.0\n",
+    "\n",
+    "# setting the Chebyshev coefficients for the basic flow \n",
+    "b[0] = 0.5  # For the Poiseuille flow (u = 1-y^2)\n",
+    "b[2] = -b[0]\n",
+    "\n",
+    "# setting the factorized vector\n",
+    "c[0] = 2.0\n",
+    "c[1:N] .= 1.0\n",
+    "\n",
+    "# setting maximum even and odd numbers from N\n",
+    "Nemax = 0\n",
+    "Nomax = 1\n",
+    "for i in 0:2:N\n",
+    "    Nemax = i\n",
+    "end\n",
+    "for i in 1:2:N\n",
+    "    Nomax = i\n",
+    "end\n",
+    "println(\"Nemax = \", Nemax, \", Nomax = \", Nomax)\n",
+    "\n",
+    "# setting matrices\n",
+    "# 支配方程式から構築される部分 (チェビシェフ次数 0<=n<=N-2 の範囲)\n",
+    "# この段階では境界条件で成り立つ係数の関係は考慮しない\n",
+    "for n in 0:N-2\n",
+    "    n2 = n^2\n",
+    "    \n",
+    "    # In v rows, a_{v,n} and a_{chi,n}\n",
+    "    for p in n+2:2:N\n",
+    "        p2 = p^2\n",
+    "        A_vv[n,p] = p * (p2 - n2)\n",
+    "    end\n",
+    "    \n",
+    "    A_vv[n,n] = - c[n] * α^2\n",
+    "    A_vchi[n,n] = - c[n]\n",
+    "    \n",
+    "    # In chi rows, a_{v,n} and a_{chi,n}\n",
+    "    for p in n+2:2:N\n",
+    "        p2 = p^2\n",
+    "        A_chichi[n,p] = p * (p2 - n2)\n",
+    "    end\n",
+    "    \n",
+    "    A_chichi[n,n] = - c[n] * α^2\n",
+    "\n",
+    "    for p in -N:N\n",
+    "        abs_n_p = abs(n - p)\n",
+    "        abs_p = abs(p)\n",
+    "        if abs_n_p <= N\n",
+    "            A_chichi[n,abs_p] = A_chichi[n,abs_p] + complex(0.0, - 0.5 * α * Re * c[abs_p] * c[abs_n_p] * b[abs_n_p])\n",
+    "        end\n",
+    "\n",
+    "        if abs_n_p+2 <= N\n",
+    "            acoef = 0.0\n",
+    "            for m in abs_n_p+2:2:N\n",
+    "                m2 = m^2\n",
+    "                acoef = acoef + m * (m2 - (n - p)^2) * b[m]\n",
+    "            end\n",
+    "            A_chiv[n,abs_p] = A_chiv[n,abs_p] + complex(0.0, 0.5 * α * Re * c[abs_p] * acoef)\n",
+    "        end\n",
+    "    end\n",
+    "    \n",
+    "    # In chi rows, lambda coefficients\n",
+    "    B_chichi[n,n] = complex(0.0, - α * Re * c[n])\n",
+    "\n",
+    "end\n",
+    "\n",
+    "# Set boundary conditions\n",
+    "for p in 0:2:N  # even mode\n",
+    "    A_vv[N-1,p] = 1.0\n",
+    "#    A_chichi[N-1,p] = 1.0\n",
+    "    A_chiv[N-1,p] = p^2\n",
+    "end\n",
+    "for p in 1:2:N  # even mode\n",
+    "    A_vv[N,p] = 1.0\n",
+    "#    A_chichi[N,p] = 1.0\n",
+    "    A_chiv[N,p] = p^2\n",
+    "end\n",
+    "\n",
+    "# Gather A_v, A_chi, B_chi -> Acal, Bcal\n",
+    "Acal[1:Ncal,1:Ncal] .= A_vv[0:N,0:N]\n",
+    "Acal[1:Ncal,Ncal+1:Ntot] .= A_vchi[0:N,0:N]\n",
+    "Acal[Ncal+1:Ntot,1:Ncal] .= A_chiv[0:N,0:N]\n",
+    "Acal[Ncal+1:Ntot,Ncal+1:Ntot] .= A_chichi[0:N,0:N]\n",
+    "Bcal[Ncal+1:Ntot,Ncal+1:Ntot] .= B_chichi[0:N,0:N]\n",
+    "\n",
+    "# Convert Transpose type to Matrix (Array) type (ドット演算して要素だけコピーしている)\n",
+    "A_J .= transpose(Acal)  # A_J = Transpose type\n",
+    "B_J .= transpose(Bcal)  # B_J = Transpose type\n",
+    "\n",
+    "F = schur(A_J,B_J)\n",
+    "Dcal = F.α ./ F.β\n",
+    "for i in 1:Ntot\n",
+    "    println(\"Schur(Eig): \", Dcal[i])\n",
+    "end\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 35,
+   "id": "1255c473",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/latex": [
+       "The analogue of IPython's \\texttt{\\%matplotlib} in Julia is to use the \\href{https://github.com/stevengj/PyPlot.jl}{PyPlot package}, which gives a Julia interface to Matplotlib including inline plots in IJulia notebooks.   (The equivalent of \\texttt{numpy} is already loaded by default in Julia.)\n",
+       "\n",
+       "Given PyPlot, the analogue of \\texttt{\\%matplotlib inline} is \\texttt{using PyPlot}, since PyPlot defaults to inline plots in IJulia.\n",
+       "\n",
+       "To enable separate GUI windows in PyPlot, analogous to \\texttt{\\%matplotlib}, do \\texttt{using PyPlot; pygui(true)}.   To specify a particular gui backend, analogous to \\texttt{\\%matplotlib gui}, you can either do \\texttt{using PyPlot; pygui(:gui); using PyPlot; pygui(true)} (where \\texttt{gui} is \\texttt{wx}, \\texttt{qt}, \\texttt{tk}, or \\texttt{gtk}), or you can do \\texttt{ENV[\"MPLBACKEND\"]=backend; using PyPlot; pygui(true)} (where \\texttt{backend} is the name of a Matplotlib backend, like \\texttt{tkagg}).\n",
+       "\n",
+       "For more options, see the PyPlot documentation.\n",
+       "\n"
+      ],
+      "text/markdown": [
+       "The analogue of IPython's `%matplotlib` in Julia is to use the [PyPlot package](https://github.com/stevengj/PyPlot.jl), which gives a Julia interface to Matplotlib including inline plots in IJulia notebooks.   (The equivalent of `numpy` is already loaded by default in Julia.)\n",
+       "\n",
+       "Given PyPlot, the analogue of `%matplotlib inline` is `using PyPlot`, since PyPlot defaults to inline plots in IJulia.\n",
+       "\n",
+       "To enable separate GUI windows in PyPlot, analogous to `%matplotlib`, do `using PyPlot; pygui(true)`.   To specify a particular gui backend, analogous to `%matplotlib gui`, you can either do `using PyPlot; pygui(:gui); using PyPlot; pygui(true)` (where `gui` is `wx`, `qt`, `tk`, or `gtk`), or you can do `ENV[\"MPLBACKEND\"]=backend; using PyPlot; pygui(true)` (where `backend` is the name of a Matplotlib backend, like `tkagg`).\n",
+       "\n",
+       "For more options, see the PyPlot documentation.\n"
+      ],
+      "text/plain": [
+       "  The analogue of IPython's \u001b[36m%matplotlib\u001b[39m in Julia is to use the PyPlot package\n",
+       "  (https://github.com/stevengj/PyPlot.jl), which gives a Julia interface to\n",
+       "  Matplotlib including inline plots in IJulia notebooks. (The equivalent of\n",
+       "  \u001b[36mnumpy\u001b[39m is already loaded by default in Julia.)\n",
+       "\n",
+       "  Given PyPlot, the analogue of \u001b[36m%matplotlib inline\u001b[39m is \u001b[36musing PyPlot\u001b[39m, since\n",
+       "  PyPlot defaults to inline plots in IJulia.\n",
+       "\n",
+       "  To enable separate GUI windows in PyPlot, analogous to \u001b[36m%matplotlib\u001b[39m, do \u001b[36musing\n",
+       "  PyPlot; pygui(true)\u001b[39m. To specify a particular gui backend, analogous to\n",
+       "  \u001b[36m%matplotlib gui\u001b[39m, you can either do \u001b[36musing PyPlot; pygui(:gui); using PyPlot;\n",
+       "  pygui(true)\u001b[39m (where \u001b[36mgui\u001b[39m is \u001b[36mwx\u001b[39m, \u001b[36mqt\u001b[39m, \u001b[36mtk\u001b[39m, or \u001b[36mgtk\u001b[39m), or you can do\n",
+       "  \u001b[36mENV[\"MPLBACKEND\"]=backend; using PyPlot; pygui(true)\u001b[39m (where \u001b[36mbackend\u001b[39m is the\n",
+       "  name of a Matplotlib backend, like \u001b[36mtkagg\u001b[39m).\n",
+       "\n",
+       "  For more options, see the PyPlot documentation."
+      ]
+     },
+     "metadata": {},
+     "output_type": "display_data"
+    }
+   ],
+   "source": [
+    "#-- Drawing\n",
+    "##########\n",
+    "#  Plot  #\n",
+    "##########\n",
+    "\n",
+    "using PyPlot\n",
+    "\n",
+    "rc(\"font\", family=\"IPAPGothic\")\n",
+    "fig = figure(\"pyplot_majorminor\",figsize=(7,5))\n",
+    "\n",
+    "x = real.(Dcal[1:Ntot])\n",
+    "y = imag.(Dcal[1:Ntot])\n",
+    "\n",
+    "#ax = gca()\n",
+    "\n",
+    "#cp = scatter(x,y)\n",
+    "scatter(x,y)\n",
+    "xlim(0.0,1.0)\n",
+    "ylim(-1.0,0.1)\n",
+    "\n",
+    "#ax = gca()\n",
+    "\n",
+    "xlabel(\"Real(λ)\")\n",
+    "ylabel(\"Imag(λ)\")\n",
+    "#if i == 1 #div(ne,2)-1\n",
+    "#   plt.colorbar(cp)\n",
+    "#end\n",
+    "grid(\"on\")\n",
+    "\n",
+    "fig.canvas.draw()\n",
+    "gcf() # Needed for IJulia to plot inline\n",
+    "\n",
+    "ttext = \"OSE spectra: Re = $(Re), α = $(α)\"\n",
+    "title(ttext)\n",
+    "\n",
+    "fname = \"test.png\"\n",
+    "savefig(fname)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "026799fc",
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Julia 1.6.1",
+   "language": "julia",
+   "name": "julia-1.6"
+  },
+  "language_info": {
+   "file_extension": ".jl",
+   "mimetype": "application/julia",
+   "name": "julia",
+   "version": "1.6.1"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/Lilly1966/OSE_solve_Don96.ipynb b/Lilly1966/OSE_solve_Don96.ipynb
new file mode 120000
index 0000000..8978b86
--- /dev/null
+++ b/Lilly1966/OSE_solve_Don96.ipynb
@@ -0,0 +1 @@
+OSE_solve.ipynb
\ No newline at end of file
diff --git a/Lilly1966/Untitled1.ipynb b/Lilly1966/OSE_solve_Orszag1971.ipynb
similarity index 100%
rename from Lilly1966/Untitled1.ipynb
rename to Lilly1966/OSE_solve_Orszag1971.ipynb
diff --git a/Lilly1966/Untitled.ipynb b/Lilly1966/OSE_solve_finite.ipynb
similarity index 100%
rename from Lilly1966/Untitled.ipynb
rename to Lilly1966/OSE_solve_finite.ipynb
diff --git a/Lilly1966/Untitled2.ipynb b/Lilly1966/Untitled2.ipynb
deleted file mode 100644
index 1fad5a2..0000000
--- a/Lilly1966/Untitled2.ipynb
+++ /dev/null
@@ -1,659 +0,0 @@
-{
- "cells": [
-  {
-   "cell_type": "markdown",
-   "id": "d7f053a4",
-   "metadata": {},
-   "source": [
-    "# Dongarra et al. (1996) D2 method\n",
-    "$n$次チェビシェフ多項式$T_n$の$0\\leq n\\leq N-4$について, \n",
-    "\\begin{equation}\n",
-    "\\begin{split}\n",
-    "&\\dfrac{1}{24} \\sum^{N}_{\\substack{p=n+4\\\\ p\\equiv n\\; (\\mathrm{mod} \\; 2)}}{p\\left[p^2(p^2-4)^2-3n^2p^4+3n^4p^2-n^2\\left(n^2-4 \\right) ^2-48\\alpha ^2(p^2-n^2) \\right] a_p} -8\\alpha ^2(n+2)(n+1)a_{n+2} +\\alpha ^4c_na_n \\\\\n",
-    "&\\quad -\\dfrac{i\\alpha R}{2} \\left[\\sum^{N}_{p=2}{a_p\\sum^{}_{\\substack{m\\equiv p\\; (\\mathrm{mod} \\; 2) \\\\ |m|\\leq p-2 \\\\ |n-m|\\leq N}}{p(p^2-m^2)c_{|n-m|}b_{|n-m|} }} -\\alpha ^2\\sum^{}_{\\substack{|p|\\leq N \\\\ |n-p|\\leq N}}{c_{|p|}a_{|p|}c_{|n-p|}b_{|n-p|}} -\\sum^{}_{\\substack{|p|\\leq N \\\\ |n-p|\\leq N}}{c_{|p|}a_{|p|}\\sum^{N}_{\\substack{m=|n-p|+2\\\\ m+n\\equiv p\\; (\\mathrm{mod} \\; 2)}}{m\\left[m^2-(n-p)^2 \\right] b_m}} \\right] \\\\\n",
-    "&\\quad +\\lambda \\left\\{i\\alpha R\\sum^{N}_{\\substack{p=n+2\\\\ p\\equiv n\\; (\\mathrm{mod} \\; 2)}}{p(p^2-n^2)a_p} -i\\alpha ^3Rc_na_n \\right\\} =0.\n",
-    "\\end{split}\n",
-    " \\label{eq:Orszag1971-17-A44}\n",
-    "\\end{equation}\n",
-    "および境界条件から得られる係数の関係式:\n",
-    "$$\n",
-    "\\sum^{N_{\\mathrm{e}}}_{\\substack{p=0\\\\ p=0\\; (\\mathrm{mod}\\; 2)}}{a_{\\upsilon ,p}} =0,\\quad \n",
-    "\\sum^{N_{\\mathrm{e}}}_{\\substack{p=0\\\\ p=0\\; (\\mathrm{mod}\\; 2)}}{a_{\\chi ,p}} =0,\\quad \n",
-    "\\sum^{N_{\\mathrm{o}}}_{\\substack{p=1\\\\ p=1\\; (\\mathrm{mod}\\; 2)}}{a_{\\upsilon ,p}} =0,\\quad \n",
-    "\\sum^{N_{\\mathrm{o}}}_{\\substack{p=1\\\\ p=1\\; (\\mathrm{mod}\\; 2)}}{a_{\\chi ,p}} =0.\n",
-    "$$\n",
-    "求める固有値方程式は\n",
-    "$$\n",
-    "(A-B\\lambda )\\textbf{x} =\\textbf{0} \\quad \\rightarrow \\quad A\\textbf{x} =\\lambda B\\textbf{x} ,\n",
-    "$$\n",
-    "$$\n",
-    "\\textbf{x} \\equiv \\left[a_{\\upsilon ,0},a_{\\upsilon ,1},\\cdots ,a_{\\upsilon ,N},\\cdots a_{\\chi ,0},a_{\\chi ,1},\\cdots ,a_{\\chi ,N} \\right] ^T.\n",
-    "$$"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 7,
-   "id": "7e7d6392",
-   "metadata": {},
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "Nemax = 200, Nomax = 199\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): 0.009959653231081049 - 48468.16070542684im\n",
-      "Schur(Eig): 0.010639510751611096 - 2607.8212490697633im\n",
-      "Schur(Eig): 0.012351539288276479 - 800.590046080902im\n",
-      "Schur(Eig): 0.015015902879698298 - 384.45342326340136im\n",
-      "Schur(Eig): 0.01865670668560613 - 225.024165049837im\n",
-      "Schur(Eig): 0.023254868059602536 - 147.71669819294004im\n",
-      "Schur(Eig): 0.02881645974978412 - 104.4202105849741im\n",
-      "Schur(Eig): 0.035328284168819016 - 77.7782209386429im\n",
-      "Schur(Eig): 0.042790145135120326 - 60.219231986395144im\n",
-      "Schur(Eig): 0.051190030989770624 - 48.043442857050216im\n",
-      "Schur(Eig): 0.06052407372056047 - 39.25389119196648im\n",
-      "Schur(Eig): 0.07078015331861712 - 32.70373530270634im\n",
-      "Schur(Eig): 0.08195179841863179 - 27.69170216772187im\n",
-      "Schur(Eig): 0.09402636094641309 - 23.772156870806974im\n",
-      "Schur(Eig): 0.10699539386551074 - 20.649194542724224im\n",
-      "Schur(Eig): 0.12084564762227566 - 18.121217196096186im\n",
-      "Schur(Eig): 0.13556718775621068 - 16.046282019972388im\n",
-      "Schur(Eig): 0.15114629288579406 - 14.322631989179538im\n",
-      "Schur(Eig): 0.16757203548479477 - 12.875408871277605im\n",
-      "Schur(Eig): 0.1848305610327682 - 11.648819539828576im\n",
-      "Schur(Eig): 0.20291060582378456 - 10.600397746194227im\n",
-      "Schur(Eig): 0.22179887197318027 - 9.69750735097785im\n",
-      "Schur(Eig): 0.24148485401338932 - 8.914624066981391im\n",
-      "Schur(Eig): 0.26195720507049675 - 8.231641208356391im\n",
-      "Schur(Eig): 0.2832082482056146 - 7.632483407168671im\n",
-      "Schur(Eig): 0.30523142737800846 - 7.104229485402162im\n",
-      "Schur(Eig): 0.32802600967215895 - 6.636361878241808im\n",
-      "Schur(Eig): 0.35159602525827593 - 6.220287241885267im\n",
-      "Schur(Eig): 0.37595602953432217 - 5.848909822285037im\n",
-      "Schur(Eig): 0.4011325641288547 - 5.516357599505128im\n",
-      "Schur(Eig): 0.427172760708608 - 5.217730231087406im\n",
-      "Schur(Eig): 0.4541508730721548 - 4.948937189436128im\n",
-      "Schur(Eig): 0.4821841601989238 - 4.706545317401889im\n",
-      "Schur(Eig): 0.5114529456517629 - 4.487684537813415im\n",
-      "Schur(Eig): 0.5422467284169451 - 4.289969833993502im\n",
-      "Schur(Eig): 0.5750855923683216 - 4.111510215453719im\n",
-      "Schur(Eig): 0.6112554918766318 - 3.9512059356024944im\n",
-      "Schur(Eig): 0.6556421923364719 - 3.8130880986298696im\n",
-      "Schur(Eig): 0.6811406872521558 - 3.7220401147493427im\n",
-      "Schur(Eig): 0.6676821488923474 - 3.608660121218405im\n",
-      "Schur(Eig): 0.6665804098495649 - 3.4880416709301314im\n",
-      "Schur(Eig): 0.6667802308133557 - 3.3711532333488363im\n",
-      "Schur(Eig): 0.6668179192758314 - 3.2564494050315274im\n",
-      "Schur(Eig): 0.666830012605901 - 3.1436975186061074im\n",
-      "Schur(Eig): 0.6668435606458105 - 3.032904279152159im\n",
-      "Schur(Eig): 0.6668585958188162 - 2.9240710578166387im\n",
-      "Schur(Eig): 0.6668751189438638 - 2.817196855044981im\n",
-      "Schur(Eig): 0.6668933011487165 - 2.712280547529687im\n",
-      "Schur(Eig): 0.6669133395517908 - 2.609320914704567im\n",
-      "Schur(Eig): 0.6669354582660756 - 2.508316621120469im\n",
-      "Schur(Eig): 0.666959912627878 - 2.409266201444541im\n",
-      "Schur(Eig): 0.6669869942827231 - 2.312168045837876im\n",
-      "Schur(Eig): 0.6670170371782566 - 2.2170203843601035im\n",
-      "Schur(Eig): 0.6670504243268697 - 2.1238212663919858im\n",
-      "Schur(Eig): 0.6670875960353477 - 2.032568540289848im\n",
-      "Schur(Eig): 0.6671290596021298 - 1.9432598271678956im\n",
-      "Schur(Eig): 0.6671754007569 - 1.8558924955316818im\n",
-      "Schur(Eig): 0.6672272971485756 - 1.770463627219107im\n",
-      "Schur(Eig): 0.14806292127708567 - 0.05163618759434957im\n",
-      "Schur(Eig): 0.6672855350586689 - 1.686969982149579im\n",
-      "Schur(Eig): 0.18943271837978434 - 0.060310465401177626im\n",
-      "Schur(Eig): 0.6673510275664712 - 1.6054079576922191im\n",
-      "Schur(Eig): 0.6674248391707343 - 1.5257735409670852im\n",
-      "Schur(Eig): 0.6675082115223181 - 1.4480622554169715im\n",
-      "Schur(Eig): 0.6676025974735393 - 1.3722691003535645im\n",
-      "Schur(Eig): 0.2706703557950911 - 0.1377834285417466im\n",
-      "Schur(Eig): 0.6677096995536513 - 1.2983884808041604im\n",
-      "Schur(Eig): 0.667831515915676 - 1.2264141288970472im\n",
-      "Schur(Eig): 0.6679703987350285 - 1.1563390138241891im\n",
-      "Schur(Eig): 0.9646425102039556 - 0.03518658306579477im\n",
-      "Schur(Eig): 0.35180487825589996 - 0.1824922602726231im\n",
-      "Schur(Eig): 0.9363517819722577 - 0.06325156250092638im\n",
-      "Schur(Eig): 0.6681291149119583 - 1.0881552375169095im\n",
-      "Schur(Eig): 0.6683109300410721 - 1.0218539239217141im\n",
-      "Schur(Eig): 0.9080563347735641 - 0.09131283496500782im\n",
-      "Schur(Eig): 0.4244641747556518 - 0.2181592421958267im\n",
-      "Schur(Eig): 0.8797556878457122 - 0.11937065543353942im\n",
-      "Schur(Eig): 0.6685196991471754 - 0.95742508243062im\n",
-      "Schur(Eig): 0.851449342351972 - 0.14742543062961255im\n",
-      "Schur(Eig): 0.490668402403938 - 0.24830312776203464im\n",
-      "Schur(Eig): 0.8231368203454562 - 0.1754777512534702im\n",
-      "Schur(Eig): 0.66875996887317 - 0.8948574571088873im\n",
-      "Schur(Eig): 0.7948177175783381 - 0.2035284187529041im\n",
-      "Schur(Eig): 0.5520070688753858 - 0.27432211719125693im\n",
-      "Schur(Eig): 0.7664917766276773 - 0.23157845822606135im\n",
-      "Schur(Eig): 0.6690371089548646 - 0.8341383859356367im\n",
-      "Schur(Eig): 0.6693574770682996 - 0.7752535724645789im\n",
-      "Schur(Eig): 0.7381588850411379 - 0.2596285881428227im\n",
-      "Schur(Eig): 0.6094786305339083 - 0.29701138828862733im\n",
-      "Schur(Eig): 0.6697285205350777 - 0.7181868981960513im\n",
-      "Schur(Eig): 0.6701590160826616 - 0.6629202654263124im\n",
-      "Schur(Eig): 0.7098628298448739 - 0.2876864831901749im\n",
-      "Schur(Eig): 0.6647784213473716 - 0.31836286892273497im\n",
-      "Schur(Eig): 0.6807217407217481 - 0.3194733972965826im\n",
-      "Schur(Eig): 0.6747420918448707 - 0.3677702817194503im\n",
-      "Schur(Eig): 0.6706593423085951 - 0.6094332030088696im\n",
-      "Schur(Eig): 0.6712414364632036 - 0.5577026125164619im\n",
-      "Schur(Eig): 0.6736288631068176 - 0.4127747913857582im\n",
-      "Schur(Eig): 0.6719191224217285 - 0.5077027995090068im\n",
-      "Schur(Eig): 0.6727091273727731 - 0.45940473289235245im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n",
-      "Schur(Eig): NaN + NaN*im\n"
-     ]
-    }
-   ],
-   "source": [
-    "#using .Lilly1966_solver\n",
-    "using OffsetArrays  # 配列要素番号を任意開始にするモジュール (行列計算時は 1 始まりの通常 matrix に copy する)\n",
-    "using LinearAlgebra\n",
-    "\n",
-    "#-- Drawing parameters    \n",
-    "N = 200  # wavenumber in Chebyshev\n",
-    "Ncal = N + 1\n",
-    "Ntot = 2 * (N + 1)\n",
-    "\n",
-    "# initialize working arrays\n",
-    "# A: non-lambda coefficients, B: lambda coefficients\n",
-    "A_J = reshape(zeros(Complex{Float64},Ntot,Ntot),Ntot,Ntot)  # Matrix for calculation\n",
-    "B_J = reshape(zeros(Complex{Float64},Ntot,Ntot),Ntot,Ntot)  # Matrix for calculation\n",
-    "### --- From here, supposing matrix element as row, col as in mathematical notation\n",
-    "Acal = reshape(zeros(Complex{Float64},Ntot,Ntot),Ntot,Ntot)  # Matrix for calculation\n",
-    "Bcal = reshape(zeros(Complex{Float64},Ntot,Ntot),Ntot,Ntot)  # Matrix for calculation\n",
-    "A_vv = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for v in v Chebyshev\n",
-    "A_vchi = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for chi in v Chebyshev\n",
-    "A_chiv = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for chi in chi Chebyshev\n",
-    "A_chichi = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for chi in chi Chebyshev\n",
-    "B_chichi = OffsetArray(Acal[1:Ncal,1:Ncal],0:N,0:N)  # coefficient matrix for chi in chi Chebyshev\n",
-    "### --- From here, supposing matrix element as col, row as in julia notation (only using these arrays for matrix calc)\n",
-    "b = OffsetArray(zeros(Float64,Ncal),0:N)  # A vector composed of Chebyshev coefficients for the basic state flow\n",
-    "c = OffsetArray(zeros(Float64,Ncal),0:N)  # factorized vector (c_0=2, c_n=1, (n>0))\n",
-    "\n",
-    "# setting parameters\n",
-    "α = 1.0\n",
-    "Re = 10000.0\n",
-    "\n",
-    "# setting the Chebyshev coefficients for the basic flow\n",
-    "b[0] = 0.5\n",
-    "b[2] = -b[0]\n",
-    "\n",
-    "# setting the factorized vector\n",
-    "c[0] = 2.0\n",
-    "c[1:N] .= 1.0\n",
-    "\n",
-    "# setting maximum even and odd numbers from N\n",
-    "Nemax = 0\n",
-    "Nomax = 1\n",
-    "for i in 0:2:N\n",
-    "    Nemax = i\n",
-    "end\n",
-    "for i in 1:2:N\n",
-    "    Nomax = i\n",
-    "end\n",
-    "println(\"Nemax = \", Nemax, \", Nomax = \", Nomax)\n",
-    "\n",
-    "# setting matrices\n",
-    "# 支配方程式から構築される部分 (チェビシェフ次数 0<=n<=N-2 の範囲)\n",
-    "# この段階では境界条件で成り立つ係数の関係は考慮しない\n",
-    "for n in 0:N-2\n",
-    "    n2 = n^2\n",
-    "    \n",
-    "    # In v rows, a_{v,n} and a_{chi,n}\n",
-    "    for p in n+2:2:N\n",
-    "        p2 = p^2\n",
-    "        A_vv[n,p] = p * (p2 - n2)\n",
-    "    end\n",
-    "    \n",
-    "    A_vv[n,n] = - c[n] * α^2\n",
-    "    A_vchi[n,n] = - c[n]\n",
-    "    \n",
-    "    # In chi rows, a_{v,n} and a_{chi,n}\n",
-    "    for p in n+2:2:N\n",
-    "        p2 = p^2\n",
-    "        A_chichi[n,p] = p * (p2 - n2)\n",
-    "    end\n",
-    "    \n",
-    "    A_chichi[n,n] = - c[n] * α^2\n",
-    "\n",
-    "    for p in -N:N\n",
-    "        abs_n_p = abs(n - p)\n",
-    "        abs_p = abs(p)\n",
-    "        if abs_n_p <= N\n",
-    "            A_chichi[n,abs_p] = A_chichi[n,abs_p] + complex(0.0, - 0.5 * α * Re * c[abs_p] * c[abs_n_p] * b[abs_n_p])\n",
-    "        end\n",
-    "\n",
-    "        if abs_n_p+2 <= N\n",
-    "            acoef = 0.0\n",
-    "            for m in abs_n_p+2:2:N\n",
-    "                m2 = m^2\n",
-    "                acoef = acoef + m * (m2 - (n - p)^2) * b[m]\n",
-    "            end\n",
-    "            A_chiv[n,abs_p] = A_chiv[n,abs_p] + complex(0.0, 0.5 * α * Re * c[abs_p] * acoef)\n",
-    "        end\n",
-    "    end\n",
-    "    \n",
-    "    # In chi rows, lambda coefficients\n",
-    "    B_chichi[n,n] = complex(0.0, - α * Re * c[n])\n",
-    "\n",
-    "end\n",
-    "\n",
-    "# Set boundary conditions\n",
-    "for p in 0:2:N  # even mode\n",
-    "    A_vv[N-1,p] = 1.0\n",
-    "    A_chichi[N-1,p] = 1.0\n",
-    "end\n",
-    "for p in 1:2:N  # even mode\n",
-    "    A_vv[N,p] = 1.0\n",
-    "    A_chichi[N,p] = 1.0\n",
-    "end\n",
-    "\n",
-    "# Gather A_v, A_chi, B_chi -> Acal, Bcal\n",
-    "Acal[1:Ncal,1:Ncal] .= A_vv[0:N,0:N]\n",
-    "Acal[1:Ncal,Ncal+1:Ntot] .= A_vchi[0:N,0:N]\n",
-    "Acal[Ncal+1:Ntot,1:Ncal] .= A_chiv[0:N,0:N]\n",
-    "Acal[Ncal+1:Ntot,Ncal+1:Ntot] .= A_chichi[0:N,0:N]\n",
-    "Bcal[Ncal+1:Ntot,Ncal+1:Ntot] .= B_chichi[0:N,0:N]\n",
-    "\n",
-    "# Convert Transpose type to Matrix (Array) type (ドット演算して要素だけコピーしている)\n",
-    "A_J .= transpose(Acal)  # A_J = Transpose type\n",
-    "B_J .= transpose(Bcal)  # B_J = Transpose type\n",
-    "\n",
-    "F = schur(A_J,B_J)\n",
-    "Dcal = F.α ./ F.β\n",
-    "for i in 1:Ncal\n",
-    "    println(\"Schur(Eig): \", Dcal[i])\n",
-    "end\n",
-    "\n",
-    "#G = eigvecs(Ccal)\n",
-    "#H = inv(G) * Ccal * G\n",
-    "#ncmax = findmax(F)[1]\n",
-    "#Ae = OffsetArray(A[0:div(N,2),0:div(N,2)],0:div(N,2),0:div(N,2))\n",
-    "#Ao = OffsetArray(A[1:div(N,2),1:div(N,2)],1:div(N,2),1:div(N,2))\n",
-    "#Be = OffsetArray(B[0:div(N,2),0:div(N,2)],0:div(N,2),0:div(N,2))\n",
-    "#Bo = OffsetArray(B[1:div(N,2),1:div(N,2)],1:div(N,2),1:div(N,2))\n",
-    "#for n in 0:div(N,2)\n",
-    "#    for p in 0:div(N,2)\n",
-    "#        Ae[n,p] = A[2*n,2*p]\n",
-    "#        Be[n,p] = B[2*n,2*p]\n",
-    "#    end\n",
-    "#end\n",
-    "#for n in 1:div(N,2)\n",
-    "#    for p in 1:div(N,2)\n",
-    "#        Ao[n,p] = A[2*n-1,2*p-1]\n",
-    "#        Bo[n,p] = B[2*n-1,2*p-1]\n",
-    "#    end\n",
-    "#end\n",
-    "#Acal[1:div(N,2)-1,1:div(N,2)-1] .= Ae[0:div(N,2)-2,0:div(N,2)-2]\n",
-    "#Bcal[1:div(N,2)-1,1:div(N,2)-1] .= Be[0:div(N,2)-2,0:div(N,2)-2]\n",
-    "#Ccal = inv(Bcal[1:div(N,2)-1,1:div(N,2)-1]) * Acal[1:div(N,2)-1,1:div(N,2)-1]\n",
-    "##Ccal = inv(Bcal) * Acal\n",
-    "#Dcal = eigvals(Ccal)\n",
-    "#F = imag.(Dcal)\n",
-    "#for i in 1:div(N,2)-1\n",
-    "#    println(F[i])\n",
-    "#end\n",
-    "\n",
-    "#Acal[1:div(N,2)-2,1:div(N,2)-2] .= Ao[1:div(N,2)-2,1:div(N,2)-2]\n",
-    "#Bcal[1:div(N,2)-2,1:div(N,2)-2] .= Bo[1:div(N,2)-2,1:div(N,2)-2]\n",
-    "#Ccal = inv(Bcal[1:div(N,2)-2,1:div(N,2)-2]) * Acal[1:div(N,2)-2,1:div(N,2)-2]\n",
-    "##Ccal = inv(Bcal) * Acal\n",
-    "#Dcal = eigvals(Ccal)\n",
-    "#F = imag.(Dcal)\n",
-    "#for i in 1:div(N,2)-2\n",
-    "#    println(F[i])\n",
-    "#end\n"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 8,
-   "id": "1255c473",
-   "metadata": {},
-   "outputs": [
-    {
-     "data": {
-      "text/plain": [
-       "Figure(PyObject <Figure size 700x500 with 1 Axes>)"
-      ]
-     },
-     "metadata": {},
-     "output_type": "display_data"
-    }
-   ],
-   "source": [
-    "#-- Drawing\n",
-    "##########\n",
-    "#  Plot  #\n",
-    "##########\n",
-    "using PyPlot\n",
-    "\n",
-    "rc(\"font\", family=\"IPAPGothic\")\n",
-    "fig = figure(\"pyplot_majorminor\",figsize=(7,5))\n",
-    "\n",
-    "x = real.(Dcal[3:Ncal-4])\n",
-    "y = imag.(Dcal[3:Ncal-4])\n",
-    "\n",
-    "ax = gca()\n",
-    "\n",
-    "cp = scatter(x,y)\n",
-    "\n",
-    "ax = gca()\n",
-    "\n",
-    "xlabel(\"X\")\n",
-    "ylabel(\"Y\")\n",
-    "#if i == 1 #div(ne,2)-1\n",
-    "#   plt.colorbar(cp)\n",
-    "#end\n",
-    "grid(\"on\")\n",
-    "\n",
-    "fig.canvas.draw()\n",
-    "gcf() # Needed for IJulia to plot inline\n",
-    "\n",
-    "#PyPlot.title(\"Re = Re[i]\")\n",
-    "\n",
-    "fname = \"test.png\"\n",
-    "#savefig(\"B-init_2d.pdf\")\n",
-    "savefig(fname)\n",
-    "#psi\n"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 6,
-   "id": "623f8c38",
-   "metadata": {},
-   "outputs": [
-    {
-     "ename": "LoadError",
-     "evalue": "UndefVarError: Ccal not defined",
-     "output_type": "error",
-     "traceback": [
-      "UndefVarError: Ccal not defined",
-      "",
-      "Stacktrace:",
-      " [1] top-level scope",
-      "   @ In[6]:27",
-      " [2] eval",
-      "   @ ./boot.jl:360 [inlined]",
-      " [3] include_string(mapexpr::IJulia.var\"#37#38\", mod::Module, code::String, filename::String)",
-      "   @ Base ./loading.jl:1094",
-      " [4] execute_code(code::String, filename::String)",
-      "   @ IJulia ~/.julia/packages/IJulia/a1SNk/src/execute_request.jl:27",
-      " [5] execute_request(socket::ZMQ.Socket, msg::IJulia.Msg)",
-      "   @ IJulia ~/.julia/packages/IJulia/a1SNk/src/execute_request.jl:86",
-      " [6] #invokelatest#2",
-      "   @ ./essentials.jl:708 [inlined]",
-      " [7] invokelatest",
-      "   @ ./essentials.jl:706 [inlined]",
-      " [8] eventloop(socket::ZMQ.Socket)",
-      "   @ IJulia ~/.julia/packages/IJulia/a1SNk/src/eventloop.jl:8",
-      " [9] (::IJulia.var\"#15#18\")()",
-      "   @ IJulia ./task.jl:411"
-     ]
-    }
-   ],
-   "source": [
-    "#-- Drawing\n",
-    "##########\n",
-    "#  Plot  #\n",
-    "##########\n",
-    "#using PyPlot\n",
-    "\n",
-    "#rc(\"font\", family=\"IPAPGothic\")\n",
-    "#fig = figure(\"pyplot_majorminor\",figsize=(7,5))\n",
-    "\n",
-    "nx = Ncal-4\n",
-    "ny = Ncal-4\n",
-    "\n",
-    "xax = reshape(zeros(nx,ny),nx,ny)\n",
-    "yax = reshape(zeros(nx,ny),nx,ny)\n",
-    "\n",
-    "for i in 1:nx\n",
-    "    xax[i,1:ny] .= i  # xmin + (xmax-xmin) * (i-1) / (nx-1)\n",
-    "end\n",
-    "for i in 1:ny\n",
-    "    yax[1:nx,i] .= i  # ymin + (ymax-ymin) * (i-1) / (ny-1)\n",
-    "#    yax[1:nx,i] .= Re[i]  # ymin + (ymax-ymin) * (i-1) / (ny-1)\n",
-    "end\n",
-    "\n",
-    "cont_val = reshape(zeros(nx,ny),nx,ny)\n",
-    "shade_val = reshape(zeros(nx,ny),nx,ny)\n",
-    "#for i in 1:nr#div(nr,2)-1:div(nr,2)+1\n",
-    "    cont_val = imag.(Ccal[1:nx,1:ny])\n",
-    "    #for j in 1:ne\n",
-    "    shade_val = imag.(A_J[1:nx,1:ny])\n",
-    "    #end\n",
-    "#draw_num2 = 2\n",
-    "#if draw_num2 == 1\n",
-    "#    cp = contourf(xax[1:nx,1:nx], xax'[1:nx,1:nx], Pftime[1:nx,1:nx,1], levels=[-15.0, -10.0, -5.0, 0.0, 5.0, 10.0, 15.0, 20.0], origin=\"image\", cmap=ColorMap(\"viridis\"), extend=\"both\")\n",
-    "    #cp = contourf(xax[1:nx,1:nx], xax'[1:nx,1:nx], Pftime[1:nx,1:nx,2], levels=[-2.0, -1.0, -0.5, -0.25, 0.25, 0.5, 1.0, 1.5, 2.0], origin=\"image\", cmap=ColorMap(\"viridis\"), extend=\"both\")\n",
-    "    #cp = contourf(xax[1:nx,1:no], xax'[1:nx,1:no], Pftime[1:nx,1:nx,2]*Hop'*inv(Hop*Pftime[1:nx,1:nx,2]*Hop'+Ro), origin=\"image\", cmap=ColorMap(\"viridis\"), extend=\"both\")\n",
-    "#elseif draw_num2 == 2\n",
-    "    ##cp = contour(xax[1:nx,1:ny], yax[1:nx,1:ny], shade_val[1:nx,1:ny], 8, levels=[-0.01, -0.005, -0.003, -0.002, -0.001, -0.0005, 0.0, 0.0005, 0.001, 0.002, 0.003, 0.005, 0.01])\n",
-    "    cp = contourf(xax[1:nx,1:ny], yax[1:nx,1:ny], cont_val[1:nx,1:ny], cmap=ColorMap(\"viridis\"), extend=\"both\")\n",
-    "    ##cp = contourf(xax[1:nx,1:ny], yax[1:nx,1:ny], cont_val[1:nx,1:ny], 8, levels=[-1.0, -0.5, -0.2, -0.1, 0.0, 0.1, 0.2, 0.5, 1.0], cmap=ColorMap(\"viridis\"), extend=\"both\")\n",
-    "    #cp = contourf(xax[1:nx,1:ny], yax[1:nx,1:ny], shade_val[1:nx,1:ny], 10, cmap=ColorMap(\"viridis\"), extend=\"both\")\n",
-    "#elseif draw_num2 == 3\n",
-    "#    cp = contourf(xoax[1:nx,1:nto-1], toax[1:nx,1:nto-1], evec_max[1:nx,1:nto-1], 10, cmap=ColorMap(\"viridis\"), extend=\"both\")\n",
-    "#elseif draw_num2 == 4\n",
-    "#    cp = contourf(xoax[1:nx,1:40], toax[1:nx,1:40], x_inc[1:nx,1:40], levels=[-2.0, -1.0, -0.5, -0.25, 0.25, 0.5, 1.0, 1.5, 2.0], cmap=ColorMap(\"viridis\"), extend=\"both\")\n",
-    "#elseif draw_num2 == 5\n",
-    "#    cp = contourf(xoax[1:no,1:40], toax[1:no,1:40], y_innov[1:no,1:40], levels=[-2.0, -1.0, -0.5, -0.25, 0.25, 0.5, 1.0, 1.5, 2.0], cmap=ColorMap(\"viridis\"), extend=\"both\")\n",
-    "#end\n",
-    "#ax.label(cp, inline=1, fontsize=10)\n",
-    "#legend()\n",
-    "ax = gca()\n",
-    "\n",
-    "xlabel(\"X\")\n",
-    "ylabel(\"Y\")\n",
-    "#if i == 1 #div(ne,2)-1\n",
-    "   plt.colorbar(cp)\n",
-    "#end\n",
-    "grid(\"on\")\n",
-    "\n",
-    "PyPlot.title(\"Re = Re[i]\")\n",
-    "\n",
-    "#########################\n",
-    "#  Set tick dimensions  #\n",
-    "#########################\n",
-    "#ax.xaxis.set_tick_params(which=\"major\",length=5,width=2,labelsize=10)\n",
-    "#ax.xaxis.set_tick_params(which=\"minor\",length=5,width=2)\n",
-    "\n",
-    "fig.canvas.draw() # Update the figure\n",
-    "gcf() # Needed for IJulia to plot inline\n",
-    "#fname = \"Re\" * string(1) * \".png\"\n",
-    "fname = \"test.png\"\n",
-    "#savefig(\"B-init_2d.pdf\")\n",
-    "#savefig(fname)\n",
-    "#psi\n",
-    "#end\n",
-    "#println(shade_val[1:na,1])"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": 47,
-   "id": "af6a72bd",
-   "metadata": {},
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "ComplexF64[0.0 + 0.0im 0.0 + 0.0im; 0.0 + 0.0im 0.0 + 0.0im]\n"
-     ]
-    },
-    {
-     "data": {
-      "text/plain": [
-       "GeneralizedSchur{ComplexF64, Matrix{ComplexF64}}\n",
-       "S factor:\n",
-       "2×2 Matrix{ComplexF64}:\n",
-       " 2.23607-0.0im      0.0+0.0im\n",
-       "     0.0+0.0im  2.23607+0.0im\n",
-       "T factor:\n",
-       "2×2 Matrix{ComplexF64}:\n",
-       " 2.23607+0.0im      0.0+0.0im\n",
-       "     0.0+0.0im  2.23607+0.0im\n",
-       "Q factor:\n",
-       "2×2 Matrix{ComplexF64}:\n",
-       " -0.447214-0.894427im  0.0+0.0im\n",
-       "       0.0+0.0im       1.0+0.0im\n",
-       "Z factor:\n",
-       "2×2 Matrix{ComplexF64}:\n",
-       " -1.0-0.0im       0.0+0.0im\n",
-       " -0.0-0.0im  0.894427+0.447214im\n",
-       "α:\n",
-       "2-element Vector{ComplexF64}:\n",
-       " 2.23606797749979 - 0.0im\n",
-       " 2.23606797749979 + 0.0im\n",
-       "β:\n",
-       "2-element Vector{ComplexF64}:\n",
-       " 2.23606797749979 + 0.0im\n",
-       " 2.23606797749979 + 0.0im"
-      ]
-     },
-     "execution_count": 47,
-     "metadata": {},
-     "output_type": "execute_result"
-    }
-   ],
-   "source": [
-    "a33 = reshape(zeros(Complex{Float64},2,2),2,2)\n",
-    "c33 = reshape(zeros(Complex{Float64},2,2),2,2)\n",
-    "println(a33)\n",
-    "a33[1,1] = 1.0 + 2.0im\n",
-    "a33[2,2] = 2.0 - 1.0im\n",
-    "b33 = transpose(a33)\n",
-    "c33 .= b33\n",
-    "schur(a33,c33)"
-   ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "id": "b28e1119",
-   "metadata": {},
-   "outputs": [],
-   "source": []
-  }
- ],
- "metadata": {
-  "kernelspec": {
-   "display_name": "Julia 1.6.1",
-   "language": "julia",
-   "name": "julia-1.6"
-  },
-  "language_info": {
-   "file_extension": ".jl",
-   "mimetype": "application/julia",
-   "name": "julia",
-   "version": "1.6.1"
-  }
- },
- "nbformat": 4,
- "nbformat_minor": 5
-}