一元非线性方程求根的算法——二分法/牛顿迭代法

爱吃肉的小春

已于 2022-04-16 16:39:03 修改

阅读量3.4k

点赞数 9

分类专栏：作业

于 2022-04-05 18:38:31 首次发布

本文链接：https://blog.csdn.net/qq_53369920/article/details/123973894

版权

算法 c语言 c++

作业专栏收录该内容

4 篇文章 1 订阅

订阅专栏

大家好，我是慢慢努力的小刘，时隔七天终于更新了。今天分享的内容是：用代码实现二分法/牛顿迭代法求一元非线性方程求根。如果对大家有帮助，请大家动动小手点个赞或关注吧。

问题描述：

选择合适的初始有根区间或者初值，分别用二分法和牛顿迭代法求出一元非线性方程 $f(x)=x^3-sinx-4x+1$ 的全部实根，其中精度ε＝10-4 。

理论基础：

拿二分法/牛顿迭代法的思想是什么呢？二分法的基本思想是：通过逐步缩小隔根区间的长度，求出根的近似值。二分法只能求实根，不能求复根及偶重根。迭代法的基本思想是：建立一个迭代序列来逼近根。迭代效果（指迭代的敛散性、收敛的速度）与迭代初值、迭代公式的收敛性有关。迭代初始值一般根据高等数学的知识或作图的方法来确定；迭代公式通过对已知方程的某种变换得到。

解决问题：

数值法求非线性方程的跟一般分三步：

判定根的存在性；
确定根的分布范围，找出隔根区间；
根的精确化，在隔根区间内找出满足精度要求的根的近似值。

确定根的范围：

利用Visual Studio 2010（2019版本都一样）中的MFC或者MATLAB来画出方程图，我这里是用 Visual Studio 2010中的MFC画图的，下面是详细步骤：

步骤①〔创建项目〕：打开Visual Studio 2010，选择“文件”－“新建”－“项目”－“Visual C++”－“MFC”－“MFC应用程序”；输入项目名称“Exp1”。

点击“确定”按钮，得到：

点击左边的“应用程序类型”，在右边的应用程序类型选择“单个文档”

点击“完成”，得到一个项目“Exp1”。

步骤②〔创建自定义函数〕：选择“类试图”面板，如下图（左）所示。右键点击“CExp1View”结点，选择“添加”，如下图（右）所示。

在“添加”窗口选择“添加函数”，输入返回类型“void”，函数名“drawCure”。在参数定义部分，输入参数类型“int”，参数名“N”，点击“添件”按钮；再输入参数类型“CDC*”，参数名“pDC”，点击“添加”，完成函数定义。

上述操作将在项目中添件函数“drawCure”的声明与定义。

注意在类的头部引入 math.h：

#include <math.h>

完整的drawCure函数内容如下：

void CEXp1View::drawCure(int N, CDC* pDC)

