Skip to content

Commit

Permalink
Add all changes made before thesis revision.
Browse files Browse the repository at this point in the history
  • Loading branch information
quantheory committed Nov 23, 2020
1 parent c9c7eae commit 18e8169
Show file tree
Hide file tree
Showing 14 changed files with 12,022 additions and 172 deletions.
501 changes: 501 additions & 0 deletions AccretionExample.ipynb

Large diffs are not rendered by default.

51 changes: 31 additions & 20 deletions complex_eigenvalues.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,9 @@
blocksize = 4096
num_files = 12
end_column = 48601
#blocksize = 4096
#num_files = 1
#end_column = 4095

splits = [(blocksize*i, blocksize*(i+1)-1) for i in range(num_files)]
splits[-1] = (splits[-1][0], end_column)
Expand All @@ -31,7 +34,7 @@
efiles.append(nc4.Dataset(name, 'r'))

plt.autoscale(tight=True)
cmap = plt.get_cmap('Reds')
cmap = plt.get_cmap('Greys')

total_evalues = 0
for efile in efiles:
Expand Down Expand Up @@ -75,7 +78,7 @@
circle_rad = np.log10(1./300.)-cutoff_exp
circle_x = circle_rad * np.cos(np.linspace(0., 2.*np.pi, 1001))
circle_y = circle_rad * np.sin(np.linspace(0., 2.*np.pi, 1001))
plt.plot(circle_x, circle_y)
plt.plot(circle_x, circle_y, color='b')
converge_rad = 1./300.
converge_x = converge_rad * (np.cos(np.linspace(0., 2.*np.pi, 1001)) - 1.)
converge_y = converge_rad * np.sin(np.linspace(0., 2.*np.pi, 1001))
Expand All @@ -84,21 +87,21 @@
plot_abs = np.maximum(np.log10(converge_abs) - cutoff_exp, 0.)
plot_x = plot_abs * np.cos(converge_angle)
plot_y = plot_abs * np.sin(converge_angle)
plt.plot(plot_x, plot_y)
plt.plot(plot_x, plot_y, color='g')
plt.axvline(x=0., color='k', linewidth=2.)
plt.axhline(y=0., color='k', linewidth=2.)
plt.axis([-max_plotval, max_plotval, -max_plotval, max_plotval])
ticks = np.arange(-max_plotval, max_plotval+1., 1)
ticks = np.arange(-max_plotval, max_plotval+1., 2)
# The "0.01" hack below is there to ensure that "0" outputs a 10 with a positive sign.
tickvals = ["${}^{{{}}}$".format(int(np.sign(i+0.01)*10),int(abs(i)+cutoff_exp)) for i in ticks]
ax = plt.gca()
ax.set_xticks(ticks)
ax.set_xticklabels(tickvals)
ax.set_xticklabels(tickvals, fontsize=16)
ax.set_yticks(ticks)
ax.set_yticklabels(tickvals)
plt.clim(vmin=0., vmax=1.e4)
ax.set_yticklabels(tickvals, fontsize=16)
plt.clim(vmin=0., vmax=8.e3)
plt.colorbar()
plt.savefig('./complex_eigenvalues.png')
plt.savefig('./complex_eigenvalues.eps')
plt.close()

midpoint = nbins // 2
Expand Down Expand Up @@ -130,43 +133,51 @@
plt.bar(bin_centers, real_zeros,
width=(bins[1]-bins[0]), color='b', label='Real zeros')
plt.legend(loc='best')
ticks = np.arange(-max_plotval, max_plotval+1., 1)
ticks = np.arange(-max_plotval, max_plotval+1., 2)
tickvals = ["${}^{{{}}}$".format(int(np.sign(i+0.01)*10),int(abs(i)+cutoff_exp)) for i in ticks]
ax = plt.gca()
ax.ticklabel_format(style='sci')
ax.set_xlim(left=bins[0], right=bins[-1])
ax.set_xticks(ticks)
ax.set_xticklabels(tickvals)
plt.savefig('./complex_vertint.png')
ax.set_xticklabels(tickvals, fontsize=16)
plt.yticks(fontsize=16)
plt.savefig('./complex_vertint.eps')
plt.close()

plt.plot(bin_centers, real_zeros / all_zeros, color='k')
ticks = np.arange(-max_plotval, max_plotval+1., 1)
ticks = np.arange(-max_plotval, max_plotval+1., 2)
tickvals = ["${}^{{{}}}$".format(int(np.sign(i+0.01)*10),int(abs(i)+cutoff_exp)) for i in ticks]
ax = plt.gca()
ax.ticklabel_format(style='sci')
ax.set_xticks(ticks)
ax.set_xticklabels(tickvals)
plt.savefig('./complex_ratio.png')
ax.set_xticklabels(tickvals, fontsize=16)
plt.yticks(fontsize=16)
plt.savefig('./complex_ratio.eps')
plt.close()

plt.bar(bin_centers, all_zeros,
width=(bins[1]-bins[0]), color='r', label='All zeros')
plt.bar(bin_centers, real_zeros + near_real_zeros,
width=(bins[1]-bins[0]), color='b', label='Real zeros')
plt.legend(loc='best')
ticks = np.arange(-max_plotval, max_plotval+1., 1)
ticks = np.arange(-max_plotval, max_plotval+1., 2)
tickvals = ["${}^{{{}}}$".format(int(np.sign(i+0.01)*10),int(abs(i)+cutoff_exp)) for i in ticks]
ax = plt.gca()
ax.ticklabel_format(style='sci')
ax.set_xlim(left=bins[0], right=bins[-1])
ax.set_xticks(ticks)
ax.set_xticklabels(tickvals)
plt.savefig('./complex_vertint_cutoff.png')
ax.set_xticklabels(tickvals, fontsize=16)
plt.yticks(fontsize=16)
plt.savefig('./complex_vertint_cutoff.eps')
plt.close()

plt.plot(bin_centers, (real_zeros + near_real_zeros) / all_zeros, color='k')
ticks = np.arange(-max_plotval, max_plotval+1., 1)
ticks = np.arange(-max_plotval, max_plotval+1., 2)
tickvals = ["${}^{{{}}}$".format(int(np.sign(i+0.01)*10),int(abs(i)+cutoff_exp)) for i in ticks]
ax = plt.gca()
ax.ticklabel_format(style='sci')
ax.set_xticks(ticks)
ax.set_xticklabels(tickvals)
plt.savefig('./complex_ratio_cutoff.png')
ax.set_xticklabels(tickvals, fontsize=16)
plt.yticks(fontsize=16)
plt.savefig('./complex_ratio_cutoff.eps')
plt.close()
24 changes: 24 additions & 0 deletions count_clusters.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,24 @@
#!/usr/bin/env python

import numpy as np
import scipy.linalg as la
import scipy.stats as stats
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import netCDF4 as nc4

CLUSTER_FILE_NAME = "/home/santos/Data/MG2_data_collection.10_cluster_labels.0001-01-06-00000.nc"

cfile = nc4.Dataset(CLUSTER_FILE_NAME, 'r')

ncluster = len(cfile.dimensions['ncluster'])

labels = cfile.variables["label"][:].flatten()

for i in range(ncluster):
counter = 0
for label in labels:
if label == i:
counter += 1
print(counter, " grid points are in cluster ", i)
Loading

0 comments on commit 18e8169

Please sign in to comment.