/*
************************************************************************
* ルンゲ・クッタ法による2次元運動方程式の解法 *
* 速度の1乗および2乗に比例する抵抗がある斜方投射 *
************************************************************************
t:時間、x[0]=x、x[1]=y:位置、x[2]=vx、x[3]=vy:速度、dt:時間の刻み幅
t0、x0、y0、vx0、vy0:初期値
v0:初速度、theta:投げ上げ角度
alpha:速度に比例する抵抗の係数
beta:速度の2乗に比例する抵抗の係数
*/
#include <stdio.h>
#include <math.h>
#define N 4 /* 方程式の数 */
double alpha, beta;
main()
{
void rungen(double t, double x[], double dt);
double f(double t, double x[], int i);
double t, x[N], v0, vx0, vy0, theta, dt, h;
double rho, rball, mass, eta;
double pi;
int j;
/* 出力ファイルの指定 */
FILE *output;
output = fopen("ball.dat","w"); /* ball.datにデータをセーブ */
/* 初速度 */
printf("Input: v0 (km/h)\n");
scanf("%lf", &v0);
printf("v0=%g\n", v0);
/* 投げ上げ角度 */
printf("Input: theta\n");
scanf("%lf", &theta);
printf("theta=%g\n", theta);
pi=3.1415926535897932385;
/* 時間の刻み幅 */
printf("Input: dt\n");
scanf("%lf", &dt);
printf("dt=%g\n", dt);
rho=1.2; /* 空気の密度 */
rball=0.0715/2.0; /* ボールの半径 */
mass=0.1445; /* ボールの質量 */
eta=1.8e-5; /* 空気の粘性係数 */
v0=v0*1000.0/3600.0;
vx0=v0*cos(pi*theta/180);
vy0=v0*sin(pi*theta/180);
alpha=6.0*pi*rball*eta/mass;
beta=0.25*pi*rho*rball*rball/mass;
/* 初期値 */
t = 0.0;
x[0] = 0.0; /* 初期位置 */
x[1] = 0.0; /* 初期位置 */
x[2] = vx0; /* 初期速度 */
x[3] = vy0; /* 初期速度 */
fprintf(output, "%f\t%f\t%f\n", t, x[0], x[1]);
while(x[1] >= 0.0)
{
/* ルンゲ・クッタ法ルーチンを呼び出す */
rungen(t, x, dt);
t += dt;
fprintf(output, "%f\t%f\t%f\n", t, x[0], x[1]);
}
printf("data stored in ball.dat\n");
fclose(output);
}
/*-----------------------end of main program--------------------------*/
void rungen(double t, double x[], double dt)
{
double f(double t, double x[], int i);
double h=dt/2.0,
t1[N], t2[N], t3[N], /* Runge-Kutta法のための */
k1[N], k2[N], k3[N],k4[N]; /* 一時変数 */
int i;
for (i=0; i<N; i++)
{
k1[i] = dt*f(t, x, i);
t1[i] = x[i]+0.5*k1[i];
}
for (i=0; i<N; i++)
{
k2[i] = dt*f(t+h, t1, i);
t2[i] = x[i]+0.5*k2[i];
}
for (i=0; i<N; i++)
{
k3[i] = dt*f(t+h, t2, i);
t3[i] = x[i]+ k3[i];
}
for (i=0; i<N; i++)
{
k4[i] = dt*f(t + dt, t3, i);
}
for (i=0; i<N; i++) x[i] += (k1[i]+2*k2[i]+2*k3[i]+k4[i])/6.0;
}
/*--------------------------------------------------------------------*/
/* 微分方程式の右辺 */
double f(double t, double x[], int i)
{
double vabs; /* 速度の絶対値 */
vabs=sqrt(x[2]*x[2]+x[3]*x[3]);
if (i == 0) return(x[2]); /* dx[0]/dx = dx/dt */
if (i == 1) return(x[3]); /* dx[1]/dx = dy/dt */
if (i == 2) return(-alpha*x[2] -beta*vabs*x[2]); /* dx[2]/dx = dvx/dt */
if (i == 3) return(-9.8-alpha*x[3] -beta*vabs*x[3]); /* dx[3]/dx = dvy/dt */
}
このプログラムを写真の式に対応させてプログラムを書き換えていただきたいです。よろしくお願いいたします。
お探しのQ&Aが見つからない時は、教えて!gooで質問しましょう!
似たような質問が見つかりました
- 数学 単振り子とルンゲ・タック法 1 2022/07/15 00:05
- その他(プログラミング・Web制作) 物理の斜方投射のシミュレーションにおける位置や速度の単位について 4 2023/05/31 09:50
- C言語・C++・C# c言語でユーザ関数を利用して複素数のべき乗と絶対値の数列を計算するプログラムが作りたいです。 3 2023/01/29 22:13
- その他(プログラミング・Web制作) 物理の斜方投射で目盛りに数値を入れたい 2 2023/05/27 06:32
- その他(プログラミング・Web制作) 物理の斜方投射の目盛り線とx軸、y軸の追加について 3 2023/05/26 21:11
- C言語・C++・C# 10個の実数に対する降順ソート結果を出力するプログラムを作りたいのですが、以下のプログラムをどう直せ 1 2022/07/09 22:16
- その他(プログラミング・Web制作) Pythonにおける物理のシミュレーションでの単位変換について 2 2023/06/02 17:11
- C言語・C++・C# プログラミングの授業の課題です 1 2023/01/17 22:15
- その他(プログラミング・Web制作) Pythonでのかんたんな物理シミュレーションについての書籍 5 2023/06/02 07:37
- その他(プログラミング・Web制作) ボールの動きがスムーズに動いてかつ目盛り線描画を維持するためには 4 2023/05/31 10:01
関連するカテゴリからQ&Aを探す
おすすめ情報
デイリーランキングこのカテゴリの人気デイリーQ&Aランキング
-
プログラムでの数字につく”f”の...
-
doubleの変数にintとintの割り...
-
float型とdouble型の変数の違い...
-
c言語で、繰り返し文の中で、0....
-
sin(x)の近似について
-
2次方程式の解を求めるプログ...
-
Cで3乗根を求める方法
-
関数におけるif文とreturn文に...
-
int とdoubleの比較
-
相互相関関数
-
C言語でポインタを用いた平均,...
-
MATLABで画像のヒストグラムを...
-
C言語で直角三角形の斜辺を求め...
-
C言語を実行すると-infが出てき...
-
C言語の問題です。
-
至急です! マクロ定義で #defi...
-
ボール同士の衝突
-
関数のプログラム
-
C言語の型による処理速度の違い
-
C#イベント中の戻り値の設定の...
マンスリーランキングこのカテゴリの人気マンスリーQ&Aランキング
-
プログラムでの数字につく”f”の...
-
doubleの変数にintとintの割り...
-
C 開放してるのにエラー(doubl...
-
Cで3乗根を求める方法
-
float型とdouble型の変数の違い...
-
至急です! マクロ定義で #defi...
-
C言語の型による処理速度の違い
-
int とdoubleの比較
-
関数におけるif文とreturn文に...
-
C言語初心者 構造体 課題について
-
c言語のコンパイルエラー canno...
-
C言語 関数プロトタイプ宣言の...
-
C言語を実行すると-infが出てき...
-
float?数字の後にLがつくもの
-
数値を指数部と仮数部に分離したい
-
difftime()について
-
浮動小数点数が表示されないん...
-
たくさんの数の平均を求める方...
-
DWORDの警告
-
-1.#IND00と出てしまうのですが...
おすすめ情報
今回は流速はないものと考えuiは0として下さい。
よろしくお願いします。