Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss

"Safe division" routine

588 views
Skip to first unread message

kat...@ibm.net

unread,
Jul 11, 1999, 3:00:00 AM7/11/99
to
Please, could anyone give me a routine or function to divide two real
(double, if possible) numbers safely? That is, one that performs the
division only if it does not produce overflow.

The type of function I have in mind is something like:

q = safe_div(a,b,altv)

where q is set at a/b in the normal case, but if a/b overflows, then the
returned value is altv. Or something of that sort.

Or perhaps, someone can suggest a better way of handling situations when a
division (or other operation) may produce overflow (or other error), without
requiring the intervention of the user.

Thank you,

Julio

Charles Crawford

unread,
Jul 12, 1999, 3:00:00 AM7/12/99
to
We use a routine similar to the one you describe. It merely returns
very large number if the denominator is too small. We have arrived at
values for "very large" and "too small" based on our particular
application. We use the routine in situations where we haven't been
able to rewrite the computation to guarantee that the denominator will
be large enough (or we haven't had time to rewrite it). This way the
program does not crash, but the numbers are out of range so the user
knows something went wrong. It is particularily useful with an
optimization algorithm where the objective function is not defined for
certain values.

-- Chuck Crawford, Toronto

Herman D. Knoble

unread,
Jul 12, 1999, 3:00:00 AM7/12/99
to
C=A/B implies C=ExpBaseb(LogBaseB(A)-LogBaseB(B))
While it may be slow, looking at the quantity Logbase(A)-LogBaseB(B)
and comparing this to a correct threshold for the platform and precision
would enable division to be performed iff no overflow would
occur. This would essentially be vendor independent (different
PC compilers offer different precisions, see:
http://www.polyhedron.com/ ).

For some compilers, it may be faster to let the overflow occur
and capture it and return "Signed Infinity" , for example, as the result.
Indeed, the (vendor dependent) program, slightly modified
from the on-line documentation:

use dflib
Double Precision A,B,C
Integer(2) Status
A=1.D200
B=1.D-200
C=A/B

C DVF dependent check for overflow during division
CALL GETSTATUSFPQQ(Status)
IF (IAND(Status, FPSW$OVERFLOW ) .NE. 0) THEN
WRITE (*,*) 'Overflow occurred.'
END IF

WRITE(*,*) 'C=',C
STOP
END

when compiled with DVF (6.0b3) and run yields the output lines:
Overflow occurred.
C= Infinity

Salford FTN95 allows one to capture this overflow and
via a vendor dependent routine, FUNCTION TRAP_EXCEPTION@,
write a routine to handle it as needed.

Similarly Lahey FTN95 offers vendor dependent routines,
Intrup and Invalop.

On Sun, 11 Jul 1999 22:07:20 -0400, <kat...@ibm.net> wrote:

-|Please, could anyone give me a routine or function to divide two real
-|(double, if possible) numbers safely? That is, one that performs the
-|division only if it does not produce overflow.
-|
-|The type of function I have in mind is something like:
-|
-|q = safe_div(a,b,altv)
-|
-|where q is set at a/b in the normal case, but if a/b overflows, then the
-|returned value is altv. Or something of that sort.
-|
-|Or perhaps, someone can suggest a better way of handling situations when a
-|division (or other operation) may produce overflow (or other error), without
-|requiring the intervention of the user.
-|
-|Thank you,
-|
-|Julio
-|


Herman D. (Skip) Knoble, Research Associate
Mailto:h...@psu.edu
Web: http://www.personal.psu.edu/hdk
Center for Academic Computing
Penn State University
214C Computer Building
University Park, PA 16802-2101
Phone:+1 814 865-0818 Fax:+1 814 863-7049

kat...@ibm.net

unread,
Jul 17, 1999, 3:00:00 AM7/17/99
to
Thank you to all that offered their suggestions for a "safe division"
routine to prevent overflow. For those that are curious about the solution
to the problem, I found useful to adopt a subroutine along the lines
suggested by William Long:

if(exponent(a) - exponent (b) >= maxexponent(a) .OR. b==0)Then
q=altv
else
q=a/b
endif

The = in the >= is to take into account the case when the fractional part of
a is 1.111... and that of b is 1.

It works very well.

Thank you again,
Julio.

0 new messages