LoginSignup
3
1

More than 3 years have passed since last update.

OpenFOAMの時系列データとsetsデータの読み込み

Last updated at Posted at 2020-01-18

概要

OpenFOAMの時系列データとsetsデータ(sampleコマンドで出力する空間分布のデータ)をPythonで読み込むための関数をつくって,matplotlibでのプロットを簡略化しました。

前提として,OpenFOAMの時系列データはリスタートの度に新しいディレクトリに保存されますが,それをcatコマンドでくっつけてひとつづきのファイルになってるとします。

Pythonのコード

時系列データの読み込み関数用関数

  • caseはケースファイル名
  • positionsには時系列データの場所の情報をリストでいれます。
  • variablesが時系列データの物理量(p, Uとか)のリスト
  • 各位置の時系列データを一括で読み込みます。
  • 「物理量_場所」の文字列をkey,実際のデータをndarrayとしてvalueにもつディクショナリデータを返します。
def readHistory(case, positions, variables):
    import numpy as np

    """
    case: case name (string data)
    positions: position name (list data)
    variables: physical quantities (list data)
    """

    data ={}

    for position in positions:
        for variable in variables:
            key = variable+"_"+position
            value = np.genfromtxt(case+"/"+position+"/"+variable, unpack=True)
            data[key] = value

    return data

setsデータの読み込み用関数

引数の説明

  • caseはケースファイル名
  • fileはsetsデータのファイル名
  • variablesが出力したい物理量(p, Uとか)のリスト
  • setsデータは複数の物理量のデータが1つのファイルに出力されるので,各物理量を個別に読み込みます。
  • 「物理量_場所」の文字列をkey,実際のデータをdataframeとしてvalueにもつディクショナリデータを返します。
def readSets(case, file, variables):
    import numpy as np
    import pandas as pd
    import glob
    import os

    """
    case: case name (string data)
    file: sets data file name (string data)
    variables: physical quantities (list data)
    """

    # list of timing
    list1 = glob.glob(case + "/sets/*")

    variable_names = file.split("_")
    # remove ".xy" 
    variable_names[-1] = variable_names[-1].split(".")[0]
    position, variable_names[0] = variable_names[0], "x"
    #print(variable_names)

    data = {}
    for variable in variables:
        key = variable+"_"+position
        data[key] = pd.DataFrame()

        for i in list1:
            path = i+"/" + file
            df = pd.read_table(path, header = None, names = variable_names)   
            data[key] = pd.concat([data[key], df[variable]], axis=1,)

        list_index=[]    
        for i in list1:
            list_index.append(str(float(os.path.basename(i))))

        data[key].columns = list_index
        data[key] = pd.concat([df[variable_names[0]], data[key].sort_index(axis=1)], axis=1)

    # return dictionary of dataframe
    return data

matplotlibでプロット

matplotlibでデータをプロットします。

データの配置は以下のようになってるとします。
test01がケース名。

test01/position-A/p
test01/position-A/U
test01/position-B/p
test01/position-B/U
test01/sets/0/center_p_k.xy
test01/sets/10/center_p_k.xy

p, Uはposition-A, Bの時系列データ
center_p_k.xyは時刻0,10の分布データ
です。

ケースが複数あって,ケース名だけ変えたプロットを簡単に作ることを想定した作りになってます。

import matplotlib.pyplot as plt
%matplotlib inline

case = "test01"
history_position = ["position-A", "position-B"]
history_variables = ["p", "U"]
sets_variables = ["p", "k"]
sets_position = "center"
dataSet_test01= [readHistory(case, history_position, history_variables), 
                  readSets(case, "center_p_k.xy", sets_variables)]

# 時系列データ一式のプロット
for i in history_position:
    for j in history_variables:
        fig = plt.figure(figsize=(8,8))
        ax = fig.add_subplot(111)

        ax.plot(dataSet_test01[0][j+"_"+i][0], dataSet_test01[0][j+"_"+i][1], label=j)
        ax.legend(bbox_to_anchor=(1.01, 1., 0., 0), loc='upper left', borderaxespad=0.,)

# 空間分布のプロット
for sets in sets_variables:
    for column in dataSet_test01[1][sets+"_"+sets_position]:
        fig = plt.figure(figsize=(8,8))
        ax = fig.add_subplot(111)

        ax.plot(dataSet_test01[1][sets+"_"+sets_position]["x"], 
            dataSet_test01[1][sets+"_"+sets_position][column], 
            label=sets+"_"+column)
        ax.legend(bbox_to_anchor=(1.01, 1., 0., 0), loc='upper left', borderaxespad=0.,)
3
1
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
3
1