3
$\begingroup$

In engineering, we usually have to solve a function numerically, which has no symbolic computation ability and may using a low percision data type like fp32 or even fp16.
Hence, we have to reform a function expression to avoid division by a small number, or remove discontinuities, to improve numerical stability on more condition.
For example, we can reform $$\frac{\frac{1}{t^2}-1}{(t+\frac{1}{t})^2}$$ to $$\frac{1-t^2}{(1+t^2)^2}$$ so we have definition when x == 0 and better numerical stability when t is small.

How can I do this work by Mathematica automatically?

For example, I have a function, which is unstable for low precision data types when f->0 and f->1, and hard to reform manually: $$g(f)=-\frac{i-ie^{-2if\pi}}{4f\pi-4f^3\pi},(f\geq0)$$ Obviously $$\lim_{f\rightarrow0}g(f)=\frac{1}{2}$$ and $$\lim_{f\rightarrow1}g(f)=-\frac{1}{4}$$ How can I reform it by Mathematica to remove the discontinuities at f == 0 and f == 1, and prevent the small number division if possible?

$\endgroup$
9
  • 3
    $\begingroup$ Instead of using the original function g[v_] := I/(4 \[Pi] v) (1 - Exp[2 \[Pi] I v])/(1 - v^2) use its Taylor series with an appropriate lenght, e.g. Series[g[v], {v, 1/2, 6}] // Normal in the whole region. $\endgroup$
    – Artes
    Commented Mar 14, 2020 at 5:20
  • 2
    $\begingroup$ Perhps you could exploit PadeApproximant, however it depends on specific aspects of your problem. If you use Series remember to evaluate with Normal in order to get an ordinary expression. $\endgroup$
    – Artes
    Commented Mar 14, 2020 at 5:32
  • 3
    $\begingroup$ E.g. PadeApproximant[g[v], {v, 0, {3, 4}}] // N yields (0.5 + (0. + 0.298574 I) v - 0.154029 v^2 - (0. + 0.0344845 I) v^3)/(1. - (0. + 2.54445 I) v - 2.72193 v^2 + (0. + 1.47873 I) v^3 + 0.354081 v^4). You have to play around a bit to decide what would be the most efficient way and simultaneously well approximated. $\endgroup$
    – Artes
    Commented Mar 14, 2020 at 5:38
  • 2
    $\begingroup$ You might also want to consider to use Taylor or Padé approximations only in small intervals of the problematic points. That will take you several branching constructs like (Piecewise is particularly useful for that within Mathematica, but in general you might have to use ifs or other switchs). You also have to fiddle around with parameters for the expansion order and the diameter of the intervals, but that can be figured out by running test in Mathematica... $\endgroup$ Commented Mar 14, 2020 at 7:38
  • 2
    $\begingroup$ Once you have found a nice polynomial, you might want to apply HornerForm to it in order to make its evaluation both more efficient and more stable. $\endgroup$ Commented Mar 14, 2020 at 7:39

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.