#include <stdio.h>
#include <stdlib.h>
#include <math.h>
double g = 9.8;
double k = 0.008;
double a_x(double vx, double vy) {
return -k * sqrt(vx * vx + vy * vy) * vx;
}
double a_y(double vx, double vy) {
return -g - k * sqrt(vx * vx + vy * vy) * vy;
}
void sim(double degree, double v0, double dt, int n_steps,
double* out_v, double* out_d, double* out_radian)
{
double radian = degree * M_PI / 180.0;
double vx0 = v0 * cos(radian);
double vy0 = v0 * sin(radian);
double* vx = (double*)malloc(sizeof(double) * n_steps);
double* vy = (double*)malloc(sizeof(double) * n_steps);
double* dx = (double*)malloc(sizeof(double) * n_steps);
double* dy = (double*)malloc(sizeof(double) * n_steps);
vx[0] = vx0;
vy[0] = vy0;
dx[0] = 0.0;
dy[0] = 0.0;
int last = 0;
for (int i = 1; i < n_steps; i++) {
if (dy[i - 1] >= 0.0) {
double k1_vx = a_x(vx[i - 1], vy[i - 1]) * dt;
double k1_vy = a_y(vx[i - 1], vy[i - 1]) * dt;
double k2_vx = a_x(vx[i - 1] + k1_vx / 2.0,
vy[i - 1] + k1_vy / 2.0) * dt;
double k2_vy = a_y(vx[i - 1] + k1_vx / 2.0,
vy[i - 1] + k1_vy / 2.0) * dt;
vx[i] = vx[i - 1] + k2_vx;
vy[i] = vy[i - 1] + k2_vy;
dx[i] = dx[i - 1] + vx[i - 1] * dt;
dy[i] = dy[i - 1] + vy[i - 1] * dt;
last = i;
}
else {
break;
}
}
*out_v = sqrt(vx[last] * vx[last] + vy[last] * vy[last]);
*out_d = dx[last];
*out_radian = radian;
free(vx);
free(vy);
free(dx);
free(dy);
}
int main() {
double v0 = 55.0;
double dt = 0.01;
int n_steps = 1000;
double degreeList[100];
double distanceList[100];
int idx = 0;
for (int degree = 10; degree <= 90; degree += 5) {
double v, d, radian;
sim(degree, v0, dt, n_steps, &v, &d, &radian);
printf("角度 %d 度 ラジアン %.5f\n", degree, radian);
printf("最終速さ %.2f m/s\n", v);
printf("飛距離 %.2f m\n\n", d);
degreeList[idx] = degree;
distanceList[idx] = d;
idx++;
}
// 最大飛距離を探す
double maxDist = distanceList[0];
double bestDeg = degreeList[0];
for (int i = 1; i < idx; i++) {
if (distanceList[i] > maxDist) {
maxDist = distanceList[i];
bestDeg = degreeList[i];
}
}
printf("最大飛距離 %.2f m は角度 %.2f 度で達成されます\n",
maxDist, bestDeg);
return 0;
}