计算机中的浮点数
计算机中以固定长度存储浮点数的方式,造成了浮点数运算过程容易产生上溢和下溢。以float32为例, 其标记位占1bit,指数位占8bit,小数部分占23bit
经典下溢场景
不满足精度导致截断误差
#include <iostream>
#include <iomanip>
using namespace std;
int main() {
float a = 1.f;
float eps = 1e-8f;
float c = a + eps;
cout << setprecision(16) << a << "+" << eps << "=" << c << endl;
return 0;
}
微小的误差很容易被放大
这里以二元一次方程的求根为例
a
x
2
+
b
x
+
c
=
0
ax^{2}+bx+c=0
ax2+bx+c=0
根据基础数学知识,你会给出这样一个解决方案
x
1
,
2
=
−
b
±
b
2
−
4
a
c
2
a
x_{1,2} = \frac{-b\pm \sqrt {b^2-4ac}}{2a}
x1,2=2a−b±b2−4ac
由此设计程序
#include<iostream>
#include<cmath>
using namespace std;
int main() {
float a, b, c;
float temp, root, r1, r2;
cout << "该程序用于求一元二次方程ax^2+bx+c=0的解" << endl;
cout << "请依次输入a b c的值(a不要为0)" << endl;
cin >> a >> b >> c;
temp = b * b - 4 * a * c;
root = sqrt(temp);
if (temp < 0) {
cout << "改方程无解\n";
return -1;
}
r1 = (-b + root) / 2.0 / a;
r2 = (-b - root) / 2.0 / a;
cout << "一元二次方程的解为:" << r1 << " , " << r2 << endl;
return 0;
}
使用这段程序输入
0.01
,
1000
,
1
0.01,1000,1
0.01,1000,1
输出
−
0.00305176
,
−
100000
-0.00305176 , -100000
−0.00305176,−100000
而保留小数点后六位有效数字应该是
−
0.001000
,
−
99999.990000
-0.001000,-99999.990000
−0.001000,−99999.990000
此时,第一项的相对误差为百分之两百,而第二项的相对误差为千万分之一。
显然,两个相近的数相减会使得运算后的有效位数变少,也就是在
a
,
c
a,c
a,c的值很小时, $ -b + \sqrt {b^2-4ac}$ 这一操作过后,使得实际结果的有效位变低了(或者说引入了较大的误差),并且这个误差会在后续的运算中被放大。
这时可以通过数学手段减少误差
既然
x
1
,
2
=
−
b
±
b
2
−
4
a
c
2
a
x_{1,2} = \frac{-b\pm \sqrt {b^2-4ac}}{2a}
x1,2=2a−b±b2−4ac
分子分母同时乘以
−
b
−
b
2
−
4
a
c
,
−
b
+
b
2
−
4
a
c
-b - \sqrt {b^2-4ac},-b + \sqrt {b^2-4ac}
−b−b2−4ac,−b+b2−4ac 可得
$x_{1} = \frac{2c}{-b - \sqrt {b^2-4ac}} , x_{2} =\frac{2c}{-b + \sqrt {b^2-4ac}} $
根据b 的正负性,两个求根公式各取一半得到新的算法
新算法的结果为
−
0.001
,
−
10000
-0.001 , -10000
−0.001,−10000
#include<iostream>
#include<cmath>
using namespace std;
int main() {
float a, b, c;
float temp, root, r1, r2;
cout << "该程序用于求一元二次方程ax^2+bx+c=0的解" << endl;
cout << "请依次输入a b c的值(a不要为0)" << endl;
cin >> a >> b >> c;
temp = b * b - 4 * a * c;
root = sqrt(temp);
if (temp < 0) {
cout << "改方程无解\n";
return -1;
}
if (b > 0) {
r1 = 2 * c / (-b - root);
r2 = (-b - root) / 2 / a;
}
else if (b < 0) {
r1 = (-b + root) / 2 / a;
r2 = 2 * c / (-b + root);
}
else {
temp = c / a;
r1 = sqrt(-temp);
r2 = -sqrt(-temp);
}
cout << "一元二次方程的解为:" << r1 << " , " << r2 << endl;
return 0;
}