Function bessel_1_asympt_beta_fast 

pub(crate) fn bessel_1_asympt_beta_fast(recip: DoubleDouble) -> DoubleDouble

Note expansion generation below: this is negative series expressed in Sage as positive, so before any real evaluation x=1/x should be applied

Generated by SageMath:

def binomial_like(n, m):
    prod = QQ(1)
    z = QQ(4)*(n**2)
    for k in range(1,m + 1):
        prod *= (z - (2*k - 1)**2)
    return prod / (QQ(2)**(2*m) * (ZZ(m).factorial()))

R = LaurentSeriesRing(RealField(300), 'x',default_prec=300)
x = R.gen()

def Pn_asymptotic(n, y, terms=10):
    # now y = 1/x
    return sum( (-1)**m * binomial_like(n, 2*m) / (QQ(2)**(2*m)) * y**(QQ(2)*m) for m in range(terms) )

def Qn_asymptotic(n, y, terms=10):
    return sum( (-1)**m * binomial_like(n, 2*m + 1) / (QQ(2)**(2*m + 1)) * y**(QQ(2)*m + 1) for m in range(terms) )

P = Pn_asymptotic(1, x, 50)
Q = Qn_asymptotic(1, x, 50)

def sqrt_series(s):
    val = S.valuation()
    lc = S[val]  # Leading coefficient
    b = lc.sqrt() * x**(val // 2)

    for _ in range(5):
        b = (b + S / b) / 2
        b = b
    return b

S = (P**2 + Q**2).truncate(50)

b_series = sqrt_series(S).truncate(30)
# see the beta series
print(b_series)

See notes/bessel_asympt.ipynb for generation