In [ ]:

    

import yt
import sys
import matplotlib
matplotlib.rcParams['figure.dpi'] = 150
import matplotlib.pyplot as plt
from matplotlib import animation
from matplotlib.animation import FuncAnimation
from matplotlib import rc_context
from tempfile import NamedTemporaryFile
from IPython.display import display, HTML
import numpy as np
import base64

yt.mylog.setLevel("ERROR")

IMG_TAG = """<img src="data:image/gif;base64,{0}" alt="some_text">"""
 
def anim_to_gif(anim):
    data="0"
    with NamedTemporaryFile(suffix='.gif') as f:
        anim.save(f.name, writer='imagemagick', fps=10);
        print(f.name)
        data = open(f.name, "rb").read()
        #print(data)
        data = base64.b64encode(data).decode()
        #print(data)
    return IMG_TAG.format(data)
 
def display_animation(anim):
    plt.close(anim._fig)
    return display(HTML(anim_to_gif(anim)))


    

In [ ]:

    

ts = yt.load('/d/d9/ychen/FLASH4/20180310_BrioWu/bw_2d_hdf5_plt_cnt_*')
#ts = yt.load('/d/d9/ychen/FLASH4/20180321_shock/oblique_shock_2d_hdf5_plt_cnt_*')


    

In [ ]:

    

field = 'density'
mi, ma = 1E-1, 1E1
log = True

plot = yt.SlicePlot(ts[0], 'z', field, width=(1,0.2), origin='native')
plot.set_zlim(field, mi, ma)
plot.set_log(field, log)

fig = plot.plots[field].figure

# animate must accept an integer frame number. We use the frame number
# to identify which dataset in the time series we want to load
def animate(i):
    ds = ts[i]
    sys.stdout.write('\r%s' % ds.basename)
    plot._switch_ds(ds)

anim = FuncAnimation(fig, animate, frames=15)


    

In [ ]:

    

display_animation(anim)


    

In [ ]:

    

field = 'pressure'
mi, ma = 1E0, 1E1
log = True

plot = yt.SlicePlot(ts[0], 'z', field, width=(1,0.2), origin='native')
plot.set_zlim(field, mi, ma)
plot.set_log(field, log)

fig = plot.plots[field].figure

# animate must accept an integer frame number. We use the frame number
# to identify which dataset in the time series we want to load
def animate(i):
    ds = ts[i]
    sys.stdout.write('\r%s' % ds.basename)
    plot._switch_ds(ds)

anim = FuncAnimation(fig, animate, frames=15)


    

In [ ]:

    

display_animation(anim)


    

In [ ]:

    

ds = ts[100]
print(ds)
ds.field_list
ds.dimensionality = 3
xl = 0.0
xr = 0.9
sl = ds.r[xl:xr,0.5,0.5]
sl['shok']


    

In [ ]:

    

fig = plt.figure(figsize=(10,4))
fields = ['density', 'pressure', 'velocity_x', 'sound_speed', 'shks', 'shok']
ln = {}
for fi in fields:
    ln[fi], = plt.plot(sl['x'], sl[fi], '.-', label=fi)
#ln['entropy'], = plt.plot(sl['x'], sl['pres']/sl['dens']**(sl['gamc']+1), '.-', label='entropy')
plt.semilogy()
plt.ylim(1E-1, 1E1)
plt.legend()
shok = sl['shok']
shok_begins, = np.where(shok[1:] - shok[:-1] == 1.0)
shok_ends, = np.where(shok[1:] - shok[:-1] == -1.0)
dx = sl['x'][1:]-sl['x'][:-1]
xx = []
widths = []
for s0, s1 in zip(shok_begins, shok_ends):
    xx.append(sl['x'][s0]+0.5*dx[s0])
    widths.append(sl['x'][s1]-sl['x'][s0])
bars = plt.bar(xx, [1E3]*len(xx), widths, color='y', alpha=0.2, bottom=1E-1)
for i, (x, width) in enumerate(zip(xx,widths)):
    bars.patches[i].set_x(x)
    bars.patches[i].set_width(width)
text = plt.text(0.03, 0.90, 't=%.4f' % ds.current_time,
    horizontalalignment='left', verticalalignment='center',
    transform=plt.axes().transAxes)
plt.legend(loc=1)


    

In [ ]:

    

def animate(i):
    ds = ts[i]
    sys.stdout.write('\r%s' % ds.basename)
    ds.field_list
    ds.dimensionality = 3
    sl = ds.r[xl:xr,0.5,0.5]
    for fi in fields:
        ln[fi].set_data(sl['x'], sl[fi])
    xx = np.zeros(len(shok_begins))
    widths = np.zeros(len(shok_begins))
    shok = sl['shok']
    shok_begins, = np.where(shok[1:] - shok[:-1] == 1.0)
    shok_ends, = np.where(shok[1:] - shok[:-1] == -1.0)
    dx = sl['x'][1:]-sl['x'][:-1]
    for i, (s0, s1) in enumerate(zip(shok_begins, shok_ends)):
        xx[i] = sl['x'][s0]+0.5*dx[s0]
        widths[i] = sl['x'][s1]-sl['x'][s0]
    for i, (x, width) in enumerate(zip(xx,widths)):
        bars.patches[i].set_x(x)
        bars.patches[i].set_width(width)
    text.set_text('t=%.4f' % ds.current_time)

anim = FuncAnimation(fig, animate, frames=100)


    

In [ ]:

    

display_animation(anim)


    

In [ ]:

    

np.sqrt(sl['shks']*(sl['gamc']+1)/2/sl['gamc']+1)


    

In [ ]:

    

sl['velocity_x']


    

In [ ]:

    

0.1/0.05


    

In [ ]:

    

sl['pres']/sl['dens']**(sl['gamc']+1)


    

In [ ]:

    

hok_begins, = np.where(shok[1:] - shok[:-1] == 1.0)
shok_ends, = np.where(shok[1:] - shok[:-1] == -1.0)