我想在 odeint C++ library 中使用 runge_kutta4
方法.我已经在 Matlab 中解决了这个问题。我在 Matlab 中使用以下代码求解 x'' = -x - g*x'
,初始值 x1 = 1
,x2 = 0
, 如下
main.m
clear all
clc
t = 0:0.1:10;
x0 = [1; 0];
[t, x] = ode45('ODESolver', t, x0);
plot(t, x(:,1));
title('Position');
xlabel('time (sec)');
ylabel('x(t)');
ODESolver.m
function dx = ODESolver(t, x)
dx = zeros(2,1);
g = 0.15;
dx(1) = x(2);
dx(2) = -x(1) - g*x(2);
end
我已经安装了 odeint 库。我使用runge_kutta4
的代码如下
#include <iostream>
#include <boost/numeric/odeint.hpp>
using namespace std;
using namespace boost::numeric::odeint;
/* The type of container used to hold the state vector */
typedef std::vector< double > state_type;
const double gam = 0.15;
/* The rhs of x' = f(x) */
void lorenz( const state_type &x , state_type &dx , double t )
{
dx[0] = x[1];
dx[1] = -x[0] - gam*x[1];
}
int main(int argc, char **argv)
{
const double dt = 0.1;
runge_kutta_dopri5<state_type> stepper;
state_type x(2);
x[0] = 1.0;
x[1] = 0.0;
double t = 0.0;
cout << x[0] << endl;
for ( size_t i(0); i <= 100; ++i){
stepper.do_step(lorenz, x , t, dt );
t += dt;
cout << x[0] << endl;
}
return 0;
}
结果如下图
我的问题是为什么结果不同?我的 C++ 代码有问题吗?
这是两种方法的第一个值
Matlab C++
-----------------
1.0000 0.9950
0.9950 0.9803
0.9803 0.9560
0.9560 0.9226
0.9226 0.8806
0.8806 0.8304
0.8304 0.7728
0.7728 0.7084
0.7083 0.6379
更新:
问题是我忘记在我的 C++ 代码中包含初始值。感谢@horchler 注意到它。在包含正确的值并按照@horchler 的建议使用 runge_kutta_dopri5
之后,结果是
Matlab C++
-----------------
1.0000 1.0000
0.9950 0.9950
0.9803 0.9803
0.9560 0.9560
0.9226 0.9226
0.8806 0.8806
0.8304 0.8304
0.7728 0.7728
0.7083 0.7084
我更新了代码以反射(reflect)这些修改。
最佳答案
odeint 中的 runge_kutta4
步进器与 Matlab 的 ode45
完全不同。 ,这是一种基于 Dormand-Prince 的自适应方案方法。要复制 Matlab 的结果,您可能应该尝试 runge_kutta_dopri5
步进器。此外,请确保您的 C++ 代码使用与 ode45
相同的绝对和相对容差(默认值分别为 1e-6
和 1e-3
).最后,看起来您可能没有在 C++ 结果中保存/打印您的初始条件。
如果您对为什么 ode45
不采取固定步骤感到困惑,即使您指定了 t = 0:0.1:10;
,请参阅 my detailed answer here .
如果您必须使用固定步长runge_kutta4
步进器,那么您需要减少 C++ 代码中的积分步长以匹配 Matlab 的输出。
关于c++ - odeint的runge_kutta4与Matlab的ode45的比较,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/26747492/