|
3 | 3 | import numpy as np |
4 | 4 | import operator |
5 | 5 |
|
| 6 | +from mpmath import mp |
| 7 | + |
6 | 8 | import numpy_quaddtype |
7 | 9 | from numpy_quaddtype import QuadPrecDType, QuadPrecision |
| 10 | +from numpy_quaddtype import pi as quad_pi |
8 | 11 |
|
9 | 12 |
|
10 | 13 | def test_create_scalar_simple(): |
@@ -1735,6 +1738,225 @@ def test_divmod_broadcasting(): |
1735 | 1738 | np.testing.assert_allclose(float(quotients[i]), expected_quotients[i], rtol=1e-14) |
1736 | 1739 | np.testing.assert_allclose(float(remainders[i]), expected_remainders[i], rtol=1e-14) |
1737 | 1740 |
|
| 1741 | +class TestTrignometricFunctions: |
| 1742 | + @pytest.mark.parametrize("op", ["sin", "cos", "tan", "atan"]) |
| 1743 | + @pytest.mark.parametrize("val", [ |
| 1744 | + # Basic cases |
| 1745 | + "0.0", "-0.0", "1.0", "-1.0", "2.0", "-2.0", |
| 1746 | + # pi multiples |
| 1747 | + str(quad_pi), str(-quad_pi), str(2*quad_pi), str(-2*quad_pi), str(quad_pi/2), str(-quad_pi/2), str(3*quad_pi/2), str(-3*quad_pi/2), |
| 1748 | + # Small values |
| 1749 | + "1e-10", "-1e-10", "1e-15", "-1e-15", |
| 1750 | + # Values near one |
| 1751 | + "0.9", "-0.9", "0.9999", "-0.9999", |
| 1752 | + "1.1", "-1.1", "1.0001", "-1.0001", |
| 1753 | + # Medium values |
| 1754 | + "10.0", "-10.0", "20.0", "-20.0", |
| 1755 | + # Large values |
| 1756 | + "100.0", "200.0", "700.0", "1000.0", "1e100", "1e308", |
| 1757 | + "-100.0", "-200.0", "-700.0", "-1000.0", "-1e100", "-1e308", |
| 1758 | + # Fractional values |
| 1759 | + "0.5", "-0.5", "1.5", "-1.5", "2.5", "-2.5", |
| 1760 | + # Special values |
| 1761 | + "inf", "-inf", "nan", |
| 1762 | + ]) |
| 1763 | + def test_sin_cos_tan(self, op, val): |
| 1764 | + mp.prec = 113 # Set precision to 113 bits (~34 decimal digits) |
| 1765 | + numpy_op = getattr(np, op) |
| 1766 | + mpmath_op = getattr(mp, op) |
| 1767 | + |
| 1768 | + quad_val = QuadPrecision(val) |
| 1769 | + mpf_val = mp.mpf(val) |
| 1770 | + |
| 1771 | + quad_result = numpy_op(quad_val) |
| 1772 | + mpmath_result = mpmath_op(mpf_val) |
| 1773 | + # convert mpmath result to quad for comparison |
| 1774 | + # Use mp.nstr to get full precision (40 digits for quad precision) |
| 1775 | + mpmath_result = QuadPrecision(mp.nstr(mpmath_result, 40)) |
| 1776 | + |
| 1777 | + # Handle NaN cases |
| 1778 | + if np.isnan(mpmath_result): |
| 1779 | + assert np.isnan(quad_result), f"Expected NaN for {op}({val}), got {quad_result}" |
| 1780 | + return |
| 1781 | + |
| 1782 | + # Handle infinity cases |
| 1783 | + if np.isinf(mpmath_result): |
| 1784 | + assert np.isinf(quad_result), f"Expected inf for {op}({val}), got {quad_result}" |
| 1785 | + assert np.sign(mpmath_result) == np.sign(quad_result), f"Infinity sign mismatch for {op}({val})" |
| 1786 | + return |
| 1787 | + |
| 1788 | + # For finite non-zero results |
| 1789 | + np.testing.assert_allclose(quad_result, mpmath_result, rtol=1e-32, atol=1e-34, |
| 1790 | + err_msg=f"Value mismatch for {op}({val}), expected {mpmath_result}, got {quad_result}") |
| 1791 | + |
| 1792 | + # their domain is [-1 , 1] |
| 1793 | + @pytest.mark.parametrize("op", ["asin", "acos"]) |
| 1794 | + @pytest.mark.parametrize("val", [ |
| 1795 | + # Basic cases (valid domain) |
| 1796 | + "0.0", "-0.0", "1.0", "-1.0", |
| 1797 | + # Small values |
| 1798 | + "1e-10", "-1e-10", "1e-15", "-1e-15", |
| 1799 | + # Values near domain boundaries |
| 1800 | + "0.9", "-0.9", "0.9999", "-0.9999", |
| 1801 | + "0.99999999", "-0.99999999", |
| 1802 | + "0.999999999999", "-0.999999999999", |
| 1803 | + # Fractional values (within domain) |
| 1804 | + "0.5", "-0.5", |
| 1805 | + # Special values |
| 1806 | + "nan" |
| 1807 | + ]) |
| 1808 | + def test_inverse_sin_cos(self, op, val): |
| 1809 | + mp.prec = 113 # Set precision to 113 bits (~34 decimal digits) |
| 1810 | + numpy_op = getattr(np, op) |
| 1811 | + mpmath_op = getattr(mp, op) |
| 1812 | + |
| 1813 | + quad_val = QuadPrecision(val) |
| 1814 | + mpf_val = mp.mpf(val) |
| 1815 | + |
| 1816 | + quad_result = numpy_op(quad_val) |
| 1817 | + mpmath_result = mpmath_op(mpf_val) |
| 1818 | + # convert mpmath result to quad for comparison |
| 1819 | + # Use mp.nstr to get full precision (40 digits for quad precision) |
| 1820 | + mpmath_result = QuadPrecision(mp.nstr(mpmath_result, 40)) |
| 1821 | + |
| 1822 | + # Handle NaN cases |
| 1823 | + if np.isnan(mpmath_result): |
| 1824 | + assert np.isnan(quad_result), f"Expected NaN for {op}({val}), got {quad_result}" |
| 1825 | + return |
| 1826 | + |
| 1827 | + # For finite non-zero results |
| 1828 | + np.testing.assert_allclose(quad_result, mpmath_result, rtol=1e-32, atol=1e-34, |
| 1829 | + err_msg=f"Value mismatch for {op}({val}), expected {mpmath_result}, got {quad_result}") |
| 1830 | + |
| 1831 | + # mpmath's atan2 does not follow IEEE standards so hardcoding the edge cases |
| 1832 | + # for special edge cases check reference here: https://en.cppreference.com/w/cpp/numeric/math/atan2.html |
| 1833 | + # atan2: [Real x Real] -> [-pi , pi] |
| 1834 | + @pytest.mark.parametrize("y", [ |
| 1835 | + # Basic cases |
| 1836 | + "0.0", "-0.0", "1.0", "-1.0", |
| 1837 | + # Small values |
| 1838 | + "1e-10", "-1e-10", "1e-15", "-1e-15", |
| 1839 | + # Medium/Large values |
| 1840 | + "10.0", "-10.0", "100.0", "-100.0", "1000.0", "-1000.0", |
| 1841 | + # Fractional |
| 1842 | + "0.5", "-0.5", "2.5", "-2.5", |
| 1843 | + # Special |
| 1844 | + "inf", "-inf", "nan", |
| 1845 | + ]) |
| 1846 | + @pytest.mark.parametrize("x", [ |
| 1847 | + "0.0", "-0.0", "1.0", "-1.0", |
| 1848 | + "1e-10", "-1e-10", |
| 1849 | + "10.0", "-10.0", "100.0", "-100.0", |
| 1850 | + "0.5", "-0.5", |
| 1851 | + "inf", "-inf", "nan", |
| 1852 | + ]) |
| 1853 | + def test_atan2(self, y, x): |
| 1854 | + mp.prec = 113 |
| 1855 | + |
| 1856 | + quad_y = QuadPrecision(y) |
| 1857 | + quad_x = QuadPrecision(x) |
| 1858 | + mpf_y = mp.mpf(y) |
| 1859 | + mpf_x = mp.mpf(x) |
| 1860 | + |
| 1861 | + quad_result = np.arctan2(quad_y, quad_x) |
| 1862 | + |
| 1863 | + # IEEE 754 special cases - hardcoded expectations |
| 1864 | + y_val = float(y) |
| 1865 | + x_val = float(x) |
| 1866 | + |
| 1867 | + # If either x is NaN or y is NaN, NaN is returned |
| 1868 | + if np.isnan(y_val) or np.isnan(x_val): |
| 1869 | + assert np.isnan(quad_result), f"Expected NaN for atan2({y}, {x}), got {quad_result}" |
| 1870 | + return |
| 1871 | + |
| 1872 | + # If y is ±0 and x is negative or -0, ±π is returned |
| 1873 | + if y_val == 0.0 and (x_val < 0.0 or (x_val == 0.0 and np.signbit(x_val))): |
| 1874 | + expected = quad_pi if not np.signbit(y_val) else -quad_pi |
| 1875 | + np.testing.assert_allclose(quad_result, expected, rtol=1e-32, atol=1e-34, |
| 1876 | + err_msg=f"Value mismatch for atan2({y}, {x}), expected {expected}, got {quad_result}") |
| 1877 | + return |
| 1878 | + |
| 1879 | + # If y is ±0 and x is positive or +0, ±0 is returned |
| 1880 | + if y_val == 0.0 and (x_val > 0.0 or (x_val == 0.0 and not np.signbit(x_val))): |
| 1881 | + assert quad_result == 0.0, f"Expected ±0 for atan2({y}, {x}), got {quad_result}" |
| 1882 | + assert np.signbit(quad_result) == np.signbit(y_val), f"Sign mismatch for atan2({y}, {x})" |
| 1883 | + return |
| 1884 | + |
| 1885 | + # If y is ±∞ and x is finite, ±π/2 is returned |
| 1886 | + if np.isinf(y_val) and np.isfinite(x_val): |
| 1887 | + expected = quad_pi / 2 if y_val > 0 else -quad_pi / 2 |
| 1888 | + np.testing.assert_allclose(quad_result, expected, rtol=1e-32, atol=1e-34, |
| 1889 | + err_msg=f"Value mismatch for atan2({y}, {x}), expected {expected}, got {quad_result}") |
| 1890 | + return |
| 1891 | + |
| 1892 | + # If y is ±∞ and x is -∞, ±3π/4 is returned |
| 1893 | + if np.isinf(y_val) and np.isinf(x_val) and x_val < 0: |
| 1894 | + expected = 3 * quad_pi / 4 if y_val > 0 else -3 * quad_pi / 4 |
| 1895 | + np.testing.assert_allclose(quad_result, expected, rtol=1e-32, atol=1e-34, |
| 1896 | + err_msg=f"Value mismatch for atan2({y}, {x}), expected {expected}, got {quad_result}") |
| 1897 | + return |
| 1898 | + |
| 1899 | + # If y is ±∞ and x is +∞, ±π/4 is returned |
| 1900 | + if np.isinf(y_val) and np.isinf(x_val) and x_val > 0: |
| 1901 | + expected = quad_pi / 4 if y_val > 0 else -quad_pi / 4 |
| 1902 | + np.testing.assert_allclose(quad_result, expected, rtol=1e-32, atol=1e-34, |
| 1903 | + err_msg=f"Value mismatch for atan2({y}, {x}), expected {expected}, got {quad_result}") |
| 1904 | + return |
| 1905 | + |
| 1906 | + # If x is ±0 and y is negative, -π/2 is returned |
| 1907 | + if x_val == 0.0 and y_val < 0.0: |
| 1908 | + expected = -quad_pi / 2 |
| 1909 | + np.testing.assert_allclose(quad_result, expected, rtol=1e-32, atol=1e-34, |
| 1910 | + err_msg=f"Value mismatch for atan2({y}, {x}), expected {expected}, got {quad_result}") |
| 1911 | + return |
| 1912 | + |
| 1913 | + # If x is ±0 and y is positive, +π/2 is returned |
| 1914 | + if x_val == 0.0 and y_val > 0.0: |
| 1915 | + expected = quad_pi / 2 |
| 1916 | + np.testing.assert_allclose(quad_result, expected, rtol=1e-32, atol=1e-34, |
| 1917 | + err_msg=f"Value mismatch for atan2({y}, {x}), expected {expected}, got {quad_result}") |
| 1918 | + return |
| 1919 | + |
| 1920 | + # If x is -∞ and y is finite and positive, +π is returned |
| 1921 | + if np.isinf(x_val) and x_val < 0 and np.isfinite(y_val) and y_val > 0.0: |
| 1922 | + expected = quad_pi |
| 1923 | + np.testing.assert_allclose(quad_result, expected, rtol=1e-32, atol=1e-34, |
| 1924 | + err_msg=f"Value mismatch for atan2({y}, {x}), expected {expected}, got {quad_result}") |
| 1925 | + return |
| 1926 | + |
| 1927 | + # If x is -∞ and y is finite and negative, -π is returned |
| 1928 | + if np.isinf(x_val) and x_val < 0 and np.isfinite(y_val) and y_val < 0.0: |
| 1929 | + expected = -quad_pi |
| 1930 | + np.testing.assert_allclose(quad_result, expected, rtol=1e-32, atol=1e-34, |
| 1931 | + err_msg=f"Value mismatch for atan2({y}, {x}), expected {expected}, got {quad_result}") |
| 1932 | + return |
| 1933 | + |
| 1934 | + # If x is +∞ and y is finite and positive, +0 is returned |
| 1935 | + if np.isinf(x_val) and x_val > 0 and np.isfinite(y_val) and y_val > 0.0: |
| 1936 | + assert quad_result == 0.0 and not np.signbit(quad_result), f"Expected +0 for atan2({y}, {x}), got {quad_result}" |
| 1937 | + return |
| 1938 | + |
| 1939 | + # If x is +∞ and y is finite and negative, -0 is returned |
| 1940 | + if np.isinf(x_val) and x_val > 0 and np.isfinite(y_val) and y_val < 0.0: |
| 1941 | + assert quad_result == 0.0 and np.signbit(quad_result), f"Expected -0 for atan2({y}, {x}), got {quad_result}" |
| 1942 | + return |
| 1943 | + |
| 1944 | + # For all other cases, compare with mpmath |
| 1945 | + mpmath_result = mp.atan2(mpf_y, mpf_x) |
| 1946 | + # Use mp.nstr to get full precision (40 digits for quad precision) |
| 1947 | + mpmath_result = QuadPrecision(mp.nstr(mpmath_result, 40)) |
| 1948 | + |
| 1949 | + if np.isnan(mpmath_result): |
| 1950 | + assert np.isnan(quad_result), f"Expected NaN for atan2({y}, {x}), got {quad_result}" |
| 1951 | + return |
| 1952 | + |
| 1953 | + if np.isinf(mpmath_result): |
| 1954 | + assert np.isinf(quad_result), f"Expected inf for atan2({y}, {x}), got {quad_result}" |
| 1955 | + assert np.sign(mpmath_result) == np.sign(quad_result), f"Infinity sign mismatch for atan2({y}, {x})" |
| 1956 | + return |
| 1957 | + |
| 1958 | + np.testing.assert_allclose(quad_result, mpmath_result, rtol=1e-32, atol=1e-34, |
| 1959 | + err_msg=f"Value mismatch for atan2({y}, {x}), expected {mpmath_result}, got {quad_result}") |
1738 | 1960 |
|
1739 | 1961 | @pytest.mark.parametrize("op", ["sinh", "cosh", "tanh", "arcsinh", "arccosh", "arctanh"]) |
1740 | 1962 | @pytest.mark.parametrize("val", [ |
|
0 commit comments