フーリエ解析：円でリャマを描く¶

フーリエ変換は、入力信号を時間領域から周波数領域へ変換する手法だ。応用範囲はきわめて広い。たとえば音響エンジニアはフーリエ変換で楽曲から不要な周波数を取り出し、逆フーリエ変換でその周波数を除いた音を得る。デジタル画像処理では、フーリエ変換を使って画像フィルタリングができる。ガウスぼかしや JPEG 圧縮がその例だ。そして最後に、これから本記事で展開するように——お絵かきである。 お絵かきとは何のことか。考え方はこうだ。何か——魚、π 記号、あるいは描けるものであればほぼ何でも——を表す一本の経路を取る。紙に鉛筆を下ろし、途中で持ち上げずに一筆で描き切る。この経路があれば、異なる周波数で円の周りを回転するベクトルを何本も端から端へつなげられ、最後のベクトルの先端が最初のスケッチを、少なくともそれにかなり近いものをなぞる。いくつか例を挙げよう：

図 1 円で直線を描く¶

直線とはできるだけかけ離れた図形で、直線が描けるのだ。これらの円を変えれば描かれる形も変わる。たとえば外側の円の半径を小さくし、回転速度を上げると正方形が描ける。

図 2 正方形¶

では三つ目の円を加えたらどうなるか。曲線はさらに細かくなる。三つの円があれば、魚に見える曲線が描ける。

図 3 魚。¶

つなげた円だけでどこまで極端なことができるか。ここにリャマがいる。描画に $1024$ 個の円を使い、それぞれが異なる周波数、半径、開始角度を持つ。一般に、（ほぼ）すべての閉曲線は円で描ける。これらの円をエピサイクルと呼ぶ。

図 4 リャマ？！¶

一般に、円を増やすほど描ける図形は複雑になる。同じ絵に対しては、円を増やすほど元の「入力」図形に近づく。このアニメーションがどう成り立っているかを理解するには、フーリエ変換（およびフーリエ級数）の背後にある数学と直観を掘り下げる必要がある。前述のとおり、円で異なる曲線を描くには三つの条件を調整できる。回転する速さ（周波数）、円の大きさ（半径）、そして開始角度（位相）だ。この三つはすべて一つの複素数で表せる。

次の例は自分で動かして、円を増やすことが最終的な図形に与える影響を体感してみてほしい。

フーリエ変換¶

だがその前に、まずフーリエ変換とは何かを簡単に説明しておくべきだろう。時間領域から周波数領域への変換だ。どういうことか。関数 $\sin (2 π \times 3 t)$ があるとして、この $\sin (2 π \times 3 t)$ をフーリエ変換すると、周波数 $\pm 3 hz$ のところに鋭いピークが現れる。

図 5 $\sin (2 π \times 3 t)$ のグラフ¶

図 6 $\sin (2 π \times 3 t)$ の周波数領域における振幅。¶

案の定、周波数 $3 hz$ と $- 3 hz$ にピークがある。要するに、扱いたい対象にどんな周波数が含まれるかを見つけ出し、それらの正弦波を組み合わせて元の対象の近似を作る、ということだ。以降の節で、フーリエ変換がどう働くかを見ていく。

複素数¶

この節は複素数とは何かの簡単な復習なので、すでに知っていれば飛ばしてよい。複素数は「実」部と「虚」部の二つからなる。すでに慣れ親しんだ数を取り、本質的にもう一本軸を足して表したものだ。複素数は $3 + 4 i$ と書き、 $3$ が実部、 $4 i$ が虚部である。この二つを足し合わせて一つの複素数になる。ここで複素数には便利な性質がある。二次元ベクトルとして扱いやすく、実部が $x$ 、虚部が $y$ で、二次元ユークリッド空間の標準的な点 $(x, y)$ に対応する。三角関数と複素指数関数の結びつきのおかげで、この表現はここでは非常に都合がよい。両者を結ぶ式をオイラーの等式と呼ぶ。本記事の後の節でさらに詳しく扱う。

線形代数¶

