科斯塔斯陣列、有限域、聲納和最難聽的音樂_風聞

**撰文 |**蔣迅

來源：和樂數學

如果我們能創作出一段沒有任何重複的音樂，沒有旋律，沒有模式，沒有比例。那麼得到的就是一首非常難聽的音樂。

從上面的視頻中，我們看到了陳奕迅把歌唱到了最難聽的境地。陳奕迅能做到這一點真的是不容易。但對數學家來説，創作出一個最難聽的歌曲並不是特別難的事情。今天我就來講講數學家是怎樣做到的。不過先要從一個看似無關的數學遊戲説起。由此我們介紹數學在通訊科學中的一個應用，然後以最難聽的音樂結束本文。

建議在閲讀這一篇之前，先閲讀筆者對哥隆尺的介紹。這樣可能對本文的思想有些幫助。不過，本文並不要求讀者具有哥隆尺的知識。

01 科斯塔斯陣列的定義

我們考慮在平面上的n×n的網格。在這些網格中，我們將放入 n 個圓點，使得在每一行和每一列上都只有一個圓點。我們可以用 0 和 1 來代替這些點：有圓點的方格上我們放入 1，否則就放入 0。滿足上述描述的最簡單的網格就是對角網格。下面是兩個

表 1

再看一個非對角網格。

表 2使用 0 和 1 表示的的好處是便於文字敍述。我們可以把一個n×n的網格看成是一個矩陣（或者稱為陣列）：從左到右分別為第 1, 2, …,n列；從上到下分別為第 1, 2, …,n行。在第i行第j列的數字就可以用ai,j或a(i,j) 來表示，即ai,j= 1 如果相應的格子裏是 1，否則ai,j= 0。注意，我們在這裏是從上往下為遞增方向，因為這個表示法與數學上的矩陣的一般表達是一致的：

我們希望使用這種矩陣的表達形式，但我們並不要求讀者有矩陣的知識。所以我們情願把它叫作陣列。

對於上面的對角網格來説，四個圓點的座標是：(1,1),(2,2),(3,3),(4,4)；另一個非對角的網格的圓點座標則是 (3,1),(4,2),(2,3),(1,4)。用矩陣的表示法，上面兩個例子分別是

在矩陣語言裏，這樣的矩陣叫作置換矩陣，因為每一個這樣的矩陣都可以通過一系列行與行之間的互換而最終變成對角矩陣。作為置換矩陣，我們也可以把它寫成數學上的“排列”。

比如上面的對角矩陣可以寫成一個平凡排列：[1,2,3,4]；非對角矩陣可以寫成 [3,4,2,1]，它是平凡排列 [1,2,3,4] 的一個排列。像我們介紹的哥隆尺，我們可以從排列的表示來構造倒三角。當初科斯塔斯就是從這個角度來構造最初的幾個科斯塔斯陣列的。我們不深入討論。

我們現在可以給出科斯塔斯陣列的定義了。數學上，科斯塔斯陣列首先是在平面上的滿足上述條件的n×n的網格上的n個點。每兩個圓點的連線形成一個線段。我們也可以把它們看作是向量，只是沒有一個明確的方向。這樣的向量一共有

個。當n = 4 時，容易看出這樣的線段一共有 6 個。科斯塔斯陣列就是要求那些向量都不相同（不平行或者不同長度）。

在上面的兩個例子中，第一個對角陣列顯然不是一個科斯塔斯陣列，比如 (1,1) 和 (2,2) 形成的向量和由 (2,2) 和 (3,3) 形成的向量的方向和長度是一樣的，它們都是與水平線構成 45 度角且長度為√2 的線段。第二個非對角陣列是一個科斯塔斯陣列。讀者可以找出它的全部

個向量並驗證它們都是不同的。如果我們把這裏的網格比作刻度尺的話，那麼科斯塔斯陣列就可以比作哥隆尺。所以我們可以把科斯塔斯陣列看作是哥隆尺在二維的推廣。

02 科斯塔斯陣列的一個等價定義

科斯塔斯陣列有一個等價的定義，而這個定義能幫助我們理解科斯塔斯陣列在應用中的意義。為此，我們對陣列A擴展如下：

