import yt
runs = ['parent-46c7166', 'parent-46c7166-run2', 'fix-7c4dd95', 'tip-5d66537']
dsa = {}
for r in runs:
dsa[r] = yt.load('%s/RD0008/RedshiftOutput0008' % r)
fields = ['density', 'pressure', 'kinetic_energy', 'thermal_energy']
extrema = {}
for f in fields:
extrema[f] = [1e20, -1e20]
for name, ds in dsa.items():
ad = ds.all_data()
this_ex = ad.quantities.extrema(fields)
for i,f in enumerate(fields):
extrema[f][0] = min(extrema[f][0], this_ex[i][0])
extrema[f][1] = max(extrema[f][1], this_ex[i][1])
nbins = 128
profiles = {}
labels = {}
for f in fields:
profiles[f] = []
labels[f] = []
for r in runs:
ad = dsa[r].all_data()
for f in fields:
profiles[f].append(yt.create_profile(ad, f, fields="cell_mass", weight_field=None))
labels[f].append(r)
for f in fields:
p = yt.ProfilePlot.from_profiles(profiles[f], labels[f])
p.show()
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
yfield = 'cell_mass'
for f in fields:
y0 = profiles[f][0][yfield]
plt.clf()
alldy = np.array([])
for i in range(1, len(profiles[f])):
dy = (profiles[f][i][yfield] - y0) / y0
alldy = np.concatenate((alldy, dy.v))
plt.semilogx(profiles[f][i].x, dy, label=labels[f][i])
ysigma = alldy[np.isfinite(alldy)].std()
ysigma2 = alldy[np.isfinite(alldy) & (np.abs(alldy) < 3*ysigma)].std()
print ysigma, ysigma2
plt.ylim(-10*ysigma2, 10*ysigma2)
plt.title('Comparing to %s' % labels[f][0])
plt.xlabel(f)
plt.ylabel('Fractional difference')
plt.legend(loc='best')
plt.show()
0.00143771131973 0.000457584887241
0.00169442646893 0.000833166144564
0.0982165365194 0.0119970738231
0.0731508505264 0.00190050347502