如何在 Maxima CAS 中计算 Brjuno 函数？

Question

我尝试为 [0, 1] 范围内的 x 中的实数计算 real 1-Brjuno function

  G(x):= block(
            [g],
        if (x=0) 
            then g : 0
            else g : float(1/x - floor(1/x)),
        return( g)
        )$
    
   bet(j,x) := block(
        [r],
        if (j=-1) 
            then r:1
            else r : product (G(x)^i, i, 0, j),
        return( r)
        )$
    
    A(i,x) := bet(i-1,x)*log(1/G(x)^i)$
    B(x):= sum(A(i,x), i, 0,10)$

存在数值错误：

                        0
    expt: undefined: 0.0

由某些 x 值的函数引起：1/2、1/5、1/10、1/3、1/9、1/4、0

我该如何解决这个问题？

Answer 1

Brjuno 函数在有理点处具有“对数奇点”，因此使用简单的“

G(x):= block(

    [g],
    if (x=0) 
        then g : 0
        else g : float(1/x - floor(1/x)),
    return( g)
    )$
    

beta(j,x) := block(

    [r],
    if (j=-1) 
        then r:1
        else r : product (G(x)^i, i, 0, j),
    return( r)
    )$
    
    
A(i,x) := beta(i-1,x)*log(1/G(x)^i)$
B(x):= sum(A(i,x), i, 0,10)$
plot2d(B(x),[x,0,1], [y, -1, 10]);

给出情节：

但该图与论文中的图1（左）不一样

Answer 2

这会产生与论文相似的情节。我使用了 wikipedia 中的一个方程和一个连分数近似值：

cflength: 10 $
B(x):= if length(x)=1 then 0 else (local(x), x[1]: 0, -log(cf2num(x))+cf2num(x)*B(rest(x))) $
cf2num(e):=ev(cfdisrep(e), numer, infeval) $
x0: makelist(i*sqrt(2)/1421, i, 1, 1000), numer $
x: map(cf, x0) $
y: map(B, x) $
draw2d(points(x0, y)) $

如何在 Maxima CAS 中计算 Brjuno 函数？

How to compute Brjuno function in Maxima CAS?

function

numerical-computing

maxima