5
$\begingroup$

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?

$\endgroup$
  • 1
    $\begingroup$ 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 $\endgroup$ – acl May 17 '12 at 22:43
  • $\begingroup$ As @acl said, try Plot[B[5, 1, rho, 1, 1], {rho, 0, 100}] ... Plank's units $\endgroup$ – Dr. belisarius May 17 '12 at 22:45
  • $\begingroup$ 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. $\endgroup$ – dws May 17 '12 at 22:51
10
$\begingroup$

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.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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