Re: a tricky limit
- To: mathgroup at smc.vnet.net
- Subject: [mg15562] Re: [mg15327] a tricky limit
- From: jtischer at col2.telecom.com.co (Cl-Jurgen Tischer)
- Date: Tue, 26 Jan 1999 13:44:42 -0500 (EST)
- Organization: Universidad del Valle
- Sender: owner-wri-mathgroup at wolfram.com
Arnold,
tricky indeed. I'm still far from having a 8 digits estimate of your
limit, but at least I think I'm closer than the other offers. My best
bet at this moment for the limit is 2.29217, but what I think I really
have is a lower bound, 2.29202448976273576, which is the value at
.9999. It took 700 seconds to calculate, and I'm afraid .99999 would
need 7000. Here is the code I used:
p[m_] := m - m/Divisors[m][[2]]
logh[x_, prec_:$MachinePrecision] :=
Module[{x1 = SetPrecision[x, Max[2*prec, prec + 10]], a, b, q,
n = 10, data, qdata}, a = -Log[1 - x1];
data = Prepend[Table[Sum[Log[1 - x1^p[m]/(m + 1)],
{m, 1000*n + 1, 1000*(n + 1)}], {n, 1, 9}],
Sum[Log[1 - x1^p[m]/(m + 1)], {m, 2, 1000}]];
While[Abs[Plus @@ Last[data]] > 10^(-(prec + 1)),
AppendTo[data, Sum[Log[1 - x1^p[m]/(m + 1)],
{m, 1000*n + 1, 1000*(n + 1)}]]; n++];
qdata = Drop[data, 1]/Drop[data, -1];
q = Fit[Rest[qdata], {1, 1/x, 1/x^2, 1/x^3}, x] /. x -> 0;
While[Abs[data[[-1]]*q/(1 - q)] > 10^(-(prec + 1)),
AppendTo[data, Sum[Log[1 - x1^p[m]/(m + 1)],
{m, 1000*n + 1, 1000*(n + 1)}]];
AppendTo[qdata, Last[data]/data[[-2]]];
q = Fit[Rest[qdata], {1, 1/x, 1/x^2, 1/x^3}, x] /. x -> 0; n++];
{N[a + Plus @@ data, prec], Abs[data[[-1]]*q/(1 - q)], n}]
h[(x_)?NumericQ, prec_:$MachinePrecision] :=
If[x <= 0, 0, Exp[logh[x, prec][[1]]]]
And here are some values (16 digits).
{{99/100, 2.263266115340667}, {991/1000, 2.266275569325346},
{124/125, 2.269289055550239}, {993/1000, 2.272303620578945},
{497/500, 2.275315070142442}, {199/200, 2.278317282731124},
{249/250, 2.281300928902045}, {997/1000, 2.284250680452356},
{499/500, 2.287137638043625}, {999/1000, 2.289893430705526},
{9991/10000, 2.290155960315155}, {9993/10000, 2.290669636995974},
{1999/2000, 2.291163421615598}, {9997/10000,
2.29162632005948},{9999/10000,2.29202448976273576`}}
I couldn't make NSum, NProduct or NLimit work, so I just fooled around
with your function and searched for a heuristic algorithm. I wouldn't
like to prove the results but I'm fairly sure they are ok. I tried to
use those values to calculate the limit with more precision but failed.
Maybe someone else will give it a try.
Jurgen
Arnold Knopfmacher wrote:
>
> I wish to obtain a numerical estimate (say 8 decimal digits) of the
> limit as x tends to 1 from below of the function
> h[x]=(Product[(1-fm[x]/(m+1)),{m,2,Infinity}])/(1-x) where
> fm[x]=x^(m-m/d) and d is the smallest divisor of m that is greater than
> 1. The problem is that when I replace Infinity by say 1000 as the upper
> limit of the product, the function blows up near 1. Visual inspection
> of the graph of h[x] for 0<x<0.9 say, suggests that the limit should
> have a value around 2.1. Can anyone help?
>
> Thanks
> Arnold Knopfmacher
> Dept of Computational and Applied Math Witwatersrand University
> South Africa