基于MATLAB的賴特函數分區算法研究及實現

收稿日期：2023-08-28

DOI：10.19850/j.cnki.2096-4706.2024.05.001

摘" 要：物理數學中的多數特殊函數可將復平面劃分區域采用不同數值技術來計算。賴特函數在分數微積分及其工程應用中有著重要作用，作為一類新型特殊函數也可使用分區算法進行數值計算。通過研究賴特函數在大參數下的漸近展開式和公式中系數的計算方法，修正積分表達式及積分半徑選擇定理的錯誤，進一步改進完善復數域賴特函數的分區算法，并利用MATLAB軟件進行編程仿真分析算法精度。實驗結果表明，分區算法的適用性廣，有良好逼近效果。

關鍵詞：分數微積分；特殊函數；漸近展開；分區算法；MATLAB

中圖分類號：TP301.6" 文獻標識碼：A" 文章編號：2096-4706（2024）05-0001-06

Research and Implementation of the Partitioning Algorithm for the Wright Function Based on MATLAB

LI Yan， YUAN Xiao

（College of Electronics and Information Engineering， Sichuan University， Chengdu" 610065， China）

Abstract： Most special functions in physical mathematics can be calculated by dividing the complex plane into different regions and using various numerical techniques. The Wright function plays an important role in fractional calculus and its engineering applications. And as a new type of special function， the Wright function can also be calculated using the partitioning algorithm. By studying the asymptotic expansion of the Wright function under large parameters and the calculation method of coefficients in the formula， correcting the errors in the integral expression and the integral radius selection theorem， the partitioning algorithm of the Wright function in the complex field has been further improved and perfected， and finally MATLAB software is used for programming simulation to analyze the accuracy of the algorithm. The experimental results show that the partitioning algorithm has wide applicability and good approximation effect.

Keywords： fractional calculus; special function; asymptotic expansion; partitioning algorithm; MATLAB

0" 引" 言

指數函數" 常出現在整數階常微分方程與偏微分方程的解析中，在應用數學及概率論中有著核心作用。與經典微積分理論相似，新型特殊函數米塔-列夫勒函數[1-3]與賴特函數等在分數微積分中扮演著類似指數函數的角色，它們都是指數函數的推廣，并繼承了指數函數的一些顯著特性。

1933年，Wright在研究分區漸近理論[4]時首次提出參數α＞0的賴特函數。1935年，Wright在研究廣義超幾何函數的漸近行為[5]時引入了福克斯-賴特函數，其定義為：

（1）

在同一論文中，他還研究了α＞0，γ＞0時的四參數賴特函數：

（2）

1940年，Wright在論文中[6]從數學角度再次研究了賴特函數Wα， β （z），并將參數α的范圍擴展到α＞-1的實數域，故賴特函數的完整定義為：

（3）

賴特函數Wα， β （z）根據α的不同取值分為兩種類型，分別為α＞0時的第一類賴特函數和-1＜α＜0時的第二類賴特函數。

隨著分數微積分理論體系的完善和拓展，賴特函數的研究涵蓋了不同領域，包括理論物理學[7，8]、概率論[9]等。在信號處理方面，賴特函數作為新型分數階微積分初等函數之一，它的導數可以處理不完全可微和高階可微信號，例如：在語音識別、圖像識別和模式識別領域常常需要處理復雜信號，利用賴特函數可以提高信號處理精度。此外，機器學習也逐漸成為該函數應用的一個熱門領域。科學家使用神經網絡等技術來預測與分析數據中存在的整體行為和趨勢，并且使用賴特函數作為一個強大工具來解決這些問題。賴特函數的應用日益增加，研究其性質及數值算法十分重要。

賴特函數的數值算法研究目前還處在一個探索階段。2008年，Luchko提出了分區算法[10]，將實數域分成三個區域利用不同數值技術來計算具有實參數的賴特函數，但這種算法的積分表達式及積分半徑選擇定理存在問題。2010年，他修正了該算法的積分半徑選擇定理中的錯誤，可該算法只被限制在負實數域[11]上，同時缺少對大參數下賴特函數的漸近級數的研究，只能作為賴特函數數值計算的參考。本文是在Luchko的分區算法基礎上，進一步研究賴特函數的漸近展開，詳細論述了其漸近展開公式中系數的計算方法，同時對積分表達式中的錯誤進行修改并推廣至復數域賴特函數的積分計算，最后利用MATLAB編寫代碼分析算法精度。