線形代数は、ベクトルを円の上で回転させて図形を描く仕組み、具体的にはあるベクトルを他のベクトルで表す仕組みを理解するうえで重要な概念だ。標準的な二次元 $(x, y)$ ユークリッド空間では、すべてのベクトルは $𝐱$ と $𝐲$ の単位ベクトルで表せる。この二つはふつう $\hat{𝐢}$ と $\hat{𝐣}$ （「i-hat」「j-hat」と読む）と書き、 $\hat{𝐢} = (1, 0)$ （ $x$ 軸）、 $\hat{𝐣} = (0, 1)$ （ $y$ 軸）である。これには二つの演算が使える。スカラー倍と加法だ。たとえばベクトル $\vec{a} = (5, - 8)$ は次のように表せる。

(1)¶

\begin{matrix} \begin{matrix} \begin{matrix} \vec{a} & = 5 \cdot \hat{𝐢} - 8 \cdot \hat{𝐣} \\ = 5 \cdot (1, 0) + - 8 \cdot (0, 1) \\ = (5, 0) + (0, - 8) \\ = (5, - 8) \end{matrix} \end{matrix} \end{matrix}

このように、スカラー倍とベクトルの和を組み合わせて別のベクトルを作る操作を、ベクトルの線形結合と呼ぶ。この例では、そのベクトルの組は単位ベクトル $\hat{𝐢}$ と $\hat{𝐣}$ である。

図 7 2つのベクトル $\vec{u}$ と $\vec{v}$ 。なす角は $90^{\circ}$ 、長さはどちらも $1$ ¶

ただし、ほかのベクトルを表すのに使えるベクトルは $\hat{𝐢}$ と $\hat{𝐣}$ だけではない。実は、互いに倍数の関係にない2つのベクトルであれば、それぞれを拡大縮小して足し合わせることで任意のベクトルを作れる。たとえば $(1, - 1)$ と $(- 2, 2)$ は使えない。 $(- 2,2) = - 2 \times (1, - 1)$ だからだ。上の図では、ベクトル $\vec{u} = (\frac{1}{\sqrt{2}}, \frac{1}{\sqrt{2}})$ と $\vec{v} = (- \frac{1}{\sqrt{2}}, \frac{1}{\sqrt{2}})$ は長さが 1 で、互いに直交している（なす角が $90^{\circ}$ ）。この2つはどちらも拡大縮小してもう一方にはなれないので線形独立であり、2次元空間全体の基底ベクトルとして使える。

内積（「ドット積」とも呼ぶ）は、選んだベクトルを使ってあるベクトルを表すとき、各ベクトルを何倍すればよいかを求めるために使う。ベクトル $\vec{𝐰} = (0.7, 1.2)$ （上の図を参照）を $\vec{u}$ と $\vec{v}$ の線形結合 $\vec{w} = c_{1} \cdot \vec{u} + c_{2} \cdot \vec{v}$ として表すには、この2つのベクトルそれぞれに掛ける数を求めればよい。具体的には、内積が与えてくれるのはこの2つの係数 $c_{1}$ と $c_{2}$ そのものだ。別の言い方をすれば、内積は「 $\vec{w}$ の中に $\vec{u}$ がどれだけ含まれているか」という問いに答える。この直感は、Fourier 変換の背後にある数学を扱うときに重要になる。ともあれ、 $\vec{w}$ と $\vec{u}$ の内積を取れば $c_{1}$ が求まり、順序は関係ない（たとえば $⟨ \vec{u}, \vec{w} ⟩ = ⟨ \vec{w}, \vec{u} ⟩$ ）。

(2)¶

\begin{matrix} \begin{matrix} \begin{matrix} c_{1} & = ⟨ \vec{u}, \vec{w} ⟩ \\ = (u_{x} \times w_{x}) + (u_{y} \times w_{y}) \\ = (\frac{1}{\sqrt{2}} \times 0.7) + (\frac{1}{\sqrt{2}} \times 1.2) \\ = 1.3435 \end{matrix} \end{matrix} \end{matrix}

つまり $\vec{u}$ を $1.3435$ 倍、すなわち元の長さの $1.3435$ 倍に引き伸ばす。したがって $\vec{w}$ のうち $\vec{u}$ に属する成分（ $c_{1}$ ）は $1.3435$ である。

(3)¶

\begin{matrix} \begin{matrix} \begin{matrix} c_{2} & = ⟨ \vec{v}, \vec{w} ⟩ \\ = (v_{x} \times w_{x}) + (v_{y} \times w_{y}) \\ = (- \frac{1}{\sqrt{2}} \times 0.7) + (\frac{1}{\sqrt{2}} \times 1.2) \\ = 0.353553 \end{matrix} \end{matrix} \end{matrix}

