C/C++实现Runge-Kutta法求解常微分方程 数值结果与解析解不符
问题根因定位
你的代码数值解和解析解偏差极大的核心错误是指数项符号计算错误,RK4迭代逻辑本身没有问题。
- 你在导数函数
f()和解析解函数ya()中,将本应为exp(-x²)的指数项错写为exp((-x)*(-x)):(-x)*(-x)的计算结果是正的x²,完全丢失了指数上的负号,直接导致导数计算、解析解计算全部错误。
正确性验证:对题目给出的解析解求导
若 $y(x) = e{-x2}(\sin x - x\cos x + 1)$,求导可得:
$y' = -2x e{-x2}(\sin x -x\cos x +1) + e{-x2} \cdot x\sin x = x e{-x2}\sin x - 2xy$
和题目给出的常微分方程完全匹配,指数项必须为-x²。
修正方案
将两个函数中的指数项修正即可,修正后的代码片段如下:
float f(float x, float y) { // 修正指数为 -x*x 对应exp(-x²) return x*exp(-x*x) * sin(x) - 2*x*y; } float ya(float x) { // 同步修正解析解的指数项 return exp(-x*x) * ( sin(x) - x*cos(x) + 1 ); }
可选优化点
- 步长计算规范:当前
n=100个离散点时,你写的h=(b-a)/(n-1)虽然能让最后一个点刚好落在x=2位置,但常规区间n等分的标准写法是h=(b-a)/n,循环内取x = a + k*h,可避免多场景下的边界计算错误。 - 精度提升:
float为单精度浮点数,有效位数仅6~7位,若需要更高计算精度可将所有浮点变量替换为double类型,对应printf输出格式符改为%lf即可。
修正后重新运行,RK4数值解和解析解的绝对误差会降到1e-6量级,不会出现运行结果中明显偏离的情况。
内容的提问来源于stack exchange,提问作者Shams
相关产品推荐
相关产品推荐

