Интерполяция экстраполирует, ой!

[kodt_rsdn]

Compiler Explorer - C++

int main() { double l = 0.00619730599999999965; double r = 0.18464830610000002031; double ratio = 1.0; double res = l + (r - l) * ratio; if (res > r) { std::cout << std::fixed << std::setprecision(20) << "l = " << l << std::endl << "r = " << r << std::endl << "(r-l) = " << (r - l) << std::endl <<...

godbolt.org

Внезапно

Comments

[fat_crocodile]

ну да, операции в числах с плавающей запятой не ассоциативны, от порядка завися.

прочитал про lerp, прикольно, не знал про такую. Видимо там просто проверка стоит if (t == 1) return b; или как-то так.

[kodt_rsdn]

Скорее всего, да, проверка.

Я сходу не смог заставить ни лерп, ни ручную формулу получить ошибку при ratio < 1, сколь угодно близкой к 1.

А не сходу - было немножко лень. Очень может быть, что ошибка денормализации D никогда не превосходит ошибку умножения M, имеющую противоположный знак?

Поэтому (x1-x0) * (1-E) ≈ (x1-x0+D) * (1-E) = (x1-x0)*(1-E) + D-DE ≈ (x1-x0)*(1-E) + D-DE-M ≈ (x1-x0)*(1-E)

[kodt_rsdn]

Операции не просто не ассоциативны, они ещё и не дистрибутивны :)

Это не поле, это даже не полукольцо...

x0 + (x1-x0)*r = x0*(1-r) + x1*r = x0+x1*r-x0*r кажется более устойчиво к денормализации... но при r близким к 0 или 1 тоже может шуметь.

Может быть, лерп как-то умеет более тщательно компенсировать шумы? Например, выбирать наиболее подходящий порядок вычислений в зависимости от знака и упорядоченности.

Э, не хочу влезать с головой в теорию ошибок.

[kodt_rsdn]

Во, смотри, что нашёл.

https://github.com/emsr/cxx_linear?tab=readme-ov-file

Там действительно выбираются разные стратегии в зависимости от знаков.

И используется полезная функция fma (floatpoint multiply add).

https://en.cppreference.com/w/cpp/numeric/math/fma

Кстати, у функции интерполяции середины - std::midpoint - тоже весёлая реализация.