龙格-库塔法解常微分方程(c++)
用龙格库塔法计算
#include
<iostream>
#include<iomanip>
#include
<cmath>
using
namespace std;
int main()
{
double a
= 1, b =
3; //a,b表示[a,b]求解区间
double
x0 = 1, y0 =
2; //x0表示初始时刻x的值,y0表示初始时刻y的值
double
x,
y; //x,y分别表示变化的时候x,y的值
double
F1, F2, F3,
F4; //F1,F2,F3,F4分别表示斜率值
double
h = 1.0f /
128; //h表示步长
cout << setiosflags(ios::left)
<< setw(25)
<< "x的值"
<< setw(25)
<< "龙格库塔计算得到的值"
<< setw(25)
<< "解析解得到的值"
<< setw(25)
<< "误差"
<< endl;
x
= x0;
y
= y0;
cout << setw(25)
<< x
<< setw(25)
<< y
<< setw(25)
<< y0
<< setw(25)
<< abs(y0-y)<< endl;
do
{
F1
= h*pow(x, -2)*(x*y - y*y);
F2
= h*pow(x + h / 2, -2)*((x + h / 2 )* (y + F1 / 2) - (y +
F1 / 2)*(y + F1 / 2));
F3
= h*pow(x + h / 2, -2)*((x + h / 2)* (y + F2 / 2) - (y +
F2 / 2)*(y + F2 / 2));
F4
= h*pow(x + h, -2)*((x + h)* (y + F3) - (y + F3)*(y +
F3));
y
+= (F1 + 2 * F2 + 2 * F3 + F4) / 6;
x
= x + h;
cout<< setw(25)
<< x
<< setw(25)
<< y
<< setw(25)
<< x / (1.0f / 2
+ log(x))
<< abs(x
/ (1.0f / 2 + log(x)) - y)
<< endl;
}
while (x<=b);
return
0;
}
最新文章
- 极简Photoshop 教程
- The first DP!
- C#小程序飞行棋地图绘制
- 在ubuntu12.04下编译android4.1.2添加JNI层出现问题
- 创建多模块maven项目
- .NET,你忘记了么?(八)—— 从dynamic到特性误用 [转]
- HTML5实现IP Camera网页输出
- 动态Order by
- TextUtils判断
- jQuery.proxy()的用法
- Dynamics CRM教程:图表的Top设置及导出修改和导入
- jwt实现
- 【Python】unittest-2-断言
- odoo开发笔记 -- wkhtmltox打印不显示中文 --ubuntu字体安装
- Django之views视图函数
- 【CF903G】Yet Another Maxflow Problem 线段树
- LeetCode刷题记录(python3)
- LeetCode--Combination Sum --ZZ
- [BZOJ4004][JLOI2015]装备购买(贪心+线性基)
- c#中类和成员的修饰符介绍