#20020: asymptotic expansion generator: singularity analysis (log-type) without
renormalization
-------------------------------------+-------------------------------------
Reporter: cheuberg | Owner:
Type: enhancement | Status: positive_review
Priority: major | Milestone: sage-7.1
Component: asymptotic | Resolution:
expansions | Merged in:
Keywords: | Reviewers: Daniel Krenn
Authors: Clemens Heuberger | Work issues:
Report Upstream: N/A | Commit:
Branch: u/dkrenn/asy | 4462b276a506a8a3439a4cad07fa1df09a4d6892
/singularity-generator-log-non- | Stopgaps:
normalized |
Dependencies: #19969 |
-------------------------------------+-------------------------------------
Comment (by cheuberg):
Replying to [comment:10 dkrenn]:
> From the docstring:
> {{{
> sage:
asymptotic_expansions._SingularityAnalysis_non_normalized_(
> ....: 'n', 1, alpha=2, beta=2, precision=4) / n
> log(n)^2 + (2*euler_gamma - 2)*log(n)
> - 2*euler_gamma + euler_gamma^2 - 1/6*pi^2 + 2
> + n^(-1)*log(n)^2 + O(n^(-1)*log(n))
>
> Be aware that the last result does *not* coincide with
[FS2009]_,
> they do have a different error term.
> }}}
> I've checked this result and it (in particular, the term
`n^(-1)*log(n)^2`) seems to be correct.
I have run numerical tests (via `compare_with_values`) for
{{{
sage: asymptotic_expansions._SingularityAnalysis_non_normalized_(
....: 'n', 1, alpha=2, beta=2, precision=6)
n*log(n)^2 + (2*euler_gamma - 2)*n*log(n) + (-2*euler_gamma +
euler_gamma^2 - 1/6*pi^2 + 2)*n + log(n)^2 + (2*euler_gamma + 1)*log(n) +
euler_gamma + euler_gamma^2 - 1/6*pi^2 + O(n^(-1)*log(n)^2)
}}}
the result is
{{{
[(5000, -0.1171037?),
(10000, -0.1069521?),
(15000, -0.101778?),
(20000, -0.098395?),
(25000, -0.095920?),
(30000, -0.093987?),
(35000, -0.09242?),
(40000, -0.09109?),
(45000, -0.08995?),
(50000, -0.08896?)]
}}}
which very much supports the additional term here (if it were not here,
the renormalization would have multiplied with `n`, which is clearly not
occurring).
I therefore submitted this as an erratum at the
[http://ac.cs.princeton.edu/errata/ book site].
For completeness, here is my code for testing.
{{{
@cached_function
def H(k):
"""
Harmonic numbers. Generating function: -log(1-z)/(1-z)
"""
if k < 1:
return 0
else:
return 1/k + H(k-1)
@cached_function
def a(k):
"""
Coefficients of (-log(1-z))^2: differentiate GF once,
obtain 2(-log(1-z))/(1-z), i.e., the harmonic numbers multiplied by 2.
Shift and divide by k in order to take differntiation into account.
"""
return 2*H(k-1)/k
@cached_function
def sum_a(k):
"""
Summatory function of a(k), thus GF (-log(1-z))^2/(1-z).
"""
if k<1:
return 0
else:
return sum_a(k-1)+a(k)
@cached_function
def sum_sum_a(k):
"""
Summatory function of sum_a(k), thus GF (-log(1-z))^2/(1-z)^2.
"""
if k<1:
return 0
else:
return sum_sum_a(k-1)+sum_a(k)
R.<z> = PowerSeriesRing(QQ, default_prec=100)
assert (log(1/(1-z)))^2 - sum(z^n*2*H(n-1)/n for n in range(1, 101)) ==
O(z^101)
assert (log(1/(1-z))/(1-z))^2 - sum(z^n*sum_sum_a(n) for n in range(1,
101)) == O(z^101)
%time result = [sum_sum_a(_) > 0 for _ in range(100, 50001, 100)]
ae = asymptotic_expansions._SingularityAnalysis_non_normalized_('n', 1,
alpha=2, beta=2, precision=6); ae
ae.compare_with_values('n', sum_sum_a, srange(5000, 50001, 5000))
}}}
--
Ticket URL: <http://trac.sagemath.org/ticket/20020#comment:11>
Sage <http://www.sagemath.org>
Sage: Creating a Viable Open Source Alternative to Magma, Maple, Mathematica,
and MATLAB
--
You received this message because you are subscribed to the Google Groups
"sage-trac" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
To post to this group, send email to [email protected].
Visit this group at https://groups.google.com/group/sage-trac.
For more options, visit https://groups.google.com/d/optout.