MathGroup Archive 2007

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: FindRoot and Bose-Einstein distribution

On Oct 30, 3:40 am, P_ter <peter_van_summe... at> wrote:
> Hello,
> I have two functions and two values:
> Clear[p, q, k, mm, mn]
> n[p_, q_] := Sum[k/(Exp[p + k q ] - 1), {k, 1, 10000}]
> m[p_, q_] := Sum[1/(Exp[p + k q ] - 1), {k, 1, Infinity}]
> mm = 0.501
> mn = 41.0959
> The first two equations are the Bose-Einstein distribution. Given mm and mn, find p and q.
> A first estimation is: p=5.086,q= 0.01226
> I know that n[p_,q_] is stable until k = 50000 for p=5.086 and q= 0.01226. My check is: n[5.086,0.01226]= 41.1969
> m[p_,q_] is also stable until 100000 and m[5.086,0.01226]= 0.502762

It is not clear to me what is meant by "stable until k=...". Do you
actually intend to do infinite summations for both m and n, and are
checking to see how far finite truncations can be computed wo some

> Everything seems ok.
> So, I tried:
> FindRoot[{m[p, q] - mm, n[p, q] - mn}, {{p, 5}, {q, 0.1}}]
> The answer from FindRoot is that the Sum does not converge.
> What did I do wrong, what went wrong, why?
> with friendly greetings,
> P_ter

There is not much hope to working with the infinite summation since
Sum cannot compute a closed form for the result. Hence anything it
does will rely on some sort of approximation via numeric summing of
some initial summands, and extrapolation for the tail. This might or
might not give reliable results. Since it is fairly easy to bound the
error of finite summations "by hand", given that {p,q} are in the
whereabouts of {5,.01}, it is probably more reliable just to truncate
the sum at some sufficiently high, but finite, bound. 10000 will do
just fine (check it).

n[p_Real, q_Real] := Sum[k/(Exp[p + k*q] - 1), {k, 1, 10000}]
m[p_Real, q_Real] := Sum[1/(Exp[p + k*q] - 1), {k, 1, 10000}]
mm = 0.501;
mn = 41.0959;

In[13]:= rt =
 FindRoot[{m[p, q] - mm, n[p, q] - mn}, {{p, 5}, {q, 0.1}}]
Out[13]= {p -> 5.09055, q -> 0.0122471}

In[14]:= m[p, q] /. rt
Out[14]= 0.501

In[15]:= n[p, q] /. rt
Out[15]= 41.0959

I'll sum more terms for n[] just to show that the result really does
not change (that is, the truncation approximation was fine).

In[26]:= Sum[1/(Exp[p + k q] - 1) /. rt, {k, 1, 10^5}] // InputForm


This is the same as I get if opnly summing to 10^4.

Daniel Lichtblau
Wolfram Research

  • Prev by Date: Re: Setting Negatives to Zero
  • Next by Date: Re: Manipulating a complex modulus expression
  • Previous by thread: Re: FindRoot and Bose-Einstein distribution
  • Next by thread: Re: Problem with element and Maximize