1" 賴特函數的漸近展開

1940年，Wright首次利用應用于積分表示的最速下降法研究賴特函數的漸近展開性質[6]，但并未考慮到斯托克斯現象。1999年，趙育求等人研究得到賴特函數Wα， β （z）的漸近公式[12，13]，并研究了斯托克斯現象，但不全面。2010年，Paris研究了福克斯-賴特函數" 的漸近展開[14]并在2014年的研究中考慮到該函數的斯托克斯現象[15]。2020年，Paris針對賴特函數的z參數取絕對值很大的實數及β取很大值的情況，建立了新的輔助函數研究其漸近展開[16]。2021年，Paris等人研究得到了兩個賴特型輔助函數Fv （z）和Mv （z）（即Mainardi函數）的漸近展開[17]。

第一類賴特函數的漸近展開可直接利用文獻[15]中福克斯-賴特函數" 的漸近展開原理得到。參數-1＜α＜0的賴特函數的漸近展開不能直接采用" 的漸近理論結果，需利用伽馬函數的反射公式來求得。表1列舉出賴特函數的漸近展開公式，其中系數Cj和Dj的計算將在第2節詳細描述。

2" 漸近展開中系數計算方法

2.1" 系數Cj的計算方法

賴特函數漸近展開式中的系數Cj可依據超幾何函數" 漸近公式的系數計算方法[14]來確定。系數Cj出現g（s） / s！的逆階乘展開式中，歸一化系cj = Cj /C0根據該展開式可改寫為：

（4）

其中 。將縮放伽馬函數[18]中z替換為αs + β得到：

（5）

將式（4）左邊式子利用式（5）通過一些常規代換可以寫為

（6）

根據式（4）與（6）可以得到：

（7）

將式（7）左側的乘積寫成" 的逆冪展開：

（8）

因為" 存在展開[19]：

（9）

故式（7）的右側式子隨著" 可化為：

（10）

對比式（8）和（10），可發現當1≤j≤M-1時有Hj = Xj，于是有：

（11）

式（11）為歸一化系數cj的計算公式，主要需要計算參數Hj與" 的值。其中" 是廣義伯努利多項式[20]。廣義伯努利數定義為 ，可以由遞推關系給出：

前3個廣義伯努利數為：

廣義伯努利多項式可以寫為：

（12）

計算參數Hj首先需將式（5）中" 與 進行逆冪展開。 有：

（13）

最后利用指數函數的泰勒展開即可求得" 的逆冪展開系數。根據縮放伽馬函數與其倒數的逆冪展開式[21]可以得到" 與" 的逆冪展開：

其中γk為斯特林系數[22，23]，可利用遞推關系精確求出。

賴特函數的漸近展開會用到p = 0，q = 1和p = 1，q = 0兩種情況，分別有：

其中σ = -α，δ = 1-β。兩種情況的算法原理及編程相似，現針對p = 0，q = 1進行具體描述：

Y0， 1（s）的系數可以利用MATLAB中conv函數算出。

y011=conv（gammaxi（alf，bet，n1，2），gammaxi（1，1，n1，2））;%gammaxi（alf，bet，n1，2）計算" 的系數；

%gammaxi（1，1，n1，2）計算" 的系數；

y01=conv（y011，gammaxi（k，q1，n1，1））;%gammaxi（k，q1，n1，1）計算" 的系數；

R0， 1（s）可以寫為：

R0， 1（s）的系數可以利用MATLAB中symsum函數與taylor函數算出：

Q=symsum（（-1）^ii/（ii*（ii+1））*（q1^ii*（-q1+（ii+1）/2）/k^ii-b^ii*（-b+（ii+1）/2）/a^ii-（ii-1）/2）*（s）^ii，ii，1，n1）;

R011=exp（Q）;

R01=taylor（R011，s，'Order'，n1）;

R01=sym2poly（R01）;

將Y0， 1（s）與R0， 1（s）的系數進行卷積，即可得到Hj的值。

2.2" 系數Dj的計算方法

賴特函數的漸近展開主要研究σ = 1/2（μ = 1）的情況，根據文獻[15]可得到：

（14）

由式（14）可看出，系數Dj主要需計算Gk， j（1/2）參數。Gk， j（1/2）出現在展開式：

（15）

其中λ j = v - j，v = m0 - X - 2β + 3/2，v選擇使m0為整數的最小數，X = x2/2。

