"""Figures for reusable divertor profile, width, total, and load-map results."""
from __future__ import annotations
import matplotlib.pyplot as plt
import numpy as np
from .eich import (DivertorEichProfile, DivertorLambdaQPoint, DivertorLoadMap,
DivertorTargetTotalsPoint, DivertorWidthScalings)
[docs]
def plot_divertor_eich_profiles(
profiles: list[DivertorEichProfile],
*,
show_particle: bool = True,
show_energy: bool = True,
show_ions: bool = True,
show_electrons: bool = True,
show_total: bool = True,
xlim: tuple[float, float] | None = None,
ylim: tuple[float, float] | None = None,
):
"""
Plot divertor particle/energy profiles with optional Eich fit overlays.
Parameters
----------
profiles : list[DivertorEichProfile]
Profiles returned by :func:`compute_divertor_eich_profiles`.
show_particle, show_energy : bool
Select which physical channels are shown. At least one must be enabled.
show_ions, show_electrons, show_total : bool
Select species/component curves. Eich fits are only overlaid for total
profiles.
xlim, ylim : tuple[float, float] | None
Optional display limits. ``xlim`` is in millimeters. ``ylim`` applies to
whichever axes are displayed.
Returns
-------
matplotlib.figure.Figure
Figure owned by the caller. GUI code can save it to PNG; scripts can
display or further customize it.
"""
if not show_particle and not show_energy:
raise ValueError("At least one of particle or energy must be enabled.")
nrows = int(show_particle) + int(show_energy)
fig, axes = plt.subplots(nrows=nrows, ncols=1, figsize=(7.4, 3.6*nrows), dpi=150, squeeze=False)
axes_flat = list(axes[:, 0])
channel_axes = []
if show_particle:
channel_axes.append(("particle", axes_flat[len(channel_axes)]))
if show_energy:
channel_axes.append(("energy", axes_flat[len(channel_axes)]))
colors = {"Outer": "#0f5ea8", "Inner": "#a04400"}
for channel, ax in channel_axes:
for profile in profiles:
_plot_channel(
ax,
profile,
channel,
show_ions=show_ions,
show_electrons=show_electrons,
show_total=show_total,
color_base=colors.get(profile.label, "#333333"),
)
ax.axvline(0.0, color="black", linewidth=0.8, alpha=0.45)
ax.set_xlabel(r"$\Delta_{sep}$ [mm]")
ax.grid(alpha=0.25, linewidth=0.4)
ax.legend(loc="best", fontsize=7, ncol=2)
if xlim is not None:
ax.set_xlim(xlim)
if ylim is not None:
ax.set_ylim(ylim)
fig.suptitle("Divertor heatdiag profile mapped to midplane radius", fontsize=11)
fig.tight_layout(pad=0.35)
return fig
[docs]
def plot_divertor_lambda_q_timeseries(
points: list[DivertorLambdaQPoint],
*,
scalings: DivertorWidthScalings | None = None,
ylim: tuple[float, float] | None = None,
show_lines: bool = True,
show_symbols: bool = False,
show_eich_uncertainty: bool = False,
show_chang_uncertainty: bool = False,
):
"""
Plot total-energy Eich ``lambda_q`` values versus heatdiag time.
Parameters
----------
points : list[DivertorLambdaQPoint]
Points returned by :func:`compute_divertor_lambda_q_timeseries`.
scalings : DivertorWidthScalings | None
Optional static upstream-midplane scaling estimates from Simulation.
ylim : tuple[float, float] | None
Optional vertical axis limits in millimeters.
show_lines, show_symbols : bool
Select line, symbol-only, or line-plus-symbol display style.
show_eich_uncertainty, show_chang_uncertainty : bool
Shade the published Eich parameter bounds and Chang fit RMS range.
Returns
-------
matplotlib.figure.Figure
Figure owned by the caller. GUI code can save it to PNG.
"""
if not show_lines and not show_symbols:
raise ValueError("Enable at least one of lines or symbols.")
if not points:
raise ValueError("No converged total-energy Eich fits found in the selected range.")
linestyle = "-" if show_lines else "None"
marker = "o" if show_symbols else None
colors = {"Outer": "#0f5ea8", "Inner": "#a04400"}
fig, ax = plt.subplots(figsize=(7.4, 4.0), dpi=150)
plotted = False
for label in ("Outer", "Inner"):
selected = [point for point in points if point.label == label]
if not selected:
continue
order = np.argsort([point.time_s for point in selected])
time_ms = np.asarray([selected[i].time_s*1.0e3 for i in order], dtype=float)
lambda_q_mm = np.asarray([selected[i].lambda_q_mm for i in order], dtype=float)
ax.plot(
time_ms,
lambda_q_mm,
linestyle=linestyle,
marker=marker,
markersize=4.0,
linewidth=1.5,
color=colors.get(label, "#333333"),
label=label,
)
plotted = True
if not plotted:
raise ValueError("No outer or inner target lambda_q points are available.")
if scalings is not None:
if show_eich_uncertainty:
bpol = scalings.bpol_midplane_t
eich_powers = (bpol**(-1.19 - 0.08), bpol**(-1.19 + 0.08))
ax.axhspan(
(0.63 - 0.08)*min(eich_powers),
(0.63 + 0.08)*max(eich_powers),
color="#246b34", alpha=0.16, linewidth=0,
label="Eich #14 parameter bounds",
)
if show_chang_uncertainty:
ax.axhspan(
scalings.chang_mm*(1.0 - 0.187),
scalings.chang_mm*(1.0 + 0.187),
color="#783c9e", alpha=0.16, linewidth=0,
label="Chang Eq. 6 RMS (18.7%)",
)
ax.axhline(scalings.eich_mm, color="#246b34", linestyle="--", linewidth=1.4, label="Eich #14")
ax.axhline(scalings.chang_mm, color="#783c9e", linestyle=":", linewidth=1.7,
label=r"Chang Eq. 6 (initial $T_i$)")
ax.set_title(r"Divertor heat-load width")
ax.set_xlabel("time [ms]")
ax.set_ylabel(r"$\lambda_q$ [mm]")
if ylim is not None:
ax.set_ylim(ylim)
ax.grid(alpha=0.25, linewidth=0.4)
ax.legend(loc="best", fontsize=8)
fig.tight_layout(pad=0.35)
return fig
[docs]
def plot_divertor_target_totals_timeseries(
points: list[DivertorTargetTotalsPoint],
*,
show_particle: bool = True,
show_power: bool = True,
show_outer: bool = True,
show_inner: bool = True,
show_total: bool = True,
show_control: bool = True,
show_ions: bool = True,
show_electrons: bool = True,
show_target_total: bool = True,
xlim_ms: tuple[float, float] | None = None,
particle_ylim: tuple[float, float] | None = None,
power_ylim_mw: tuple[float, float] | None = None,
):
"""
Plot integrated target particle flux and power versus heatdiag time.
Parameters
----------
points : list[DivertorTargetTotalsPoint]
Points returned by :func:`compute_divertor_target_totals_timeseries`.
show_particle, show_power : bool
Select particle-flux and power panels.
show_outer, show_inner, show_total, show_control : bool
Select target/control trace groups. ``Total`` is the inner-plus-outer
divertor sum, while ``Control`` is the full-wall control quantity after
excluding the heatdiag garbage bin.
show_ions, show_electrons, show_target_total : bool
Select ion, electron, and ion-plus-electron component lines for the
inner and outer target groups.
xlim_ms : tuple[float, float] | None
Optional horizontal axis range in milliseconds.
particle_ylim : tuple[float, float] | None
Optional particle-flux vertical range in particles per second.
power_ylim_mw : tuple[float, float] | None
Optional power vertical range in megawatts.
Returns
-------
matplotlib.figure.Figure
Figure owned by the caller. GUI code can save it to PNG.
"""
if not show_particle and not show_power:
raise ValueError("At least one of particle flux or power must be enabled.")
if not show_outer and not show_inner and not show_total and not show_control:
raise ValueError("At least one target/control trace must be enabled.")
if not points:
raise ValueError("No divertor target total points found in the selected range.")
nrows = int(show_particle) + int(show_power)
fig, axes = plt.subplots(nrows=nrows, ncols=1, figsize=(7.4, 3.4*nrows), dpi=150, squeeze=False)
axes_flat = list(axes[:, 0])
panel_axes = []
if show_particle:
panel_axes.append(("particle", axes_flat[len(panel_axes)]))
if show_power:
panel_axes.append(("power", axes_flat[len(panel_axes)]))
trace_specs = _target_totals_trace_specs(
show_outer=show_outer,
show_inner=show_inner,
show_total=show_total,
show_control=show_control,
show_ions=show_ions,
show_electrons=show_electrons,
show_target_total=show_target_total,
)
if not trace_specs:
raise ValueError("Enable at least one target/control or ion/electron/target-total line.")
for panel, ax in panel_axes:
plotted = False
for spec in trace_specs:
selected = [point for point in points if point.label == spec["label"]]
if not selected:
continue
order = np.argsort([point.time_s for point in selected])
time_ms = np.asarray([selected[i].time_s*1.0e3 for i in order], dtype=float)
if panel == "particle":
y = np.asarray([_target_totals_particle_value(selected[i], spec["component"]) for i in order], dtype=float)
else:
y = np.asarray([_target_totals_power_value(selected[i], spec["component"])/1.0e6 for i in order], dtype=float)
ax.plot(
time_ms,
y,
linestyle=spec["linestyle"],
linewidth=1.5,
color=spec["color"],
label=spec["plot_label"],
)
plotted = True
if not plotted:
raise ValueError("No enabled target/control traces are available.")
ax.set_xlabel("time [ms]")
if panel == "particle":
ax.set_ylabel("Particle flux [ptl/s]")
if particle_ylim is not None:
ax.set_ylim(particle_ylim)
else:
ax.set_ylabel("Power [MW]")
if power_ylim_mw is not None:
ax.set_ylim(power_ylim_mw)
if xlim_ms is not None:
ax.set_xlim(xlim_ms)
ax.grid(alpha=0.25, linewidth=0.4)
ax.legend(loc="best", fontsize=8)
fig.suptitle("Divertor heatdiag target totals", fontsize=11)
fig.tight_layout(pad=0.35)
return fig
[docs]
def plot_divertor_load_map(
load_map: DivertorLoadMap,
*,
xlim: tuple[float, float] | None = None,
ylim: tuple[float, float] | None = None,
clim: tuple[float, float] | None = None,
contour_count: int = 256,
cmap: str = "seismic",
):
"""
Plot a mapped divertor load map as filled contours without contour lines.
Parameters
----------
load_map : DivertorLoadMap
Map returned by :func:`compute_divertor_load_map`.
xlim, ylim : tuple[float, float] | None
Optional display ranges. ``xlim`` is in millimeters. ``ylim`` is in
milliseconds for the time view and degrees for the toroidal view.
clim : tuple[float, float] | None
Optional color scale in displayed units. Energy maps are displayed in
MW m^-2; particle maps are displayed in ptl m^-2 s^-1.
contour_count : int
Number of filled contour levels. Values below two are promoted to two.
cmap : str
Matplotlib colormap name. The default is ``"seismic"``.
Returns
-------
matplotlib.figure.Figure
Figure owned by the caller. GUI code can save it to PNG.
"""
contour_count = max(2, int(contour_count))
x = np.asarray(load_map.x_mm, dtype=float)
y = np.asarray(load_map.y_values, dtype=float)
z = np.asarray(load_map.values, dtype=float)
if z.ndim != 2:
raise ValueError(f"Expected a 2D load map, got shape={z.shape}")
if load_map.channel == "energy":
z_plot = z/1.0e6
cbar_label = r"Energy flux density [MW m$^{-2}$]"
else:
z_plot = z
cbar_label = r"Particle flux density [ptl m$^{-2}$ s$^{-1}$]"
if load_map.view == "time":
y_plot = y*1.0e3
y_label = "time [ms]"
else:
y_plot = y
y_label = "toroidal angle [deg]"
finite = z_plot[np.isfinite(z_plot)]
if finite.size == 0:
raise ValueError("Load map contains no finite values.")
if clim is None:
zmin = float(np.nanmin(finite))
zmax = float(np.nanmax(finite))
else:
zmin, zmax = clim
if not np.isfinite(zmin) or not np.isfinite(zmax):
raise ValueError("Color scale limits must be finite.")
if zmin == zmax:
pad = max(abs(zmin)*1.0e-6, 1.0e-12)
zmin -= pad
zmax += pad
if zmin > zmax:
zmin, zmax = zmax, zmin
levels = np.linspace(zmin, zmax, contour_count)
fig, ax = plt.subplots(figsize=(7.4, 4.4), dpi=150)
contour = ax.contourf(x, y_plot, z_plot, levels=levels, cmap=cmap, extend="both")
ax.set_xlabel(r"$\Delta_{sep}$ [mm]")
ax.set_ylabel(y_label)
if xlim is not None:
ax.set_xlim(xlim)
if ylim is not None:
ax.set_ylim(ylim)
title_channel = "energy" if load_map.channel == "energy" else "particle"
ax.set_title(f"{load_map.target_label} divertor {load_map.component} {title_channel} load")
cbar = fig.colorbar(contour, ax=ax)
cbar.set_label(cbar_label)
fig.tight_layout(pad=0.35)
return fig
def _target_totals_trace_specs(
*,
show_outer: bool,
show_inner: bool,
show_total: bool,
show_control: bool,
show_ions: bool,
show_electrons: bool,
show_target_total: bool,
) -> list[dict[str, str]]:
"""Return line specifications for target-total plots."""
specs: list[dict[str, str]] = []
def add_target(label: str, color: str) -> None:
"""Append enabled component traces for one divertor target."""
if show_ions:
specs.append({"label": label, "component": "ion", "plot_label": f"{label} ions", "color": color, "linestyle": "--"})
if show_electrons:
specs.append({"label": label, "component": "electron", "plot_label": f"{label} electrons", "color": color, "linestyle": ":"})
if show_target_total:
specs.append({"label": label, "component": "total", "plot_label": f"{label} total", "color": color, "linestyle": "-"})
if show_outer:
add_target("Outer", "#0f5ea8")
if show_inner:
add_target("Inner", "#a04400")
if show_total:
specs.append({"label": "Total", "component": "total", "plot_label": "Total divertor", "color": "#20845a", "linestyle": "-"})
if show_control:
specs.append({"label": "Control", "component": "total", "plot_label": "Wall control", "color": "#333333", "linestyle": "--"})
return specs
def _target_totals_particle_value(point: DivertorTargetTotalsPoint, component: str) -> float:
"""Return one particle-flux component from an integrated totals point."""
if component == "ion":
return point.particle_i_rate_per_s
if component == "electron":
return point.particle_e_rate_per_s
return point.particle_rate_per_s
def _target_totals_power_value(point: DivertorTargetTotalsPoint, component: str) -> float:
"""Return one power component from an integrated totals point."""
if component == "ion":
return point.power_i_w
if component == "electron":
return point.power_e_w
return point.power_w
def _plot_channel(
ax,
profile: DivertorEichProfile,
channel: str,
*,
show_ions: bool,
show_electrons: bool,
show_total: bool,
color_base: str,
) -> None:
"""Draw one particle or energy channel for a target profile.
Species curves are drawn with dashed/dotted styles, the total is drawn as a
solid curve, and the Eich fit overlay plus fitted ``lambda``/``S`` text are
shown only when ``show_total`` is enabled and a fit is available.
"""
if channel == "particle":
yi = profile.particle_i
ye = profile.particle_e
fit = profile.fit_particle
ylabel = r"Particle flux density [m$^{-2}$ s$^{-1}$]"
else:
yi = profile.energy_i/1.0e6
ye = profile.energy_e/1.0e6
fit = profile.fit_energy
ylabel = r"Energy flux density [MW m$^{-2}$]"
if show_ions:
ax.plot(profile.x_mm, yi, linestyle="--", linewidth=1.0, label=f"{profile.label} ions")
if show_electrons:
ax.plot(profile.x_mm, ye, linestyle=":", linewidth=1.1, label=f"{profile.label} electrons")
if show_total:
total = yi + ye
ax.plot(profile.x_mm, total, linewidth=1.6, label=f"{profile.label} total")
if fit is not None:
params, fit_xy = fit
y_fit = fit_xy[:, 1]/1.0e6 if channel == "energy" else fit_xy[:, 1]
ax.plot(fit_xy[:, 0], y_fit, color=color_base, linewidth=1.2, alpha=0.85)
ax.text(
0.02,
0.92 if profile.label == "Outer" else 0.82,
f"{profile.label}: lambda={params[2]*1.0e3:.2f} mm, S={params[3]*1.0e3:.2f} mm",
transform=ax.transAxes,
fontsize=8,
color=color_base,
)
ax.set_ylabel(ylabel)