2

In a programme I'm running, at a certain point there's a multipplication of variables that gives underflow...

For example $c=c_1\times c_2$. Is there anyway to check if that multiplication gives underflow (error message General::munfl), before it gives a warning, and by how much it would be necessary to multiply $c$ for it not to give an underflow warning? Or if not possible, how could I just consider $c$ as zero, without that specific warning being showed?

Michael E2
  • 235,386
  • 17
  • 334
  • 747
An old man in the sea.
  • 2,527
  • 1
  • 18
  • 26
  • 1
    There are some pointers here: https://mathematica.stackexchange.com/questions/170416/new-generalmunfl-error-and-loss-of-precision – Roman May 03 '19 at 08:27

2 Answers2

7

An example:

c1 = c2 = 2^-1022 // N
(* 2.22507*10^-308 *)
c1*c2
(* General::munfl warning *)
(* 0. *)

To suppress the underflow warning and replace the result by zero:

Quiet[c1*c2, General::munfl]
(* 0. *)

Or switch this message off globally:

Off[General::munfl];
c1*c2
(* 0. *)

Multiplying each number by $2^{512}$ fixes the underflow for sure. But then you may run into overflow issues on the other end. Here's a diagnostic multiplication function:

On[General::munfl];
mymult::stretch = 
  "underflow detected - please multiply both factors by at least `1`.";
mymult[a_?MachineNumberQ, b_?MachineNumberQ] := 
  Quiet[
    Check[a*b, 
      Message[mymult::stretch, 4*Exp[-512 Log[2] - (Log[a] + Log[b])/2]];
      $Failed, 
      General::munfl],
    General::munfl]

mymult[c1, c2]
(* mymult::stretch: underflow detected - please multiply both factors by at least 1.34`*^154. *)
(* $Failed *)

Maybe you could work with the logarithms of these numbers instead? This way you're much less likely to run into over-/underflow issues.

lc1 = Log[c1];
lc2 = Log[c2];
lc1 + lc2
(* -1416.79 *)
Roman
  • 47,322
  • 2
  • 55
  • 121
0

One way to avoid "Machine underflow" issues it to not use machine numbers. Instead, use exact or extended precision numbers. For instance, with exact numbers:

c1 = c2 = 2^-1022;
c1 c2

1/2019812879456937956294679793041871997527756416857217752008146589220290946179243180824825088220182091480544872557618626218382472446905682568443009524153017695039429835456312255734387359399353256674753602399004223017299513665163734760114880896154760654411352865752269065180473493221316613037972024945245649095119645836854271401292810924160285593428511002207895128629862853708189137044278769634391162054011069795371475232403866084849896947018852869025231100827080451951695355294426263107822318857933207716854908911291043620940374829272062414888470322899339833471475133464576850290332404912809509262097240850691224764416

Or, with extended precision numbers:

e1 = e2 = N[2^-1022, 10];
e1 e2

4.950953676*10^-616

Most functions in Mathematica will work with both exact and extended precision numbers, although speed of computation will be lowered.

A typical example where extended precision numbers are useful:

c1 = N[Exp[700]]
c2 = N[Exp[-380]]

c1 c2^2

1.01423*10^304

9.29174*10^-166

General::munfl: 9.29174*10^-166^2 is too small to represent as a normalized machine number; precision may be lost.

0.

Using extended precision numbers instead:

Activate @ SetPrecision[Inactivate[c1 c2^2], 10]

8.756511*10^-27

Carl Woll
  • 130,679
  • 6
  • 243
  • 355