Question

Sequence

Sequence，在英文当中是“一连串”的意思，在数学当中则是“一连串数字”的意思，中文译作“数列”。

像这样就是一个数列：

4 -1 6 0 9

我们习惯用一维阵列储存一个数列：

数学家谈数列，是谈等差数列、等比数列、数列求和、递迴数列。计算学家谈数列，则是谈积分、微分、排序、搜寻。

积分（前缀和）

从头开始、连续数字和。

原本数列  4 -1  6  0  9
积分之后  4  3  9  9 18

UVa 10324 10994 983

微分（差分）

相邻数字差。

原本数列  4 -1  6  0  9
微分之后  4 -5  7 -6  9

Sequence（Under Construction!）

数列与多项式互相转换（z-transform）

(2 -5 1 0 4)  <--->  2x⁰ - 5x¹ + x² + 0x³ + 4x⁴
(a₀ a₁ a₂ a₃ a₄)  <--->  a₀x⁰ + a₁x¹ + a₂x² + a₃x³ + a₄x⁴

最后来个趣味题目，两步猜出多项式的係数：

http://www.matrix67.com/blog/archives/4064

加减乘除模

点积

摺积

Number Theoretic Transform

数论转换是一对一函数，输入和输出都是一串馀数。

e 的纯虚数次方会不断绕圈，馀数的次方也会不断绕圈。于是有人把 e^i⋅2π/N改成原根的次方。

http://www.apfloat.org/prim.html

输出入每个数值都是馀数，皆大于等于零、小于模数。当输入数值很大，那麽模数也必须足够大。

数论转换的特色是完全没有浮点数误差！

Polynomial Multiplication（Under Construction!）

多项式加法与减法 = 数列加法与减法

多项式加法
        a₀ x⁰ +       a₁ x¹ +       a₂ x²
+)      b₀ x⁰ +       b₁ x¹ +       b₂ x²
—————————————————————————————————————————
   (a₀+b₀) x⁰ +  (a₁+b₁) x¹ +  (a₂+b₂) x²

数列加法
(a₀ a₁ a₂) + (b₀ b₁ b₂) = (a₀+b₀ a₁+b₁ a₂+b₂)

多项式乘法与除法 = 数列摺积与反摺积

多项式乘法，是数列的什麽呢？数学家以不断位移的点积，兜出多项式乘法的结果，然后命名为摺积。

多项式乘法
                         a₀ x⁰ +   a₁ x¹ +   a₂ x²
×)                       b₀ x⁰ +   b₁ x¹ +   b₂ x²
——————————————————————————————————————————————————
                       a₀b₂ x² + a₁b₂ x³ + a₂b₂ x⁴
             a₀b₁ x¹ + a₁b₁ x² + a₂b₁ x³
   a₀b₀ x⁰ + a₁b₀ x¹ + a₂b₀ x²
——————————————————————————————————————————————————
     c₀ x⁰ +   c₁ x¹ +   c₂ x² +   c₃ x³ +   c₄ x⁴

其中
c₀ = a₀b₀
c₁ = a₀b₁ + a₁b₀
c₂ = a₀b₂ + a₁b₁ + a₂b₀
c₃ = a₁b₂ + a₂b₁
c₄ = a₂b₂

数列摺积
(a₀ a₁ a₂) convolution (b₀ b₁ b₂) = (c₀ c₁ c₂ c₃ c₄)

其中
c₀ = (-  -  a₀ a₁ a₂ -  - ) dot (b₂ b₁ b₀ -  -  -  - )
c₁ = (-  -  a₀ a₁ a₂ -  - ) dot (-  b₂ b₁ b₀ -  -  - )
c₂ = (-  -  a₀ a₁ a₂ -  - ) dot (-  -  b₂ b₁ b₀ -  - )
c₃ = (-  -  a₀ a₁ a₂ -  - ) dot (-  -  -  b₂ b₁ b₀ - )
c₄ = (-  -  a₀ a₁ a₂ -  - ) dot (-  -  -  -  b₂ b₁ b₀)