{

    pDC->TextOutW(250,10,_T("非线性函数 f(x)=x^3-sin(x)-4x+1"));//文本输出

    //设置画笔，将影响画线的样式

    CPen *pnewPen,*poldPen;

    pnewPen=new CPen();

    pnewPen->CreatePen(PS_SOLID,3,RGB(0,0,0));

    poldPen=pDC->SelectObject(pnewPen);

       //画坐标系

    pDC->MoveTo(0,384);    pDC->LineTo(1366,384);//x轴

    pDC->MoveTo(683,0);    pDC->LineTo(683,768);//Y轴

    //画刻度

1   int Ex =50, Ey =50;//坐标放大倍数

    CString str; //输出刻度值用

    for (int i=-8;i<=8;i++){

       int xx = i*Ex+683;

       pDC->MoveTo(xx,384-4);

       pDC->LineTo(xx,384);

       str.Format(_T("%d"),i);

       pDC->TextOutW(xx,384+4,str);

    }

    for (int i=-6;i<=6;i++){

       if (i==0){

               continue;     //原点已经标定

       }

       int yy = -i*Ey+384;

       pDC->MoveTo(683,yy);

       pDC->LineTo(683+4,yy);

       str.Format(_T("%d"),i);

       pDC->TextOutW(683-20,yy,str);

    }

    //画曲线,先移动到本曲线的第一个点

    //原始坐标

    float x1 = -4;

    float y1 = x1*x1*x1-sin(x1)-4*x1+1;

    //放大

    x1 = x1*Ex;

    y1 = -y1*Ey;

    //取整

    int xx1 = (int)(x1+683);

    int yy1 = (int)(y1+384);

    //移动到第一个点

    pDC->MoveTo(xx1,yy1);

    for(int i=0;i<=N;i++){

                  //原始坐标

       float x0=-4+8.0/N*i;

       float y0= x0*x0*x0-sin(x0)-4*x0+1 ;

               //坐标放大

               x0=x0*Ex;

       y0=-y0*Ey;

               //坐标转化为整型

       int x1=(int)(x0+683);

       int y1=(int)(y0+384);

       pDC->LineTo(x1,y1);

    }

    //重置画笔

    pDC->SelectObject(poldPen);

}

在“CExp1View”类的“OnDraw(CDC* pDC)”函数中加入代码〔注意修改输入参数，去掉注释〕：

点击下图中的三角形，运行项目

运行结果:

由此得到三个隔根区间：（－3,-2）、（0,0.5)和（1.5，2.5）

二份法的代码：

头文件的内容(VS2019)：

#pragma once
#include<iostream>
#include<math.h>
#include<iomanip>
using namespace std;
class bnary
{
public:
	bnary(double a,double b);
	~bnary();
	int sign(double s);//判断函数值的正负;
	double reustbnary(double left, double right,double(*p)(double));//方法代码块
	void show(double left, double right, int c, double d, int count);//c是传入sign函数的返回值，d是区间的中点值;
	void show2(double a, double b);

private:
	double a;
	double b;
};

cpp文件：

#include"ErFen.h"
double e=pow(10,-4);
bnary::bnary(double a,double b)
{
	this->a = a;
	this->b = b;
}
bnary::~bnary()
{

}
double Fx(double x)
{
	return x * x * x - sin(x) - 4 * x + 1;
}
int bnary::sign(double s)
{
	if (s > 0)
	{
		return 1;//代表符号为'+'
	}
	else if (s < 0)
	{
		return -1;//符号维‘-’;
	}
	else
	{
		return 0;
	}
}

