iRIC v4対応ソルバー fortran編

はじめに

本記事ではiRIC v4に対応したソルバー記述方法について紹介する。言語はfortran、２次元構造格子の計算プログラムを前提とする。サンプルプログラムやコンパイルは以下を参考にして欲しい。

>>> こちら

1. ファイルオープン～初期化

iRICソフトウェア上でプログラム実行ボタンがクリックされると、”Case1.cgn"というファイルを引数としてソルバープログラム(***.exe)が実行される。そのため、プログラムではその引数からファイル名を取得し、ファイルをオープンする。既に計算結果が格納されている場合は削除してから計算を実行したいので、既存計算結果をクリアするためのサブルーチンをcallする。計算途中に計算を停止したい場合に備えた準備をする。

    write(*,'(a)') "---------- start program ---------"

! get argument
    if(iargc() == 1) then
        call getarg(1, cgnsfile)
        
    else
        write(*,"(a)") "Error: no argument."
        write(*,"(a)") "usage: "
        write(*,"(a)") "***.exe Case1.cgn"
        stop
    end if

! open file
    call cg_iric_open(trim(cgnsfile), IRIC_MODE_MODIFY, FID, IER)

! clear previous result data
	call cg_iric_clear_sol(fid, ier) 

! set cancel option
	call iric_initoption(IRIC_OPTION_CANCEL, ier)

2. 計算条件読込

計算条件の読込み関数は以下の種類が準備されている。"definition.xml"に定義した計算条件に応じて読込部分を作成する。

! integer
    call cg_iric_read_integer(fid, 'ival', ival, ier)

! real
    call cg_iric_read_real(fid, 'rval', rval, ier)

! string
	call cg_iric_read_string(fid, 'strval', strval, ier)

! function
	call cg_iric_read_functionalsize(fid, 'func_a', size, ier)
    allocate(xval(1:size), yval(1:size))
    call cg_iric_read_functionalwithname(fid, 'func_a', 'xval', xval, ier)
    call cg_iric_read_functionalwithname(fid, 'func_a', 'yval', yval, ier)

! grid data
    call cg_iric_read_grid2d_str_size(fid, isize, jsize, ier)
    allocate(xx(1:isize, 1:jsize), yy(1:isize, 1:jsize))
    call cg_iric_read_grid2d_coords(fid, xx, yy, ier)

! node attributes
    allocate(zz(1:isize, 1:jsize))
    call cg_iric_read_grid_real_node(fid, 'elevation', zz, ier)

! cell attributes
    allocate(zz_c(1:isize-1, 1:jsize-1), sn_c(1:isize-1, 1:jsize-1), obst_c(1:isize-1, 1:jsize-1)) 
    call cg_iric_read_grid_real_cell(fid, 'elevation_c', zz_c, ier)
    call cg_iric_read_grid_real_cell(fid, 'sn_c', sn_c, ier)
    call cg_iric_read_grid_real_cell(fid, 'obst_c', obst_c, ier)

3. タイムループと計算結果出力

計算のタイムループ＝時間更新は大凡以下の形となる。
iRICソフトウェア上で「更新」ボタンがクリックされたときや、キャンセルボタンがクリックされたときのアクションは毎タイムステップ確認するのが良い。計算結果出力は、cg_iric_write_sol_startとcg_iric_write_sol_endで挟んだ外部サブルーチンに設定するのが良いだろう。

    do icount = istart, iend

! set time
        ctime = icount * dt
    
! check update button clicked
        call cg_iric_check_update(fid, ier)

! check cancel button clicked
        call iric_check_cancel(ier)
        if(ier /= 0)then
            write(*,*) ""
            write(*,*) "solver stops because the cancel button was clicked."
            call cg_iric_close(fid, ier)
            stop
        end if
    
! ouput start
        if(mod(icount, ioutput) == 0)then
            write(*,"(a7, f20.3, a)") "time = ", ctime, "[sec]"
            call cg_iric_write_sol_start(fid, ier)
            call cg_iric_write_sol_time(fid, time, ier)
            call write_results()
            call cg_iric_write_sol_end(fid, ier)
        end if
    
!--------------------------------------------------!--------------------------------------------------
! start time integration calculation
    
    
    
    
    
! end time integration calculation
!--------------------------------------------------!--------------------------------------------------
    end do

    call cg_iric_close(fid, ier)
    write(*,*) ""
    write(*,'(a)') "---------- program normal end ---------"
    stop

end program

4. 計算結果出力

１タイムステップで１つの値(上流端流量、下流端水位など）、計算格子（格子形状が時間に依存して変化しない場合は不要）の出力は以下。
また、計算結果は"node"、"cell"、"iface"、"jface"の４箇所でそれぞれで出力することができる。いずれもは配列のインデックスを指定して出力するのが良いだろう。

       
subroutine call write_results()

! write baseiterative
    call cg_iric_write_sol_baseiterative_real(fid, "Q[m3_s]", qp, ier)

! write grid
    call cg_iric_write_sol_grid2d_coords(fid, xx(1:isize, 1:jsize), yy(1:isize, 1:jsize), ier)

! write node values
    call cg_iric_write_sol_node_real(fid, "Depth[m]", dep(1:isize, 1:jsize), ier)

! write cell values
    call cg_iric_write_sol_cell_real(fid, "Elevation_c[m]", eta_c(1:isize-1, 1:jsize-1), ier)
    call cg_iric_write_sol_cell_real(fid, "WaterSurf_c[m]", wl_c(1:isize-1, 1:jsize-1), ier)
    call cg_iric_write_sol_cell_real(fid, "Depth_c[m]", dep_c(1:isize-1, 1:jsize-1), ier)
    call cg_iric_write_sol_cell_real(fid, "Velocity[m_s]", vv_c(1:isize-1, 1:jsize-1), ier)
    call cg_iric_write_sol_cell_real(fid, "Usta[m_s]", usta_c(1:isize-1, 1:jsize-1), ier)
    call cg_iric_write_sol_cell_real(fid, "Tausta[-]", tausta_c(1:isize-1, 1:jsize-1), ier)
    call cg_iric_write_sol_cell_real(fid, "Qb[m2_s]", qb_c(1:isize-1, 1:jsize-1), ier)
    call cg_iric_write_sol_cell_real(fid, "Cc[-]", cc_c(1:isize-1, 1:jsize-1), ier)
    call cg_iric_write_sol_cell_real(fid, "Cb[-]", cb_c(1:isize-1, 1:jsize-1), ier)

! write edgei values
    call cg_iric_write_sol_iface_real(fid, "Vx[m_s]", vx(1:isize, 1:jsize-1), ier)
    call cg_iric_write_sol_iface_real(fid, "Qbx[m2_s]", qbx(1:isize, 1:jsize-1), ier)

! write edgej values
    call cg_iric_write_sol_jface_real(fid, "Vy[m_s]", vy(1:isize-1, 1:jsize), ier)
    call cg_iric_write_sol_jface_real(fid, "Qby[m2_s]", qby(1:isize, 1:jsize-1), ier)

end subroutine

【参考】構造格子と位置名

  |                |
-node----jface---node-------
  |                |
  |                |
iface   cell      iface
  |                |
  |                |
-node----jface----node-------
  |                |
j
↑---> i

おわりに

fortranで２次元の構造格子プログラムを記述するときの参考になれば幸いである。