同様に、 $\vec{v}$ を $0.353553$ 倍、つまり元のサイズの $0.353553$ 倍に縮める。よって $\vec{w}$ のうち $\vec{v}$ の成分（ $c_{2}$ ）は $0.353553$ である。

最終的に得られる線形結合は

(4)¶

\begin{matrix} \begin{matrix} \begin{matrix} \vec{w} & = 1.3435 \cdot \vec{u} + 0.353553 \cdot \vec{v} \\ = (0.95, 0.95) + (- 0.25, 0.25) \\ = (0.7, 1.2) \end{matrix} \end{matrix} \end{matrix}

これは確かに最初のベクトル $\vec{w}$ に等しい。

図 8 $\vec{u}$ と $\vec{v}$ の線形結合で $\vec{w}$ を得る¶

任意の2つのベクトルの内積は次のように定義される。

(5)¶

\begin{matrix} \begin{matrix} \begin{matrix} ⟨ \vec{v_{1}}, \vec{v_{2}} ⟩ & = {v_{1}}_{x} \times {v_{2}}_{x} + {v_{1}}_{y} \times {v_{2}}_{y} \\ = ‖ \vec{v_{1}} ‖ \times ‖ \vec{v_{2}} ‖ \times \cos (θ) \end{matrix} \end{matrix} \end{matrix}

ここで $θ$ は $v_{1}$ と $v_{2}$ のなす角、 $‖ \vec{v} ‖$ はベクトルの大きさである。大きさはノルムとも呼ばれ、 $‖ \vec{v} ‖ = \sqrt{v_{x}^{2} + v_{y}^{2}}$ と定義される。この式には見覚えがあるはずだ。小中学校の数学で習った三平方の定理そのものである。

図 9 直角三角形¶

内積には、ちょうど我々の役に立つ重要な性質がいくつかある。まず、任意のベクトルと自身の内積はそのベクトルの大きさの2乗に等しい（ $\vec{v_{1}} \cdot \vec{v_{1}} = ‖ \vec{v_{1}} ‖^{2}$ ）。例えば、

(6)¶

\begin{matrix} ⟨ [5, 7], [5, 7] ⟩ = {\sqrt{5^{2} + 7^{2}}}^{2} = {\sqrt{74}}^{2} = 74 \end{matrix}

次に、互いに直交する（直角に交わる）2つのベクトルの内積は $0$ になる。 $\cos (\frac{π}{2}) = 0$ だからである。したがって、 $\vec{v_{1}}, \vec{v_{2}}$ が直交するとき、

(7)¶

\begin{matrix} ⟨ \vec{v_{1}}, \vec{v_{2}} ⟩ = ‖ \vec{v_{1}} ‖ \times ‖ \vec{v_{2}} ‖ \times 0 = 0 \end{matrix}

先ほどの「 $\vec{v_{1}}$ の中に $\vec{v_{2}}$ がいくつあるか」という直感的な理解に照らせば、この結果が $0$ になるのは妥当である。この2つのベクトルは直交しているからだ。

連続関数の内積¶

前節では2つのベクトルの内積を扱ったが、この定義は連続関数に拡張できる。鍵となる違いは和の取り方にある。前述のとおり、2次元の内積の式は

(8)¶

\begin{matrix} ⟨ \vec{a}, \vec{b} ⟩ = a_{x} \times b_{x} + a_{y} \times b_{y} \end{matrix}

である。より正確には、2つのベクトルの対応する項を足し合わせている。このときベクトルの長さは有限（ここでは $2$ ）なので、和は離散的である。この式は任意の次元について次のようにも書ける。

(9)¶

\begin{matrix} ⟨ \vec{a}, \vec{b} ⟩ = \sum_{n = 1}^{N} a_{n} \times b_{n} \end{matrix}

ここで $N$ は次元である。このシグマ記号による和は離散区間上で取られる。例えば3次元空間なら、和の区間は $[1, 3]$ である。ただし、この和が使えるのは $ℝ^{n}$ のベクトルに限る。長さ $n$ の複素ベクトル $ℂ^{n}$ の式もよく似ているが、2番目のベクトルの複素共役を取る点が異なる。

(10)¶

