現在(2013-07-02 (火) 18:18:57)作成中です。 既に書いてある内容も大幅に変わる可能性が高いので注意。
神戸大学 大学院システム情報学研究科 計算科学専攻 陰山 聡
【目次】
次にアニメーションによって収束の様子を確認しよう。 そのためのデータ(連番つきファイル群)を書き出すためのプログラム heat_animation_x_prof.f90 を以下に示す。 このプログラムをpgf90コマンドでコンパイルし、 実行すればdataディレクトリに連番ファイルが出力される。
module constants implicit none integer, parameter :: m=61, nmax=5000 integer, parameter :: SP = kind(1.0) integer, parameter :: DP = selected_real_kind(2*precision(1.0_SP)) end module constants module print use constants implicit none contains subroutine print__profile(f) real(DP), dimension(0:m+1), intent(in) :: f integer :: i integer :: counter = 0 ! automatic save attribute character(len=*), parameter :: base = "../data/temp.j=middle." character(len=4) :: serial_num write(serial_num,'(i4.4)') counter open(10,file=base//serial_num) do i = 0 , m+1 write(10,*) i, f(i) end do close(10) counter = counter + 1 end subroutine print__profile end module print program heat1_animation_x_prof use constants use print implicit none integer :: i,j,n real(DP), dimension(:,:), allocatable :: u, un real(DP) :: h, heat=1.0_DP real(DP) :: random allocate(u(0:m+1,0:m+1)) allocate(un(m,m)) h = 1.0_DP/(m+1) u(:,:) = 0.0_DP call print__profile(u(:,m/2+1)) do n=1, nmax do j=1, m do i=1, m un(i,j)=(u(i-1,j)+u(i+1,j)+u(i,j-1)+u(i,j+1))/4.0_DP+heat*h*h end do end do u(1:m,1:m)=un(1:m,1:m) if (mod(n,100)==0) call print__profile(u(:,m/2+1)) end do deallocate(u,un) end program heat1_animation_x_prof
gnuplotのスクリプトを使えばアニメーションも簡単にできる。 ただし、コマ数が多くなると手でスクリプトを書くのは難しくなるので、 「gnuplotスクリプト自動生成プログラム」を作ったほうがよい。 以下はその例である。
以下のコードを heat1_animation_x_prof_gp_generator.f90 というファイルに保存し、ソースコードの冒頭にあるUsageに従ってgnuplotによる可視化アニメーションを作れ。
! ! heat1_animation_x_prof_gp_generator.f90 ! ! by Akira Kageyama ! at Kobe University ! on 2010.07.07 ! for data visualization of heat1_animation_x_prof.f90 ! in the Lecture "Computational Science" (2010) ! ! Usage: ! (1) check the mesh size m and the counter_end. ! (2) pgf95 this_code ! (3) ./a.out > anyname ! (4) gnuplot anyname ! program heat1_animation_x_prof_gp_generator implicit none integer, parameter :: m = 61 integer :: counter, counter_end = 50 character(len=*), parameter :: base='../data/temp.j=middle.' character(len=4) :: serial_num print *, "#" print *, "# gnuplot script generated by heat1_animation_x_prof_gp_generator.f90" print *, "#" print *, " " print *, "set xlabel 'j' # x-axis" print *, "set ylabel 'temperature' # y-axis" print *, "set xrange[0:", m+1, "] # j-grid min & max" print *, "set yrange[0.0:0.5] # temperature min & max" do counter = 0 , counter_end write(serial_num,'(i4.4)') counter print *, "plot '"//base//serial_num//"' w lp" if ( counter==0) then print *, "pause 5" else print *, "pause 1" end if end do print *, "pause -1" end program heat1_animation_x_prof_gp_generator
・・・未定・・・
今回の演習内容はどうでしたか?(どれか一つ、一度だけ押してください。)
「お名前」欄は空欄で可。
コメントはありません。 コメント/6.1 1次元可視化?
as of 2023-12-11 (月) 17:27:07 (1431)