void bnary::show(double left, double right, int c, double d,int count)
{
	cout << setiosflags(ios::fixed) << setprecision(6) << endl;
	cout<<setw(5) <<count<<setw(23) << left << setw(15) << right << setw(12) << d << setw(6) << c << setw(18) << fabs(left - right)/2<<endl;
	//cout << setw(5)<< count <<setw(23)<<left<<setw(14)<<right<<setw(6)<<d<<setw(6)<<c<<setw(18)<< fabs(left - right)/2<< endl;
}
double bnary::reustbnary(double left, double right,double(*p)(double))
{
	int c;//储存函数符号的返回值;
	double mid;//区间的中点值;
	mid = (left + right) / 2;
	int count = 0;
	c = sign(Fx(mid));//计算中点的正负;
	while (1)//非递归
	{
		mid = (left + right) / 2;
		c = sign(Fx(mid));//计算中点的正负;
		if (c == sign(Fx(left)))
		{
			count++;
			this->show(left, right, c, mid,count);
			left = mid;
		}
		else
		{
			count++;
			this->show(left, right, c, mid,count);
			right = mid;
		}
		if (fabs(right - left)< e)
		{
			return mid;
		}
	}
	//if (fabs(left - right)< e)//递归 
	//{
	//	//this->show(left, right, c, mid);
	//	return mid;
	//}
	//if (c == sign(result(left)))//符号与区间左端点一样;
	//{
	//	this->show(left,right,c,mid);
	//	return reustbnary(mid, right);
	//}
	//if (c == sign(result(right)))
	//{
	//	this->show(left, right, c, mid);
	//	return reustbnary(left,mid);
	//}
}
void bnary::show2(double a, double b)
{
	cout << "初始区间为：[" << a << "," << b << "]" << endl;
	cout << setiosflags(ios::fixed) << setprecision(6) << endl;
	cout << "f(a)=" << Fx(a) << setw(10) << "f(b)=" << Fx(b)<< endl << endl;
	cout << "二分的次数" << setw(16) << "[a,b]-a" << setw(16) << "[a,b]-b" << setw(10) << "mid" << setw(12) << "f(x)+/-" << setw(15) << "1/2|b-a|" << endl;
}
int main()
{
	bnary A(-3.0, -2.0), B(0.0, 1.0), C(1.0, 2.0);
	cout << "二分法求根结果如下：" << endl;
	A.show2(-3.0, -2.0);
	double s;
	s= A.reustbnary(-3.0, -2.0, Fx);
	cout << "the root is:" << setw(12) << s << endl << endl << endl;
	B.show2(0.0, 1.0);
	s = B.reustbnary(0.0, 1.0,Fx);
	cout << "the root is:" << setw(12) << s << endl << endl << endl;
	C.show2(1.0, 2.0);
	s=C.reustbnary(1.0, 2.0,Fx);
	cout << "the root is:" << setw(12) << s << endl << endl << endl;
	return 0;
}

牛顿迭代法的代码:

头文件：

#pragma once
#include<iostream>
#include<math.h>
#include<iomanip>
using namespace std;
double e = pow(10, -4);
class NiuDon
{
public:
	NiuDon(double a);
	double result(double x, double(*p)(double));
	~NiuDon();

private:
	double a;
};

cpp代码：

#include"NewTon.h"
NiuDon::NiuDon(double x)
{
	this->a = x;
}

NiuDon::~NiuDon()
{

}
double Fx(double x)
{
	return (pow(x, 3) - sin(x) - 4 * x + 1) / (3 * pow(x, 2) - cos(x) - 4);
}
double NiuDon::result(double x,double(*p)(double))
{
	int count = 0;
	double x1, x2 = 0.0;
	x1 = x;
	while (1)
	{
		count++;
		double c = 3 * pow(x1, 2) - cos(x) - 4;
		if (c != 0)
		{
			double d = Fx(x1);
			x2 = x1 - d;
			if (fabs(x2 - x1) < e)
			{
				cout << setiosflags(ios::fixed) << setprecision(6) << endl;
				cout<<"第"<<count<<"次迭代"<<setw(12) << x2 << setw(16) << "|x2-x1|=" << fabs(x2 - x1) << endl;
				return x2;
				break;
			}
			else
			{
				cout << setiosflags(ios::fixed) << setprecision(6) << endl;
				cout << "第" << count << "次迭代" << setw(12) << x2 << setw(16) << "|x2-x1|=" << fabs(x2 - x1) << endl;
				x1 = x2;
				x2 = 0.0;
			}
		}
		else
		{
			cout << "x1=" << x1 << "的值歧异" << endl;
			break;
		}
	}
}
int main()
{
	NiuDon A(-2.5), B(0.5), C(1.5);
	double s = 0.0;
	cout <<"迭代次数"<<setw(15) << "迭代初值:"<< setw(16) << "|x2-x1|" << endl;
	s = A.result(-2.5,Fx);
	cout << "root=" << s << endl << endl;
	cout << "迭代次数" << setw(15) << "迭代初值:"<< setw(16) << "|x2-x1|" << endl;
	s = B.result(0.5, Fx);
	cout << "root=" << s << endl << endl;
	cout << "迭代次数"<< setw(15) << "迭代初值:"<< setw(16) << "|x2-x1|" << endl;
	s = C.result(1.5, Fx);
	cout << "root" << s << endl;
	return 0;
}