0
$\begingroup$

When asking Mathematica to integrate a function with Integrate, I'm getting drastically different results depending on whether I have a proportionality constant inside or outside of the Integrate command.

To summarize: I'm attempting to define a probability distribution from the function ds:

ds[e_, th_]:= a^2r^2(P[e, th]^2)(P[e, th]+P[e, th]^(-1)-1+Cos[th]^2)/2

with

P[e_, th_]:=1/(1+(e/(m*c^2))(1-Cos[th]))

and a, r, c, and m real constants, as follows:

c = 3*10^(8);
h = 1.05*10^(-34);
m = 9.11*10^(-31);
r = h/(m*c);
a = 7.30*10^(-3);

I'd like the distribution to give the probability density function for values of th from 0 to 2pi. However, for this to work properly with ProbabilityDistribution, ds needs to be normalized. I should be able to do this by integrating ds over the relevant interval (for a specific e, in this case 10^(-15)) and dividing ds by the resulting constant:

s=Integrate[ds[10^(-15), th],{th, 0, 2 Pi}]

1.66118*10^(-25)

ds2[th_]:=ds[10^(-15), th]/s

However, something seems to be broken in this solution. If I define a probability distribution in terms of ds2, it doesn't act as though the distribution is normalized, and instead returns values far outside the desired range. To diagnose the problem, I've tried integrating ds divided by s. Here's what I'm getting: When s is outside the integral, I have

Integrate[ds[10^(-15), th], {th, 0, 2 Pi}]/s

1.

as expected, but when s is inside the integral, I get

Integrate[ds[10^(-15), th]/s, {th, 0, 2 Pi}]

0.000217878

Because s is a constant, moving it inside the integral shouldn't affect the result, so clearly something is going wrong. I've retried the computation with common functions instead of something complicated like ds, and don't have any problems there. What might the problem be?

$\endgroup$
3
  • $\begingroup$ Welcome to Mathematica.SE! I suggest that: 1) You take the introductory Tour now! 2) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! 3) As you receive help, try to give it too, by answering questions in your area of expertise. $\endgroup$ – bbgodfrey Mar 29 '15 at 4:23
  • $\begingroup$ I do not obtain 1.66118*10^(-25) for s but rather an expression dependent on your various constants. Please provide the values of those constants, if you wish readers to try to reproduce your result for s. Of course, Integrate[ds[10^(-15), th], {th, 0, 2 Pi}]/s still yields 1. $\endgroup$ – bbgodfrey Mar 29 '15 at 4:49
  • $\begingroup$ Values have been added $\endgroup$ – TheMac Mar 29 '15 at 4:55
2
$\begingroup$

You are having round-off problems. One way to eliminate them is by replacing your approximate real numbers by rational numbers, because Mathematica computes results with exact numbers to arbitrary precison.

c = 3*10^(8); h = 105*10^(-36); m = 911*10^(-33); r = h/(m*c); a = 730*10^(-5);

Then, both normalized integrals yield precisely 1.

$\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.