\begin{matrix} ⟨ \vec{a}, \vec{b} ⟩ = \sum_{n = 1}^{N} a_{n} \times b_{n} \end{matrix}

ここで複素共役を取ることで、この式は「2番目の引数について反線形」になる。つまり $\vec{a} \cdot (\vec{b} + \vec{c}) = \vec{a} \cdot \vec{b} + \vec{a} \cdot \vec{c}$ である。ただしスカラー倍について線形ではなく、スカラーの共役を外に出す必要がある。 $\vec{a} \cdot (b \times \vec{c}) = b \times (\vec{a} \cdot \vec{c})$ となる。一方、1番目の引数については線形である。つまり、内積の1番目の引数からスカラー倍をくくり出すとき、そのスカラーの共役を取る必要はない。一般に、複素共役を取る側の引数について反線形、もう一方の引数について線形である。

連続区間上で和を取る必要があり、しかもその区間はいくらでも小さくできる。まさに積分の出番だ。

(11)¶

\begin{matrix} ⟨ f (t), g (t) ⟩ = \frac{1}{T} \int_{0}^{T} f (t) g (t) d t \end{matrix}

唯一の条件は、関数 $f (t)$ と $g (t)$ が区間 $[0, T]$ 上で二乗可積分であること、つまり絶対値の二乗の積分が $[0, T]$ 上で有限であることだ。

(12)¶

\begin{matrix} f (t) \in L^{2} (0, T) if \int_{0}^{T} | f (t) |^{2} d t is finite. \end{matrix}

注意 $L^{2}$ は二乗可積分関数の集合、 $L^{2} (0, T)$ は区間 $[0, T]$ 上の二乗可積分関数の集合である。

2D 空間の内積の定義に寄せて、見た目を似せるなら、この積分をリーマン和（シグマ記号）で書ける。

(13)¶

\begin{matrix} \frac{1}{T} \int_{0}^{T} f (t) g (t) d t = \lim_{N \to \infty} \frac{1}{T} \sum_{k = 0}^{N} f (\frac{T \times k}{N}) \times g (\frac{T \times k}{N}) \times \frac{T}{N} \end{matrix}

ここで $\frac{T}{N}$ が $d t$ にあたり、各ステップで掛ける無限小量だ。関数 $f$ と $g$ の引数は、積分の上端 $T$ に、現在の和のステップ数と極限で近づく値を掛けたもので、これにより $f$ と $g$ の引数は常に積分区間内に収まる。ステップ数が小さいうち（ $k = 0, 1, 2$ など）は $\frac{T \times k}{N}$ は $0$ に非常に近く、引数も $0$ に非常に近い。 $k$ が $N$ に近づくと $\frac{T \times k}{N}$ は $1$ に近づき、 $T$ に近づく。

複素指数で正弦と余弦を表す¶

複素指数を使うと正弦と余弦を都合よく扱える。三角法と複素平面の間には意外な関係があり、それを表す公式がオイラーの等式だ。以下がその式で、複素数の節で約束したとおりである。

(14)¶

\begin{matrix} e^{i x} = \cos (x) + i \times \sin (x) \end{matrix}

ここで $x$ は円をどれだけ進むかを表す。たとえば $e^{\frac{π}{4} i}$ は円を $\frac{π}{4}$ ラジアン進み、 $\frac{\sqrt{2}}{2} + \frac{\sqrt{2}}{2} i$ に落ち着く。また $e^{π i}$ は円を $π$ ラジアン進み、 $- 1 + 0 i$ に落ち着く。これがオイラーのもう一つの有名な等式だ。

(15)¶

\begin{matrix} e^{π i} = - 1 \end{matrix}

図 10 $e^{i x}$ と三角関数の関係¶

なぜこうなるのか。 $e$ を虚数乗すると、なぜ単位円を回るのか。答えは物理にある。時刻 $t$ における位置 $p$ が

(16)¶

\begin{matrix} p (t) = e^{i t} \end{matrix}

なら、その点の速度は $p (t)$ の一階微分だ。

(17)¶

\begin{matrix} v (t) = \frac{d p}{d t} = i \times e^{i t} \end{matrix}

