From 9a07326d47a0baa8b607111506147523d7943054 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Fri, 23 Feb 2024 14:11:23 -0500 Subject: [PATCH] some updates to the massive star slice plots (#2759) now we offset t based on the driving time also add a time-sequence plot --- .../analysis/massive_star_multi.py | 18 +- .../analysis/massive_star_sequence.py | 178 ++++++++++++++++++ 2 files changed, 191 insertions(+), 5 deletions(-) create mode 100755 Exec/science/massive_star/analysis/massive_star_sequence.py diff --git a/Exec/science/massive_star/analysis/massive_star_multi.py b/Exec/science/massive_star/analysis/massive_star_multi.py index 29f6bf3829..4722093e33 100755 --- a/Exec/science/massive_star/analysis/massive_star_multi.py +++ b/Exec/science/massive_star/analysis/massive_star_multi.py @@ -22,6 +22,13 @@ def make_plot(plotfile, prefix="plot", width_frac=0.1, ds = CastroDataset(plotfile) + t_drive = 0.0 + if "[*] castro.drive_initial_convection_tmax" in ds.parameters: + t_drive = ds.parameters["[*] castro.drive_initial_convection_tmax"] + elif "castro.drive_initial_convection_tmax" in ds.parameters: + t_drive = ds.parameters["castro.drive_initial_convection_tmax"] + print(t_drive) + xmin = ds.domain_left_edge[0] xmax = ds.domain_right_edge[0] @@ -70,8 +77,9 @@ def make_plot(plotfile, prefix="plot", width_frac=0.1, if f == "enuc": # now do a contour of NSE - sp.annotate_contour("in_nse", ncont=1, clim=(0.5, 0.5), take_log=False, - plot_args={"colors": "k", "linewidths": 2}) + if ("boxlib", "in_nse") in ds.derived_field_list: + sp.annotate_contour("in_nse", ncont=1, clim=(0.5, 0.5), take_log=False, + plot_args={"colors": "k", "linewidths": 2}) sp.set_axes_unit("cm") @@ -79,7 +87,7 @@ def make_plot(plotfile, prefix="plot", width_frac=0.1, cbar_location=cbar_location, cbar_pad="2%") fig.subplots_adjust(left=0.05, right=0.95, bottom=0.025, top=0.975) - fig.text(0.02, 0.02, f"time = {float(ds.current_time):8.2f} s", + fig.text(0.02, 0.02, f"$t - \\tau_\\mathrm{{drive}}$ = {float(ds.current_time) - t_drive:6.1f} s", transform=fig.transFigure) fig.set_size_inches(size) fig.tight_layout() @@ -87,7 +95,7 @@ def make_plot(plotfile, prefix="plot", width_frac=0.1, if layout[0] >= layout[1]: extra = "_vertical" - fig.savefig(f"{prefix}_{os.path.basename(plotfile)}_w{width_frac:04.2f}{extra}.pdf", pad_inches=0) + fig.savefig(f"{prefix}_{os.path.basename(plotfile)}_w{width_frac:04.2f}{extra}.pdf", pad_inches=0.1, bbox_inches="tight") if __name__ == "__main__": @@ -95,7 +103,7 @@ def make_plot(plotfile, prefix="plot", width_frac=0.1, p = argparse.ArgumentParser() p.add_argument("--vertical", action="store_true", help="plot 2x2 or 1x4") - p.add_argument("--width_fraction", type=float, default=0.1, + p.add_argument("--width-fraction", type=float, default=0.1, help="fraction of domain to show") p.add_argument("plotfile", type=str, nargs=1, help="plotfile to plot") diff --git a/Exec/science/massive_star/analysis/massive_star_sequence.py b/Exec/science/massive_star/analysis/massive_star_sequence.py new file mode 100755 index 0000000000..96b4218085 --- /dev/null +++ b/Exec/science/massive_star/analysis/massive_star_sequence.py @@ -0,0 +1,178 @@ +#!/usr/bin/env python3 + +import argparse +import os + +import matplotlib +import matplotlib.pyplot as plt +import numpy as np +from mpl_toolkits.axes_grid1 import ImageGrid + +import yt +from yt.frontends.boxlib.api import CastroDataset +# assume that our data is in CGS +from yt.units import cm + +matplotlib.use('agg') + +times = [50, 100, 150, 200, 250] + +def find_files(plist): + + mask = np.zeros(len(times)) + files_to_plot = [] + for pfile in plist: + for k, t in enumerate(times): + if mask[k]: + continue + print(pfile) + ds = CastroDataset(pfile) + if ds.current_time >= t: + files_to_plot.append(pfile) + mask[k] = 1.0 + + return files_to_plot + + +def make_plot(field, pfiles, *, + width_frac=0.1, + size=(19.2, 10.8)): + + fig = plt.figure() + + if len(pfiles) > 8: + nrows = 3 + ncols = (len(pfiles) + 1)//3 + elif len(pfiles) >= 4: + nrows = 2 + ncols = (len(pfiles) + 1)//2 + else: + nrows = 1 + ncols = len(pfiles) + + grid = ImageGrid(fig, 111, nrows_ncols=(nrows, ncols), + axes_pad=0.2, cbar_pad=0.05, + label_mode="L", cbar_mode="single") + + # get some metadata + ds = CastroDataset(pfiles[0]) + + t_drive = 0.0 + if "[*] castro.drive_initial_convection_tmax" in ds.parameters: + t_drive = ds.parameters["[*] castro.drive_initial_convection_tmax"] + elif "castro.drive_initial_convection_tmax" in ds.parameters: + t_drive = ds.parameters["castro.drive_initial_convection_tmax"] + print(t_drive) + + xmin = ds.domain_left_edge[0] + xmax = ds.domain_right_edge[0] + + ymin = ds.domain_left_edge[1] + ymax = ds.domain_right_edge[1] + + xctr = 0.0 * xmin + L_x = xmax - xmin + + yctr = 0.5 * (ymin + ymax) + L_y = ymax - ymin + + f = field + + if f == "Temp" or f == "abar": + text_color = "white" + else: + text_color = "black" + + for i in range(nrows * ncols): + + if i < len(pfiles): + pf = pfiles[i] + else: + grid[i].remove() + continue + + ds = CastroDataset(pf) + + sp = yt.SlicePlot(ds, "theta", f, + center=[xmin + 0.5*width_frac*L_x, yctr, 0.0*cm], + width=[width_frac*L_x, width_frac*L_y, 0.0*cm], + fontsize="12") + + sp.set_buff_size((2400,2400)) + + sp.annotate_text((0.05, 0.05), + f"$t - \\tau_\\mathrm{{drive}}$ = {float(ds.current_time) - t_drive:8.2f} s", + coord_system="axis", + text_args={"color": text_color}) + + + if f == "Ye": + sp.set_zlim(f, 0.46, 0.5) + sp.set_log(f, False) + sp.set_cmap(f, "magma_r") + elif f == "abar": + sp.set_log(f, False) + sp.set_cmap(f, "viridis") + elif f == "enuc": + sp.set_log(f, True, linthresh=1.e12) + sp.set_zlim(f, -1.e20, 1.e20) + sp.set_cmap(f, "bwr") + elif f == "MachNumber": + sp.set_zlim(f, 1.e-4, 0.3) + sp.set_cmap(f, "plasma") + elif f == "magvel": + sp.set_zlim(f, 100.0, 2.e7) + sp.set_cmap(f, "viridis") + elif f == "magvort": + sp.set_cmap(f, "magma") + sp.set_zlim(f, 1.e-2, 5) + + if f == "enuc": + # now do a contour of NSE + if ("boxlib", "in_nse") in ds.derived_field_list: + sp.annotate_contour("in_nse", levels=1, clim=(0.5, 0.5), take_log=False, + plot_args={"colors": "k", "linewidths": 2}) + + sp.set_axes_unit("cm") + + plot = sp.plots[field] + plot.figure = fig + plot.axes = grid[i].axes + plot.cax = grid.cbar_axes[i] + if i < len(pfiles)-1: + grid[i].axes.xaxis.offsetText.set_visible(False) + if i > 0: + grid[i].axes.yaxis.offsetText.set_visible(False) + + sp._setup_plots() + + fig.set_size_inches(size) + fig.tight_layout() + extra = "" + + fig.savefig(f"massive_star_{f}_w{width_frac:04.2f}.pdf", pad_inches=0.1, bbox_inches="tight") + + +if __name__ == "__main__": + + p = argparse.ArgumentParser() + p.add_argument("--var", type=str, default="Temp", + help="variable to plot") + p.add_argument("--vertical", action="store_true", + help="plot 2x2 or 1x4") + p.add_argument("--width-fraction", type=float, default=0.1, + help="fraction of domain to show") + p.add_argument("plotfiles", type=str, nargs="+", + help="list of plotfiles to plot") + + args = p.parse_args() + + if args.vertical: + size = (7.5, 8.5) + else: + size = (19.2, 8.5) + + plist = find_files(args.plotfiles) + + make_plot(args.var, plist, + width_frac=args.width_fraction, size=size)