也就是將陣列A往四個方向用 0 無限延伸。我們稱之為A的擴展陣列。我們定義陣列A的非週期自相關函數（autocorrelation function）為

對擴展陣列A’=(ai,j)，我們把一個陣列向右移動u單位，同時向下移動v單位，得到的是A’u,v=(ai+v,j+u). 顯然，CA就是把陣列A’和它的 (u,v) 平移陣列A’u,v在它們的相同座標上的值兩兩相乘，然後相加這些乘積。這個數越大，那麼它的自相關就越高。如果 (u,v)=(0,0)，那麼A’與A’u,vCA(u,v)=n. 注意到當 |u| ≥n或 |v| ≥n時，A’u,v已經移到A’可能取 1 的範圍之外，因而這時必有CA(u,v)= 0. 所以我們只需考慮 |u| <n或 |v| <n的情況。也就是説，我們可以把CA看作是一個 (2n-1)×(2n-1) 的陣列。回到我們前面的一個對角陣列和非對角陣列的例子上，為説話方便，我們把它們分別記為A和B。這時，n = 4. 經簡單計算，我們可以得到下面兩個 7×7 陣列：

注意座標 (0,0) 在這兩個陣列的中心，而且在座標 (0,0) 上的值都是 4。比較兩個例子，我們還可以看到，CA中對角線上的值是中心向兩邊逐漸遞減的，而CB上除了在中心有一個 4 以外，其他的值都是 0 和 1。事實上，我們可以證明滿足這個性質的置換陣列就一定是科斯塔斯陣列，反之亦然。容易看出，當n充分大時，科斯塔斯陣列只有在原點 (0,0) 的自相關函數值達到n。對其他平移後的點，其自相關函數值都非常小；而對角陣列則不具備這個性質。

在本節的最後，我們順便定義兩個n×n陣列A和B的互相關函數（cross-correlation function）：

03 科斯塔斯和已知的一些科斯塔斯陣列

圖 1. 約翰・科斯塔斯

約翰・科斯塔斯（JohnP. Costas）是美國電子工程師。1947 年，他從普渡大學畢業。這時正值第二次世界大戰。他加入了美國海軍，成了一名雷達軍官。後來他進入麻省理工學院，研究干擾濾波和線性系統編碼。在學校裏，他得到了美國應用數學家諾伯特・維納（NorbertWiener）、意裔美籍計算機科學家羅伯特・法諾（Robert Mario Fano）、美國電子工程學家傑羅姆・威斯納（Jerome Bert Wiesner）和中國現代早期電機工程學家李鬱榮。從 1951 年開始，他受僱於通用電氣公司。1980 年代初，他轉到 Cogent Systems 公司直到退休。科斯塔斯最著名的成就還不是科斯塔斯陣列，而是他在 1950 年代發明的對現代數字通信產生深遠影響的科斯塔斯循環（Costasloop）。進入 1960 年代後，他解決了聲納系統效果不佳的問題。他的解就是本文的科斯塔斯陣列。我們將在稍後做詳細介紹。

找到科斯塔斯陣列比找到哥隆尺相對容易一些，因為在二維上自由度比一維大一些。科斯塔斯在 1975 年用手在一張紙上湊出了一個 12×12 的科斯塔斯陣列。由於他無法找到更大的例子，他懷疑對n > 12 可能不存在這樣的陣列。但後來人們發現了一些算法，可以得到任意大的科斯塔斯陣列。下面的表格給出前 36 階的科斯塔斯陣列的數量表。

表3. 科斯塔斯陣列一覽表

目前，從n = 1 到 29 的所有科斯塔斯陣列都已經找到。在 29 之後還沒有一個n得到了全部科斯塔斯陣列。我們用斜體字表示我們得到的是科斯塔斯陣列的數目的下限。特別讓人們意外的是，至今人們還沒有找到n= 32 和 33 時哪怕一個科斯塔斯陣列。人們也無法證明不存在n = 32 和 33 時的科斯塔斯陣列。甚至有人估算説，當n = 32 時，用現在已知的算法和現有的設備，找到全部科斯塔斯陣列需要 45000 年的計算機時間！所以至今是否對任意的正整數n來説都存在科斯塔斯陣列還是一個未解之謎。人們仍然在努力尋找新的算法。我們不打算囊括全部這些算法，而只是介紹一下最早的兩個算法。這兩個構造法都有很強的數學背景。

