Return sqrt(x * y) computed with improved accuracy and without overflow/underflow.
(x: float, y: float)
| 1769 | |
| 1770 | |
| 1771 | def _sqrtprod(x: float, y: float) -> float: |
| 1772 | "Return sqrt(x * y) computed with improved accuracy and without overflow/underflow." |
| 1773 | |
| 1774 | h = sqrt(x * y) |
| 1775 | |
| 1776 | if not isfinite(h): |
| 1777 | if isinf(h) and not isinf(x) and not isinf(y): |
| 1778 | # Finite inputs overflowed, so scale down, and recompute. |
| 1779 | scale = 2.0 ** -512 # sqrt(1 / sys.float_info.max) |
| 1780 | return _sqrtprod(scale * x, scale * y) / scale |
| 1781 | return h |
| 1782 | |
| 1783 | if not h: |
| 1784 | if x and y: |
| 1785 | # Non-zero inputs underflowed, so scale up, and recompute. |
| 1786 | # Scale: 1 / sqrt(sys.float_info.min * sys.float_info.epsilon) |
| 1787 | scale = 2.0 ** 537 |
| 1788 | return _sqrtprod(scale * x, scale * y) / scale |
| 1789 | return h |
| 1790 | |
| 1791 | # Improve accuracy with a differential correction. |
| 1792 | # https://www.wolframalpha.com/input/?i=Maclaurin+series+sqrt%28h**2+%2B+x%29+at+x%3D0 |
| 1793 | d = sumprod((x, h), (y, -h)) |
| 1794 | return h + d / (2.0 * h) |
| 1795 | |
| 1796 | |
| 1797 | def _normal_dist_inv_cdf(p, mu, sigma): |
no outgoing calls
no test coverage detected
searching dependent graphs…