对齐一下
c₀:                     c₁:                     c₂:
(-  -  a₀ a₁ a₂ -  - )  (-  -  a₀ a₁ a₂ -  - )  (-  -  a₀ a₁ a₂ -  - )
(b₂ b₁ b₀ -  -  -  - )  (-  b₂ b₁ b₀ -  -  - )  (-  -  b₂ b₁ b₀ -  - )

c₃:                     c₄:
(-  -  a₀ a₁ a₂ -  - )  (-  -  a₀ a₁ a₂ -  - )
(-  -  -  b₂ b₁ b₀ - )  (-  -  -  -  b₂ b₁ b₀)

Convolution

等长的两串数列，做“对应项乘法”运算，得到一串数列。时间複杂度 O(N)。

Pointwise Product:
                                   a₀  a₁  a₂  
(a₀ a₁ a₂) multiply (b₀ b₁ b₂) = ( ×   ×   ×  ) = (a₀b₀ a₁b₁ a₂b₂)
                                   b₀  b₁  b₂  

等长的两串数列，做“点积”运算，得到一个值：对应项相乘后求和。时间複杂度 O(N)。

Dot Product:
                            a₀   a₁   a₂
(a₀ a₁ a₂) dot (b₀ b₁ b₂) = ×  + ×  + ×  = a₀b₀ + a₁b₁ + a₂b₂
                            b₀   b₁   b₂

两串数列，做“摺积”运算，得到一个数列：第二串数列头尾颠倒，穷举各种对齐方式，各做一次点积。至于第二串数列为何要头尾颠倒？正是为了迎合多项式乘法！时间複杂度 O(AB)。

Convolution:

(a₀ a₁ a₂ a₃ a₄) convolution (b₀ b₁ b₂) = (c₀ c₁ c₂ c₃ c₄ c₅ c₆)

c₀ = (-  -  a₀ a₁ a₂ a₃ a₄ -  - ) dot (b₂ b₁ b₀ -  -  -  -  -  - )
c₁ = (-  -  a₀ a₁ a₂ a₃ a₄ -  - ) dot (-  b₂ b₁ b₀ -  -  -  -  - )
c₂ = (-  -  a₀ a₁ a₂ a₃ a₄ -  - ) dot (-  -  b₂ b₁ b₀ -  -  -  - )
c₃ = (-  -  a₀ a₁ a₂ a₃ a₄ -  - ) dot (-  -  -  b₂ b₁ b₀ -  -  - )
c₄ = (-  -  a₀ a₁ a₂ a₃ a₄ -  - ) dot (-  -  -  -  b₂ b₁ b₀ -  - )
c₅ = (-  -  a₀ a₁ a₂ a₃ a₄ -  - ) dot (-  -  -  -  -  b₂ b₁ b₀ - )
c₆ = (-  -  a₀ a₁ a₂ a₃ a₄ -  - ) dot (-  -  -  -  -  -  b₂ b₁ b₀)

摺积是很多次点积，第二串头尾颠倒，各种位移。

多项式乘法有交换率、结合率、分配律，当然摺积也有。

Deconvolution

“反摺积”就是摺积的反运算，解摺积式子。

☐ convolution b = c   --->   c deconvoution b = ☐

化作矩阵乘法，求反矩阵即可。

http://dsp.stackexchange.com/questions/15096/ 。

但是我不清楚如何求得摺积的反元素。

多项式循环乘法 = 数列循环摺积

紧接著引入“循环”的概念！多项式相乘，刻意让次方循环。

多项式循环乘法
                         a₀ x⁰ +   a₁ x¹ +   a₂ x²
×)                       b₀ x⁰ +   b₁ x¹ +   b₂ x²
—————————————————————————————————————————————————————
                       a₀b₂ x² + a₁b₂ x³ + a₂b₂ x⁴
             a₀b₁ x¹ + a₁b₁ x² + a₂b₁ x³
   a₀b₀ x⁰ + a₁b₀ x¹ + a₂b₀ x²
—————————————————————————————————————————————————————
   a₁b₂ x⁰ + a₂b₂ x¹ + a₀b₂ x²
   a₂b₁ x⁰ + a₀b₁ x¹ + a₁b₁ x²
   a₀b₀ x⁰ + a₁b₀ x¹ + a₂b₀ x²
