關于位移水仙花數的研究

楊健，韋佳玉，蔣劍軍

（銅陵學院數學與計算機學院，銅陵244061）

0 引言

整數世界的神秘讓人自小就會產生無窮的求知欲望。水仙花數是不少人初中時代播下的一顆種子，到大學時代對它的認知徹底清晰起來。所謂水仙花數是指一個3位數，它各位數字的3次方之和等于該數自身[1]。利用枚舉法或計算機編程，知水仙花數共有四個：153、370、371 和 407。

尋找水仙花數，用數學語言表達就是求解下述關于x,y,z的不定方程：

其中x,y,z都是介于0到9的整數，且x≠0。

有了（1）式的數學語言，初中時代播下的求知種子，到大學時代發了芽。方程（1）是三位整數的情況，它在n位的情形是怎么樣的呢？即求下述關于(n;x1,x2,…,xn)的不定方程的解（或全部解）：

其中xi(1,2…,n)都是介于 0到 9的整數，且x1≠0。

滿足方程（2）的n位正整數也稱為水仙花數，它是經典的三位水仙花數往高位整數上的拓展。

眾所周知，方程（2）至少有解（3;1，5，3）、（3;3，7，0）、（3;3，7，1）和（3;4，0，7）分別對應著上面所列的 4 個經典的水仙花數。當n＞3呢？查閱文獻知，2002年林宣治[2]及2015年衛洪春[3]分別求出了方程（2）的所有解，即得到了所有的水仙花數，其中最大的水仙花數有39 位：115132219018763992565095597973971522401。

因為水仙花數的趣味性，以及在計算機程序設計教學中的特殊地位，近年來依然有學者將它作為各種計算機語言程序設計教學中的范例[4-8]，或做著水仙花數的科普工作[9]。

文獻[2]及文獻[3]對水仙花數的徹底解決并沒有阻擋我們對水仙花數未知世界的探索。本文進一步將水仙花數推廣為更為廣闊的概念——位移水仙花數，研究該類水仙花數的性質，并設計基于MATLAB矩陣運算的快速算法求解滿足一定條件的所有的位移水仙花數。

1 位移水仙花數

1.1 定義

一個n位正整數若等于它各個位上數字的n+d次方之和，則稱該n位正整數為指數位移d的n位水仙花數，簡稱位移d水仙花數。此時稱d為位移。設是一個位移d水仙花數，則下述（3）式成立：

其中xi(1,2…n)都是介于0到9的整數且x1≠0。

例如，4150=45+15+55+05，故4150是位移d=1水仙花數；再如，194979=15+95+45+95+75+95，故194979是位移d=-1的水仙花數。

當d＞0,d=0,d＜0時分別稱為正位移、0位移、負位移的水仙花數。易知，0位移水仙花數就是一般意義下的水仙花數。

1.2 性質

因為位移水仙花數比水仙花數多一個位移參數，所以具有豐富的性質。

性質1【有限性】設d為一給定整數，則位移d水仙花數的個數是有限的。

易得：

上式兩邊取以10為底的對數，得：

將（4）式整理為如下形式：

時位移d的n位水仙花數是不存在的，故位移d的水仙花數的個數是有限的。

性質2【n和d的關系性質】位移水仙花數中，位移d與位數n滿足如下關系式：

證明（6）式的左邊不等式已由（5）式給出。下面證明（6）式的右邊不等式。

對上式右邊放大后易得：

按上述放大所得的數n·9n+d，其位數至多增加1，即n·9n+d至多是n+1位整數，即：

上式兩邊取以10為底的對數，得：

這就完成了性質2的證明。

性質2意義非常豐富。首先，（6）式可視化如下；然后，由（6）式得到d的最小值為-1；最后，由（6）式，應用MATLAB編程計算給定d對應的n值的范圍如表1所示。

表1 給定d對應的n值的范圍

注意，應用（6）式計算得到的n的上界和下界，并非n的上確界和下確界。例如，當d=0時，n的下確界是 1，上確界是 39（見文獻[1-2]）。

圖1

在描述性質3之前，先提出并證明下面的引理。為方便描述下述引理和性質3，我們引入如下記號：

引理 1【7整除1(n)的判別】（1）；（2）當6|n時，7||1(n)。

證明：

（1）簡單計算知，當n=1,2,3,4,5時，7?1(n)；當n=6時，7|1(n)。

充分性：設6|n，可寫為n=6k。則：

上式右端的k個加項都能被7整除，故能被7整除。

必要性：設7|1(n)。

設n=6k+r,0＜r≤6，下面證明r=6。

將1(n)改寫為：

上式右端能被7整除，故左端中必有7|1(r)，從而得r=6。

