1
$\begingroup$

I have a complex but analytical expression that should always evaluate to zero for certain values, but I'm getting strange numerical noise I think I can avoid. The problem seems to occur in a subsection involving logs. Here's the offending section in Mathematica code:

  rn = 250*10^-6; e = 0.75;
  f = e*rn;

  test1[r_, t_] := Log[(r^2 - (f^2) (Cos[t]^2) )/(rn^2  - f^2)];
  t1out = Plot[test1[rn*Sqrt[1 - (e*Sin[t])^2], t], {t, 0, 2*Pi}]

Now, I know when $r = r_{n}\sqrt{1 - (e\sin(t))^2}$, then for any value of $t$,

$\frac{r^2 - f^2\cos^2 (t)}{r_{n}^2 - f^2} = \frac{r_{n}^2 - f^2(\sin^2(t) + \cos^2 (t))}{r_{n}^2 - f^2} = 1$

and so the log should always be zero. And for some values of $t$this holds - but for others, I get weird oscillating errors...

errors

These are very small (order $10^{-16}$) but I think they're avoidable, especially as uncorrected they might mess up some other outputs I need. Is there a way to eradicate this numerical noise, given there's an analytic answer?

$\endgroup$
  • 1
    $\begingroup$ Since you are working with machine numbers, this numerical noise is unavoidable. But you can wrap test1 in Chop for suppressing the noise.. $\endgroup$ – Fred Simons May 17 '18 at 11:56
  • $\begingroup$ Though this question is a "simple mistake" as per the close reasons, machine numbers are lately one of my favorite topics. +1 $\endgroup$ – LLlAMnYP May 17 '18 at 19:19
  • $\begingroup$ Thanks, very interesting! $\endgroup$ – DRG May 21 '18 at 13:57
3
$\begingroup$

It is a precision issue.

rn = 250*10^-6;
e = 3/4; (* use exact numbers for parameters *)
f = e*rn;

test1[r_, t_] := Log[(r^2 - (f^2) (Cos[t]^2))/(rn^2 - f^2)];

Use arbitrary precision rather than machine precision.

t1out = Plot[test1[rn*Sqrt[1 - (e*Sin[t])^2], t], {t, 0, 2*Pi}, 
  WorkingPrecision -> 10]

enter image description here

EDIT: Alternatively, use Simplify and Evaluate

t1out = Plot[
  test1[rn*Sqrt[1 - (e*Sin[t])^2], t] // Simplify // Evaluate, {t, 0, 2*Pi}]
$\endgroup$
  • $\begingroup$ Cheers Bob - makes sense now! I figured it'd evaluate first but doing this forces it to do so! $\endgroup$ – DRG May 21 '18 at 13:57
  • $\begingroup$ @DRG - If you evaluate Attributes[Plot] you will see that Plot has the attribute HoldAll so evaluation does not occur initially unless forced by use of Evaluate. $\endgroup$ – Bob Hanlon May 21 '18 at 14:25

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.