最常用:四阶龙格库塔方法
例:
#include <iostream> #include <fstream> #include <cstring> double func(double x) { return x * (1 - x); } void RK(double x0,double x[20]) { double dt{ 0.1 }; double k1, k2, k3, k4; int count{ 0 }; int flag{ 0 }; for (int count{ 0 };count < 100;++count) { if (count % 5 == 0) { flag = count / 5; x[flag] = x0; } k1 = func(x0) * dt; k2 = func(x0 + 0.5 * k1) * dt; k3 = func(x0 + 0.5 * k2) * dt; k4 = func(x0 + 0.5 * k3) * dt; x0 += (1.0 / 6.0) * (k1 + 2 * k2 + 2 * k3 + k4); } } void show(double* x, int num) { for (int i = 0;i < num;i++) { std::cout << x[i] << '\n'; } } void savetxt(double* x,int num,std::string name) { std::ofstream write; write.open(name,std::ostream::app); //以添加模式打开文件 for (int i = 0;i < num;i++) { write << x[i] << '\n'; } write.close(); } int main() { double x[20]{}; RK(0.5,x); //show(x, 20); savetxt(x, 20, "k2.txt"); RK(1.0, x); savetxt(x, 20, "K.txt"); RK(0.2, x); savetxt(x, 20, "less than k2.txt"); RK(0.8, x); savetxt(x, 20, "more than k2.txt"); RK(1.5, x); savetxt(x, 20, "more than k.txt"); }
使用vector:
/* x0:初值 dt:步长 D:总时长 p:控制输出的数量 */ vecRow RK2(double x0,double dt,int D,int p) { double k1, k2, k3, k4; int count{ 0 }; int flag{ 0 }; int t{ static_cast<int>(D / dt) }; int num{ static_cast<int>(t / p) }; vecRow x(num); for (int count{ 0 };count < t;++count) { if (count % 5 == 0) { flag = count / p; x[flag] = x0; } k1 = func(x0) * dt; k2 = func(x0 + 0.5 * k1) * dt; k3 = func(x0 + 0.5 * k2) * dt; k4 = func(x0 + 0.5 * k3) * dt; x0 += (1.0 / 6.0) * (k1 + 2 * k2 + 2 * k3 + k4); } return x; }
标签:count,求解,int,double,c++,func,微分方程,x0,dt From: https://www.cnblogs.com/zhimingyiji/p/17080783.html