一、问题描述

已知一元二次方程 a x 2 + b x + c = 0 ( a ≠ 0 ) ax^2 + bx + c = 0(a\neq0) ax2+bx+c=0(a=0)，求方程的实数根。在数学上，我们很容易可以得到方程的实数根为：
x 1 = − b + b 2 − 4 a c 2 a (1) x_1=\frac{-b+\sqrt{b^2-4ac}}{2a} \tag 1 x1=2a−b+b2−4ac (1)
x 2 = − b − b 2 − 4 a c 2 a (2) x_2=\frac{-b-\sqrt{b^2-4ac}}{2a} \tag 2 x2=2a−b−b2−4ac (2)
其中， b 2 − 4 a c ≥ 0 b^2-4ac\geq0 b2−4ac≥0
然而，在计算机上，由于计算机字长的限制，有些浮点数无法绝对精确地存储。并且，一些计算过程可能会导致精度损失。例如接下来这两个例子。
( A ) (A) (A)设有两个数 p = 3.1415926536 p=3.1415926536 p=3.1415926536和 q = 3.1415957341 q=3.1415957341 q=3.1415957341，两者几乎相等，而且都有11位有效数字精度。它们的差为 p − q = − 0.0000030805 p-q=-0.0000030805 p−q=−0.0000030805，只有5位有效数字精度！
( B ) (B) (B)设 f ( x ) = x ( x + 1 − x ) , g ( x ) = x x + 1 + x f(x)=x(\sqrt{x+1}-\sqrt{x}),g(x)=\frac{x}{\sqrt{x+1}+\sqrt{x}} f(x)=x(x+1 −x ),g(x)=x+1 +x x。显然，当 x ≥ 0 x\geq0 x≥0时， f ( x ) f(x) f(x)可以通过分母有理化得到 g ( x ) g(x) g(x)，在数学上两者是等价的。现在，我们采用6位有效数字精度和舍入法比较 f ( 500 ) f(500) f(500)与 g ( 500 ) g(500) g(500)的计算结果。 f ( 500 ) = 500 × ( 501 − 500 ) = 500 × ( 22.3830 − 22.3607 ) = 500 × 0.0223 = 11.1500 f(500)=500\times(\sqrt{501}-\sqrt{500})=500\times(22.3830-22.3607)=500\times0.0223=11.1500 f(500)=500×(501 −500 )=500×(22.3830−22.3607)=500×0.0223=11.1500。 g ( 500 ) = 500 501 + 500 = 500 22.3830 + 22.3607 = 500 44.7437 = 11.1748 g(500)=\frac{500}{\sqrt{501}+\sqrt{500}}=\frac{500}{22.3830+22.3607}=\frac{500}{44.7437}=11.1748 g(500)=501 +500 500=22.3830+22.3607500=44.7437500=11.1748。 g ( 500 ) g(500) g(500)比 f ( 500 ) f(500) f(500)计算精度更高！(真实值为11.174755300747198…，舍入到6位有效数字，与 g ( 500 ) g(500) g(500)的计算值相同）
精度损失的现象在我们编程过程中很有可能会出现，有时可能会使计算结果与正确结果大相径庭，产生灾难性的后果，我们必须引起足够重视。精度损失的现象有时候是没办法避免的，有时候则可以通过一些数值计算算法或技巧来避免。显然，式 ( 1 ) (1) (1)和式 ( 2 ) (2) (2)的分子在一定条件下，可能是两个相近的数相减，影响计算精度。那么，如何编写一个计算精度较高的程序来求一元二次方程的实数根呢？

二、推导步骤

编写一个计算精度较高的程序来求一元二次方程 a x 2 + b x + c = 0 ( a ≠ 0 ) ax^2 + bx + c = 0(a\neq0) ax2+bx+c=0(a=0)的实数根，至少需要考虑两种情形： ( A ) ∣ b ∣ ≈ b 2 − 4 a c ( B ) c = 0 (A)\left|b\right| \approx\sqrt{b^2-4ac}\ \ (B)c=0 (A)∣b∣≈b2−4ac (B)c=0
( A ) (A) (A)若 ∣ b ∣ ≈ b 2 − 4 a c \left|b\right| \approx\sqrt{b^2-4ac} ∣b∣≈b2−4ac ，当 b > 0 b>0 b>0时，式 ( 1 ) (1) (1)中， − b + b 2 − 4 a c -b+\sqrt{b^2-4ac} −b+b2−4ac 会导致计算精度损失，对式 ( 1 ) (1) (1)分子有理化得：
x 1 = − 2 c b + b 2 − 4 a c (3) x_1=\frac{-2c}{b+\sqrt{b^2-4ac}} \tag 3 x1=b+b2−4ac −2c(3)
当 b ≤ 0 b\leq0 b≤0时，式 ( 2 ) (2) (2)中， − b − b 2 − 4 a c -b-\sqrt{b^2-4ac} −b−b2−4ac 会导致计算精度损失，对式 ( 2 ) (2) (2)分子有理化得：
x 2 = − 2 c b − b 2 − 4 a c (4) x_2=\frac{-2c}{b-\sqrt{b^2-4ac}} \tag 4 x2=b−b2−4ac −2c(4)
( B ) (B) (B)若 c = 0 c=0 c=0，在数学上， b 2 − 4 a c = ∣ b ∣ \sqrt{b^2-4ac}=\left|b\right| b2−4ac =∣b∣，但计算机在数值计算过程中，对 b b b进行了开方和平方根运算，会引入计算误差(或截断误差)和舍入误差，导致此时 b 2 − 4 a c ≠ ∣ b ∣ \sqrt{b^2-4ac}\neq\left|b\right| b2−4ac =∣b∣，使方程根的计算误差变大。因而，当 c = 0 c=0 c=0时，需特殊处理，此时方程的两实数根为：

{ x 1 = 0 x 2 = − b a (5) \left\{ \begin{array}{c} x_1=0 \\ \tag 5 x_2=\frac{-b}{a}\end{array}\right. {x1=0x2=a−b(5)

三、 C C C代码

#include <math.h>
#include <float.h>
#include <stdio.h>#define UINT32 unsigned int
#define ERR_NO_ERROR 0x00000000
#define ERR_NAN 0x00000001
#define ERR_INF 0x00000002
#define MAX(a, b) ((a) > (b)) ? (a) : (b)
#define MIN(a, b) ((a) < (b)) ? (a) : (b)/*************************************************
Function: is_number
Description: 判断浮点数是否为nan
Input: 浮点数x
Output: 无
Return: 若浮点数x为nan返回0，否则返回1
Author: Marc Pony(marc_pony@163.com)
*************************************************/
int is_number(float x)
{ return (x == x);
}/*************************************************
Function: is_finite_number
Description: 判断浮点数是否为inf
Input: 浮点数x
Output: 无
Return: 若浮点数x为inf返回0，否则返回1
Author: Marc Pony(marc_pony@163.com)
*************************************************/
int is_finite_number(float x)
{ return (x >= -FLT_MAX && x <= FLT_MAX);
}/*************************************************
Function: solve_quadratic_equation
Description: 求一元二次方程(a*x^2 + b*x + c = 0)的所有实数根
Input: 方程的系数 p = {c, b, a}
Output: 方程的所有实数根x, 实数根的个数rootCount
Return: 错误号
Author: Marc Pony(marc_pony@163.com)
*************************************************/
UINT32 solve_quadratic_equation(float p[], float x[], int *rootCount)
{int i;float a, b, c, delta, sqrtDelta;const float ZERO = FLT_MIN;  // min normalized positive value（1.175494351e-38F）const float EPS = FLT_MIN;UINT32 errNo = ERR_NO_ERROR;*rootCount = 0;for (i = 0; i < 3; i++){if (!is_number(p[i])){errNo = ERR_NAN;return errNo;}if (!is_finite_number(p[i])){errNo = ERR_INF;return errNo;}}a = p[2];b = p[1];c = p[0];if (fabs(a - 0.0) < EPS){if (fabs(b - 0.0) > EPS){x[0] = -c / b;*rootCount = 1;}}else{b /= a;c /= a;a = 1.0;delta = b * b - 4.0 * a * c;if (delta > ZERO){if (fabs(c - 0.0) < EPS) //若c = 0,由于计算误差,sqrt(b*b - 4*a*c）不等于|b|{x[0] = 0.0;x[1] = -b / a;}else{sqrtDelta = sqrt(delta);if (b > 0.0){x[0] = (-2.0 * c) / (b + sqrtDelta);    //避免两个很接近的数相减,导致精度丢失x[1] = (-b - sqrtDelta) / (2.0 * a);}else{x[0] = (-b + sqrtDelta) / (2.0 * a);x[1] = (-2.0 * c) / (b - sqrtDelta);  //避免两个很接近的数相减,导致精度丢失}}*rootCount = 2;}else if (fabs(delta - 0.0) < EPS){x[0] = x[1] = -b / (2.0 * a);*rootCount = 2;}else{*rootCount = 0;}}return errNo;
}void main(void)
{float x[2], p[3];int rootCount;float x1, x2, a, b, c;UINT32 errNo = ERR_NO_ERROR;//一元二次方程测试//（1）(x - 1000)*(x - 0.001) = x^2 - 1000.001*x + 1 = 0a = 1;b = -1000.001;c = 1;//(2) 3*x ^ 2 - 1000000*x + 0 = 0//a = 3;//b = -1000000;//c = 0;//(3) 1.0e-10*x^2 - 2.0e-10*x + 1.0e-10 = 0//a = 1.0e-10;//b = -2.0e-10;//c = 1.0e-10;p[0] = c;p[1] = b;p[2] = a;errNo = solve_quadratic_equation(p, x, &rootCount);x1 = (-b - sqrt(b * b - 4.0 * a * c)) / (2.0 * a);x2 = (-b + sqrt(b * b - 4.0 * a * c)) / (2.0 * a);printf("原始方法求得方程的根:x1=%f, x2=%f\n本文方法求得方程的根:x1=%f, x2=%f\n",MIN(x1, x2), MAX(x1, x2), MIN(x[0], x[1]), MAX(x[0], x[1]));
}

(1)当 a = 1 , b = − 1000.001 , c = 1 a=1,b=-1000.001,c=1 a=1,b=−1000.001,c=1时，方程根的准确值为:0.001和1000，程序运行结果：

(2)当 a = 3 , b = − 1000000 , c = 0 a=3,b=-1000000,c=0 a=3,b=−1000000,c=0时，方程根的准确值为:0和333333.33…，程序运行结果：

(3)当 a = 1.0 e − 10 , b = − 2.0 e − 10 , c = 1.0 e − 10 a=1.0e^{-10},b=-2.0e^{-10},c=1.0e^{-10} a=1.0e−10,b=−2.0e−10,c=1.0e−10时，方程根的准确值为:1和1，程序运行结果：

四、总结

一元二次方程实数根的c语言实现，在很多人看来非常简单，根本不值一提，然而就是这么简单的问题，要写出计算精度高、简洁优雅、可读性强、健壮的代码也不是一件简单的事。

五、参考文献/资料

Numerical Methods Using MATLAB Fourth Edition. John Mathews, Kurtis D.Fink
数值方法(MATLAB版)(第四版) 周璐，陈渝，钱方等译

一元二次方程高精度实数根(C语言)相关推荐

一元三次方程高精度实数根(C语言)
文章目录一.问题描述二.CCC代码三.参考文献/资料一.问题描述已知一元三次方程ax3+bx2+cx+d=0(a≠0)ax^3 + bx^2 + cx+d = 0(a\neq0)ax3+ ...
java求一元二次方程“ax2＋bx＋c＝0，（a≠0）”的根。要求：通过命令行参数或者Scanner类输入a、b和c的值，借助Math类计算并输出该一元二次方程的实数根。
编程题:求一元二次方程"ax2+bx+c＝0,(a≠0)"的根.要求:通过命令行参数或者Scanner类输入a.b和c的值,借助Math类计算并输出该一元二次方程的实数根. 提示: ...
1205: 求一元二次方程的实数根（C）
一.题目 acm.wust.edu.cn/problem.php?id=1205&soj=0 二.分析一元二次方程有三个系数a.b.c,两个根x1.x2,以及d(德尔塔): a.b.c均为实 ...
C语言例题——求一元二次方程的实数根
#include "stdio.h" #include "math.h"void main() {//方程的系数int a, b, c;//△和方程的解doub ...
python中一元二次方程的判别式_一元二次方程难？根的判别式常考这四题型，若只知其一难考高分...
学习有方法,解题有技巧.对于九年级学生来说,现在的学习可以说是进入"白热化"阶段,如何在紧张的学习中脱颖而出.在同样的努力前提下,对重要考点进行全面总结,不留知识死角也是要提上日程 ...
一元三次方程的实数根
一元三次方程的实数根一元三次方程的一般形式一元三次方程的一般形式为 ax3+bx2+cx+d=0,a≠0. ax^3+bx^2+cx+d=0,\quad a\neq 0. 令 x=y−b3ax=y ...
js条件语句，用if...else if....else方程ax2+bx+c=0一元二次方程。求根
if 语句 - 只有当指定条件为 true 时,使用该语句来执行代码 if...else 语句 - 当条件为 true 时执行代码,当条件为 false 时执行其他代码 if...else if... ...
编写lisp程序解一元二次方程_怎样用C语言编一个解一元二次方程的程序？
展开全部用三个变量,来接一元二次方程的三个参数, 然后62616964757a686964616fe78988e69d8331333365633839根据这三个参数计算△=b*b-4ac,再根据△判 ...
2021.11.17 关于一元二次方程组求根
题目:定义一个名为equation的函数,实现求解一元二次方程组(标准形式)的根.equation接收四个参数,第一个参数是二次项的系数(a)第二个参数是一次项的系数(b)第三个参数是常数项(c)第四 ...

一元二次方程高精度实数根(C语言)

文章目录