有没有更有效的方法来嵌套对数？

Question

这是这里发布的两个问题的继续，

Declaring a functional recursive sequence in Matlab

Nesting a specific recursion in Pari-GP

长话短说，我构造了一个函数族，用于求解四次泛函方程。我已经证明了这些东西是全纯的。现在是时候制作图表了，或者至少是一些可以通过的代码来评估这些东西。我已经设法达到了大约13个有效数字的精度，但如果我试图获得更多，我会遇到一个特定的错误。这个错误实际上只是一个溢出错误。但这是一个特殊的溢出错误；Pari-GP似乎不喜欢嵌套对数。

我的特殊数学函数是通过取一些大的东西(考虑顺序e^e)来产生小的东西(顺序e^(-n))来近似的。数学本质上需要较大值的样本来产生这些较小的值。奇怪的是，随着我们越来越接近数值近似(大约13位有效数字)，我们也越来越接近溢出，因为我们需要如此大的值来获得这13位有效数字。我是一个糟糕透顶的程序员；我想知道是否有一些我看不到的工作。

/*
This function constructs the approximate Abel function
The variable z is the main variable we care about; values of z where real(z)>3 almost surely produces overflow errors
The variable l is the multiplier of the approximate Abel function
The variable n is the depth of iteration required
n can be set to 100, but produces enough accuracy for about 15 
The functional equation this satisfies is exp(beta_function(z,l,n))/(1+exp(-l*z)) = beta_function(z+1,l,n); and this program approaches the solution for n to infinity
*/

beta_function(z,l,n) =
{
    my(out = 0);
    for(i=0,n-1,
        out = exp(out)/(exp(l*(n-i-z)) +1));
    out;
}

/*
This function is the error term between the approximate Abel function and the actual Abel function
The variable z is the main variable we care about
The variable l is the multiplier
The variable n is the depth of iteration inherited from beta_function
The variable k is the new depth of iteration for this function
n can be set about 100, still; but 15 or 20 is more optimal. 
Setting the variable k above 10 will usually produce overflow errors unless the complex arguments of l and z are large. 
Precision of about 10 digits is acquired at k = 5 or 6 for real z, for complex z less precision is acquired. k should be set to large values for complex z and l with large imaginary arguments.
*/

tau_K(z,l,n,k)={
    if(k == 1, 
      -log(1+exp(-l*z)), 
      log(1 + tau_K(z+1,l,n,k-1)/beta_function(z+1,l,n)) - log(1+exp(-l*z))
    )
}

/*
This is the actual Abel function
The variable z is the main variable we care about
The variable l is the multiplier
The variable n is the depth of iteration inherited from beta_function
The variable k is the depth of iteration inherited from tau_K
The functional equation this satisfies is exp(Abl_L(z,l,n,k)) = Abl_L(z+1,l,n,k); and this function approaches that solution for n,k to infinity
*/

Abl_L(z,l,n,k) ={
    beta_function(z,l,n) + tau_K(z,l,n,k);
}

这是我已经证明是全纯的近似函数的代码；但不幸的是，我的代码太可怕了。在这里，附加了一些预期输出，其中您可以看到函数方程满足大约10 - 13个有效数字。

Abl_L(1,log(2),100,5)

%52 = 0.1520155156321416705967746811

exp(Abl_L(0,log(2),100,5))

%53 = 0.1520155156321485241351294757

Abl_L(1+I,0.3 + 0.3*I,100,14)

%59 = 0.3353395055605129001249035662 + 1.113155080425616717814647305*I

exp(Abl_L(0+I,0.3 + 0.3*I,100,14))

%61 = 0.3353395055605136611147422467 + 1.113155080425614418399986325*I

Abl_L(0.5+5*I, 0.2+3*I,100,60)

%68 = -0.2622549204469267170737985296 + 1.453935357725113433325798650*I

exp(Abl_L(-0.5+5*I, 0.2+3*I,100,60))

%69 = -0.2622549205108654273925182635 + 1.453935357685525635276573253*I 

现在，你会注意到，我必须为不同的值改变k值。当自变量z，l离实轴更远时，我们可以使k非常大(我们必须这样做才能获得良好的精度)，但它最终仍会溢出；通常情况下，一旦我们获得了大约13-15个有效数字，函数就会开始爆炸。你会注意到，设置k = 60，意味着我们取60个对数。这听起来已经是个坏主意了，哈哈。但是从数学上讲，值Abl_L(z，l，无穷大，无穷大)正是我想要的函数。我知道这一定很奇怪；嵌套的无限for循环听起来像是胡说八道，哈哈。

我想知道是否有人能想出一种方法来避免这些溢出错误，并获得更高的准确性。在一个完美的世界中，这个对象几乎肯定会收敛，并且这个代码是完美的(尽管它可能有点慢)；但是我们可能需要无限地增加堆栈大小。理论上这是非常好的；但在现实中，它不仅仅是不切实际的。作为一个程序员，有没有办法解决这个问题呢？

