Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Planck's law dependent on frequency rho is as follows:

   B[T_, h_, rho_, k_, c_] := (2 h rho^3)/c^2 1/(Exp[h rho/(T k)] - 1)

As you can see, the denominator can only be zero for rho = 0. All the rest are nonzero constants. Curiously if we try to plot it

Plot[B[5800, 6.626 10^-34, rho, 1.38 10^-23, 299792458], {rho, 
  0.1 10^-6, 3 10^-6}]

It will return "Infinite expression 1/0. encountered" and the plot will be blank. In fact, substituting any value for rho appears to result in this error (out = ComplexInfinity). I have no experience with precision, but I don't see any other possible reason. Anyway to solve this?

share|improve this question
1  
when doing numerical calculations, it's better to work in dimensionless units, like here: mathematica.stackexchange.com/a/2057/16 so as to avoid underflow/overflow problems. so if I were you I'd set $h=k_B=c=1$ and go from there –  acl May 17 '12 at 22:43
    
As @acl said, try Plot[B[5, 1, rho, 1, 1], {rho, 0, 100}] ... Plank's units –  belisarius May 17 '12 at 22:45
    
Also, not sure what you're doing but there is not much 5800K black-body radiation at 10^-6 Hertz. A better place to look might be around 10^15 Hz. –  dws May 17 '12 at 22:51
add comment

1 Answer 1

up vote 10 down vote accepted

To understand what is going on, look at this:

N[b[5800, 6626 10^-37, rho, 138 10^-25, 299792458]]

Mathematica graphics

look at the exponent in the second term in the denominator. If you put $\rho\approx 10^{-6}$, you're trying to calculate the difference between two numbers that are extremely close to each other, namely, between $1$ and $\exp(\alpha)$ with $\alpha\approx 10^{-21}$. So you run into numerical accuracy problems with machine numbers; see here and at the end of this answer.

In general, when doing numerical calculations for physical problems, it's better to work in dimensionless units, so as to avoid underflow/overflow problems. so if I were you I'd set $h=k_B=c=1$ (this amounts to a choice of units) and go from there.

share|improve this answer
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.