フーコーの振り子実験

フーコーの振り子実験の数値解析によるシミュレーションは数多くの方法が提案され,Web上でも様々に紹介されて 来ましたが,地球の自転を示す正確な方法はなかなか見当たりません。教科書の「Cで理解を深める基礎力学」に 単純で正確なフーコーの振り子実験のシミュレーション法の解説を記載しました。

シミュレーションの手順

長さ10 m の振り子を北緯30度の地点で@静止させた状態から北極に向けて速度を与えた場合と,A真南側に持ち上げて 静かに離す場合の2つの4時間のシミュレーション(時間刻み0.01 ms,1.44×10^9ステップ)を行った。
double fi=pi/6, df=0;
f0="Foucau35e"; fg="Foucau35e.emf"; f1="Foucault35e"; fg1="Foucault35e.emf";
の部分で/* */を付け替えて上と下を取り替えれば@とAの条件が入れ替わる。条件@では正確に時を刻むが,条件Aでは 円錐振り子の成分が加わってしまう。図はそれぞれ最初の5振動と最後の5振動を示している。
〇 C プログラム
#include
#include
#define GP "gnuplot"
#define N 1440000000
#define N1 3000000
#define M 10000
int main(){
int n, n0, nn, nm;
double rn, tn, so, co, x, y, z, xn, xn1, yn, yn1, zn, zn1, uf, ue;
double uff, g=9.823, o=7.27E-5, pi=3.1415, c=pi/3, R0=6.37E+6;
double sc=sin(c), cc=cos(c), dt=1.0E-5, r=10.0, Ro2=R0*o*o, Rdt=1.0E-8;
/* double fi=0, df=5.0E-6; */
double fi=pi/6, df=0;
char *f0, *fg, *f1, *fg1;
FILE *PGP, *FF0, *FF1;
/* f0="Foucau35x"; fg="Foucau35x.emf"; f1="Foucault35x"; fg1="Foucault35x.emf"; */
f0="Foucau35e"; fg="Foucau35e.emf"; f1="Foucault35e"; fg1="Foucault35e.emf";
FF0=fopen(f0,"w"); FF1=fopen(f1,"w");
xn=-sin(c-fi); xn1=-sin(c-fi-df); yn=0; yn1=0; zn=-cos(c-fi); zn1=-cos(c-fi-df);
for(n=1;n so=sin(o*tn); co=cos(o*tn);
x=2*r*xn-r*xn1+(Ro2-g)*sc*co*dt*dt;
y=2*r*yn-r*yn1+(Ro2-g)*sc*so*dt*dt;
z=2*r*zn-r*zn1-g*cc*dt*dt;
xn1=xn; yn1=yn; zn1=zn; rn=sqrt(x*x+y*y+z*z); xn=x/rn; yn=y/rn; zn=z/rn;
if(n < N1){
if(n0==0){ uf=-so*xn+co*yn; ue=cc*(co*xn+so*yn)-sc*zn;
fprintf(FF0,"%12.3f %15.8f %15.8f %15.8f %15.8f %15.8f \n",tn, xn, yn, zn, ue, uf);
}}
if( n > N-N1) {
if(n0==0){ uf=-so*xn+co*yn; ue=cc*(co*xn+so*yn)-sc*zn;
fprintf(FF1,"%12.3f %15.8f %15.8f %15.8f %15.8f %15.8f \n",tn, xn, yn, zn, ue, uf);
}}
}
fclose(FF0); fclose(FF1);
PGP=popen(GP,"w");
fprintf(PGP,"set term emf size 424,300 \n");
fprintf(PGP,"set output '%s' \n unset key \n",fg);
fprintf(PGP,"set xlabel 'ue ' \n set ylabel 'uf ' \n");
fprintf(PGP," plot '%s' using 5:6 with linespoints pt 7 \n",f0);
fflush(PGP);
fprintf(PGP,"set term emf size 424,300 \n");
fprintf(PGP,"set output '%s' \n unset key \n set size ratio -1 \n",fg1);
fprintf(PGP,"set xlabel 'ue ' \n set ylabel 'uf ' \n");
fprintf(PGP," plot '%s' using 5:6 with linespoints pt 7 \n",f0);
fprintf(PGP,"replot '%s' using 5:6 with linespoints pt 8 \n",f1);
fflush(PGP);
pclose(PGP);
return 0;
}