/*----------------------------------------*/ /* downward     */ /*---------------------------------------*/ #include #include #define n 100 /* 一辺あたりのメッシュの数 */ #define Length 1.0e1 /* 矩形の一辺の長さ */ #define cs 1e0 /* 波の伝播速度 */ #define pi 3.141592 /* π */ #define tcal 10. main(void) { int i,j,k,num; double factor; double delta_L,dt; double f[n+2][2], t; double tout,dtout; /* f の最後の添え字0,1はテキストのk,k+1に対応する */ FILE *output_file; /* 時間ステップなど */ delta_L = (float)Length/((float)n+1); dt = 0.1*delta_L/cs; factor=cs*dt/delta_L; output_file=fopen("output2.dat","w"); /* 出力ファイルオープン*/ /* ---------------- 初期化 及び境界条件 ----------------- */ t=0; tout=0.1; dtout=0.1; /* 境界はメッシュ番号 0 とメッシュ番号 n+1 にとる */ for(i=0; i < n+2 ; i++){ f[i][0]=0.0; /* 初期化及び境界の値を0に */ } for(i=0; i < n/2 ; i++){ f[i][0]=1.; /* 初期条件設定 */ } num=1; /* この段階でt=0はk=1に対応することに注意 */ /* -------------- 初期化 及び境界条件終わり ------------- */ /* -- 初期条件の出力 -- */ for(i=0; i < n+2 ; i++){ fprintf(output_file,"%f %f\n",(float)delta_L*i,f[i][0]); } fprintf(output_file,"\n \n"); fflush(output_file); /* --------------------*/ /* -------- 時間ステップの繰り返しループ --------- */ while(t<=tcal){ for(i=1; i < n+1 ; i++){ f[i][1]=f[i][0]-factor*(f[i][0]-f[i-1][0]); } t=t+dt; num+=1; for(i=1; i < n+1 ; i++){ f[i][0]=f[i][1]; } /* -------- 出力ルーチン -----------*/ if(t>=tout){ /* t>tout出力 */ printf("number of iteration = %d %f\n",num,t); for(i=0; i < n+2 ; i++){ fprintf(output_file,"%f %f %f\n",(float)delta_L*i,(float)f[i][0]); } fprintf(output_file,"\n \n"); fflush(output_file); tout+=dtout; } /* -------- 出力ルーチン終わり -----------*/ } /* -------- 時間ステップ繰り返しループ 終わり----------- */ }