\(\require{physics}\)

クエット流れ可視化プログラム

import numpy as np
import matplotlib.pyplot as plt
import glob
from matplotlib import animation

# extract data from a file
def get_data(file_name):
  columns = []
   for line in open(file_name).readlines():
    items = line.rstrip().split()
    if items[0] == "time":
      time = float(items[1])
    elif items[0] == "dz":
      dzs = [float(x) for x in items[1:]]
    elif items[0] == "limit_depth":
      pass
    elif items[0] == "items":
      pass
    elif items[0] == "column":
      columns.append({"number"     :int(items[1]),    "k_offset"   :int(items[2]),
                      "center_x"   :float(items[3]),  "center_y"   :float(items[4]),
                      "delta_x"    : float(items[5]), "delta_y"    :float(items[6]),
                      "btm_k"      :int(items[7]),    "sfc_k"      :int(items[8]),
                      "btm_z"      :float(items[9]),  "sfc_z"      :float(items[10]),
                      "container_i":int(items[11]),   "container_j":int(items[12]),
                      "container_i":int(items[13]),   "container_j":int(items[14]),
                      "status"     :items[15],        "variables"  :{}})
    else:
      columns[int(items[0])]["variables"][items[1]]=[float(x) for x in items[2:]]
  return time, dzs, columns

if __name__ == '__main__':

  files = sorted(glob.glob("couette*.dat"))

  # analytical solution
  h       = 1.0
  U       = 1
  z_level = np.arange(0.01, 1.0,0.02)
  fig = plt.figure()
  img0, = plt.plot([ U*z/h for z in z_level],z_level, label="Exact")
  imgs  = [[img0, img0]]
  for n,file in enumerate(files):
    t, dzs, columns = get_data(file)
    u = columns[0]["variables"]["u"]
    img, = plt.plot(u,z_level,"x",label="numerical")
    if n==0: plt.legend()
    plt.title("Couette flow")
    plt.xlabel("U (m/s)")
    plt.ylabel("z (m)")
    imgs.append([img0, img])

  ani = animation.ArtistAnimation(fig, imgs, interval=200, blit=True,repeat_delay=1000)
  ani.save('couette.gif', writer="imagemagick")