（2）因為1(6)=3×7×11×13×17，即 7||1(6)，所以當6|n時，有

從而得7||1(n)。

注：引理1說明，7至多是1(n)的單因子。

同理可得下述兩個引理。

引理2【3恰好整除1(n)的判別】3||1(n)?3||n。

注：引理2說明，當且僅當3||n時，3是1(n)的單因子。

引理 3【9 整除1(n)的判別】（1）；（2）當9|n時，9||1(n)。

注：引理3說明，9也至多是1(n)的單因子。

性質3【位異性】位數大于1的位移型水仙花數，至少存在兩個位上的數字是互異的。

反證法。假設存在整數d(≥-1)，使m(n)得是位移d的n位水仙花數，即：

注意到，m(n)=m·1(n)，結合上式得：

顯然，（8）式中m≠1,5且m不能為偶數。于是m只可能是 3、7、9。

經檢驗，（8）式右端n+d-1≥2，故若m是 3、7、9之一，則必然有m是1(n)的至少2重因子，這矛盾于上面的三個引理。

性質3說明，位數大于1的各位數字相同的正整數不會是位移型水仙花數，從而不會是水仙花數。這是下述性質4的基礎。

性質4【惟一性】設某個n位的位移d水仙花數從高位到低位的數字依次是x1,x2,…,xn，則這n個數字的其他任意異于x1x2…xn的排列所成的n位正整數都不是位移d水仙花數。

證明：

設x1,x2,…,xn是位移d的n位水仙花數從高位到低位的數字，xi1,xi2,…,xin（其中xi1≠0）是x1,x2,…,xn的異于x1,x2,…,xn的一個排列，令：

因為xi1,xi2,…,xin和x1,x2,…,xn互異的兩個排列，必然有A≠B。

因為A是位移d的水仙花數，故：

故B不是位移d的水仙花數。

3 算法實現

3.1 算法設計

計算泛水仙花數，計算量超大是其特點。所以在設計基于MATLAB快速計算泛水仙花數的算法方面，須充分考慮MATLAB的計算特點。眾所周知，MAT?LAB計算特點是數據以矩陣為單元，循環效率相對較低。所以在設計算法時必須將輸入、輸出和存儲數據都矩陣化，降低對循環的依賴，以提高效率。本文正是基于MATLAB擅長矩陣運算的特點設計算法以計算位移型水仙花數。算法描述如下。

算法：基于MATLAB的計算位數不超過len、位移為d的所有水仙花數的矩陣算法

毋庸諱言，本算法有很大的局限，并不能應對任意位數的水仙花數的求解。

3.2 計算結果

由d與n的關系圖知，d=-1的位移水仙花數最高位不超過34位，實際最大值為8，一共有2個：194979，14459929；d=1的位移水仙花數最高位不超過84位，因為個人電腦算力的原因，本文計算了10位以內的所有d=1的位移水仙花數，獲得兩個該類水仙花數：4150，4151；d=2的位移水仙花數位數最小為58，最大為108，本文沒能計算得到該類水仙花數；更大的位移已沒有計算能力應對。

表2 計算結果

4 水仙花數的分類

現在我們將位移型水仙花數簡稱為水仙花數，它包括了一般意義上的水仙花數（即3位經典水仙花數和高位水仙花數），那么對水仙花數的描述有2個特征：位數n和位移d；進而對水仙花數的分類就可按位數n分類或按位移d分類。那么，按哪個特征分類更好呢？我們認為按位移分類更好，因為，簡單的枚舉計算就知道n=2時沒有的任何位移的水仙花數。目前計算得到的結果按d分類如表3。

即位移-1的水仙花數有2個，位移0的水仙花數有88個，位移1的水仙花數最高位不超過84，計算得到10位以內有2個。需要說明的是，由于1的特殊性，我們不將1視為任何位移的水仙花數。

表3 水仙花數的分類

5 結語

本文首先將水仙花數推廣為更一般的位移水仙花數，使得水仙花數成為我們研究范疇的一個特例。為了實際上討論和描述的方便，本文可以從一開始就稱呼位移水仙花數為水仙花數。然后對水仙花數的性質進行了比較細致的研究，得到的性質既是水仙花數的特征，也十分有趣。最后對水仙花數的分類進行了一些討論。在對水仙花數分類的過程中，不禁產生了一些疑惑，其中最大的疑惑就是：

問題1是否存在任意位移的水仙花數？

顯然，無論計算機技術發展到什么程度都無法從計算的角度解決上述問題。于是，我們把上述問題等價地轉化為如下數學問題：

問題2下述關于(n;d,x1,x2,…,xn)的不定方程：

的解是否是有限的？其中0≤xk≤9,k=1,2,…,n，且x1≠0。

期待數學工作者給出問題2的解答。