Edited at

[C言語] 斜方投射をルンゲクッタ法で計算

More than 3 years have passed since last update.

ルンゲクッタでやる必要を特に感じないけど、まあ勉強ということで。


サンプル

地上5mから速度15m/sで鉛直上方へ投げ上げた質点の運動をルンゲクッタ法で計算する。


sample.c

#include <stdio.h>

#include <math.h>
#define g 9.80665

int main(void) {

double x0, x1, x2, v1, v2, v0, t, dt, x, k1, k2, k3, k4, k, h1, h2, h3, h4, h;

x0 = 5;
v0 = 15;

x1 = x0;
v1 = v0;
dt = 0.01;

for (int i=0; i<=100; ++i) {
t = i*dt;

k1 = dt*v1;
h1 = dt*(-g);

k2 = dt*(v1 + h1/2.);
h2 = dt*(-g);

k3 = dt*(v1 + h2/2.);
h3 = dt*(-g);

k4 = dt*(v1 + h3/2.);
h4 = dt*(-g);

k = (k1+ 2*k2 + 2*k3 + k4)/6.;
x2 = x1 + k;

h = (h1 + 2*h2 + 2*h3 +h4)/6.;
v2 = v1 + h;

x = x0 + v0*t - g*t*t/2.;

printf("t = %f, 数値計算 -> x = %f, 理論値 -> x = %f\n", t, x1, x);

x1 = x2;
v1 = v2;
}
return 0;
}



問題

速度10m/sで地面から角度60°で投げ上げた質点の運動をルンゲクッタ法で計算し、軌道をグラフにする。


problem.c

#include <stdio.h>

#include <math.h>
#define g 9.80665

int main(void) {

double y1, y2, v1, v2, v0, t, dt, x, y, rad, deg, k1, k2, k3, k4, k, h1, h2, h3, h4, h;

printf("初速度 [m/s] v0 = ");
scanf("%lf", &v0);

printf("角度(°) deg = ");
scanf("%lf", &deg);

rad = deg*(M_PI/180.);

y1 = 0;
v1 = v0*sin(rad);
dt = 0.01;

FILE * data;
data = fopen("problem_data.csv", "w");

for (int i=0; i<=200; ++i) {
t = i*dt;

k1 = dt*v1;
h1 = dt*(-g);

k2 = dt*(v1 + h1/2.);
h2 = dt*(-g);

k3 = dt*(v1 + h2/2.);
h3 = dt*(-g);

k4 = dt*(v1 + h3/2.);
h4 = dt*(-g);

k = (k1+ 2*k2 + 2*k3 + k4)/6.;
y2 = y1 + k;

h = (h1 + 2*h2 + 2*h3 + h4)/6.;
v2 = v1 + h;

x = v0*cos(rad)*t;
y = v0*sin(rad)*t - (1./2.)*g*t*t;

fprintf(data,"%f, %f\n", x, y1);
printf("t = %f, 数値計算 -> x = %f, y = %f, 理論値 -> x = %f, y = %f\n", t, x, y1, x, y);

y1 = y2;
v1 = v2;
}
return 0;
}



実行結果