@samitouri / QOSamiQemu / commits / d2a727e347

fpu: Drop parts_log2

Use parts64_log2 at each call site. That leaves parts128_log2 unused, so move the whole function back to softfloat.c and specialize for FloatParts64. Reviewed-by: Philippe Mathieu-Daudé <philmd@linaro.org> Signed-off-by: Richard Henderson <richard.henderson@linaro.org>

Richard Henderson committed Apr 26, 2026 at 20:40 UTC d2a727e34732e57fb40d9f42d0c09a6d3a6c6ddf

2 files changed +125 -139

fpu/softfloat-parts.c.inc

-131

@@ -1671,134 +1671,3 @@ static void partsN(scalbn)(FloatPartsN *a, int n, float_status *s)

1671 g_assert_not_reached();

1672 }

1673 }

1674 -

1675 -/*

1676 - * Return log2(A)

1677 - */

1678 -static void partsN(log2)(FloatPartsN *a, float_status *s, const FloatFmt *fmt)

1679 -{

1680 - uint64_t a0, a1, r, t, ign;

1681 - FloatPartsN f;

1682 - int i, n, a_exp, f_exp;

1683 -

1684 - if (unlikely(a->cls != float_class_normal)) {

1685 - switch (a->cls) {

1686 - case float_class_denormal:

1687 - if (!a->sign) {

1688 - /* -ve denormal will be InvalidOperation */

1689 - float_raise(float_flag_input_denormal_used, s);

1690 - }

1691 - break;

1692 - case float_class_snan:

1693 - case float_class_qnan:

1694 - partsN(return_nan)(a, s);

1695 - return;

1696 - case float_class_zero:

1697 - float_raise(float_flag_divbyzero, s);

1698 - /* log2(0) = -inf */

1699 - a->cls = float_class_inf;

1700 - a->sign = 1;

1701 - return;

1702 - case float_class_inf:

1703 - if (unlikely(a->sign)) {

1704 - goto d_nan;

1705 - }

1706 - return;

1707 - default:

1708 - g_assert_not_reached();

1709 - }

1710 - }

1711 - if (unlikely(a->sign)) {

1712 - goto d_nan;

1713 - }

1714 -

1715 - /* TODO: This algorithm looses bits too quickly for float128. */

1716 - g_assert(N == 64);

1717 -

1718 - a_exp = a->exp;

1719 - f_exp = -1;

1720 -

1721 - r = 0;

1722 - t = DECOMPOSED_IMPLICIT_BIT;

1723 - a0 = a->frac_hi;

1724 - a1 = 0;

1725 -

1726 - n = fmt->frac_size + 2;

1727 - if (unlikely(a_exp == -1)) {

1728 - /*

1729 - * When a_exp == -1, we're computing the log2 of a value [0.5,1.0).

1730 - * When the value is very close to 1.0, there are lots of 1's in

1731 - * the msb parts of the fraction. At the end, when we subtract

1732 - * this value from -1.0, we can see a catastrophic loss of precision,

1733 - * as 0x800..000 - 0x7ff..ffx becomes 0x000..00y, leaving only the

1734 - * bits of y in the final result. To minimize this, compute as many

1735 - * digits as we can.

1736 - * ??? This case needs another algorithm to avoid this.

1737 - */

1738 - n = fmt->frac_size * 2 + 2;

1739 - /* Don't compute a value overlapping the sticky bit */

1740 - n = MIN(n, 62);

1741 - }

1742 -

1743 - for (i = 0; i < n; i++) {

1744 - if (a1) {

1745 - mul128To256(a0, a1, a0, a1, &a0, &a1, &ign, &ign);

1746 - } else if (a0 & 0xffffffffull) {

1747 - mul64To128(a0, a0, &a0, &a1);

1748 - } else if (a0 & ~DECOMPOSED_IMPLICIT_BIT) {

1749 - a0 >>= 32;

1750 - a0 *= a0;

1751 - } else {

1752 - goto exact;

1753 - }

1754 -

1755 - if (a0 & DECOMPOSED_IMPLICIT_BIT) {

1756 - if (unlikely(a_exp == 0 && r == 0)) {

1757 - /*

1758 - * When a_exp == 0, we're computing the log2 of a value

1759 - * [1.0,2.0). When the value is very close to 1.0, there

1760 - * are lots of 0's in the msb parts of the fraction.

1761 - * We need to compute more digits to produce a correct

1762 - * result -- restart at the top of the fraction.

1763 - * ??? This is likely to lose precision quickly, as for

1764 - * float128; we may need another method.

1765 - */

1766 - f_exp -= i;

1767 - t = r = DECOMPOSED_IMPLICIT_BIT;

1768 - i = 0;

1769 - } else {

1770 - r |= t;

1771 - }

1772 - } else {

1773 - add128(a0, a1, a0, a1, &a0, &a1);

1774 - }

1775 - t >>= 1;

1776 - }

1777 -

1778 - /* Set sticky for inexact. */

1779 - r |= (a1 || a0 & ~DECOMPOSED_IMPLICIT_BIT);

1780 -

1781 - exact:

1782 - partsN(sint_to_float)(a, a_exp, 0, s);

1783 - if (r == 0) {

1784 - return;

1785 - }

1786 -

1787 - memset(&f, 0, sizeof(f));

1788 - f.cls = float_class_normal;

1789 - f.frac_hi = r;

1790 - f.exp = f_exp - frac_normalize(&f);

1791 -

1792 - if (a_exp < 0) {

1793 - partsN(sub_normal)(a, &f);

1794 - } else if (a_exp > 0) {

1795 - partsN(add_normal)(a, &f);

1796 - } else {

1797 - *a = f;

1798 - }

1799 - return;

1800 -

1801 - d_nan:

1802 - float_raise(float_flag_invalid, s);

1803 - partsN(default_nan)(a, s);

1804 -}