04 有限域和原根

這裏我們要説到的數學背景有一段悲催的故事。這個故事始於 200 年前的法國。一位年輕人伽羅瓦（Evariste Galois）為了解決困擾五次多項式的根式解問題另闢奇徑，搞出來一個所有當時的大數學家都無法理解的思路。政治上，他強烈支持共和，受到保皇勢力的打壓。更奇特的是，在他 21 歲的時候為一位女子與人決鬥。決鬥前，他匆忙寫下了他的數學成果，然後委託他的朋友務必找到出版的地方。他的稿子沒有得到高斯（JohannKarl Friedrich Gauss）、雅可比（Carl Gustav Jacob Jacobi）的賞識。所幸的是，這個稿子後來被劉維爾（Joseph Liouville）的肯定，終於發展成為了數學界的豐碑“伽羅瓦理論”。他的故事已經出現許多數學科普作品中。我們後面會看到卡斯塔斯陣列在通訊中的應用。估計伽羅瓦沒有想到的是，他的數學成果能在二百多年後應用到通訊領域。

在數學上，實數的全體構成一個“域”。所謂域，就是一個代數結構，在它裏面可以進行加、減、乘和除（除數不為零）運算。比如説兩個實數相加還是實數，兩個實數相除也還是實數，只要除數不是零。運算結果仍然保留在這個代數結構裏。我們看到，域的概念是數域以及四則運算的推廣。

實數域是一個無限域。但並不是所有的域都是無限的。有限域也被稱為伽羅瓦域（Galois field），很顯然是為了紀念這位偉大的法國數學家伽羅瓦。有限域就是具有加減乘除運算的包含有限個元素的域。有限域最常見的例子是當p為素數時，整數對p取模。我們把它記為 GF(p)。它的元素可以用 0,…,p-1 表示。我們假定對這些元素做四則運算時在取模的意義下進行。也就是説，一旦一個運算結果達到了p，就將這個數歸零；運算結果超過了p時就把這個數減去p，直到其結果落入到 0 到 p-1 的範圍之內。這種運算叫作模運算（mod），一般用“≡”表示，但是在代數學裏也會用“=”表示，只要不會引起混淆。以 GF(7) 為例，它的元素為 0,…, 6。我們有 1 + 4 ≡ 5(mod 7)，4 + 5 = 9 ≡ 2(mod 7)。有了模運算後，我們就可以引入歐拉函數的概念。定義歐拉函數φ(p) 為與p互素並小於或等於p的那些正整數的個數。在我們考慮的例子中，p為素數，所以總是有φ(p)=p-1。注意歐拉函數並不只是對素數定義的。在後面的科斯塔斯陣列的算法中也會用到更一般的整數的歐拉函數。

我們還需要引入數論中原根的概念。對於兩個互素的正整數g和p，數論中有歐拉定理保證必存在正整數d≤p-1，使得gd≡ 1(modp)。我們不妨假設d是滿足上述條件中的最小的那個正整數。如果這個d滿足φ(p)=d，那麼基數g叫作模p的原根。我們以p= 7 為例，顯然有φ(p)= 6。我們説g= 3 是一個原根，因為

而當基數為 2 時，23= 8 ≡ 1(mod 7)，但是指數 3 < 6 =φ(p)。從而 2 不是模 7 的一個原根。我們看到，3k模 7 的週期為 6。在這個週期裏，模 7 後的餘數是 3, 2, 6, 4, 5, 1。這正好是 1, 2, 3, 4, 5, 6 的一個置換。這樣做出的置換是一個科斯塔斯陣列。現在我們可以介紹科斯塔斯陣列的構造法了。至今科斯塔斯陣列的構造法仍然是一個活躍的科研領域。但限於篇幅，我們只介紹兩個最早出現的方法：衞曲構造法和藍波-哥倫布構造法。

