c - C中弹丸(二次阻力)的四阶RK方法的递归值

Closed. This question is off-topic. It is not currently accepting answers. Learn more。
想改进这个问题吗？Update the question所以堆栈溢出的值小于aa>。
去年关门了。
我试图用四阶RK方法对弹壳的抛射运动进行编码。我在下面的代码中遇到了一个问题。

for (i=1; i<=n; i=i+1)
{
    ax[i]= (-kk[i]*v[i]*vx[i])/m;
    ay[i]= ((-kk[i]*v[i]*vy[i])-(m*g))/m;
    k1= ax[i];
    m1= ay[i];
    k2= (-kk[i]*v[i]*(vx[i]+(k1*delt/2)))/m;
    m2=  (-kk[i]*v[i]*(vy[i]+(m1*delt/2)))/m;
    k3= (-kk[i]*v[i]*(vx[i]+(k2*delt/2)))/m;
    m3=  (-kk[i]*v[i]*(vy[i]+(m2*delt/2)))/m;
    k4= (-kk[i]*v[i]*(vx[i]+(k3*delt)))/m;
    m4=  (-kk[i]*v[i]*(vy[i]+(m3*delt)))/m;
    vx[i+1]= vx[i]+((1/6)*(k1+2*k2+2*k3+k4)*delt);
    vy[i+1]= vy[i]+((1/6)*(m1+2*m2+2*m3+m4)*delt);
    x[i+1]= vx[i+1]*delt;
    y[i+1]= vy[i+1]*delt;
    xx[i+1]= x[i+1] + xx[i];
    yy[i+1]= y[i+1] + yy[i];
    t[i+1]= t[i]+delt;
    v[i+1]= sqrt((vx[i+1]*vx[i+1])+(vy[i+1]*vy[i+1]));
    if ((yy[i+1]<0)&&(i!=1))
        //if (i==100)
    {
        nn=i;
        i=n;
    }
    printf ("%f\t%f\t%f\t%f\t%f\t%f\t\n", t[i], vx[i], vy[i], v[i], xx[i], yy[i]);

在上述代码中，ax、ay、vx和vy是x和y方向上的加速度和速度，delt表示时间步长。xx和yy是相应时间步的位置。
在使用四阶RK方法之后，我尝试在下一个时间步中找到速度，以获得回路中的相应值，直到弹丸击中地面。但每次迭代和重复，我只能得到与vx[1]和vy[1]相同的值，我是C编程的初学者，还处于使用循环的初始阶段。

最佳答案

问题在于：

vx[i+1]= vx[i]+((1/6)*(k1+2*k2+2*k3+k4)*delt);
vy[i+1]= vy[i]+((1/6)*(m1+2*m2+2*m3+m4)*delt);

1/6-两个数字都是整数，因此这将导致0生成语句：

vx[i+1]= vx[i];
vy[i+1]= vy[i];

使它们浮点数。既然你经常计算，就把它设为常数。

const double one_sixth = 1.0 / 6.0;

那么

vx[i+1]= vx[i]+(one_sixth * (k1+2*k2+2*k3+k4)*delt);
vy[i+1]= vy[i]+(one_sixth * (m1+2*m2+2*m3+m4)*delt);