C++で常微分方程式(ホイン法)
初速 150 km/h で, 45°の角度で投げたボールの軌跡をホイン法で計算する
(空気抵抗を考慮しない)
#include <iostream> #include <iomanip> #include <math.h> using namespace std; // 重力加速度 const double g = 9.8; // 水平方向の速度 double fx(double t, double vx) { // 水平方向の速度は, 初速のまま return vx; } // 鉛直方向の速度 double fy(double t, double vy) { // 初速から落下速度を引く return vy - (t * g); } int main(void) { // 角度 double degree = 45; double radian = degree * M_PI / 180.0; // 初速 150 km/h -> 秒速に変換 double v = 150 * 1000 / 3600; // 水平方向の速度 double vx = v * cos(radian); // 鉛直方向の速度 double vy = v * sin(radian); // 時間間隔(秒) double dt = 0.01; // 経過秒数 double t = 0.0; // 位置 double x = 0.0; double y = 0.0; // Heun法 for (int i = 1; y >= 0.0; i++) { // 位置 // Euler法 // x = x + dt * fx(t, vx); // y = y + dt * fy(t, vy); // Heun法 x = x + (dt * fx(t, vx) + dt * fx(t + dt, vx)) / 2; y = y + (dt * fy(t, vy) + dt * fy(t + dt, vy)) / 2; // 経過秒数 t = i * dt; // 経過秒数 cout << setw(8) << fixed << setprecision(5) << t << "\t"; // 位置 (Heun法) cout << setw(8) << fixed << setprecision(5) << x << "\t"; cout << setw(8) << fixed << setprecision(5) << y << "\t"; // 位置 (正解) cout << setw(8) << fixed << setprecision(5) << vx * t << "\t"; cout << setw(8) << fixed << setprecision(5) << vy * t - (g * t * t) / 2 << endl; } return 0; }
0.01000 0.28991 0.28942 0.28991 0.28942 0.02000 0.57983 0.57787 0.57983 0.57787 0.03000 0.86974 0.86533 0.86974 0.86533 0.04000 1.15966 1.15182 1.15966 1.15182 0.05000 1.44957 1.43732 1.44957 1.43732 0.06000 1.73948 1.72184 1.73948 1.72184 0.07000 2.02940 2.00539 2.02940 2.00539 0.08000 2.31931 2.28795 2.31931 2.28795 0.09000 2.60922 2.56953 2.60922 2.56953 0.10000 2.89914 2.85014 2.89914 2.85014 省略 5.80000 168.14999 3.31399 168.14999 3.31399 5.81000 168.43991 3.03502 168.43991 3.03502 5.82000 168.72982 2.75506 168.72982 2.75506 5.83000 169.01973 2.47412 169.01973 2.47412 5.84000 169.30965 2.19221 169.30965 2.19221 5.85000 169.59956 1.90931 169.59956 1.90931 5.86000 169.88948 1.62544 169.88948 1.62544 5.87000 170.17939 1.34058 170.17939 1.34058 5.88000 170.46930 1.05474 170.46930 1.05474 5.89000 170.75922 0.76793 170.75922 0.76793 5.90000 171.04913 0.48013 171.04913 0.48013 5.91000 171.33904 0.19135 171.33904 0.19135 5.92000 171.62896 -0.09840 171.62896 -0.09840
参考文献