blob: aabefdb5bec5e5c8e0f134621fe2fda5c342f757 (
plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
|
# Copyright (c) 2024 Cisco and/or its affiliates.
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at:
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
"""Module holding functions for avoiding rounding underflows.
Some applications wish to manipulate non-negative numbers
which are many orders of magnitude apart.
In those circumstances, it is useful to store
the logarithm of the intended number.
As math.log(0.0) raises an exception (instead returning -inf or nan),
and 0.0 might be a result of computation caused only by rounding error,
functions of this module use None as -inf.
TODO: Figure out a more performant way of handling -inf.
The functions handle the common task of adding or subtracting
two numbers where both operands and the result is given in logarithm form.
There are conditionals to make sure overflow does not happen (if possible)
during the computation."""
import math
def log_plus(first, second):
"""Return logarithm of the sum of two exponents.
Basically math.log(math.exp(first) + math.exp(second))
which avoids overflow and uses None as math.log(0.0).
TODO: replace with scipy.special.logsumexp ? Test it.
:param first: Logarithm of the first number to add (or None if zero).
:param second: Logarithm of the second number to add (or None if zero).
:type first: float
:type second: float
:returns: Logarithm of the sum (or None if zero).
:rtype: float
"""
if first is None:
return second
if second is None:
return first
if second > first:
retval = second + math.log(1.0 + math.exp(first - second))
else:
retval = first + math.log(1.0 + math.exp(second - first))
return retval
def log_minus(first, second):
"""Return logarithm of the difference of two exponents.
Basically math.log(math.exp(first) - math.exp(second))
which avoids overflow and uses None as math.log(0.0).
TODO: Support zero difference?
TODO: replace with scipy.special.logsumexp ? Test it.
:param first: Logarithm of the number to subtract from (or None if zero).
:param second: Logarithm of the number to subtract (or None if zero).
:type first: float
:type second: float
:returns: Logarithm of the difference.
:rtype: float
:raises RuntimeError: If the difference would be non-positive.
"""
if first is None:
raise RuntimeError("log_minus: does not support None first")
if second is None:
return first
if second >= first:
raise RuntimeError("log_minus: first has to be bigger than second")
factor = -math.expm1(second - first)
if factor <= 0.0:
msg = "log_minus: non-positive number to log"
else:
return first + math.log(factor)
raise RuntimeError(msg)
def safe_exp(log_value):
"""Return exponential of the argument, or zero if the argument is None.
:param log_value: The value to exponentiate.
:type log_value: NoneType or float
:returns: The exponentiated value.
:rtype: float
"""
if log_value is None:
return 0.0
return math.exp(log_value)
|