【工程数学】若干种解常微分方程的算法

来源：互联网发布：阿里云客服薪资计算编辑：程序博客网时间：2024/05/16 23:01

// ConsoleAppDifferential_Euler_Solu.cpp : 定义控制台应用程序的入口点。///**函数功能：欧拉法解常微分方程*函数原形：void Differential_Euler_Solu(double a,double b) *参数：double a,double b：所要计算的区间*返回值：无*时间复杂度：O(n)*备注：此法解微分方程精度较低*日期：2014/12/12*原创：是*作者：EbowTang*Email：tangyibiao520@163.com*/#include "stdafx.h"#include "iostream"#include <iomanip>using namespace std;double Func(double x,double y);double RealFunc(double x,double y) ;void Differential_Euler_Solu(double ,double);int _tmain(int argc, _TCHAR* argv[]){cout<<"此微分方程为："<<endl;cout<<"式子一：y'=y/x-2*y*y"<<endl;cout<<"式子二：y(0)=0"<<endl;Differential_Euler_Solu(0,3);system("pause");return 0;}double Func(double x,double y)  {  if (x!=0.0)return (y/x-2*y*y);  elsereturn 1;}  double RealFunc(double x)//此式子为该微分方程的真实解形式  {  return x/(1+x*x);} void Differential_Euler_Solu(double a,double b)  {   double oldX=0.0,oldY=0.0; double newX=0.0,newY=0.0;  double h=0.2;//步长 for ( int k = 0 ; k <=(b-a)/h ; k++ )       {    oldY = newY; oldX = newX;                           newY = oldY+h*Func(oldX,oldY);//欧拉法的迭代形式 newX = oldX+h;//x递增 cout<<"........................"<<endl; cout<<"第"<<k+1<<"次迭代计算的结果为"<<endl; cout<<"x="<<setprecision(5)<<oldX; cout<<"  真实y="<<setprecision(5)<<RealFunc(oldX); cout<<"  近似y="<<setprecision(5)<<oldY; cout<<"  此次误差为="<<setprecision(5)<<oldY-RealFunc(oldX)<<endl; }}

<pre name="code" class="cpp">// ConsoleAppDifferential_Improve_Euler_Solu.cpp : 定义控制台应用程序的入口点。///**函数功能：改进的欧拉法解常微分方程*函数原形：void Differential_Improve_Euler_Solu(double a,double b) *参数：double a,double b：所要计算的区间*返回值：无*时间复杂度：O(n)*备注：此法解微分方程精度相比欧拉法精度明显提高*日期：2014/12/12*原创：是*作者：EbowTang*Email：tangyibiao520@163.com*/#include "stdafx.h"#include "iostream"#include <iomanip>using namespace std;double Func(double x,double y);double RealFunc(double x,double y) ;void Differential_Improve_Euler_Solu(double ,double);int _tmain(int argc, _TCHAR* argv[]){cout<<"此微分方程为："<<endl;cout<<"式子一：y'=y/x-2*y*y"<<endl;cout<<"式子二：y(0)=0"<<endl;Differential_Trapezoid_Solu(0,3);system("pause");return 0;}double Func(double x,double y)  {  if (x!=0.0)return (y/x-2*y*y);  elsereturn 1;}  double RealFunc(double x)//此式子为该微分方程的真实解形式  {  return x/(1+x*x);} void Differential_Improve_Euler_Solu(double a,double b)  {   double oldX=0.0,oldY=0.0; double newX=0.0,newY=0.0;  double h=0.2;//步长 for ( int k = 0 ; k <=(b-a)/h ; k++ )       {    oldY = newY; oldX = newX;     newY = oldY+h*Func(oldX,oldY); newX = oldX+h;//x递增 newY = oldY+(h/2)*(Func(oldX,oldY)+Func(newX,newY));//梯形法的迭代形式 cout<<"........................"<<endl; cout<<"第"<<k+1<<"次迭代计算的结果为"<<endl; cout<<"x="<<setprecision(5)<<oldX; cout<<"  真实y="<<setprecision(5)<<RealFunc(oldX); cout<<"  近似y="<<setprecision(5)<<oldY; cout<<"  此次误差为="<<setprecision(5)<<oldY-RealFunc(oldX)<<endl; }}

<pre name="code" class="cpp">// ConsoleAppRungeKuttaSolu.cpp : 定义控制台应用程序的入口点。///**函数功能：三阶以及四阶RungeKutta法解常微分方程*函数原形：void Differential_RungeKutta_Solu(double a,double b) *参数：double a,double b：所要计算的区间*返回值：无*时间复杂度：O(n)*备注：此法解微分方程精度相比欧拉法精度提高*日期：2014/12/12*原创：是*作者：EbowTang*Email：tangyibiao520@163.com*/#include "stdafx.h"#include "iostream"#include <iomanip>using namespace std;double Func(double x,double y);double RealFunc(double x,double y) ;void Differential_RungeKutta_Solu(double ,double);int _tmain(int argc, _TCHAR* argv[]){cout<<"此微分方程为："<<endl;cout<<"式子一：y'=y/x-2*y*y"<<endl;cout<<"式子二：y(0)=0"<<endl;cout<<"/*********************四阶经典龙格法的华丽分割线*******************/"<<endl;Differential_RungeKutta_Solu(0,3);system("pause");return 0;}double Func(double x,double y)  {  if (x!=0.0)return (y/x-2*y*y);  elsereturn 1;}  double RealFunc(double x)//此式子为该微分方程的真实解形式  {  return x/(1+x*x);} void Differential_RungeKutta_Solu(double a,double b)  {  double oldX=0.0,oldY=0.0;double newX=0.0,newY=0.0;double K1=0.0,K2=0.0,K3=0.0,K4;double h=0.2;//步长for ( int k = 0 ; k <=(b-a)/h ; k++ )      {   oldY = newY;oldX = newX;  /************经典四阶龙格法迭代形式*********************/K1=Func(oldX,oldY);K2=Func(oldX+h/2,oldY+h/2*K1);K3=Func(oldX+h/2,oldY+h/2*K2);K4=Func(oldX+h,oldY+h*K3);newY=oldY+h/6*(K1+2*K2+2*K3+K4);newX = oldX+h;//x递增/******经典三阶龙格法迭代形式K1=Func(oldX,oldY);K2=Func(oldX+h/2,oldY+h/2*K1);K3=Func(oldX+h,oldY-h*K1+2*h*K2);newY=oldY+h/6*(K1+4*K2+K3);newX = oldX+h;//x递增*********************************/cout<<"第"<<k+1<<"次迭代计算的结果为"<<endl;cout<<"x="<<setprecision(5)<<oldX;cout<<"  真实y="<<setprecision(5)<<RealFunc(oldX);cout<<"  近似y="<<setprecision(5)<<oldY;cout<<"  此次误差为="<<setprecision(5)<<oldY-RealFunc(oldX)<<endl;}}

0 0