05 勞埃德・衞曲和衞曲構造法

圖2. 勞埃德・衞曲

衞曲陣列是最早的一個構造方法。其實，這個方法是由美國數學家和編碼專家埃德加・吉爾伯特（Edgar Gilbert）在 1965 年，即科斯塔斯發現科斯塔斯陣列的同時發現的。當然這個時候他並不知道科斯塔斯的工作。所以他的出發點是不同的：拉丁方陣 (Latin square)。可惜的是，科斯塔斯發現了科斯塔斯陣列但是沒有開發一套計算方法，而吉爾伯特設計了一個巧妙的方法卻不知道科斯塔斯陣列。由於他們兩人沒有出現交集而使得吉爾伯特的工作被埋沒了。一直到 1982 年，吉爾伯特的構造法才重新被美國應用數學家和信息科學家勞埃德・衞曲（Lloyd R. Welch）發現。

衞曲於 1951 年畢業於伊利諾伊大學厄巴納-香檳分校數學系，又在 1958 年從加州理工學院獲得數學博士學位。他的博士論文的題目是：卷積積分的排序和最大化。比較自相關函數和卷積，我們可以感覺到他在讀書的時候就已經為他以後的工作鋪墊了道路。畢業後他在噴氣推進實驗室、國防分析研究所和南加州大學工作。衞曲的主要貢獻是尋找隱馬爾可夫模型的未知參數的鮑姆-衞曲算法（Baum-Welch algorithm）和一種用於高效地解碼 BCH 碼與裏德-所羅門碼的伯利坎普－衞曲算法（Berlekamp-Welch algorithm）。

衞曲並沒有投入到科斯塔斯陣列的研究。原來科斯塔斯在用紙筆反覆湊答案失敗之後，他轉向哥倫布求救。這是 1970 年代後期的事情。哥倫布一方面告訴科斯塔斯，這個問題以前還沒有人研究過（其實他是不知道吉爾伯特的工作）；他同時向他在南加州大學的同事和合作者衞曲詢問。這兩個人早就在算法學上長期合作。早在 1968 年，他們就共同提出了在編碼理論裏至今未解決的“哥倫布-衞曲”猜想，而且這方面的工作就與伽羅瓦域緊密相關。衞曲意識到了科斯塔斯的問題可以用伽羅瓦域的結果來做。這就是衞曲構造法。1982 年，哥倫布在和赫伯特・泰勒合寫的一篇關於科斯塔斯陣列的評述中首次介紹了這個算法。隨後哥倫布給出了嚴格的證明。

衞曲構造法：設p是一個素數，g是其原根。令A= [ai,j] 為一個p-1 階的置換陣列（矩陣），滿足下列條件：

上面我們舉了一個p= 7 例子。我們已經看到 GF(7) 的原根g= 3，而且我們得到了一個由模餘數形成的置換 3, 2, 6, 4, 5, 1。可以驗證相應的 6 階陣列A就是：

這個科斯塔斯陣列就是一個用衞曲構造法生成的。我們可以把它記作 [1 3 2 6 4 5]。細心的讀者會注意到，這個陣列不是想象中的 [3 2 6 4 5 1]。其實，衞曲構造法還可以稍微推廣一點。假定c是一個整數。令A= [ai,j] 為一個p-1 階的陣列（矩陣），滿足下列條件：

顯然，前面的定義是當c= 0 時的特例。如果c= 1，那麼可以看到

而這正是我們預期的陣列 [3 2 6 4 5 1]。這兩陣列的區別僅僅是在水平方向的一個平移。所以從本質上説它們是等價的。

06 藍波-哥倫布構造法

第二個早期科斯塔斯陣列構造法也是從有限域出發的。這就是哥倫布介紹的藍波-哥倫布構造法（Lempel-Golomb construction）。