對式（15）中" 映射進行級數反演時，有：

（16）

更高系數g （k）可利用Mathematic軟件中的InverseSeries函數獲得。根據式（16）有：

（17）

利用式（15）與（17）進行系數卷積即可求得Gk， j（1/2）的值，更高系數可利用MATLAB求得。

GA1=（1+RT）^RJ;%GA1為 ；

GA=taylor（GA1，m，'Order'，n1）;%對" 進行泰勒展開；

GA1=sym2poly（GA）;%獲得多項式系數；

GG=conv（GA1，CRT1）;%CRT1為" 的冪展開系數，兩者卷積即可求得" 的相反數；

3" 分區算法

將復平面分為三個區域：A）| z |≤q1，0＜q1＜1，B）q1＜| z |≤q2，C）| z |＞q2，分別使用泰勒級數、積分表達式及漸近展開式來計算賴特函數。其中漸近展開在第1～2節已做具體論述，但無法從理論上預測評估精度。

3.1" 泰勒級數

首先考慮A區域| z |≤q1，采用泰勒級數進行計算。

定理3.1，當| z |≤q1，0＜q1＜1，對于規定精度ε＞0有：

（18）

其中截斷點滿足：

3.2" 積分表達

在區域B中，采用賴特函數的積分表達式定理3.1進行計算。被積函數" 在積分區間上有界且積分區間也有界，故積分是有限的，可以直接通過MATLAB中integral函數以規定的精度ε＞0來計算。被積函數K涉及無窮上限，需對積分半徑進行截斷，得到積分半徑選擇定理3.2。

定理3.1一般情況下，賴特函數Wα， β （z）的積分表達式為：

（19）

當-1＜α＜0且β＜1，有：

（20）

當-1＜α＜0且β = 1，有：

（21）

其中：

（22）

定理3.2當" 時，設z = x + iy，積分表達式

1）當" 時：

2）當" 時：

4" 實驗仿真

本節利用MATLAB對賴特函數的分區算法進行編程仿真，主要考察算法近似值與函數真實值的相對誤差η = | （W - W真） / W真|×100%。取特例：

（23）

利用分區算法分別繪制 、 與真實值之間的相對誤差，如圖1、2所示。其中，分區算法在| z |≤0.95內采用泰勒級數計算，在0.95＜

| z |≤6內采用積分表達式計算，當| z |＞6時采用截斷指數M = 10的漸近展開式計算。由圖可觀察：分區算法計算" 時在| z |較小時使用積分表達式能取得高達1×10-13數量級的精度，但隨著| z |不斷增大后使用漸近展開式精度有所下降，但仍能保持1×10-11數量級的精度；對于 ，分區算法在整個矩形區域內至少取得1×10-12數量級的精度，對于大部分取值其精度能達到1×10-14數量級。

圖1" "與函數真實值的相對誤差

圖2" "與函數真實值的相對誤差

5" 結" 論

分數階微積分由于具有時間記憶性和全局相關性被廣泛應用于信號與信息處理、分數階圖像處理、控制系統等工程領域中來描述復雜現象。賴特函數作為一類新型特殊函數在分數微積分中有著重要地位，研究其性質與數值算法對于工程應用很有必要。本文將復平面劃分三個區域，采用泰勒級數、積分表達和漸近展開數值技術來計算賴特函數，并基于MATLAB進行編程仿真。實驗結果表明：分區算法的適用性廣，能計算" 范圍內的參數；在計算精度方面，分區算法會受到參數α，β取值的影響而呈現不一樣的精度，但均能達到較高精度；分區算法的積分表達式近似值在| z |較小時能有較高精度，但當賴特函數的值小于設置的精度時，積分表達式結果的誤差會不斷增加，漸近展開式不會受此限制，能夠較精確地得到無限靠近0的賴特函數值。

分區算法研究仍存在許多值得探索的問題，比如：

1）算法中的區域劃分參數q2及漸近展開式中的截斷參數M的最佳取值還需要繼續研究。

2）算法中積分表達式的準確性跟設置的誤差容限與數值積分方式有很大聯系，故還需研究最佳數值積分算法及誤差容限參數設置。

3）由于數值積分及漸近展開需計算系數，該算法的運算速度較慢，如何提高算法運算效率還需進一步研究。

4）該算法僅實現賴特函數， 范圍內的計算如何將范圍擴展到" 仍待研究。

參考文獻：

