15

Python は、既存の浮動小数点値の最下位ビットをインクリメントした結果の浮動小数点値を取得する関数を提供していますか?

std::nextafterC++11で追加された機能に似たものを探しています。

4

3 に答える 3

26

ここに 5 つの (実際には 4.5 分の 1) 可能な解決策があります。

解決策 1: Python 3.9 以降を使用する

2020 年 10 月にリリースされた Python 3.9 には、math.nextafterこの機能を直接提供する新しい標準ライブラリ関数が含まれていますmath.nextafter(x, math.inf)。正の無限大に向かって次の浮動小数点数を取得するために使用します。例えば：

>>> from math import nextafter, inf
>>> nextafter(100.0, inf)
100.00000000000001

メソッドによって提供される 16 進数表現を見ると、この関数が実際に次のフロートアップを生成していることを確認するのが少し簡単になります。float.hex

>>> 100.0.hex()
'0x1.9000000000000p+6'
>>> nextafter(100.0, inf).hex()
'0x1.9000000000001p+6'

math.ulpPython 3.9 では、値とゼロから離れた次の値の差を与える、密接に関連し、頻繁に役立つコンパニオン関数も導入されています。

>>> from math import ulp
>>> nextafter(100.0, inf) - 100.0
1.4210854715202004e-14
>>> ulp(100.0)
1.4210854715202004e-14

解決策 2: NumPy を使用する

Python 3.9 以降を持っていなくても、NumPy にアクセスできる場合は、numpy.nextafter. 通常の Pythonfloatの場合、セマンティクスは のセマンティクスと一致します (ただし、NumPy は Python よりずっとmath.nextafter前にこの機能を利用できるため、Python のセマンティクスは NumPy のセマンティクスと一致すると言った方が適切でしょう)。

>>> from numpy import nextafter, inf
>>> nextafter(100.0, inf)
100.00000000000001

解決策 3: C をnextafter自分でラップする

C はnextafter関数を指定しますmath.h(たとえば、C99 のセクション 7.12.11.3 を参照)。これはまさに、Python >= 3.9 がそのmathモジュールでラップして公開する関数です。Python 3.9 以降を持っていない場合は、ctypesまたはを使用cffiして C を動的に呼び出すnextafterか、単純なCythonラッパーまたは C を公開する Python C 拡張機能を記述できnextafterます。これを行う方法の詳細は、別の場所ですでに十分に説明されています。この質問に対する @ Endophageの回答、および同様の StackOverflow の質問に対するこの回答(この質問は重複として閉じられています)。

struct解決策 4:モジュールによるビット操作

Python が IEEE 754 浮動小数点を使用しているという (実際にはほとんど常に安全な) 仮定を作成する意思がある場合、提供する Python 関数を作成するのは非常に簡単nextafterです。すべてのコーナーケースを正しく行うには、少し注意が必要です。

IEEE 754 バイナリ浮動小数点形式は巧妙に設計されているため、ある浮動小数点数から「次の」浮動小数点数への移動は、ビット表現をインクリメントするのと同じくらい簡単です。これは、範囲内の任意の数値に対して機能し、[0, infinity)指数境界とサブノーマルを越えて機能します。完全な浮動小数点範囲をカバーするバージョンを生成するにはnextUp、負の数、無限大、nan、および負のゼロを含む 1 つの特殊なケースも処理する必要があります。nextUp以下は、 Pythonでの IEEE 754 の関数の標準準拠バージョンです。すべてのコーナーケースをカバーしています。

import math
import struct

def nextup(x):
    # NaNs and positive infinity map to themselves.
    if math.isnan(x) or (math.isinf(x) and x > 0):
        return x

    # 0.0 and -0.0 both map to the smallest +ve float.
    if x == 0.0:
        x = 0.0

    n = struct.unpack('<q', struct.pack('<d', x))[0]
    if n >= 0:
        n += 1
    else:
        n -= 1
    return struct.unpack('<d', struct.pack('<q', n))[0]

nextDownとthenの実装は次のnextAfterようになります。(これnextAfterは IEEE 754 で指定された関数ではないことに注意してください。そのため、IEEE の特殊値で何が起こるかについては、少し当て推量があります。ここでは、Python のdecimal.Decimalクラスが基づいている IBM Decimal Arithmetic 標準に従っています。)

def nextdown(x):
    return -nextup(-x)

def nextafter(x, y):
    # If either argument is a NaN, return that argument.
    # This matches the implementation in decimal.Decimal
    if math.isnan(x):
        return x
    if math.isnan(y):
        return y

    if y == x:
        return y
    elif y > x:
        return nextup(x)
    else:
        return nextdown(x)

(部分的な) 解決策 5: 浮動小数点演算

xが小さすぎない正の値であり、IEEE float754 binary64 形式とセマンティクスを想定する場合は、驚くほど簡単な解決策があります。次の float up からxis x / (1 - 2**-53)、次の float down からxisx * (1 - 2**-53)です。

さらに詳しく説明すると、次のすべてが当てはまるとします。

  • IEEE 754 コーナー ケース (ゼロ、無限大、サブノーマル、ナン) は気にしません。
  • IEEE 754 binary64 浮動小数点形式だけでなく、IEEE 754 binary64セマンティクスも想定できます。つまり、すべての基本的な算術演算は、現在の丸めモードに従って正しく丸められます。
  • さらに、現在の丸めモードが IEEE 754 のデフォルトの round-tie-to-even モードであると仮定できます。

