初始值为负时Runge-Kutta(RK4)算法行为异常问题排查
RK4算法求解一阶微分方程的异常问题
我实现了用于求解一阶微分方程的C++版RK4算法,代码如下:
using ODE_Function = std::function<double/*dy/dx*/(double/*x*/, double/*y*/)>; template<int length> //"length" is the number of data points double* rk4(ODE_Function fxn, double y0, double x0, double h) { double* y = new double[length]; y[0] = y0; //initialize output array double x = x0; //n is our discretized x variable //main loop for (int n = 1; n < length; n++) {//y0 is already known, so we start from the second data point double k1 = fxn(x, y[n-1]); double k2 = fxn(x + h/2, y[n-1] + h*k1/2); double k3 = fxn(x + h/2, y[n-1] + h*k2/2); double k4 = fxn(x + h, y[n-1] + h*k3); y[n] = y[n-1] + h/6*(k1 + 2*k2 + 2*k3 + k4); x += h; } return y; }
但当初始值x0为负时,算法行为异常:理论上以下三个测试案例的结果都应符合(-x²)-2(x+1)的趋势,仅观察窗口不同,但只有最后一个案例符合预期。第一个案例结果呈负指数趋势,第二个案例结果呈正指数趋势,测试代码如下:
测试案例1:x0 = -1
double fxn(double x, double y) { return x*x + y; } int main() { const int length = 100; double y0 = -2; double x0 = -1; //LOOK HERE double h = 0.1; double* result = rk4<length>(fxn, y0, x0, h); }
测试案例2:x0 = -5
double fxn(double x, double y) { return x*x + y; } int main() { const int length = 100; double y0 = -2; double x0 = -5; //LOOK HERE double h = 0.1; double* result = rk4<length>(fxn, y0, x0, h); }
测试案例3:x0 = 0
double fxn(double x, double y) { return x*x + y; } int main() { const int length = 100; double y0 = -2; double x0 = 0; //LOOK HERE double h = 0.1; double* result = rk4<length>(fxn, y0, x0, h); }
内容的提问来源于stack exchange,提问作者Aryan MP
相关产品推荐
相关产品推荐