[1] GORENFLO R，KILBAS A A，MAINARDI F，et al. Mittag-Leffler functions，related topics and applications [M].Berlin：Springer，2020.

[2] PODLUBNY I.分數微積分：理論基礎與應用導論 [M].袁曉，譯.北京：電子工業出版社，2021：106-113.

[3] 方宇孟.米塔-列夫勒函數的高精度快速算法和顯示 [D].成都：四川大學，2022.

[4] WRIGHT E M. The asymptotic expansion of the generalized Bessel function [J].Proceedings of the London Mathematical Society，1935，38（1）：257-270.

[5] WRIGHT E M.The asymptotic expansion of the generalized hyper geometric function [J].Journal of the London Mathematical Society，1935，10（4）：287–293.

[6] WRIGHT E M. The generalized Bessel function of order greater than one [J].The Quarterly Journal of Mathematics，1940，11（1）：36-48.

[7] MAINARDI F，CONSIGLIO A. The Wright functions of the second kind in Mathematical Physics [J/OL].Mathematics，2020，8（6）：（2020-06-01）.https：//api.semanticscholar.org/CorpusID：219762246.

[8] POVSTENKO Y. Some applications of the wright function in continuum physics：A survey [J/OL].Mathematics，2021，9（2）：（2021-01-19）.https：//api.semanticscholar.org/CorpusID：234148001.

[9] CONSIGLIO A，LUCHKO Y，MAINARDI F. Some notes on the Wright functions in probability theory [J].WSEAS Transactions on Mathematics，2019，18：389-393.

[10] LUCHKO Y. Algorithms for evaluation of the Wright function for the real arguments’values [J].Fractional calculus and applied analysis，2008，11（1）：57-75.

[11] LUCHKO Y，Trujillo J，Velasco M P. The Wright function and its numerical evaluation [J].International Journal of Pure and Applied Mathematics，2010，64（4）：567-575.

[12] WONG R，ZHAO Y Q. Smoothing of Stokes' discontinuity for the generalized Bessel function [J].Proceedings of the Royal Society of London. Series A：Mathematical，Physical and Engineering Sciences，1999，455（1984）：1381-1400.

[13] WONG R，ZHAO Y Q. Smoothing of Stokes's discontinuity for the generalized Bessel function. II [J].Proceedings of the Royal Society of London. Series A：Mathematical，Physical and Engineering Sciences，1999，455（1988）：3065-3084.

[14] PARIS R B. Exponentially small expansions in the asymptotics of the Wright function [J].Journal of computational and applied mathematics，2010，234（2）：488-504.

[15] PARIS R B. Exponentially small expansions of the Wright function on the Stokes lines [J].Lithuanian Mathematical Journal，2014，54：82-105.

[16] PARIS R B. Asymptotic expansion of the Wright function for large variable and parameter [J/OL].arXiv：2110.06690 [math.CA].（2021-10-13）.https：//arxiv.org/abs/2110.06690.

[17] PARIS R B，Consiglio A，Mainardi F. On the asymptotics of Wright functions of the second kind [J].Fractional Calculus and Applied Analysis，2021，24（1）：54-72.

[18] NEMES G. Error bounds and exponential improvements for the asymptotic expansions of the gamma function and its reciprocal [J].Proceedings of the Royal Society of Edinburgh Section A：Mathematics，2015，145（3）：571-596.

[19] OLVER F W J. NIST handbook of mathematical functions hardback and CD-ROM [M].Cambridge：Cambridge university press，2010.

[20] ELEZOVI? N. Generalized Bernoulli polynomials and numbers，revisited [J].Mediterranean Journal of Mathematics，2016，13（1）：141-151.

[21] PARIS R B，KAMINSKI D. Asymptotics and mellin-barnes integrals [M].Cambridge：Cambridge University Press，2001.

[22] WRENCH J W. Concerning two series for the gamma function [J].Mathematics of Computation，1968，22（103）：617-626.

[23] PARIS R B. On the Asymptotic Expansion of Γ（x），Lagrange's Inversion Theorem and the Stirling Coefficients [J/OL].arXiv：1405.3423 [math.CA].（2014-05-14）.https：//arxiv.org/abs/1405.3423.

作者簡介：李燕（1995—），女，漢族，四川南充人，碩士研究生，研究方向：分數微積分理論與應用；袁曉（1964—），男，漢族，四川中江人，副教授，博士，研究方向：現代信息信號處理、分數微積分理論與應用、現代電路與系統理論與技術。