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

最新推荐文章于 2023-08-28 10:07:14 发布

原创

最新推荐文章于 2023-08-28 10:07:14 发布 · 2.6k 阅读

13 ·

CC 4.0 BY-SA版权

本文介绍了一种使用C语言实现的一元三次方程实数根求解方法，采用盛金公式进行高精度计算。文章详细阐述了求解过程，并提供了完整的代码示例，包括对一元二次方程实数根的求解函数。

文章目录

一、问题描述
二、 $C$ 代码
三、参考文献/资料

一、问题描述

已知一元三次方程 $ax^3 + bx^2 + cx+d = 0(a\neq0)$ ，求方程的实数根。笔者用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[