在这一点上，我唯一的另一个选择就是尝试创建一个强力算法来发现这个函数的泰勒级数；但我在这方面的运气并不好。这个过程非常独特，试图用泰勒级数来解决这个问题，会让我们回到起点。除非，这里有人能想出一个奇妙的方法，从这个表达式中恢复泰勒级数。

老实说，我对所有的建议，任何评论都持开放态度。我束手无策；我想知道这是不是唯一的解决方案之一就是无限期地增加堆栈大小(这绝对有效)。这不仅仅是因为我在处理大量的数字。我需要越来越大的值来计算一个小的值。出于这个原因，我想知道是否有一些我看不到的快速工作。Pari-GP输出的错误总是与tau_K一起出现的，所以我想知道这是不是进行了次优化编码；我应该在它迭代时添加一些东西来减小堆栈大小。或者，如果有可能的话。再说一次，我是一个糟糕的程序员。我需要有人向我解释，就像我在幼儿园一样。

任何需要澄清的帮助、评论和问题，都非常受欢迎。在这一点上，我就像一只追逐自己尾巴的狗；想知道为什么他不能取1000个对数，哈哈。

致以问候。

编辑：

我想我应该补充一下，我可以产生任意的精度，但是我们必须在左半平面上保持z的论点。如果变量为n,k = -real(z)，那么我们可以通过将n设置为我们想要的大小来产生任意精度。这里有一些输出来解释这一点，其中我使用了\p 200，并且我们在这个级别上几乎是相等的(减去一些数字)。

Abl_L(-1000,1+I,1000,1000)

%16 = -0.29532276871494189936534470547577975723321944770194434340228137221059739121428422475938130544369331383702421911689967920679087535009910425871326862226131457477211238400580694414163545689138863426335946 + 1.5986481048938885384507658431034702033660039263036525275298731995537068062017849201570422126715147679264813047746465919488794895784667843154275008585688490133825421586142532469402244721785671947462053*I

exp(Abl_L(-1001,1+I,1000,1000))

%17 = -0.29532276871494189936534470547577975723321944770194434340228137221059739121428422475938130544369331383702421911689967920679087535009910425871326862226131457477211238400580694414163545689138863426335945 + 1.5986481048938885384507658431034702033660039263036525275298731995537068062017849201570422126715147679264813047746465919488794895784667843154275008585688490133825421586142532469402244721785671947462053*I

Abl_L(-900 + 2*I, log(2) + 3*I,900,900)

%18 = 0.20353875452777667678084511743583613390002687634123569448354843781494362200997943624836883436552749978073278597542986537166527005507457802227019178454911106220050245899257485038491446550396897420145640 - 5.0331931122239257925629364016676903584393129868620886431850253696250415005420068629776255235599535892051199267683839967636562292529054669236477082528566454129529102224074017515566663538666679347982267*I

exp(Abl_L(-901+2*I,log(2) + 3*I,900,900))

%19 = 0.20353875452777667678084511743583613390002687634123569448354843781494362200997943624836883436552749978073278597542986537166527005507457802227019178454911106220050245980468697844651953381258310669530583 - 5.0331931122239257925629364016676903584393129868620886431850253696250415005420068629776255235599535892051199267683839967636562292529054669236477082528566454129529102221938340371793896394856865112060084*I

Abl_L(-967 -200*I,12 + 5*I,600,600)

%20 = -0.27654907399026253909314469851908124578844308887705076177457491260312326399816915518145788812138543930757803667195961206089367474489771076618495231437711085298551748942104123736438439579713006923910623 - 1.6112686617153127854042520499848670075221756090591592745779176831161238110695974282839335636124974589920150876805977093815716044137123254329208112200116893459086654166069454464903158662028146092983832*I

exp(Abl_L(-968 -200*I,12 + 5*I,600,600))

%21 = -0.27654907399026253909314469851908124578844308887705076177457491260312326399816915518145788812138543930757803667195961206089367474489771076618495231437711085298551748942104123731995533634133194224880928 - 1.6112686617153127854042520499848670075221756090591592745779176831161238110695974282839335636124974589920150876805977093815716044137123254329208112200116893459086654166069454464833417170799085356582884*I

问题是，我们不能一遍又一遍地应用exp，然后期望保持相同的精度。问题在于exp，当你在复杂的平面上迭代它时，它会显示出如此多的混乱行为，这是注定要工作的。

Answer 1

好吧，我回答了我自己的问题。@user207421发布了一条评论，我不确定它是否意味着我所认为的意思，但我认为它让我到达了我想要的地方。我假设exp不会继承其参数的精确度，但显然这是真的。所以我所需要的就是定义，

Abl_L(z,l,n,k) ={
    if(real(z) <= -max(n,k),
    beta_function(z,l,n) + tau_K(z,l,n,k),
    exp(Abl_L(z-1,l,n,k)));
}

从这里开始，一切都工作得很好；当然，对于我需要它的地方。所以，我回答了我自己的问题，它很简单。我只需要一条if语句。

无论如何，感谢阅读这篇文章的人。