—————————————————————————————————————————————————————
     c₀ x⁰ +   c₁ x¹ +   c₂ x²

其中
c₀ = a₁b₂ + a₂b₁ + a₀b₀
c₁ = a₂b₁ + a₀b₁ + a₁b₁
c₂ = a₀b₀ + a₁b₀ + a₂b₀

也就是
(a₀ a₁ a₂)      (a₀ a₁ a₂)      (a₀ a₁ a₂)
    dot    x⁰ +     dot    x¹ +     dot    x²
(b₀ b₂ b₁)      (b₁ b₀ b₂)      (b₂ b₁ b₀)

数列循环摺积
(a₀ a₁ a₂) circular convolution (b₀ b₁ b₂) = (c₀ c₁ c₂)

其中
c₀ = (a₀ a₁ a₂) dot (b₀ b₂ b₁)
c₁ = (a₀ a₁ a₂) dot (b₁ b₀ b₂)
c₂ = (a₀ a₁ a₂) dot (b₂ b₁ b₀)

Circular Convolution

等长的两串数列，做“循环摺积”运算，得到一串数列：与摺积大同小异，超出数列的部分，改成头尾循环。时间複杂度 O(N^2)。

第二串数列由第 0 项到第 N-1 项轮流作为首项、头尾颠倒、头尾循环。

                  circular
(a₀ a₁ a₂ a₃ a₄) convolution (b₀ b₁ b₂ b₃ b₄) = (c₀ c₁ c₂ c₃ c₄)

c₀ = (a₀ a₁ a₂ a₃ a₄) dot (b₀ b₄ b₃ b₂ b₁)
c₁ = (a₀ a₁ a₂ a₃ a₄) dot (b₁ b₀ b₄ b₃ b₂)
c₂ = (a₀ a₁ a₂ a₃ a₄) dot (b₂ b₁ b₀ b₄ b₃)
c₃ = (a₀ a₁ a₂ a₃ a₄) dot (b₃ b₂ b₁ b₀ b₄)
c₄ = (a₀ a₁ a₂ a₃ a₄) dot (b₄ b₃ b₂ b₁ b₀)

c0:               c1:               c2:
(a₀ a₁ a₂ a₃ a₄)  (a₀ a₁ a₂ a₃ a₄)  (a₀ a₁ a₂ a₃ a₄)
(b₀ b₄ b₃ b₂ b₁)  (b₁ b₀ b₄ b₃ b₂)  (b₂ b₁ b₀ b₄ b₃)

c3:               c4:
(a₀ a₁ a₂ a₃ a₄)  (a₀ a₁ a₂ a₃ a₄)
(b₃ b₂ b₁ b₀ b₄)  (b₄ b₃ b₂ b₁ b₀)

Convolution Theorem【尚无正式名称】

多项式与数列互相转换。

a₀x⁰ + a₁x¹ + a₂x²  <--->  (a₀ a₁ a₂)

数列转换成多项式。

(a₀ a₁ a₂)  ---->  a₀x⁰ + a₁x¹ + a₂x²

数列转换成多项式，可以看成点积，可以看成线性变换。

线性变换矩阵
A = [ x⁰ x¹ x² ]

数列
    [ a₀ ]
a = [ a₁ ]
    [ a₂ ]

数列转换成多项式
                  [ a₀ ]
Aa = [ x⁰ x¹ x² ] [ a₁ ] = a₀x⁰ + a₁x¹ + a₂x²
                  [ a₂ ]

这种线性变换有个特性：“变换前的摺积 = 变换后的乘法”。

A = [ x⁰ x¹ x² ]

    [ a₀ ]       [ b₀ ]
a = [ a₁ ]   b = [ b₁ ]
    [ a₂ ]       [ b₂ ]

p = Aa = (a₀x⁰ + a₁x¹ + a₂x²)
q = Ab = (b₀x⁰ + b₁x¹ + b₂x²)

pq = Aa Ab = (a₀x⁰ + a₁x¹ + a₂x²) (b₀x⁰ + b₁x¹ + b₂x²)
           = (a convolution b) 然后添上 [ x⁰ x¹ x² x³ x⁴ ]
           = A(a convolution b)    a 与 b 末端补零，A 末端补项。