これで十分だ。速度は位置と同じ形をしていて、実数が虚数に変わっただけである。 $i$ を掛けることは反時計回りに $90^{\circ}$ 回転させることに等しい。たとえば $a \in ℝ$ なら、 $i \times a$ は元の $a$ が虚軸上に移ったものだ。しかし $a$ がもともと虚数ならどうか。このとき $a = i \times k$ （ $k$ はある実数）である。さらに $i$ を掛けると $i \times a = i^{2} \times k$ となる。 $i$ の定義より $i^{2}$ は $- 1$ なので、 $i \times a = - k$ 。これは実数の符号を反転させただけ、つまり $180^{\circ}$ の回転だ。 $i$ を掛けることが $90^{\circ}$ の回転に等しいのだから、速度関数は常に現在位置の方向と直交する。こうして描かれる曲線は円であり、 $f (t) = e^{i t} where t \in [0, 2 π]$ が円を描くのはそのためである。

補足

(18)¶

\begin{matrix} \frac{1}{2 π} \int_{0}^{2 π} \cos (t) \sin (t) d t = 0 \end{matrix}

$\cos t$ と $\sin t$ は直交関数であり、これによって $\cos$ と $\sin$ がフーリエ変換の基底を成す。

図 11 は $\cos$ と $\sin$ が直交関数であることを示している。¶

本質的にここで行っているのは、2 次元複素空間 $(1,0), (0, 1) i$ から、関数を基底 $(. . ., e^{i \frac{- 4 π}{T} t}, e^{i \frac{- 2 π}{T} t}, 1, e^{i \frac{2 π}{T} t}, e^{i \frac{4 π}{T} t}, . . .)$ とする無限次元空間への基底変換である。つまり、複素数を生成する関数から、無限次元の複素ベクトルを生成する新しい関数へと写像している。オイラーの公式を使って基底を $(. . ., \cos (\frac{- 2 π}{T}), \sin (\frac{- 2 π}{T}), 1, \cos (\frac{2 π}{T}), \sin (\frac{2 π}{T}), . . .)$ に変換すれば、これは $ℝ$ 上の無限次元基底としても表せる。線形演算であり、ある関数空間から別の関数空間への写像であることから、これを線形作用素と呼ぶ。

方形波の近似¶

図 12 方形波。¶

方形波は信号処理における「hello world」のようなものだ。方形波は $[- π, π]$ 上の周期波として定義され、半周期が「オン」すなわち $1$ 、残りの半周期が「オフ」すなわち $0$ となる。

(19)¶

\begin{matrix} square (t) = sgn (\sin (t)) \end{matrix}

ここで sgn は符号関数であり、次のように定義される。

(20)¶

\begin{split} \begin{equation}\text{sgn}(x)=\begin{cases}1 &\text{if}\; x>0 \\ 0 &\text{if}\; x=0 \\ -1 &\text{if}\; x ## リャマを描く¶ 描きたいものを十分になめらかな閉曲線として表せるなら、前節で示したのと同じ原理でほとんど何でも描ける。リャマを表す関数を得る方法は 2 通りある。ベクターグラフィックスプログラムでリャマを描いて SVG として書き出し、その SVG をそのまま使う方法だ。SVG はもともと数学的に定義された曲線（直線、ベジェ曲線など）でできている。もう一つは、その SVG を等間隔にサンプリングして点列を得て、その点列に離散フーリエ変換（DFT）を適用してフーリエ係数を求める方法である。本記事では後者を採る。こちらのほうが直接的だからだ。 ``` const numSamples = 400; let pts = []; function entry() { const svgUrl = './images/llama.svg'; const response = await fetch(svgUrl); const text = await response.text(); pts = await sampleSvgPoints(text, numSamples); } entry(); // Might not sample exactly numPts number of sample points, but it will be // close. Not a concern using the O(N^2) DFT. Using FFT, that can be // concerning since the ideal number of point samples is an exponent of 2. async function sampleSvgPoints(xml: string, numPts: number): Promise { const parser = new DOMParser(); const doc = parser.parseFromString(xml, 'application/xml'); let samplePts = []; const paths = Array.from(doc.querySelectorAll('path')); const totalLength = paths.reduce((acc, v) => acc + v.getTotalLength(), 0); const step = (totalLength / numPts); paths.forEach(path => { const sample = []; const n = M.floor(path.getTotalLength()); for (let i = 0; i numPts) samplePts.pop(); return samplePts; } function toPt(pt: SVGPoint): [number, number] { return [pt.x, pt.y]; } ``` SVG をサンプリングして点を取る次に、これらのサンプリング点に離散フーリエ変換（DFT）を適用し、必要なフーリエ係数を求める。上で方形波を扱ったのと同じ要領だ。離散フーリエ変換の式は次のとおり。 (33)¶\[ \begin{equation}c_k=\frac{1}{N}\sum_{j=0}^{N}{p(j)e^{-i\frac{2\pi}{N}kj}}\end{equation} ParseError: Expected & or \\ or \cr or \end at position 121: …\text{if}\; x #̲# リャマを描く¶ 描きたい…

