### runge-kutta with 2d kinematics

I got to make a runge-kutta program that can solve where an object is at an x position and a y position. Where dx/dt = vx, dy/dt = vy, dvx/dt = 0, and dvy/dt = -g. The initial launch is at 20 m/s with a 30 degrees angle. X and y are suppose to start a 0. The equation for x is x = vx0 * t and y is y = vy0 * t - (1/2) * g * t^2. I am having a problem with how to use the initial condition of vx and vy with the rate of change it has. Here is the code I got so far:

 ``123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103`` ``````#include #include double ddt[4]; void f(double t, double v[4]) { ddt[0] = v[2] * t; ddt[1] = v[3] * t - (0.5) * 9.81 * pow(t,2); ddt[2] = 0 * t; ddt[3] = -9.81 * t; } int main() { using namespace std; void f(double t, double v[4]); double h, t0, v1_0, v2_0, v3_0, v4_0; int i, j; const int nos = 10; h = 0.1; t0 = 0.0; v1_0 = 0.0; v2_0 = 0.0; v3_0 = 20 * cos(30); v4_0 = 20 * sin(30); double t[nos + 1]; double v[4]; t[0] = t0; v[0] = v1_0; v[1] = v2_0; v[2] = v3_0; v[3] = v4_0; double v_n[4]; cout << "x is: " << v[0] << endl; cout << "y is: " << v[1] << endl << endl; for(i = 0; i < nos; i++) { t[i + 1] = t[i] + h; v_n[0] = v[0]; v_n[1] = v[1]; v_n[2] = v[2]; v_n[3] = v[3]; double added[4]; double k1[4], k2[4], k3[4], k4[4]; double v_npl[4]; f(t[i], v_n); for(j = 0; j < 4; j++) { k1[j] = h * ddt[j]; added[j] = v_n[j] + (0.5 * k1[j]); } f(t[i] = 0.5 * h, added); for(j = 0; j < 4; j++) { k2[j] = h * ddt[j]; added[j] = v_n[j] + (0.5 * k2[j]); } f(t[i] + 0.5 * h, added); for(j = 0; j < 4; j++) { k3[j] = h * ddt[j]; added[j] = v_n[j] + k3[j]; } f(t[i] + h, added); for(j = 0; j < 4; j++) { k4[j] = h * ddt[j]; } for(j = 0; j < 4; j++) { v_npl[j] = v_n[j] + ((k1[j] + 2 * k2[j] + 2 * k3[j] + k4[j]) / 6.0); } v[0] = v_npl[0]; v[1] = v_npl[1]; v[2] = v_npl[2]; v[3] = v_npl[3]; cout << "x is: " << v[0] << endl; cout << "y is: " << v[1] << endl << endl; } return 0; } /* doesn't seem to give "+" to y and y never returns to 0 */``````

And here is the output:

 ``` x is: 0 y is: 0 x is: 0.0231377 y is: -0.152415 x is: 0.0514171 y is: -0.339338 x is: 0.0848383 y is: -0.562826 x is: 0.123401 y is: -0.825019 x is: 0.167106 y is: -1.12814 x is: 0.215952 y is: -1.47448 x is: 0.26994 y is: -1.86644 x is: 0.32907 y is: -2.30648 x is: 0.393341 y is: -2.79714 x is: 0.462754 y is: -3.34106 ```

The problem I'm having is that y is negative and it won't return to 0. Anyone know what I am doing wrong? Thank you helping me if you do.
It was a simple radian-degree error. Here is the code:

 ``123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102`` ``````#include #include #define PI 3.1415926 double ddt[4]; void f(double t, double v[4]) { ddt[0] = v[2] * t; ddt[1] = v[3] * t - (0.5) * 9.81 * pow(t,2); ddt[2] = 0; ddt[3] = -9.81; } int main() { using namespace std; void f(double t, double v[4]); double h, t0, v1_0, v2_0, v3_0, v4_0; int i, j; const int nos = 10; h = 0.1; t0 = 0.0; v1_0 = 0.0; v2_0 = 0.0; v3_0 = 20 * cos(30 * PI/180); v4_0 = 20 * sin(30 * PI/180); double t[nos + 1]; double v[4]; t[0] = t0; v[0] = v1_0; v[1] = v2_0; v[2] = v3_0; v[3] = v4_0; double v_n[4]; cout << "x is: " << v[0] << endl; cout << "y is: " << v[1] << endl << endl; for(i = 0; i < nos; i++) { t[i + 1] = t[i] + h; v_n[0] = v[0]; v_n[1] = v[1]; v_n[2] = v[2]; v_n[3] = v[3]; double added[4]; double k1[4], k2[4], k3[4], k4[4]; double v_npl[4]; f(t[i], v_n); for(j = 0; j < 4; j++) { k1[j] = h * ddt[j]; added[j] = v_n[j] + (0.5 * k1[j]); } f(t[i] = 0.5 * h, added); for(j = 0; j < 4; j++) { k2[j] = h * ddt[j]; added[j] = v_n[j] + (0.5 * k2[j]); } f(t[i] + 0.5 * h, added); for(j = 0; j < 4; j++) { k3[j] = h * ddt[j]; added[j] = v_n[j] + k3[j]; } f(t[i] + h, added); for(j = 0; j < 4; j++) { k4[j] = h * ddt[j]; } for(j = 0; j < 4; j++) { v_npl[j] = v_n[j] + ((k1[j] + 2 * k2[j] + 2 * k3[j] + k4[j]) / 6.0); } v[0] = v_npl[0]; v[1] = v_npl[1]; v[2] = v_npl[2]; v[3] = v_npl[3]; cout << "x is: " << v[0] << endl; cout << "y is: " << v[1] << endl << endl; } return 0; } ``````
Topic archived. No new replies allowed.