次方值乘上任意倍率，此特性一样成立。

A = [ x⁰ x¹ x² ]

A = [ x⁰ x⁵ x¹⁰ ]

A = [ x⁰ x⁻⁷ x⁻²¹ ]

A = [ x⁰ x⁰ x⁰ ]

所有元素一齐乘上任意倍率，此特性一样成立。

A = [ 7x⁰ 7x¹ 7x² ]

A = [ -x⁰ -x⁵ -x¹⁰ ]

A = [ 0 0 0 ]

从一个横列推广到一个矩阵，此特性一样成立。

    [ 7x⁰ 7x¹  7x²   ]
A = [ -x⁰ -x⁵  -x¹⁰  ]
    [  x⁰  x⁻⁷  x⁻²¹ ]

Aa multiply Ab = A(a convolution b)

结论就是：

Cicular Convolution Theorem【尚无正式名称】

接下来继续补强矩阵，除了满足“变换前的摺积 = 变换后的乘法”，也要满足“变换前的乘法 = 变换后的摺积”！

从数学来看，补强性质，达成了对称之美；从计算学来看，补强限制，极可能产生特殊演算法。

那麽就得让 A 拥有反矩阵 A⁻¹，而且 A 和 A⁻¹都具备上个段落提到的特性。一种尝试是正规正交矩阵：A⁻¹ = Aᵀ，前述特性变成了行与列同时成立，容易观察。

实数系统下，x 次方渐增，数值越来越大，导致基底不可能垂直（内积不可能为零）。很不幸的，这种矩阵不存在。

             [ x⁰ x⁰ x⁰ x⁰ .. ]
     -1   T  [ x⁰ x¹ x² x³ .. ]
A = A  = A = [ x⁰ x² x⁴ x⁶ .. ]
             [ x⁰ x³ x⁶ x⁹ .. ]
             [ :  :  :  :     ]

Aa multiply Ab = A(a convolution b)
Aa convolution Ab = A(a multiply b)

折衷的方式是令 x 的次方产生循环，让数值能够变小，使得基底互相垂直、内积是零。（用数学术语来说：从 open set 改成 closed set。）

複数系统下，把 x 设定成 e^i⋅2π/N（或其倒数 e^-i⋅2π/N），令 x 的次方产生循环。此即“ 傅立叶转换 ”。

Fourier Transform: x = e-i⋅2π/N; k = 1/sqrt(N); N is size of matrix

      [ x⁰ x⁰ x⁰ x⁰ .. ]          [ x⁻⁰ x⁻⁰ x⁻⁰ x⁻⁰ .. ]
      [ x⁰ x¹ x² x³ .. ]    -1    [ x⁻⁰ x⁻¹ x⁻² x⁻³ .. ]
A = k [ x⁰ x² x⁴ x⁶ .. ]   A  = k [ x⁻⁰ x⁻² x⁻⁴ x⁻⁶ .. ]
      [ x⁰ x³ x⁶ x⁹ .. ]          [ x⁻⁰ x⁻³ x⁻⁶ x⁻⁹ .. ]
      [ :  :  :  :     ]          [ :   :   :   :      ]

Aa multiply Ab = A(a circular convolution b)
Aa circular convolution Ab = A(a multiply b)

馀数系统下，则是把 x 设定成任意一个原根（或其倒数），令 x 的次方产生循环。此即“数论转换”。

Number Theoretic Transform: x = primitive root (mod n)

    [ x⁰ x⁰ x⁰ x⁰ .. ]        [ x⁻⁰ x⁻⁰ x⁻⁰ x⁻⁰ .. ]
    [ x⁰ x¹ x² x³ .. ]    -1  [ x⁻⁰ x⁻¹ x⁻² x⁻³ .. ]
A = [ x⁰ x² x⁴ x⁶ .. ]   A  = [ x⁻⁰ x⁻² x⁻⁴ x⁻⁶ .. ]
    [ x⁰ x³ x⁶ x⁹ .. ]        [ x⁻⁰ x⁻³ x⁻⁶ x⁻⁹ .. ]
    [ :  :  :  :     ]        [ :   :   :   :      ]