ここで $N$ は変換対象のサンプリング点リストに含まれる点の総数である。今回の場面では、 $p (k)$ は要するにインデックスで配列の要素を取り出すことにあたる。たとえば $p$ がサンプリング点のリストなら、Javascript では p[k] となる。DFT とフーリエ級数は実は同じもので、違いは、等間隔の有限個のサンプリング点について和を取るか、刻み幅を $0$ に近づける極限、言い換えれば和の項数を $\infty$ に近づける極限を取るかにある。DFT はまた、時間領域の $N$ 個の等間隔な複素サンプリング点を周波数領域の表現に変換する方法でもあり、つまり各整数倍の正弦周波数の強度を求める。これはまさに、 $N$ 次元基底 $(. . ., e^{i \frac{- 4 π}{T} t}, e^{i \frac{- 2 π}{T} t}, 1, e^{i \frac{2 π}{T} t}, e^{i \frac{4 π}{T} t}, . . .)$ に基底変換を直接適用することであり、数節前の Complex Exponentials as Sines and Cosines の節で説明したとおりである。

DFT は 2 重ループでかなり簡単に実装できる。外側のループで各係数 $c_{k}$ を計算し、内側のループで数値積分の和を取り、 $N$ で割って各周波数 $k$ の平均強度を求める。この方法は非常に遅いが、点の数がそれほど多くないうちは使える。計算量は $𝒪 (k j)$ 、最悪の場合で $𝒪 (N^{2})$ である。より効率的なアルゴリズムを使えばはるかに良くできる。たとえば Radix-2 高速フーリエ変換（FFT）なら平均計算量は $𝒪 (N \log (N))$ だが、この記事はすでに長すぎるので、それはまた別の機会の話題にとっておく。

import * as M from 'mathjs';

export type DFTData = M.Complex & {
 freq: number;
 radius: number;
 phase: number;
};

export async function dft(cps: M.Complex[]): Promise {
 const X = [];
 const N = cps.length;
 for (let i = 0; i

$𝒪 (N^{2})$ DFT implementation in TypeScript.

係数が手に入ったので、曲線をパラメータ化して描画できるようにする。

function parameterize(cn: M.Complex[]) {
 return async (t: number): Promise =>
 await sigma(k =>
 M.multiply(M.complex(cn[k - 1].re, cn[k - 1].im), 
 M.exp(M.multiply(M.i, (k - 1) * t))), [1, numSamples]) as M.Complex;
}

この係数でラマを近似したとき、正弦項を増やすほどラマに近づいていくはずだ。次のアニメーションで確認できる。ここでいう「モード」は使う正の周波数の個数を指し、 $0$ の周波数と対応する負の周波数も併せて使っている。

Fig. 14 $1 \to 500$ modes of a llama approximation.¶

このアニメーションを見ると、作業の大部分は低い周波数で済んでいることが分かる。すぐにラマらしい形になるからだ。高周波を足していくにつれて、細かい部分が少しずつ補われる。

In order to recreate the image at the beginning of the article, with the circles moving around in organized chaos drawing the outline of the llama, we need some code to go through and draw each frequency’s circle. To do that we need the radius and phase which we have included in the DFTData type in the first code block. To do the actual graphics representation, we’ll be using the p5.js library. Note that the code depends on variables defined in the above code blocks. I’ll also provide the full source at the end.

function draw() {
 for (let i = 0; i M.pi * 2)
 points.pop();
}

time += ((2 * M.pi) / pts.length);
}

これで記事冒頭のアニメーションが描ける。

Fourier アニメーションの完全なソースコードはこちら：

adekau/fourier — フーリエ変換で描く。

](https://github.com/adekau/fourier)

Index Notes on Classical & Intuistionistic Logic