关于对浮点数执行的操作的准确性。
我正在尝试编写一个返回第 n 个勒让德多项式的系数的函数。 // 存储以前的多项式以使递归过程更快 无序_地图 我正在尝试编写一个返回第 n 个勒让德多项式的系数的函数。 // store previous polynomials to make the recursive process faster unordered_map<int, vector<double>> legendre_cache; // Recursive function to compute the coefficients of the nth Legendre polynomial vector<double> legendre_polynomial(int n) { // uses Bonnet's recursion formula P_{n+1} = ( (2n+1) x P_n - n P_{n-1} )/(n+1) // Check if the result is already in the cache if (legendre_cache.find(n) != legendre_cache.end()) { return legendre_cache[n]; } // Base cases if (n == 0) { return { 1.0 }; // P0(x) = 1 } if (n == 1) { return { 1.0, 0.0 }; // P1(x) = x } // Get coefficients for P_(n-1) and P_(n-2) vector<double> Pn_minus_1 = legendre_polynomial(n - 1); vector<double> Pn_minus_2 = legendre_polynomial(n - 2); // Get the sizes of the vectors to n+1 Pn_minus_1.push_back(0.0); // P_{n-1} gets multiplied by x Pn_minus_2.insert(Pn_minus_2.begin(), 0.0); // P_{n-2} needs two leading digits Pn_minus_2.insert(Pn_minus_2.begin(), 0.0); // Calculate coefficients for P_n vector<double> Pn(n + 1, 0.0); // Initialize Pn with n+1 zeros // Using the recurrence relation to fill coefficients for (int k = 0; k <= n; k++) { Pn[k] += ((2 * n - 1) * Pn_minus_1[k] - (n - 1) * Pn_minus_2[k]) / n; } // Store the result in the cache before returning legendre_cache[n] = Pn; return Pn; } 该算法有效,因为我将其结果与 n = 10 以内的其他来源进行了比较。 然而,这个函数很快就会变得不准确,而且它的递归性质并没有给它带来好处: 在这里,我评估我的多项式并将其与“boost/math”库中的类似函数进行比较: int main() { cout << fixed << setprecision(15); // Example inputs int n = 48; double x = 0.9; cout << "Result of my function:" << endl; vector<double> coeffs = legendre_polynomial(n); double result = 0; int l = coeffs.size(); for (int i = l - 1; i >= 0; i--) { result += coeffs[l - i - 1] * pow(x, i); } cout << result << endl; cout << "Result of boost library:" << endl; cout << boost::math::legendre_p(n, x) << endl; return 0; } n 越大,我的函数就越不准确,直到 n = 45 左右时完全崩溃: 测试代码的输出: Result of my function: 0.151876441197771 Result of boost library: -0.106892037065314 我从 boost 网站 知道 boost 函数使用与我相同的递归方法。 因此,假设他们的功能有效,那么也应该可以改进我的功能。 感谢您的帮助。 编辑1: 需要明确的是,我评估多项式只是为了证明我的系数是错误的;我主要需要多项式的系数(勒让德高斯求积)。 据我了解库中的函数,我认为那里从未提到过系数;递归公式仅应用于特定的 x 值。 坦率地说,我还没有达到重写代码来提取系数的水平。 编辑2: 我所说的“该算法是有效的,因为我将它的结果与其他来源的结果进行了比较,最多为n = 10”,我怀疑如果没有准确性错误,该函数将会起作用。我认为由于舍入、除法、浮点导致的错误会因为递归而像滚雪球一样不断增加,直到函数变得不可用。 (我的问题是如何处理这个问题) 我认为您无法像您一样计算勒让德多项式在给定点的值,即直接使用多项式的系数。事实上,您将面临数值错误,因为当多项式的阶数很大时,double无法处理x与“巨大”系数的乘法。 相反,您可以直接依赖这些多项式的三项递推关系来计算给定点的值: #include <iostream> #include <iomanip> #include <boost/math/special_functions/legendre.hpp> int main(int argc, char* argv[]) { std::cout << std::fixed << std::setprecision(15); int n = argc>=2 ? atoi(argv[1]) : 48; double x = 0.9; double p0=1; double p1=x; for (size_t k=1; k<n; k++) { double tmp = ( (2*k+1)*x*p1 - k*p0) / (k+1); p0=p1; p1=tmp; } std::cout << "Result of my function:" << std::endl; std::cout << p1 << std::endl; std::cout << "Result of boost library:" << std::endl; std::cout << boost::math::legendre_p(n, x) << std::endl; return 0; } 可能的输出: Result of my function: -0.106892037065314 Result of boost library: -0.106892037065314
我用 0 或 1 填充 numpy 数组,具体取决于计算的数组索引值是大于还是小于给定值。 这会产生一个由 0 和 1 组成的数组。 Matplotlib 不这么认为 比较...
在一系列浮点算术运算之后与浮点数进行小于等于比较是否有“最佳实践”? 我在 R 中有以下示例(尽管问题
C++ 是否有比 float 或 double 更准确的数据类型,或者我是否只需要满足于我的计算将会失败的事实? 编辑:正如李斯特先生所指出的,我的问题是……
比 float 或 double 更准确的数据类型? C++
C++ 是否有比 float 或 double 更准确的数据类型,或者我是否只需要满足于我的计算将会失败的事实? 编辑:正如李斯特先生所指出的,我的问题是尊重......
#包括 无效 get_eps(浮点 *eps) { *每股收益=1.0f; while ((1.0f + *eps / 2.0f) != 1.0f) { *每股收益/= 2.0f; } } int main(){ 浮动 eps; get_eps(&eps);
所以我知道浮点精度(以及像 1.1 这样的东西如何不能用二进制精确表达)等等,但我想知道:那么,数学相关的库如何实现无限精度...
这是一个自我回答的问题,鼓励在 Stackoverflow 上分享知识。因此,这是我之前关于精确计算
这是一个自我回答的问题,鼓励在 Stackoverflow 上分享知识。因此,这是我之前关于精确计算原理的自我回答问题的后续......
(指数)缩放互补误差函数,通常由 erfcx 指定,在数学上定义为 erfcx(x) := ex2 erfc(x)。它经常出现在物理学中的扩散问题中,例如......
float 0.7f 的准确值为 0.69999... 所以我认为 0.7f * 100f 的结果会低于 70,比如 69.99999... 但结果恰好是 70。 浮点乘法涉及这样的吗
在java中将浮点值四舍五入到最接近的整数的最有效方法是什么?
我已经看到很多关于与舍入浮点值相关的SO的讨论,但没有考虑到效率方面的可靠问答。所以这里是: 舍入最有效(但正确)的方法是什么
Neumaier 求和是 Kahan 求和的改进,用于对浮点数数组进行精确求和。 将 numba 导入为 nb @nb.njit def neumaier_sum(arr): s = arr[0] c = 0.0 对于我在范围(1,...
SQL Server 如何在 FLOAT(53) 列中存储超过 2^24 的整数值?
我有一个 SQL Server 数据库,我可以从数据类型为 float 的列中读取超过 2^24 的整数值。这让我很困惑,因为这些值应该太大/太长才能准确......
在PHP中,我愿意在一些操作后比较浮点数,但它不能正确显示。 例如: $a = 0.2; 如果(($a - 0.2)=== 0) 返回真; 别的 返回假; 它回来了
12.7 - 20 + 7.3 = -8.8817841970013e-016 我在读《Programming in Lua》这本书时遇到了这个问题,我不知道如何让计算结果为零。 当我格式化资源时...
如何纠正 Floor(t / dt) * dt != t 中的数值不稳定?
我有两个浮点数 t 和 dt。典型值为 t = 1 和 dt = .01f。我想在从 0 到 t 的每个“时间”执行一个函数 foo,步长为 dt。如果碰巧 t 不可整除...
Python 内置 round() 函数在 2.4 和 2.7 之间的变化
Python 中内置的 round() 函数在 2.4 和 2.7 之间有变化吗? Python 2.4: Python 2.4.6(#1,2009 年 2 月 12 日,14:52:44) >>> f = 1480.39499999999998181010596454143524169921875 >&g...
ceil(double(int)*double) 或 ceil(float(int*double))?
我有这个C++代码: 整数x; 双 y; // 做某事 int z = ceil(x*y); 如果 x = 410,y = 1.1,则 z = 452,而它应该是 451。我更改为 z = ceil(float(x*y)) 并且它有效。我朋友建议...
[虽然这是一个自我回答的问题,但我会很乐意投票并接受任何替代答案,该答案要么为相同数量的计算工作提供卓越的准确性,要么减少金额...