Aa multiply Ab = A(a circular convolution b)
Aa circular convolution Ab = A(a multiply b)

原本的特性，也变得循环：“变换前的循环摺积 = 变换后的乘法”、“变换前的乘法 = 变换后的循环摺积”！

引入时域与频域的观念：“频域的乘法 = 时域的循环摺积”、“频域的循环摺积 = 时域的乘法”！

Circular Convolution 的快速演算法
Polynomial Circular Multiplication 的快速演算法

正向傅立叶转换、数论转换需时 O(NlogN)，对应项相乘需时 O(N)，逆向傅立叶转换、数论转换需时 O(NlogN)。总时间複杂度为 O(NlogN)。

傅立叶转换的弱点是浮点数误差。数论转换的弱点是数值大小不得超过模数大小。

注意到，快速傅立叶转换、快速数论转换，数列长度必须是 2 的次方。当数列长度不是 2 的次方，千万不能直接补零到 2 的次方。

正确的循环摺积计算结果：
            circular
(a₀ a₁ a₂) convolution (b₀ b₁ b₂) = (c₀ c₁ c₂)

c₀ = a₀b₀ + a₁b₂ + a₂b₁
c₁ = a₀b₁ + a₁b₀ + a₂b₂
c₂ = a₀b₂ + a₁b₁ + a₂b₀

长度补到 2 的次方，计算结果完全不同：
              circular
(a₀ a₁ a₂ 0) convolution (b₀ b₁ b₂ 0) = (d₀ d₁ d₂ d₃)

d₀ = a₀b₀ + a₂b₂
d₁ = a₀b₁ + a₁b₀
d₂ = a₀b₂ + a₁b₁ + a₂b₀
d₃ = a₁b₂ + a₂b₁

正确的方式：先补零直到不受循环影响，再补零直到长度是 2 的次方，最后让输出数列循环。

想要计算      circular
  (a₀ a₁ a₂) convolution (b₀ b₁ b₂) = (c₀ c₁ c₂)

改为计算          circular
　(a₀ a₁ a₂ 0 0) convolution (b₀ b₁ b₂ 0 0) = (d₀ d₁ d₂ d₃ d₄)

长度补到 2 的次方         circular
　(a₀ a₁ a₂ 0 0 0 0 0) convolution (b₀ b₁ b₂ 0 0 0 0 0)
= (d₀ d₁ d₂ d₃ d₄ 0 0 0)

最后让输出数列循环
  c₀ = d₀ + d₃
  c₁ = d₁ + d₄
  c₂ = d₂

Convolution 的快速演算法
Polynomial Multiplication 的快速演算法

运用循环摺积，计算摺积。

想要计算
  (a₀ a₁ a₂ a₃) convolution (b₀ b₁ b₂)
= (c₀ c₁ c₂ c₃ c₄ c₅)

改为计算             circular
  (a₀ a₁ a₂ a₃ 0 0) convolution (b₀ b₁ b₂ 0 0 0)
= (c₀ c₁ c₂ c₃ c₄ c₅ 0)

长度补到 2 的次方           circular
  (a₀ a₁ a₂ a₃ a₄ 0 0 0) convolution (b₀ b₁ b₂ 0 0 0 0 0)
= (c₀ c₁ c₂ c₃ c₄ c₅ 0 0)

截断输出数列至正确长度
  (c₀ c₁ c₂ c₃ c₄ c₅ 0 0) -> (c₀ c₁ c₂ c₃ c₄ c₅)

范例：大数乘法

大数乘法即是多项式乘法！

傅立叶转换、数论转换得以计算大数乘法，时间複杂度从 O(N^2) 降为 O(NlogN)。

Polynomial Factorization（Under Construction!）

多项式与函数点互相转换（point-value presentation）

根据“ 多项式内插 ”的“唯一解定理”，N 项多项式，等同 N 个相异函数点。

N 项多项式的加减法，等同 N 个函数点的加减法。

N 项与 M 项多项式的乘法（次方变成 N+M-1），等同 N+M-1 个函数点的乘法。