亞伯拉罕・藍波（Abraham Lempel）出生于波蘭。他分別於 1963 年、1965 年和 1967 年從以色列理工學院獲得學士、碩士和博士學位。然後，他作為研究助理前往南加州大學。在那裏開始了與哥倫布的長期合作。1969 年，他加入了位於馬薩諸塞州薩德伯裏的斯佩裏蘭德研究中心任研究員。1971 年，他回到以色列理工學院，在那裏擔任計算機科學教授直到退休。他同時還擔任 IBM 沃森研究中心的職務。他最著名的工作是在無損數據壓縮算法的兩個算法“LZ77 與 LZ78”，而且這個工作也是基於有限域的性質。在那段時間裏，他在有限域方面有許多研究。順便提一句，另有一位叫作泰瑞・衞曲（Terry Welch）的美國計算機學家把無損數據壓縮算法做了改進，這個新算法稱為“藍波-立夫-衞曲（Lempel-Ziv-Welch）編碼法”。

藍波的算法也是哥倫布與泰勒 1982 年的同一篇論文中首先披露的。稍後哥倫布完成了證明。哥倫布和泰勒介紹的是一個推廣了的藍波構造法。我們在這裏也先介紹藍波-哥倫布構造法，然後作為一個特例給出藍波構造法。

藍波-哥倫布構造法：設p為一個素數，n為一個正整數，q=pn> 2 為p的一個冪。再設φ和ρ為有限域 GF(q) 的兩個（可能相同的）原根。令A= [ai,j] 為一個q-2 階的置換陣列（矩陣），滿足下列條件：

當φ=ρ時就是藍波提出的原始情形。讓我們看一個例子：令p= 11,n= 1，從而q=pn= 11。計算可知φ= 2，ρ= 7 是 GF(11) 的兩個原根：

於是，使用藍波-哥倫布構造法，我們可以得到下面的置換陣列：

表 4

我們把計算留給讀者。

07 科斯塔斯陣列在聲納的應用

前面我們説過，科斯塔斯是在研究聲納時發現的科斯塔斯陣列的。現在我們就來談談聲納與科斯塔斯陣列的關係。

故事還是從科斯塔斯開始。從麻省理工大學博士畢業後，他受命為海軍解決用聲納偵查潛艇效率不高的問題。由於電磁波在水中衰減的速率非常的高，無法做為偵測的訊號來源，以聲波探測水面下的人造物體成為運用最廣泛的手段。無論是潛艇或者是水面船隻，都利用這項技術的衍生系統，探測水底下的物體，或者是以其作為導航的依據。聲納的工作原理是它發出音響訊號，藉由這個訊號接觸物體後反射回來聲波的變化，做為計算這個物體的相對方位與距離的資料。這種方法就是利用了多普勒效應。

多普勒效應（Doppler effect）是波源和觀察者有相對運動時，觀察者接受到波的頻率與波源發出的頻率並不相同的現象。可能你早就注意到，遠方急駛過來的火車鳴笛聲變得尖細（即頻率變高，波長變短），而離我們而去的火車鳴笛聲變得低沉（即頻率變低，波長變長）。這就是多普勒效應的現象。這一現象最初是由奧地利物理學家多普勒（Christian Doppler）1842 年發現的。

在實際應用中，人們一般是在連續的幾個相同的時間段裏發出不同頻率的聲波系列，然後收集反射回來的聲波。當收到的一個從移動物體反射回來的迴音時，這個系統會與它發出的音頻系列各個時間和頻率上的平移做比較，從中找到與這個反射回來的信號高度吻合的那個時間和頻率的平移。從時間的平移人們可以計算出物體的距離範圍；從頻率的變化，人們可以計算出物體移動的速度。

這個系列可以用n×n的置換陣列 [ai,j] 來表示，其中 j 對應於時間區間t1,t2, …,tn，i對應於頻率f1,f2, …,fn，滿足：ai,j= 1 當且僅當頻率fi在時間段tj裏發射。當這個設定好的頻率系列發出去後，接收器稍後收到一組迴音。這組迴音的頻率會因物體的運動速度而有些變化；而發射時間與接受時間之間的時間差則決定於發射器與反射回音的物體之間的距離。

