digital.cpp

高阶微分方程分解成为两个方程后

//沈硕 自动化5班 3004203132

#include <iostream.h>
#include <math.h>
#include <iomanip.h>
#define Pi 3.1415926535
double fun1(double x,double y,double z)//高阶微分方程分解后的第一个方程
{
	return -98*sin(y);
}
double fun2(double x,double y,double z)//高阶微分方程分解后的第二个方程
{
	return z;
}
class Ola//改进欧拉法&龙格－库塔 解 高阶微分方程
{
	friend double fun1(double x,double y,double z);
	friend double fun2(double x,double y,double z);
private:double x,y1,y2,h;//y1改进欧拉法解，y2龙格－库塔解
		int num;//迭代次数
		double z1,z2;//z是中间变量
public:Ola(double x1,double yi,double zi)//赋初值
	   {
		   x=x1;y1=y2=yi;z2=z1=zi;
		   cout<<"初值："<<endl<<setw(20)<<"X="<<x<<setw(20)<<"Y="<<y1<<endl;
	   }
	   void num_get(int n,double x0,double xn)//迭代次数,步长
	   {
		   num=n;h=(xn-x0)/num;
	   }
	   void diedai()
	   {
		   cout<<"迭代次数"<<setw(10)<<"X"<<setw(20)<<"改进欧拉";
		   cout<<setw(20)<<"龙格－库塔"<<endl;
		   double yp,yc;
		   double zp,zc;//改进欧拉法中间变量
		   double k1,k2,k3,k4,l1,l2,l3,l4;//龙格－库塔法中间变量
		   
		   for(int i=0;i<num;i++)
		   {
			   //改进欧拉法迭代格式
			   zp=z1+h*fun1(x,y1,z1);yp=y1+h*fun2(x,y1,z1);
			   zc=z1+h*fun1(x+h,yp,zp);yc=y1+h*fun2(x+h,yp,zp);
			   z1=0.5*(zp+zc);y1=0.5*(yp+yc);
			   //龙格－库塔法迭代格式
			   k1=fun1(x,y2,z2);  l1=fun2(x,y2,z2);
			   k2=fun1(x+h/2,y2+h/2*k1,z2+h/2*l1);  l2=fun2(x+h/2,y2+h/2*k1,z2+h/2*l1);
			   k3=fun1(x+h/2,y2+h/2*k2,z2+h/2*l2);  l3=fun2(x+h/2,y2+h/2*k2,z2+h/2*l2);
			   k4=fun1(x+h,y2+h*k3,z2+h*l3);  l4=fun2(x+h,y2+h*k3,z2+h*l3);
			   y2=y2+h/6*(l1+2*l2+2*l3+l4);  z2=z2+h/6*(k1+2*k2+2*k3+k4);
			   
			   x=x+h;
			   cout<<i+1<<setw(20)<<x<<setw(20)<<y1<<setw(20)<<y2<<endl;
		   }
	   }

};
		   
void fun_print()
{
	cout<<"微分方程：(d^2 x)/(d t^2)=-g/L*sin(x)"<<endl;
	cout<<"其中初值："<<endl<<"\tx(0)=Pi/4"<<endl;
	cout<<"\tg=9.8\t(d x)/(d t)=0\tL=0.1"<<endl<<endl;
}
void main()
{//解高阶微分方程
	fun_print();
	Ola a(0,Pi/4,0);
	a.num_get(50,0,0.7);//计算一个周期内的值
	a.diedai();


	
}