N 项与 M 项多项式的除法（次方变成 N-M+1），若整除，则 N-M+1 个函数点皆整除。若不整除，则无法处理。

Cicular Convolution Theorem【尚无正式名称】

这个版本比较容易记诵。教科书採用这个版本。

傅立叶矩阵碰巧是 Vandermonde Matrix！可援引多项式内插！

正向傅立叶转换，可以看成多项式求值：已知 N 个多项式係数，x 代入下述数值，求得 N 个函数值。

逆向傅立叶转换，可以看成多项式内插：已知 N 个函数值，x 等于下述数值，求得 N 个多项式係数。

傅立叶转换：x 代入 e-i⋅(2π/N)⋅0、e-i⋅(2π/N)⋅1、……、e-i⋅(2π/N)⋅(N-1)。
数论转换：x 代入 r⁰、r¹、……、rN-1。

傅立叶矩阵是正规正交矩阵、一对一函数，唯一解定理成立！因此，N 项多项式的循环乘法，等同 N 个函数点的乘法。

Polynomial Factorization

N 项多项式恰有 N-1 个根，可能有重根，可能是虚数。

取 N-1 个根当作函数点！但是少了一个函数点，况且无法处理重根。

“因式分解”。

“质式 Irreducible Polynomial”

“生成器 Primitive Polynomial”

演算法（Kronecker's Algorithm）

原理是“ Lagrange Interpolation ”。

      N-1         x  - rj
p(x) = ∑  di  ∏  ---------
      i=0    j≠i  ri - rj

两个多项式整除，所有对应点也会整除。

【待补文字】

Zero 与 Pole

多项式摺积与反摺积 = 数列？？与？？

似乎没人研究。有兴趣可参考“ Normalized B-splines ”。

二项式定理、多项式定理

巴斯卡三角形、排列组合、摺积。

Toeplitz Matrix（Under Construction!）

Circulant Matrix

“循环矩阵”。每行每列皆循环的矩阵。

[ 1 5 3 7 2 ]
[ 2 1 5 3 7 ]
[ 7 2 1 5 3 ]
[ 3 7 2 1 5 ]
[ 5 3 7 2 1 ]

“循环摺积”与“循环矩阵”完全等价！

[ b₀ b₄ b₃ b₂ b₁ ] [ a₀ ]   [ c₀ ]
[ b₁ b₀ b₄ b₃ b₂ ] [ a₁ ]   [ c₁ ]
[ b₂ b₁ b₀ b₄ b₃ ] [ a₂ ] = [ c₂ ]
[ b₃ b₂ b₁ b₀ b₄ ] [ a₃ ]   [ c₃ ]
[ b₄ b₃ b₂ b₁ b₀ ] [ a₄ ]   [ c₄ ]

求值、求解、乘法、反矩阵，运用傅立叶转换或数论转换，时间複杂度为 O(NlogN)。

Toeplitz Matrix（Diagonal-constant Matrix）

“常对角矩阵”。每一条左上右下斜线，都是同一个元素的矩阵。以下只讨论方阵。

“循环矩阵”是“常对角矩阵”的特例。

[ 1 5 3 2 1 ]
[ 2 1 5 3 2 ]
[ 7 2 1 5 3 ]
[ 4 7 2 1 5 ]
[ 8 4 7 2 1 ]

“常对角矩阵”与“摺积”不等价。“摺积”可以化作“常对角矩阵”，反之则不然，除非矩阵无限大。

                               [  0 ]
[ b₂ b₁ b₀  0  0  0  0  0  0 ] [  0 ]   [ c₀ ]
[  0 b₂ b₁ b₀  0  0  0  0  0 ] [ a₀ ]   [ c₁ ]
[  0  0 b₂ b₁ b₀  0  0  0  0 ] [ a₁ ]   [ c₂ ]
[  0  0  0 b₂ b₁ b₀  0  0  0 ] [ a₂ ] = [ c₃ ]
[  0  0  0  0 b₂ b₁ b₀  0  0 ] [ a₃ ]   [ c₄ ]
[  0  0  0  0  0 b₂ b₁ b₀  0 ] [ a₄ ]   [ c₅ ]
[  0  0  0  0  0  0 b₂ b₁ b₀ ] [  0 ]   [ c₆ ]
                               [  0 ]