次に、量1 - 2**-53は として正確に表現できfloat、正の非正規化 Python float が与えられるとx、x / (1 - 2**-53)に一致しnextafter(x, inf)ます。同様に、x * (1 - 2**-53)は に一致しますが、 が最小の正の正常値 でnextafter(x, -inf)あるコーナー ケースを除きます。x2**-1022

これを使用する際に注意すべき点が 1 つあります。式はシステムの数学ライブラリから your2**-53を呼び出すため、正しく丸められることを期待するのは一般的に安全ではありません。この定数を計算するより安全な方法はたくさんありますが、そのうちの 1 つは を使用することです。次に例を示します。powpowfloat.fromhex

>>> d = float.fromhex('0x1.fffffffffffffp-1')  # 1 - 2**-53, safely
>>> d
0.9999999999999999
>>> x = 100.0
>>> x / d  # nextup(x), or nextafter(x, inf)
100.00000000000001
>>> x * d  # nextdown(x), or nextafter(x, -inf)
99.99999999999999

これらのトリックは、正確な 2 の累乗のような厄介なケースを含め、浮動小数点数の通常の範囲で正しく機能します。

証明のスケッチ: が正の法線にx / d一致することを示すために、正しさに影響を与えずに 2 のべき乗でスケーリングできるため、証明では一般性を失うことなく を仮定できます。の正確な数学的値(浮動小数点数ではなく実数と見なす)を書くと、は に等しくなります。不等式 と組み合わせると、float は間隔内で正確に離れているため、 に最も近い float が であることを保証するのに十分であると結論付けることができます。の証明も同様です。nextafter(x, inf)x0.5 <= x < 1.0zx / dz - xx * 2**-53 / (1 - 2**-53)0.5 <= x <= 1 - 2**-532**-54 < z - x <= 2**-532**-53[0.5, 1.0]znextafter(x, inf)x * d

于 2012-05-03T06:10:52.267 に答える
2

アップデート：

これは重複した質問であることが判明しました（検索「c++ nextafter python」の結果＃2としてGoogleに表示されます）：Increment a python float value by the small possible amount

受け入れられた答えは、いくつかの確かな解決策を提供します。

元の答え:

確かにこれは完璧な解決策ではありませんが、cython を数行使用するだけで、既存の C++ 関数をラップして Python で使用することができます。以下のコードをコンパイルしましたが、ubuntu 11.10 ボックスで動作します。

まず、.pyx ファイル (私は nextafter.pyx と呼びます) は、C++ へのインターフェイスを定義します。

cdef extern from "cmath":
    float nextafter(float start, float to)

def pynextafter(start, to):
    cdef float float_start = float(start)
    cdef float float_to = float(to)
    result = nextafter(start, to)
    return result

次に、setup.py で拡張機能のビルド方法を定義します。

from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext 

ext_modules=[
    Extension("nextafter",
        ["nextafter.pyx"],
        libraries=[],
        library_dirs=[],
        include_dirs=[],
        language="c++",
    )
]

setup(
    name = "nextafter",
    cmdclass = {"build_ext": build_ext},
    ext_modules = ext_modules
)

それらが同じディレクトリにあることを確認してから、でビルドしpython setup.py build_ext --inplaceます。nextafter の他のバリエーション (ダブルなど) を拡張機能に追加する方法を理解していただければ幸いです。ビルドしたら、nextafter.so が必要です。同じディレクトリで python を起動する (またはパスのどこかに nextafter.so を配置する) と、 を呼び出すことができるはずですfrom nextafter import pynextafter。

楽しみ！

于 2012-05-02T20:19:23.543 に答える
0

http://docs.python.org/library/stdtypes.html#float.hexをチェックしてください

次の次のことをあまり知らない実装を試してみましょう。

まず、16 進文字列から 16 進部分と指数を抽出する必要があります。

def extract_parts(hex_val):
    if not hex_val.startswith('0x1.'):
        return None
    relevant_chars = hex_val[4:]
    if not len(relevant_chars) > 14 and relevant_chars[13] == 'p':
        return None
    hex_portion = int(relevant_chars[:13], 16)
    if relevant_chars[14] == '+':
        p_val = int(relevant_chars[15:])
    elif relevant_chars[14] == '-':
        p_val = -int(relevant_chars[15:])
    else:
        return None
    return (hex_portion, p_val)

次に、正または負の方向にインクリメントする方法が必要です (16 進文字列は既に整数に変換されていると仮定しますhex_portion)。

def increment_hex(hex_portion, p_val, direction):
    if hex_portion == 0 and direction == -1:
        new_hex = 'f' * 13
        p_val -= 1
    elif hex_portion == int('f' * 13, 16) and direction == 1:
        new_hex = '0' * 13
        p_val += 1
    else:
        new_hex = hex(hex_portion + direction)[2:].rstrip('L').zfill(13)

    if len(new_hex) != 13:
        return None
    return format_hex(new_hex, p_val)

上記で使用した、受け入れ可能な 16 進文字列と指数を再フォーマットするためのヘルパー関数が必要です。

def format_hex(hex_as_str, p_val):
    sign = '-' if p_val < 0 else '+'
    return '0x1.%sp%s%d' % (hex_as_str, sign, p_val)

最後に、実装するにはnextafter:

def nextafter(float_val):
    hex_equivalent = float_val.hex()
    hex_portion, p_val = extract_parts(hex_equivalent)

    direction = 1
    new_hex_equiv = increment_hex(hex_portion, p_val, direction)
    return float.fromhex(new_hex_equiv)
于 2012-05-02T20:09:22.213 に答える