楊健，韋佳玉，蔣劍軍

（銅陵學(xué)院數(shù)學(xué)與計(jì)算機(jī)學(xué)院，銅陵244061）

0 引言

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

尋找水仙花數(shù)，用數(shù)學(xué)語言表達(dá)就是求解下述關(guān)于x,y,z的不定方程：

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

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

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

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

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

因?yàn)樗苫〝?shù)的趣味性，以及在計(jì)算機(jī)程序設(shè)計(jì)教學(xué)中的特殊地位，近年來依然有學(xué)者將它作為各種計(jì)算機(jī)語言程序設(shè)計(jì)教學(xué)中的范例[4-8]，或做著水仙花數(shù)的科普工作[9]。

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

1 位移水仙花數(shù)

1.1 定義

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

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

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

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

1.2 性質(zhì)

因?yàn)槲灰扑苫〝?shù)比水仙花數(shù)多一個(gè)位移參數(shù)，所以具有豐富的性質(zhì)。

性質(zhì)1【有限性】設(shè)d為一給定整數(shù)，則位移d水仙花數(shù)的個(gè)數(shù)是有限的。

易得：

上式兩邊取以10為底的對(duì)數(shù)，得：

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

時(shí)位移d的n位水仙花數(shù)是不存在的，故位移d的水仙花數(shù)的個(gè)數(shù)是有限的。

性質(zhì)2【n和d的關(guān)系性質(zhì)】位移水仙花數(shù)中，位移d與位數(shù)n滿足如下關(guān)系式：

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

對(duì)上式右邊放大后易得：

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

上式兩邊取以10為底的對(duì)數(shù)，得：

這就完成了性質(zhì)2的證明。

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

表1 給定d對(duì)應(yīng)的n值的范圍

注意，應(yīng)用（6）式計(jì)算得到的n的上界和下界，并非n的上確界和下確界。例如，當(dāng)d=0時(shí)，n的下確界是 1，上確界是 39（見文獻(xiàn)[1-2]）。

圖1

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

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

證明：

（1）簡單計(jì)算知，當(dāng)n=1,2,3,4,5時(shí)，7?1(n)；當(dāng)n=6時(shí)，7|1(n)。

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

上式右端的k個(gè)加項(xiàng)都能被7整除，故能被7整除。

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

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

將1(n)改寫為：

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

（2）因?yàn)?(6)=3×7×11×13×17，即 7||1(6)，所以當(dāng)6|n時(shí)，有

從而得7||1(n)。

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

同理可得下述兩個(gè)引理。

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

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

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

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

性質(zhì)3【位異性】位數(shù)大于1的位移型水仙花數(shù)，至少存在兩個(gè)位上的數(shù)字是互異的。

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

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

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

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

性質(zhì)3說明，位數(shù)大于1的各位數(shù)字相同的正整數(shù)不會(huì)是位移型水仙花數(shù)，從而不會(huì)是水仙花數(shù)。這是下述性質(zhì)4的基礎(chǔ)。

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

證明：

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

因?yàn)閤i1,xi2,…,xin和x1,x2,…,xn互異的兩個(gè)排列，必然有A≠B。

因?yàn)锳是位移d的水仙花數(shù)，故：

故B不是位移d的水仙花數(shù)。

3 算法實(shí)現(xiàn)

3.1 算法設(shè)計(jì)

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

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

毋庸諱言，本算法有很大的局限，并不能應(yīng)對(duì)任意位數(shù)的水仙花數(shù)的求解。

3.2 計(jì)算結(jié)果

由d與n的關(guān)系圖知，d=-1的位移水仙花數(shù)最高位不超過34位，實(shí)際最大值為8，一共有2個(gè)：194979，14459929；d=1的位移水仙花數(shù)最高位不超過84位，因?yàn)閭€(gè)人電腦算力的原因，本文計(jì)算了10位以內(nèi)的所有d=1的位移水仙花數(shù)，獲得兩個(gè)該類水仙花數(shù)：4150，4151；d=2的位移水仙花數(shù)位數(shù)最小為58，最大為108，本文沒能計(jì)算得到該類水仙花數(shù)；更大的位移已沒有計(jì)算能力應(yīng)對(duì)。

表2 計(jì)算結(jié)果

4 水仙花數(shù)的分類

現(xiàn)在我們將位移型水仙花數(shù)簡稱為水仙花數(shù)，它包括了一般意義上的水仙花數(shù)（即3位經(jīng)典水仙花數(shù)和高位水仙花數(shù)），那么對(duì)水仙花數(shù)的描述有2個(gè)特征：位數(shù)n和位移d；進(jìn)而對(duì)水仙花數(shù)的分類就可按位數(shù)n分類或按位移d分類。那么，按哪個(gè)特征分類更好呢？我們認(rèn)為按位移分類更好，因?yàn)椋唵蔚拿杜e計(jì)算就知道n=2時(shí)沒有的任何位移的水仙花數(shù)。目前計(jì)算得到的結(jié)果按d分類如表3。

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

表3 水仙花數(shù)的分類

5 結(jié)語

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

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

顯然，無論計(jì)算機(jī)技術(shù)發(fā)展到什么程度都無法從計(jì)算的角度解決上述問題。于是，我們把上述問題等價(jià)地轉(zhuǎn)化為如下數(shù)學(xué)問題：

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

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

期待數(shù)學(xué)工作者給出問題2的解答。