[ b₀  0  0  0  0 ]          [ c₀ ]
[ b₁ b₀  0  0  0 ] [ a₀ ]   [ c₁ ]
[ b₂ b₁ b₀  0  0 ] [ a₁ ]   [ c₂ ]
[  0 b₂ b₁ b₀  0 ] [ a₂ ] = [ c₃ ]
[  0  0 b₂ b₁ b₀ ] [ a₃ ]   [ c₄ ]
[  0  0  0 b₂ b₁ ] [ a₄ ]   [ c₅ ]
[  0  0  0  0 b₂ ]          [ c₆ ]

常对角矩阵的演算法较为複杂，以下开始介绍。

加法

两个 Toeplitz Matrix 相加，时间複杂度 O(N)。

2N-1 个数字就能记录一个 Toeplitz Matrix。

乘法

两个 Toeplitz Matrix 相乘，时间複杂度 O(N^2)。

http://stackoverflow.com/questions/15889521/

乘法的反运算：反矩阵

时间複杂度 O(N^2)。比高斯乔登消去法 O(N^3) 快。

Toeplitz Matrix 的反矩阵，通常不是 Toeplitz Matrix，但是可以表示成 Circulant Matrix 与上三角 Toeplitz Matrix 的乘积。

http://www.sciencedirect.com/science/article/pii/S0893965907000535

求值

时间複杂度 O(NlogN)。比普通矩阵求值 O(N^2) 快。

填充元素成为循环矩阵，化为循环摺积，套用傅立叶转换、数论转换。最后从中抽取正确答案。

                     [ 1 5 3 7 2 ] [ 1 ]
[ 1 5 3 ] [ 1 ]      [ 2 1 5 3 7 ] [ 2 ]      (1 2 3 0 0)
[ 2 1 5 ] [ 2 ] ---> [ 7 2 1 5 3 ] [ 3 ] --->      ⍟
[ 7 2 1 ] [ 3 ]      [ 3 7 2 1 5 ] [ 0 ]      (1 2 7 3 5)
                     [ 5 3 7 2 1 ] [ 0 ]
Toeplitz matrix      circulant matrix      circular convolution

求值的反运算：求解

时间複杂度 O(N^2)。比高斯消去法 O(N^3) 快。

https://en.wikipedia.org/wiki/Levinson_recursion

另外还有时间複杂度 O(N(logN)^2) 的演算法。我没有涉猎。

Toeplitz Decompoisition

任意方阵都可以化为一连串 Toeplitz Matrix 的乘积。

http://www.stat.uchicago.edu/~lekheng/work/toeplitz.pdf

Eigenvalue 与 Fourier Transform

循环矩阵，特徵向量是傅立叶矩阵，特徵值是第一个横条的傅立叶转换。

http://math.stackexchange.com/questions/25126/

1. 两个循环矩阵相乘  = 总是可以用 fourier matrix 作对角化 = 对角线矩阵相乘
　 (还是循环矩阵)        (C = F D F⁻¹) (F⁻¹ = Fᵀ)　      (对应的 eigenvalue 相乘)

2. 两个数列的循环摺积     = fft                         = 对应项相乘

3. 两个多项式的循环乘法   = 多项式求值，以 e^-itn 取样　　= 对应的点座标相乘

4. 两个分解式的循环乘法   = 这是甚麽东西？
   (2n 个根融合成 n 个根)

tridiagonal and Toeplitz matrix's eigenvalues:
a + 2 sqrt(bc) cos(k pi / (n+1))   for k = 1~n

eigenvectors of fourier matrix is gaussian function
eigenvalues of fourier matrix is +1 -1 +i -i
F^4 = I

Convolution（Under Construction!）

Dirichlet Convolution

convolution: 加法分解
Dirichlet convolution: 乘法分解

Dyadic Convolution

http://book.huihoo.com/pdf/algorithms-for-programmers/afp.pdf   p362 $9.6
http://iris.elf.stuba.sk/JEEEC/data/pdf/09-10_102-10.pdf

Hadamard matrix

http://en.wikipedia.org/wiki/Hadamard_matrix