概要
以前に『Fortranでもデータベースしたい!』という題で記事を投稿しました。本稿ではそれに続いて、数値計算プログラムのデータ管理をLibpq-FortranとPostgreSQLで行う方法について例示したいと思います。
Libpq-Fortranとは
データベースサーバーPostgreSQLへのFortranインターフェースです。
Fortran2003以降のC相互運用機能を利用してPostgreSQLのCライブラリlibpqの関数をFortranから呼び出せるようにしました。libpqのFortranラッパーとも言えるでしょう。
これを利用するとFortranプログラムからPostgreSQLのサーバーに接続してクエリを実行し、データを読み書きすることができます。
数値計算のデータ管理
Fortranプログラムを扱う人のその目的は、多くの場合、数値計算プログラムを書いて結果を評価することでしょう。数値シミュレーションを例に、このプロセスをデータに注目して分けると以下のようになります(各サブプロセスで扱うデータファイルを括弧内に書いています)。
- ソースコードを書く(ソースファイル)
- シミュレーションを実行する(入力ファイルと出力ファイル)
- 結果の図をプロット(描画スクリプトと画像ファイル)
- 1に戻る
さて、これらのファイルはどのように取り扱えば良いでしょうか?
1番目のソースコードや3番目のスクリプトの管理には、GitやSVNなどのバージョン管理システムを使用することができます。
しかし一方で、2番のシミュレーションの初期値(入力ファイル)や計算結果(出力ファイル)を管理するためのベストプラクティスは自明ではありません。単純にファイルシステムで管理すると仮定すれば、以下のような方式が考えられるでしょう。
- 識別できるファイル名をつける
- 例:
hogefuga_202312022400_A.dat
- 実行した日時などを入れる。
- できるだけキーワードを入れて長い名前にする。
- 例:
- ディレクトリによるカテゴリ分けをする。
- 例:
result/hoge-scheme/fugapiyo_20231202.dat
- 例:
しかしこれらの方法も、並列処理をすると衝突してしまったり、ディレクトリパスが長くなってしまったりして、データを取り出すのに手間取ることになります。また、そもそも長い名前を考えるのは多少の面倒でもあります。
このような手間は、データ数が小さいうちは手動で間に合いますが、数が膨大になると管理しきれない場面が出てくることでしょう。
入出力データのデーターベース管理
前のセクションで、数値計算の結果を管理するのは悩みの種になることについて述べました。このセクションではデータベース管理システムの一つであるPostgreSQLを利用して、それらのデータを管理する方法について詳しく見ていきます。
Libpq-Fortranを使用すると、Fortranプログラムからデータベースサーバーにアクセスすることが可能になります。以下では実際にPostgreSQLデータベースにアクセスして、どのように数値計算の結果を登録するかについて例示します。
例:一次元移流方程式の数値シミュレーション
一次元移流方程式の数値シミュレーションを例にとり、実際にパラメーターを色々変えて実行し、結果をデータベースに登録してみましょう。
一次元移流方程式
$$\frac{\partial u}{\partial t} + \frac{\partial f}{\partial x} = 0, f= cu,$$
を陽的オイラー法と、Lax-FriedrichスキームおよびLax-Wendroffスキームでそれぞれ離散化すると、次のような差分方程式が得られます。
$$
u^{n+1}_{j}= u^{n}_{j} - \frac{\Delta t}{\Delta x}\left(\tilde f^{n+\frac{1}{2}}_{j+\frac{1}{2}}-\tilde f^{n-\frac{1}{2}}_{j-\frac{1}{2}}\right)
$$
ここで$\tilde f$は数値流速であり、各スキームで以下のように取ります。
Lax-Friedrichスキーム:
$$
\tilde f^{n+\frac{1}{2}}_{j+\frac{1}{2}}=\frac{1}{2}\left[ \left(f^n_{j+1}+f^n_j\right)-\frac{\Delta x}{\Delta t}\left(u^n_{j+1}-u^n_j\right)\right],
$$
Lax-Wendroffスキーム:
$$
\tilde f^{n+\frac{1}{2}}_{j+\frac{1}{2}} = \frac{1}{2}\left[\left(1 - c\frac{\Delta t}{\Delta x}\right)f^n_{j+1}+ \left(1+ c\frac{\Delta t}{\Delta x}\right)f^n_j\right],
$$
なお、ここで$$f^n_{j} = f(u^n_j) = cu^n_j$$です。
また、$c$は特性速度、$\Delta t$は時間ステップ幅、$\Delta x$は空間刻み幅と呼ばれる量になり、この方程式のシミュレーション1回を特徴づけるパラメーターになります。
ユーザーとデータベースを作る
以下のFortranプログラムでは、PostgreSQLサーバーに接続して、2つのクエリを実行しています。
program create_database
use, intrinsic :: iso_c_binding
use :: libpq
implicit none
type(c_ptr) :: conn, res
character(:), allocatable :: conninfo, query
conninfo = "host=localhost user=postgres" ! 接続情報
!(必要なら管理者ユーザーpostgresのパスワードをpassword=以下に記述する)
conn = PQconnectdb(conninfo) ! 接続の実行
! エラーハンドリング
if (PQstatus(conn) /= CONNECTION_OK) then
print *, PQerrorMessage(conn)
error stop
end if
! ユーザー作成のクエリ
query = "create role shinobu with createdb login password 'hogehoge';"
res = PQexec(conn, query) ! クエリ実行
! エラーハンドリング
if (PQresultStatus(res) /= PGRES_COMMAND_OK ) then
print *, PQresultErrorMessage(res)
call PQfinish(conn)
error stop
end if
! 結果オブジェクトの解放
call PQclear(res)
! データベース作成のクエリ
query = "create database simulation with owner shinobu;"
res = PQexec(conn, query) ! クエリの実行
! エラーハンドリング
if (PQresultStatus(res) /= PGRES_COMMAND_OK ) then
print *, PQresultErrorMessage(res)
end if
call PQclear(res) ! 結果オブジェクトの解放
call PQfinish(conn) ! 接続オブジェクトの解放
end program create_database
データベースへ接続
最初にconninfo
変数に接続情報(ホスト名localhost
、ユーザー名postgres
)を記述し、関数PQconnectdb
を呼び出して、戻り値をtype(c_ptr)
の変数conn
に代入しています。
その後のif
文はPQstatus(conn) /= CONNECTION_OK
という条件文を評価して、エラーハンドリングを行っています。
ユーザー作成のクエリ
次にユーザー作成のクエリの部分では、create role shinobu with createdb login password 'hogehoge';
という文字列をquery
変数に格納した後、PQexec
関数を呼び出して、その結果をtype(c_ptr)
型変数res
に代入しています。
ここでquery
変数の中身にはデータベースサーバーに対するコマンドを含んでいて、以下のような意味を持ちます。
-
create role shinobu
:shinobu
という名前でロール(ユーザー)を作成するコマンド -
with createdb
:当該ロールにデータベースを作成する権限を付与する。 -
login password 'hogehoge'
:当該ロールのログインパスワードをhogehoge
に設定する。
ここで、ユーザー名とパスワードは例示したものではなく各自で置き換えてください。
PQexec
関数は、第一引数にPQconnectdb
で返された接続オブジェクトへのポインタを、第二引数に上記のクエリを取り、結果オブジェクトへのポインタtype(c_ptr)
型の値を返します。ここではtype(c_ptr)
型の変数res
に代入しています。
コマンドが成功した場合はPQresultStatus(res)
がPGRES_COMMAND_OK
の値を返します。失敗した場合にはこれ以外の値を返し関数PQresultErrorMessage(res)
がエラーメッセージを返します。
成功した場合は、もうres
オブジェクトは不要なのでPQclear(res)
関数を呼び出して、このオブジェクトを解放します。
データベース作成のクエリ
次にデータベース作成のクエリの文字列create database simulation with owner shinobu;
をquery
変数に格納しています。
ここでクエリに含まれるコマンドの意味は次のようなものです:
-
create database simulation
:simulation
という名のデータベース(テーブルの集合の名前)を作成する。 -
with owner shinobu
:当該データベース(simulation
)の所有者をshinobu
にする。
その後、この文字列をPQexec
関数に渡して呼び出しています。エラーハンドリングについては上と同様です。
成功した場合、結果のオブジェクトres
と接続オブジェクトconn
を解放して、このプログラムは終了です。
この時点で、ターミナルからpsql -h localhost -U shinobu -d simulaiton
を実行して、ここで作成したデータベースsimulation
に接続することができるはずです。
テーブルを作る
次に、実際にシミュレーションの情報を書き込むためのテーブルを用意します。
conninfo
に記述する接続情報には、前のセクションで作成したロール・データベース・パスワードを使用します。
program create_table
use :: libpq
use, intrinsic :: iso_c_binding
implicit none
type(c_ptr) :: conn, res
character(:), allocatable :: conninfo, query
integer :: i
conninfo = "host=localhost user=shinobu dbname=simulation password=hogehoge"
conn = PQconnectdb(conninfo)
if (PQstatus(conn) /= 0) then
print *, PQerrorMessage(conn)
error stop
end if
query = "create table advection &
& (sim_id decimal(4,0) not null, &
& start_date date , &
& sim_name varchar(100) not null, &
& cs double precision not null, &
& dt double precision not null, &
& dx double precision not null, &
& courant double precision , &
& length double precision , &
& nx integer not null, &
& nstep integer not null, &
& dat_file varchar(256) not null, &
& primary key(sim_id));"
res = PQexec(conn, query)
if (PQresultStatus(res) /= PGRES_COMMAND_OK) then
print *, PQresultErrorMessage(res)
end if
call PQclear(res)
call PQfinish(conn)
end program create_table
変数query
には最初にcreate table advection
とあり、ここでadvection
という名前のテーブルを定義することを宣言しています。その後、括弧の中に属性名・型・オプションを1つの属性(テーブルの列)としてカンマ区切りで定義を並べています。各属性の記述の意味は次の通りです。
-
sim_id
: 実験ID、整数型(4桁を小数点以下0桁精度)、非NULL制約1あり -
start_date
:実行した日付、日付型 -
sim_name
:シミュレーションの名前、最大100桁の可変長文字列型、非NULL制約あり -
cs
:特性速度の数値、倍精度実数型、非NULL制約あり -
dt
:時間ステップ幅、倍精度実数型、非NULL制約あり -
dx
:空間刻み幅、倍精度実数型、非NULL制約あり -
courant
:クーラン数=$c\Delta t/\Delta x$、倍精度実数型 -
length
:計算領域長さ、倍精度実数型 -
nx
:計算領域の分割数、整数型、非NULL制約あり -
nstep
:計算した時間ステップの総数、非NULL制約あり -
dat_file
:データファイルのパス名、最大256桁の可変長文字列型、非NULL制約あり
最後にこのテーブルでは一意にレコードを特定できる属性として属性sim_id
を選ぶことにするので、一意性制約primary key(sim_id)
を記述しています。
これらを一つの変数query
に含ませた状態で、PQexec
を実行します。
このプログラムを実行した後に、ターミナルからpsql -h localhost -U shinobu -d simulation
を実行した上で、psql
のプロンプトからselect * from advection;
を問い合わせて、以下のような属性名の一覧を得られれば成功です。
$ psql -h localhost -U shinobu -d simulation
psql (15.4)
Type "help" for help.
simulation=> select * from advection;
sim_id | start_date | sim_name | cs | dt | dx | courant | l | nx | nstep | dat_file
--------+------------+---------------+----+--------+-------+---------+---+-----+-------+--------------------
以下のセクションではこのテーブルに対してデータの登録や取り出しを試していきます。
シミュレーションを実行する。
データを保持しておくデータベースとテーブルが用意できたので、シミュレーションを実行していきます。
ここで実行するプログラムは以下の要素から構成されています。
- 主プログラム:1個の実験をパラメータを変更して全部で32個実行するループを持つ。
- 内部サブルーチン
-
main_loop
サブルーチン:引数にシミュレーションを特徴づけるパラメータを受け取り、1個のシミュレーションを実行する。 -
data_output
サブルーチン:出力ファイルの第1列にx、第2列にuの値を書き出し、区切りを2つの空行とするデータセットを書き込む。 -
data_registration
サブルーチン:main_loop
の後に、実験のパラメータと出力ファイルへのパスをデータベースに登録する
-
時間積分するメインのループ
数値計算のループはmain_loop
というサブルーチンに含まれているものとします。このサブルーチンの引数に、1回の実験を特徴づけるパラメータを渡して、変数filename
にはファイル名を作成して格納し、その名前の.dat
ファイルに結果を書き込みます。
subroutine main_loop(cs, dt, dx, nstep, filename, offset_k, k, prefix)
implicit none
real(real64), intent(in) :: cs, dt, dx ! cs is short for characteristic speed
integer(int32), intent(in) :: k
character(:), allocatable, intent(inout) :: filename
integer :: offset_k
character(*), intent(in) :: prefix
integer, intent(out) :: nstep
character(4) :: code
integer :: nth_out
real(real64) :: tout, time
real(real64) :: Courant
integer :: i
integer :: uni
real(real64), allocatable :: x(:)
real(real64), allocatable :: u(:), f(:), u_new(:), f0(:)
write(code, '(1i4.4)') k + (offset_k - 1)
Courant = cs*dt/dx
write(error_unit, '(a, i0, a)') "### START ", k, "-th Simulation ###"
write(error_unit, '(a, f10.3)') " Courant number: ", Courant
! Initialize
uni = 100
nth_out = 0
nstep = 1
time = 0d0
tout = 0
allocate(x(0:nx), source=0d0)
allocate(u(0:nx), source=0d0)
allocate(u_new(0:nx), source=0d0)
allocate(f0(0:nx), source=0d0)
allocate(f(0:nx), source=0d0)
! 空間のグリッドをセットする
x(0) = 0d0
do i = 1, nx
x(i) = x(0) + dx*dble(i)
end do
! 初期条件をセットする
do i = 0, nx
if (x(i) < 0.1d0) then
u(i) = 1d0
else
u(i) = 0d0
end if
end do
filename = 'result/'//trim(adjustl(prefix))//code//'.dat'
open(uni, file=filename, form='formatted', status='replace')
call data_output(uni, x, u, tout, time, 0, nth_out)
! 時間積分
do
nstep = nstep + 1
time = time + dt
f0(:) = cs*u(:) ! 数値流束
! Lax-Friedrich スキーム
if (trim(prefix) == 'LF_') then
do i = 1, nx-1
f(i) = 0.5d0*(f0(i+1)+f0(i)) - 0.5d0*dx/dt *(u(i+1) - u(i))
end do
end if
! Lax-Wendroff スキーム
if (trim(prefix) == 'LW_') then
do i = 1, nx-1
f(i) = 0.5d0*( (1d0-Courant)*f0(i+1) + (1d0+Courant)*f0(i) )
end do
end if
! 陽的オイラー法
do i = 1, nx-1
u_new(i) = u(i) - Courant*(f(i) - f(i-1))
end do
! uのアップデート
u(2:nx-1) = u_new(2:nx-1)
! 境界条件
u(1) = u(2)
u(nx) = u(nx-1)
if (time >= tout) call data_output(uni, x, u, tout, time, nstep, nth_out)
if (time > tend) exit
end do
close(uni)
write(error_unit, fmtStop) nstep, time
write(error_unit, '(a)') "### Normal STOP ###"
end subroutine main_loop
main_loop
サブルーチンに含まれるdoループ内部に、セクションの冒頭で式で示したLax-FriedrichスキームとLax-Wendroffスキームの数値流束の計算を引数prefix
の値に応じて切り替えて実行するようにしています。
ファイル出力
データ書き込み用のサブルーチンは以下の通りです。
subroutine data_output (uni, x, u, tout, time, nstep, nth_out)
implicit none
integer, intent(in) :: uni
real(real64), intent(in) :: x(:), u(:)
real(real64), intent(inout) :: tout
real(real64), intent(in) :: time
integer, intent(in) :: nstep
integer, intent(inout) :: nth_out
integer :: i
do i = 1, nx
write(uni, fmtOut) x(i), u(i)
end do
write(uni, '(/)') !空行を2つ書き込む(gnuplot用のデータセットの区切り)
tout = tout + tout_interval
write(error_unit, fmtErrOut) nstep, time, nth_out
nth_out = nth_out + 1
end subroutine data_output
データを登録する
データベースにデータを登録するためのサブルーチンは次の通りです。ここでのクエリはinsert文と呼ばれるもので、テーブルにレコードを挿入します。
subroutine data_registration(conn, simu_name, filename, offset_k, k, cs, dt, dx, nstep, date)
use :: iso_c_binding
use :: libpq
implicit none
type(c_ptr), intent(in) :: conn
character(*), intent(in) :: simu_name
character(*), intent(in) :: filename
integer, intent(in) :: k, offset_k
real(real64), intent(in) :: cs, dt, dx
integer, intent(in) :: nstep
character(*), intent(in) :: date
character(:), allocatable :: query
character(256) :: values(11)
type(c_ptr) :: res
real(real64) :: Courant
character(4) :: code
integer :: i
Courant = cs*dt/dx
! 文字列配列に値を格納する
write(code, '(1i4.4)') k + (offset_k-1)
values(1) = code
values(2) = trim(adjustl(date))
values(3) = trim(simu_name)
write(values(4), '(f16.8)') cs
write(values(5), '(f16.8)') dt
write(values(6), '(f16.8)') dx
write(values(7), '(f16.8)') Courant
write(values(8), '(f16.8)') L
write(values(9), '(i0)') nx
write(values(10), '(i0)') nstep
values(11) = trim(adjustl(filename))
do i = 1, 9
values(i) = trim(adjustl(values(i)))
end do
! プレースホルダーには文字列配列valuesの要素の値が展開される
query = 'insert into advection &
& (sim_id, start_date, sim_name, cs, dt, dx, courant, l, nx, nstep, dat_file) &
& values ($1, $2, $3, $4, $5, $6, $7, $8, $9, $10, $11)'
res = PQexecParams(conn, query, 11, [(0, i=1,11)], values(:))
if (PQresultStatus(res) /= PGRES_COMMAND_OK) then
print *, PQresultErrorMessage(res)
end if
call PQclear(res)
end subroutine data_registration
データベースへの登録(データの挿入)を行うには、PQexecParams
関数を使用してクエリを実行します。このPQexecParams
関数は、クエリ文字列にはドル記号と数字からなるプレースホルダー(コード例では$1
から$11
)を記述して、そこの値は文字列の配列として別に渡すことができます。PQexecParams
関数に渡す引数の内容は次の通りです。
PQexecParams 関数 |
引数の中身 |
---|---|
第1引数 | 接続オブジェクトへのポインタ(conn ) |
第2引数 | パラメータのプレースホルダーを含むクエリ文字列(query ) |
第3引数 | 第2引数に含まれるパラメータの個数(整数型) |
第4引数 | パラメータの個数だけの長さを持つゼロ値配列 (ここでのゼロはパラメータの型の解釈をサーバーに委ねるの意味) |
第5引数 | パラメータの個数だけの長さを持つ文字列型の配列 |
data_registration
サブルーチンの24行目から34行目の間に、第5引数にわたす予定の文字列の配列value(1:11)
にデータの値を書き込んでいます。
メインプログラムを実行する
ここまでで必要な材料はすべて揃いました。メインプログラムを書いてシミュレーションを実行します。メインプログラムの全体はGitHubのこちらのリンクから入手することができます。
program main
use, intrinsic :: iso_fortran_env
use, intrinsic :: iso_c_binding
use :: libpq
implicit none
real(real64), parameter :: tend = 0.6d0
real(real64), parameter :: tout_interval = 0.1d0
integer, parameter :: NX = 200 ! 分割格子の数
real(real64), parameter :: L = 1d0 ! 計算領域長さ
character(*), parameter :: fmtOut = '(e10.3, 1x, e10.3, 1x, i0)'
character(*), parameter :: fmtErrOut= "(' output ', 'step=', i8, ' time=', e10.3, ' nth_out =', i3)"
character(*), parameter :: fmtStop = "(' stop ', 'step=', i8, ' time=', e10.3)"
character(:), allocatable :: filename
integer :: k
real(real64) :: cs, dt, dx
integer :: nstep, offset_k
type(c_ptr) :: conn
character(:), allocatable :: conninfo
character(10) :: date = '2023-12-02'
conninfo = 'host=localhost dbname=simulation user=shinobu'
conn = PQconnectdb(conninfo)
if (PQstatus(conn) /= 0) then
print *, PQerrorMessage(conn)
error stop
end if
! Lax-Friedrichスキームの実験を16件実行する(id = 1:16)。
offset_k = 1
do k = 1, 16
cs = 1d0
dt = 5d-4 + 3d-4*dble(k-1)
dx = L/dble(nx)
call main_loop(cs, dt, dx, nstep, filename, offset_k, k, 'LF_')
call data_registration(conn, 'Lax-Friedrich', filename, offset_k, k, cs, dt, dx, nstep, date)
end do
! Lax-Wendroffスキームの実験を16件実行する(id = 17:32)。
offset_k = 17
do k = 1, 16
cs = 1d0
dt = 5d-4 + 3d-4*dble(k-1)
dx = L/dble(nx)
call main_loop(cs, dt, dx, nstep, filename, offset_k, k, 'LW_')
call data_registration(conn, 'Lax-Wendroff', filename, offset_k, k, cs, dt, dx, nstep, date)
end do
contains
...
end program main
データベースに登録されているのを確認する
ターミナルからpsql -h localhost -U shinobu -d simulation
を実行してpsqlを起動し、select * from advection;
を実行すると32行のレコードが挿入されていることが確認できます。
以上で、数値計算の結果のファイルをデータベースサーバーに登録することができました。次のセクションでは登録されたデータを取り出してFortranから図を描画する方法を述べたいと思います。
データベースからデータを取り出して図にプロットする
PostgreSQLにアクセスできるプログラミング言語ならFortranでなくとも、セクションタイトルのことをできますが、ここではLibpq-Fortranでデータをサーバーから読み込んでgnuplotで図を描画する方法について説明したいと思います2。
ここでの方針は、simulation
データベースのadvection
テーブルから一部のレコードを抽出して、パス名だけを出力し、文字列型変数にスペース区切りで連結してプロットスクリプトのコマンドライン引数として渡します。
以下がplot.f90
の例です。
program plot
use :: libpq
use, intrinsic :: iso_c_binding
implicit none
character(:), allocatable :: conninfo, query, cmd
type(c_ptr) :: conn, res
integer :: i
conninfo = "host=localhost dbname=simulation user=shinobu"
conn = PQconnectdb(conninfo)
if (PQstatus(conn) /= 0) then
print *, PQerrorMessage(conn)
error stop
end if
! 実行するコマンドの先頭を記述する
cmd = "gnuplot -c plot.gp "
query = "select dat_file from advection where sim_name = 'Lax-Wendroff' and courant > 0.55;"
res = PQexec(conn, query)
if (PQresultStatus(res) /= PGRES_COMMAND_OK) then
print *, PQresultErrorMessage(res)
call PQfinish(conn)
error stop
end if
do i = 0, PQntuples(res) - 1
! クエリで取り出したデータをcmdに追記していく
cmd = cmd // ' '// PQgetvalue(res, i, 0)
end do
! コマンドを実行していく
print *, cmd
call system(cmd)
call PQclear(res)
call PQfinish(conn)
end program plot
ここで、"select dat_file from advection where sim_name = 'Lax-Wendroff' and courant > 0.55;"
というクエリは、以下のような意味を持ちます。
query の一部分 |
意味 |
---|---|
from advection |
advection テーブルから、 |
where sim_name ='Lax-Wendroff' and courant > 0.55 |
属性sim_name がLax-Wendroff であり、かつcourant が0.55より大きいレコードを抽出し、 |
select dat_file |
dat_file 属性の値(すなわちパス名)を出力する。 |
このクエリの結果はc_ptr
型変数res
の示すオブジェクトから、libpqの関数を使用して取り出します。PQexec
の後のdoループでは、0からPQntuples(res) - 1
までのループを指示しています。そしてそのループの中で、PQgetvalue(res, i, 0)
で第i
番目のレコードの第0列のデータを取り出し、cmd
変数に連結しています。
その後、system
サブルーチンを呼び出し、cmd
に含まれるコマンドを実行しています。
実際に実行されるコマンドは次のようになるでしょう。
gnuplot -c plot.gp result/LW_0024.dat result/LW_0025.dat result/LW_0026.dat result/LW_0027.dat result/LW_0028.dat result/LW_0029.dat result/LW_0030.dat result/LW_0031.dat result/LW_0032.dat
ここで、plot.gp
は次のようなGnuplotスクリプトになります。
set term pngcairo size 1600,1200
set termoption noenhanced
set output "plot9.png"
set xlabel "x"
set ylabel "u"
set xrange [0:1]
set yrange [-0.5:2.0]
set multiplot layout 3,3
do for [i=1:9] {
path = ARGV[i]
set title path
unset xlabel
unset ylabel
if ( i % 3 == 1) {
set ylabel "u"
}
if ( i > 6 ) {
set xlabel "x"
}
plot path index 0 using 1:2 with lines title 't = 0.0', \
path index 1 using 1:2 with lines title 't = 0.1', \
path index 2 using 1:2 with lines title 't = 0.2', \
path index 3 using 1:2 with lines title 't = 0.3', \
path index 4 using 1:2 with lines title 't = 0.4', \
path index 5 using 1:2 with lines title 't = 0.5', \
path index 6 using 1:2 with lines title 't = 0.6
}
unset multiplot
このプログラム(plot.f90
)を実行すると次のような画像を出力することができます。
以上で、Libpq-FortranでデータをPostgreSQLサーバーから読み込んで、gnuplotスクリプトで図を描画することができました。
おわりに
ライブラリLibpq-Fortranを使用して、PostgreSQLで数値シミュレーションのデータを管理し、そのデータをサーバーから取り出して使用する方法について述べました。
この方法は同じプログラムをパラメーターを変えて大量に実行するような場面で活躍すると思われます。筆者はある差分法スキームの安定性を評価するためにパラメーターの全探索をしたときに(手で管理するのを諦めて)データベースを活用する方法を選択したことがありました。そのときはLibpq-Fortranがなかったので別の方法(Pythonだったかな?)で行っていましたが、今では本稿で説明した通りFortranプログラムから直接登録することができるようになりました。
今回扱った例では、1個のシミュレーションが1つのレコードに対応していたので比較的簡単に管理することができました。一方、1つのシミュレーションに対して複数のファイルが出力されるような場合を想定すると、データ管理の難易度はやや上がります。これについてはまたいずれ書きたいと思います。
コードのテスト環境
- Windows 10, GFortran MinGw-W64 13.1.0, Gnuplot 6.0-rc1
- Gentoo Linux, GFortran 13.2.1, Gnuplot 5.4.6