如果沒有任何雜音的話，這個結果應該就很好了。但在實際應用中雜音是避免不了的。所以我們必須能夠區分背景雜音和目標物體反射的迴音。為此，我們必須將收到的信號與發射信號的所有 (2n-1)(2n-1) 個組合一一比較，並希望在這些組合中只有那個對應於目標物體的位置和速度的平移與其有較高的互相關性。因此，在設計這組頻率信號時，我們必須讓我們置換陣列使得其所有的非平凡位移（即零位移）都與其本身具有低相關性。

我們希望以上的討論已經把科斯塔斯的思想展現出來了。科斯塔斯陣列在通訊上的研究至今活躍，在中國也是同樣。我們不打算更深入地從電子工程的角度談科斯塔斯陣列是如何實施的。但是想指出中國學者在這方面也有一些工作。我們在參考文獻裏引用了幾篇，比如周彥鵬和景曉軍對 FH-OFDMA通信系統的簡介。

08 最難聽的音樂

現在讓我們回到文章開頭最難聽的音樂。我不知道陳奕迅是如果即興創作出一首那麼難聽的音樂，但我猜測他的目標是打破觀眾的常規期待─優美的樂曲。那麼什麼是優美的音樂呢？我想重複是一個關鍵。畢達哥拉斯早在 2500 年前就發現音調之間的比例。一首樂曲有旋律，有主題。通過旋律和主題的重複。比如在貝多芬的第五交響曲中的著名的“噠吶吶吶”四音符動機音型主題在交響樂中僅在第一樂章裏就有數百次。所以這種重複的設置對於美麗來説非常重要。那麼，如果我們能創作出一段沒有任何重複的音樂，沒有旋律，沒有模式，沒有比例。那麼得到的就是一首非常難聽的音樂。

這個思想其實早在 20 世紀 30-50 年代開始就有人嘗試過。著名的奧地利裔美國作曲家和音樂理論家阿諾德・勳伯格（Arnold Schoenberg，1874－1951）提出“解放不和諧”的思想，希望能從音調結構中釋放音樂。可惜他在科斯塔斯解決了如何在數學上創造這些結構的問題之前 10 年就去世了。

經過了上面對科斯塔斯陣列的討論，我們應該不難想到切入點就是科斯塔斯陣列。鋼琴上有 88 個鍵，從最左邊的低音 A（啦）到最右邊的高音 C（哆）。我們可以把科斯塔斯在聲納中的頻率換成鋼琴的琴鍵。正好p= 89 是一個素數。上面的討論讓我們知道可以構造一個 (p-1)×(p-1)= 88×88 的科斯塔斯陣列。又g= 3 是 GF(89) 的一個原根。所以我們可以運用衞曲構造法來做。簡單計算得到：

從這個計算，我們得到下面的科斯塔斯陣列：

把它換成樂譜就是世界上首個無模式鋼琴奏鳴曲（Costas Golomb No. 1: The Perfect Ping）：

在這個曲子中，88 個琴鍵的每一個都只彈奏一次，並且是按上述科斯塔斯陣列的順序。

這段音樂是美國數學家、愛爾蘭都柏林大學學院工學院高級講師斯科特・裏卡德（Scott Rickard）為我們創作的。裏卡德在麻省理工學院獲得數學和計算機與工程學的兩個本科學位（1992 年和 1993 年）和電子工程和計算機科學的碩士學位，又在普林斯頓大學獲得應用數學和計算數學的碩士和博士學位（2000 年和 2003 年）。他在 2011 年的 TED 大會上介紹了這首樂曲。如果你能看到他的演講視頻，那麼就可以欣賞由新世界交響樂團室內音樂系主任邁克爾・林維（MichaelLinville）演奏的鋼琴獨奏。我不知道讀者會怎樣比較陳奕迅和林維的表演，但我們知道林維演奏的樂曲只有數學家可以創作出來。

本文經授權轉載自微信公眾號“和樂數學”，發表在《中國工業與應用數學學會通訊》2019 年第 3 期。

特 別 提 示

1. 進入『返樸』微信公眾號底部菜單“精品專欄“，可查閲不同主題系列科普文章。

2. 『返樸』提供按月檢索文章功能。關注公眾號，回覆四位數組成的年份+月份，如“1903”，可獲取2019年3月的文章索引，以此類推。