fpu/softfloat.c

+125 -8

@@ -779,12 +779,6 @@ static float128 QEMU_FLATTEN float128_pack_raw(const FloatParts128 *p)

779 FloatParts128 *: parts128_##NAME, \

780 FloatParts256 *: parts256_##NAME)

781

782 -static void parts64_log2(FloatParts64 *a, float_status *s, const FloatFmt *f);

783 -static void parts128_log2(FloatParts128 *a, float_status *s, const FloatFmt *f);

784 -

785 -#define parts_log2(A, S, F) \

786 - PARTS_GENERIC_64_128(log2, A)(A, S, F)

787 -

782 /*

783 * Helper functions for softfloat-parts.c.inc, per-size operations.

784 */

@@ -4851,12 +4845,135 @@ floatx80 floatx80_sqrt(floatx80 a, float_status *s)

4845 /*

4846 * log2

4847 */

4848 +

4849 +static void parts64_log2(FloatParts64 *a, float_status *s, const FloatFmt *fmt)

4850 +{

4851 + uint64_t a0, a1, r, t, ign;

4852 + int i, n, a_exp, f_exp;

4853 +

4854 + if (unlikely(a->cls != float_class_normal)) {

4855 + switch (a->cls) {

4856 + case float_class_denormal:

4857 + if (!a->sign) {

4858 + /* -ve denormal will be InvalidOperation */

4859 + float_raise(float_flag_input_denormal_used, s);

4860 + }

4861 + break;

4862 + case float_class_snan:

4863 + case float_class_qnan:

4864 + parts64_return_nan(a, s);

4865 + return;

4866 + case float_class_zero:

4867 + float_raise(float_flag_divbyzero, s);

4868 + /* log2(0) = -inf */

4869 + a->cls = float_class_inf;

4870 + a->sign = 1;

4871 + return;

4872 + case float_class_inf:

4873 + if (unlikely(a->sign)) {

4874 + goto d_nan;

4875 + }

4876 + return;

4877 + default:

4878 + g_assert_not_reached();

4879 + }

4880 + }

4881 + if (unlikely(a->sign)) {

4882 + goto d_nan;

4883 + }

4884 +

4885 + a_exp = a->exp;

4886 + f_exp = -1;

4887 +

4888 + r = 0;

4889 + t = DECOMPOSED_IMPLICIT_BIT;

4890 + a0 = a->frac_hi;

4891 + a1 = 0;

4892 +

4893 + n = fmt->frac_size + 2;

4894 + if (unlikely(a_exp == -1)) {

4895 + /*

4896 + * When a_exp == -1, we're computing the log2 of a value [0.5,1.0).

4897 + * When the value is very close to 1.0, there are lots of 1's in

4898 + * the msb parts of the fraction. At the end, when we subtract

4899 + * this value from -1.0, we can see a catastrophic loss of precision,

4900 + * as 0x800..000 - 0x7ff..ffx becomes 0x000..00y, leaving only the

4901 + * bits of y in the final result. To minimize this, compute as many

4902 + * digits as we can.

4903 + * ??? This case needs another algorithm to avoid this.

4904 + */

4905 + n = fmt->frac_size * 2 + 2;

4906 + /* Don't compute a value overlapping the sticky bit */

4907 + n = MIN(n, 62);

4908 + }

4909 +

4910 + for (i = 0; i < n; i++) {

4911 + if (a1) {

4912 + mul128To256(a0, a1, a0, a1, &a0, &a1, &ign, &ign);

4913 + } else if (a0 & 0xffffffffull) {

4914 + mul64To128(a0, a0, &a0, &a1);

4915 + } else if (a0 & ~DECOMPOSED_IMPLICIT_BIT) {

4916 + a0 >>= 32;

4917 + a0 *= a0;

4918 + } else {

4919 + goto exact;

4920 + }

4921 +

4922 + if (a0 & DECOMPOSED_IMPLICIT_BIT) {

4923 + if (unlikely(a_exp == 0 && r == 0)) {

4924 + /*

4925 + * When a_exp == 0, we're computing the log2 of a value

4926 + * [1.0,2.0). When the value is very close to 1.0, there

4927 + * are lots of 0's in the msb parts of the fraction.

4928 + * We need to compute more digits to produce a correct

4929 + * result -- restart at the top of the fraction.

4930 + * ??? This is likely to lose precision quickly, as for

4931 + * float128; we may need another method.

4932 + */

4933 + f_exp -= i;

4934 + t = r = DECOMPOSED_IMPLICIT_BIT;

4935 + i = 0;

4936 + } else {

4937 + r |= t;

4938 + }

4939 + } else {

4940 + add128(a0, a1, a0, a1, &a0, &a1);

4941 + }

4942 + t >>= 1;

4943 + }

4944 +

4945 + /* Set sticky for inexact. */

4946 + r |= (a1 || a0 & ~DECOMPOSED_IMPLICIT_BIT);

4947 +

4948 + exact:

4949 + parts64_sint_to_float(a, a_exp, 0, s);

4950 + if (r != 0) {

4951 + FloatParts64 f = {

4952 + .cls = float_class_normal, .frac = r

4953 + };

4954 + f.exp = f_exp - frac_normalize(&f);

4955 +

4956 + if (a_exp < 0) {

4957 + parts64_sub_normal(a, &f);

4958 + } else if (a_exp > 0) {

4959 + parts64_add_normal(a, &f);

4960 + } else {

4961 + *a = f;

4962 + }

4963 + }

4964 + return;

4965 +

4966 + d_nan:

4967 + float_raise(float_flag_invalid, s);

4968 + parts64_default_nan(a, s);

4969 +}

4970 +

4971 float32 float32_log2(float32 a, float_status *status)

4972 {

4973 FloatParts64 p;

4974

4975 float32_unpack_canonical(&p, a, status);

4859 - parts_log2(&p, status, &float32_params);

4976 + parts64_log2(&p, status, &float32_params);

4977 return float32_round_pack_canonical(&p, status);

4978 }

4979

@@ -4865,7 +4982,7 @@ float64 float64_log2(float64 a, float_status *status)

4982 FloatParts64 p;

4983

4984 float64_unpack_canonical(&p, a, status);

4868 - parts_log2(&p, status, &float64_params);

4985 + parts64_log2(&p, status, &float64_params);

4986 return float64_round_pack_canonical(&p, status);

4987 }

4988