Re: FindRoot and Bose-Einstein distribution
- To: mathgroup at smc.vnet.net
- Subject: [mg82792] Re: FindRoot and Bose-Einstein distribution
- From: Jean-Marc Gulliet <jeanmarc.gulliet at gmail.com>
- Date: Wed, 31 Oct 2007 06:09:37 -0500 (EST)
- Organization: The Open University, Milton Keynes, UK
- References: <fg6ql3$dmc$1@smc.vnet.net>
P_ter wrote:
> 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
> 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?
Mixing symbolic and numeric result might not be the best approach. If,
say, Sum[1/(Exp[p + k q] - 1), {k, 1, Infinity}] has no closed form (as
Mathematica claims), you are doomed to fail with this approach.
You may be more successful by using purely numerical methods and
controlling the precision as in the following example (though I do not
know whether the result {p -> 3.74867, q -> 0.0464572} makes sense).
In[1]:= Block[{$MinPrecision = 200, $MaxPrecision = 200},
n[p_?NumericQ, q_?NumericQ] :=
NSum[k/(Exp[p + k q] - 1), {k, 1, 10000}];
m[p_?NumericQ, q_?NumericQ] :=
NSum[1/(Exp[p + k q] - 1), {k, 1, 10000000}];
mm = 0.501;
mn = 41.0959;
FindRoot[{m[p, q] - mm, n[p, q] - mn}, {{p, 5}, {q, 0.1}}]
]
During evaluation of In[1]:= NIntegrate::slwcon: Numerical \
integration converging too slowly; suspect one of the following: \
singularity, value of the integration is 0, highly oscillatory \
integrand, or WorkingPrecision too small. >>
During evaluation of In[1]:= NIntegrate::ncvb: NIntegrate failed to \
converge to prescribed accuracy after 9 recursive bisections in k \
near {k} = {120.855}. NIntegrate obtained -4.98592*10^7 and \
136833.21646408358` for the integral and error estimates. >>
During evaluation of In[1]:= SequenceLimit::seqlim: The general form \
of the sequence could not be determined, and the result may be \
incorrect. >>
During evaluation of In[1]:= NIntegrate::slwcon: Numerical \
integration converging too slowly; suspect one of the following: \
singularity, value of the integration is 0, highly oscillatory \
integrand, or WorkingPrecision too small. >>
During evaluation of In[1]:= NIntegrate::ncvb: NIntegrate failed to \
converge to prescribed accuracy after 9 recursive bisections in k \
near {k} = {74.3418}. NIntegrate obtained -4.99961*10^7 and \
2600.0866205544694` for the integral and error estimates. >>
During evaluation of In[1]:= NIntegrate::slwcon: Numerical \
integration converging too slowly; suspect one of the following: \
singularity, value of the integration is 0, highly oscillatory \
integrand, or WorkingPrecision too small. >>
During evaluation of In[1]:= General::stop: Further output of \
NIntegrate::slwcon will be suppressed during this calculation. >>
During evaluation of In[1]:= NIntegrate::ncvb: NIntegrate failed to \
converge to prescribed accuracy after 9 recursive bisections in k \
near {k} = {231.1}. NIntegrate obtained -4.99785*10^7 and \
24827.157779411496` for the integral and error estimates. >>
During evaluation of In[1]:= General::stop: Further output of \
NIntegrate::ncvb will be suppressed during this calculation. >>
During evaluation of In[1]:= SequenceLimit::seqlim: The general form \
of the sequence could not be determined, and the result may be \
incorrect. >>
During evaluation of In[1]:= SequenceLimit::seqlim: The general form \
of the sequence could not be determined, and the result may be \
incorrect. >>
During evaluation of In[1]:= General::stop: Further output of \
SequenceLimit::seqlim will be suppressed during this calculation. >>
Out[1]= {p -> 3.74867, q -> 0.0464572}
Regards,
--
Jean-Marc