diff --git a/src/py4vasp/_calculation/electronic_minimization.py b/src/py4vasp/_calculation/electronic_minimization.py index ba79ff8b..88cc6bca 100644 --- a/src/py4vasp/_calculation/electronic_minimization.py +++ b/src/py4vasp/_calculation/electronic_minimization.py @@ -19,7 +19,37 @@ ) from py4vasp._raw.data_db import ElectronicMinimization_DB from py4vasp._third_party import graph -from py4vasp._util import check +from py4vasp._util import check, index, select + +# VASP labels that clash with the selection grammar (parentheses mean nesting) are +# exposed under a grammar-safe alias. The raw label is kept for display / printing. +_LABEL_TOKEN_OVERRIDES = {"rms(c)": "rms_c"} + +_SELECTION_ERROR_MESSAGE = """\ +Please choose a selection including at least one of the following keywords: +N, E, dE, deps, ncg, rms, rms_c""" + +# Energy-change series shown on the left axis of the convergence overview, beside the +# "E" distance to the converged energy. Given as (column token, display label) pairs. +_ENERGY_CHANGE_TOKENS = [("dE", "dE"), ("deps", "d eps")] +# Columns not plotted on the secondary residual axis: the iteration index, the energies, +# and the number of Hamiltonian evaluations. Every remaining column (rms and the density +# / orthonormalization residual, which VASP labels rms(c) or ort) is a residual. +_NON_RESIDUAL_TOKENS = {"N", "E", "dE", "deps", "ncg"} + +# Convergence entries whose magnitude is below this threshold are treated as +# not-yet-computed and reported as NaN. This most notably affects rms(c), which VASP +# reports as zero until density updates begin after the NELMDL delay. +_SANITY_THRESHOLD = 1e-16 + +# Overlay flagging the points whose sign is atypical for their series (a log axis only +# shows the magnitude, so the sign is otherwise invisible). Drawn as small crosses in a +# neutral colour. "Atypical" means the minority sign, so a cleanly converging quantity +# (e.g. an always-negative dE) is not flagged at all. +_UNUSUAL_SIGN_LABEL = "unusual sign" +_UNUSUAL_MARKER = "x" +_UNUSUAL_MARKER_SIZE = 8 +_UNUSUAL_COLOR = "#4d4d4d" _TO_DATABASE_SUPPRESSED_EXCEPTIONS = ( exception.Py4VaspError, @@ -48,7 +78,8 @@ def __str__(self) -> str: label_rep = "{}\t\t{}\t\t{}\t\t{}\t\t{}\t{}\t\t{}\n" string = "" labels = [label.decode("utf-8") for label in getattr(self._raw_data, "label")] - data = self.to_dict() + # print the raw OSZICAR data verbatim, without the read() sanity filter + data = self._extract(selection=None, sanitize=False) electronic_iterations = data["N"] if not self._more_than_one_ionic_step(electronic_iterations): electronic_iterations = [electronic_iterations] @@ -59,7 +90,10 @@ def __str__(self) -> str: for electronic_step in range(electronic_steps): _data = [] for label in self._raw_data.label: - _values_electronic = data[label.decode("utf-8")] + token = _LABEL_TOKEN_OVERRIDES.get( + label.decode("utf-8"), label.decode("utf-8") + ) + _values_electronic = data[token] if not self._more_than_one_ionic_step(_values_electronic): _values_electronic = [_values_electronic] _value = _values_electronic[ionic_step][electronic_step] @@ -69,35 +103,152 @@ def __str__(self) -> str: return string def to_dict(self, selection=None) -> dict: - """Extract convergence data and return as a dict.""" - return_data = {} + """Extract convergence data and return as a dict. + + The selection is parsed with the standard :mod:`py4vasp._util.select` grammar and + evaluated with :class:`py4vasp._util.index.Selector`, so multiple columns + (``"E, dE"``) and compositions (``"E + dE"``) are supported. Entries whose + magnitude is below a small threshold are treated as not-yet-computed and reported + as NaN (see :data:`_SANITY_THRESHOLD`). + """ + return self._extract(selection, sanitize=True) + + def _extract(self, selection, sanitize) -> dict: + tokens = self._tokens() + step_arrays = self._step_arrays() + column_map = {1: {token: column for column, token in enumerate(tokens)}} + selectors = [index.Selector(column_map, array) for array in step_arrays] + is_none = np.all([array.is_none() for array in step_arrays]) if selection is None: - keys_to_include = self._from_bytes_to_utf(self._raw_data.label) + selections = [(token,) for token in tokens] else: - labels_as_str = self._from_bytes_to_utf(self._raw_data.label) - if selection not in labels_as_str: - message = """\ -Please choose a selection including at least one of the following keywords: -N, E, dE, deps, ncg, rms, rms(c)""" - raise exception.RefinementError(message) - keys_to_include = [selection] - for key in keys_to_include: - return_data[key] = self._read(key) + selections = list(select.Tree.from_selection(selection).selections()) + return_data = {} + for sel in selections: + try: + key = selectors[0].label(sel) + columns = [selector[sel] for selector in selectors] + except exception.IncorrectUsage as error: + raise exception.RefinementError(_SELECTION_ERROR_MESSAGE) from error + if sanitize: + columns = [self._sanitize(column) for column in columns] + values = [list(column) for column in columns] + if len(values) == 1: + values = values[0] + return_data[key] = values if not is_none else {} return return_data - def to_graph(self, selection="E") -> graph.Graph: - """Graph the change in parameter with iteration number.""" + def _sanitize(self, column): + """Report not-yet-computed (near-zero) entries as NaN.""" + column = np.array(column, dtype=float) + column[np.abs(column) < _SANITY_THRESHOLD] = np.nan + return column + + def to_graph(self, selection=None) -> graph.Graph: + """Graph the convergence data against the iteration number. + + Without a selection this produces a convergence overview: the energy changes + on a logarithmic left axis and the residuals on a logarithmic secondary axis. + With a selection, the chosen columns are plotted as-is on a logarithmic axis. + """ + if selection is None: + return self._overview_graph() + return self._selection_graph(selection) + + def _overview_graph(self) -> graph.Graph: data = self.to_dict() - series = graph.Series(data["N"], data[selection], selection) - from py4vasp._util import select as sel_util + iterations = np.array(data["N"], dtype=float) + # (signed values, label, secondary axis) for each series of the overview; the + # distance to the converged energy is positive while E is above its final value + specs = [(self._energy_distance(data), "E - E_final", False)] + specs += [ + (np.array(data[token], dtype=float), label, False) + for token, label in _ENERGY_CHANGE_TOKENS + ] + specs += [ + (np.array(data[token], dtype=float), token, True) + for token in self._tokens() + if token not in _NON_RESIDUAL_TOKENS + ] + series = [self._make_series(iterations, *spec) for spec in specs] + series += self._unusual_sign_overlays(iterations, specs) + return graph.Graph( + series=series, + xlabel="Iteration number", + ylabel="Energy change (eV)", + y2label="Residual", + yscale="log", + y2scale="log", + ) - ylabel = " ".join(s.capitalize() for s in selection.split("_")) + def _unusual_sign_overlays(self, iterations, specs): + """One markers-only series per axis flagging the atypical-sign points (whose + magnitude is plotted but whose sign the log axis hides); the overlays on the two + axes share a single legend entry.""" + overlays = [] + for y2 in (False, True): + points = [] + for values, _, axis in specs: + if axis != y2: + continue + mask = self._unusual_sign_mask(values) + if mask.any(): + points.append((iterations[mask], np.abs(values[mask]))) + if not points: + continue + xs, ys = zip(*points) + overlays.append( + graph.Series( + np.concatenate(xs), + np.concatenate(ys), + _UNUSUAL_SIGN_LABEL, + marker=graph.Marker(_UNUSUAL_MARKER, _UNUSUAL_MARKER_SIZE), + color=_UNUSUAL_COLOR, + y2=y2, + show_legend=not overlays, + ) + ) + return overlays + + @staticmethod + def _unusual_sign_mask(values): + """Boolean mask of the points whose sign is the minority for this series; empty + when the finite values do not change sign.""" + num_positive = np.sum(values > 0) + num_negative = np.sum(values < 0) + if num_positive == 0 or num_negative == 0: + return np.zeros(len(values), dtype=bool) + return values > 0 if num_positive < num_negative else values < 0 + + def _energy_distance(self, data): + # E - E_final is positive while the energy is above its converged value; the + # final point is zero and therefore dropped so it does not vanish on the axis + values = np.array(data["E"], dtype=float) + values = values - values[-1] + values[-1] = np.nan + return values + + def _selection_graph(self, selection) -> graph.Graph: + iterations = np.array(self.to_dict("N")["N"], dtype=float) + series = [ + self._make_series(iterations, np.array(values, dtype=float), label) + for label, values in self.to_dict(selection).items() + ] return graph.Graph( - series=[series], + series=series, xlabel="Iteration number", - ylabel=ylabel, + ylabel="Convergence data", + yscale="log", ) + def _make_series(self, x, values, label, y2=False) -> graph.Series: + """Plot |values| labelled "|label|" if any value is negative (a log axis cannot + show negative numbers); otherwise plot the values as they are.""" + if np.any(values < 0): + values = np.abs(values) + label = f"|{label}|" + return graph.Series(x, values, label, y2=y2) + def is_converged(self) -> np.ndarray: is_elmin_converged = self._raw_data.is_elmin_converged[self._steps_or_last] converged = is_elmin_converged == 0 @@ -144,22 +295,20 @@ def _steps_or_last(self): def _more_than_one_ionic_step(self, data): return any(isinstance(_data, list) for _data in data) - def _read(self, key): + def _tokens(self): + """Selectable tokens: the raw labels with grammar-clashing ones aliased.""" + labels = self._from_bytes_to_utf(self._raw_data.label) + return [_LABEL_TOKEN_OVERRIDES.get(label, label) for label in labels] + + def _step_arrays(self): + """Split the convergence data into per-ionic-step arrays for the chosen steps.""" data = getattr(self._raw_data, "convergence_data") iteration_number = data[:, 0] split_index = np.where(iteration_number == 1)[0] data = np.vsplit(data, split_index)[1:][self._steps_or_last] if isinstance(self._steps, slice): - data = [raw.VaspData(_data) for _data in data] - else: - data = [raw.VaspData(data)] - labels = [label.decode("utf-8") for label in self._raw_data.label] - data_index = labels.index(key) - return_data = [list(_data[:, data_index]) for _data in data] - is_none = [_data.is_none() for _data in data] - if len(return_data) == 1: - return_data = return_data[0] - return return_data if not np.all(is_none) else {} + return [raw.VaspData(_data) for _data in data] + return [raw.VaspData(data)] def _get_electronic_steps_info(self) -> tuple: if check.is_none(self._raw_data.convergence_data): @@ -251,13 +400,17 @@ def to_dict(self, selection=None) -> dict: """Convenient alias for :py:meth:`read`. Please read the documentation there.""" return self.read(selection=selection) - def to_graph(self, selection="E") -> graph.Graph: - """Graph the change in parameter with iteration number. + def to_graph(self, selection=None) -> graph.Graph: + """Graph the convergence data against the iteration number. Parameters ---------- selection: str - Choose strings consistent with the OSZICAR format + Choose strings consistent with the OSZICAR format (N, E, dE, deps, ncg, rms, + rms_c), optionally composed with the standard selection grammar (e.g. + "E + dE"). Without a selection a convergence overview is produced with the + energy changes on a logarithmic left axis and the residuals on a logarithmic + secondary axis. Returns ------- diff --git a/src/py4vasp/_demo/electronic_minimization.py b/src/py4vasp/_demo/electronic_minimization.py index aafba339..097e0e28 100644 --- a/src/py4vasp/_demo/electronic_minimization.py +++ b/src/py4vasp/_demo/electronic_minimization.py @@ -7,9 +7,15 @@ def electronic_minimization(): random_convergence_data = np.random.rand(9, 3) + # the total energy converges downward toward its final value, so the distance + # E - E_final stays positive throughout the minimization + random_convergence_data[:, 0] = -8.0 - np.linspace(0.0, 1.0, 9) iteration_number = np.arange(1, 10)[:, np.newaxis] ncg = np.random.randint(4, 10, (9, 1)) random_rms = np.random.rand(9, 2) + # density updates (and hence rms(c)) only start after the NELMDL delay, so the + # first few electronic steps report a value of zero for this column + random_rms[:5, 1] = 0.0 convergence_data = np.hstack( [iteration_number, random_convergence_data, ncg, random_rms] ) diff --git a/src/py4vasp/_third_party/graph/__init__.py b/src/py4vasp/_third_party/graph/__init__.py index bd4a218c..1724841e 100644 --- a/src/py4vasp/_third_party/graph/__init__.py +++ b/src/py4vasp/_third_party/graph/__init__.py @@ -9,7 +9,7 @@ from .graph import Graph from .mixin import Mixin from .plot import plot -from .series import Series +from .series import Marker, Series go = import_.optional("plotly.graph_objects") pio = import_.optional("plotly.io") diff --git a/src/py4vasp/_third_party/graph/graph.py b/src/py4vasp/_third_party/graph/graph.py index d8862a91..ff9ea8a3 100644 --- a/src/py4vasp/_third_party/graph/graph.py +++ b/src/py4vasp/_third_party/graph/graph.py @@ -124,14 +124,20 @@ class Graph(Sequence): "Tuple specifying the visible range of the x-axis as (min, max)." xticks: Optional[dict] = None "Dictionary mapping tick positions (keys) to their labels (values) for the x-axis." + xscale: Optional[str] = None + 'Scale of the x-axis, e.g. "log". Defaults to None (plotly\'s automatic linear scale).' xsize: Optional[int] = 720 "Width of the figure in pixels." ylabel: Optional[str] = None "Label for the y-axis. For subplots, provide a list of labels corresponding to each subplot." yrange: Optional[tuple] = None "Tuple specifying the visible range of the y-axis as (min, max)." + yscale: Optional[str] = None + 'Scale of the y-axis, e.g. "log". Defaults to None (plotly\'s automatic linear scale).' y2label: Optional[str] = None "Label for the secondary y-axis (only applicable when series use the secondary y-axis)." + y2scale: Optional[str] = None + 'Scale of the secondary y-axis, e.g. "log". Defaults to None (automatic linear scale).' ysize: Optional[int] = 540 "Height of the figure in pixels." title: Optional[str] = None @@ -316,17 +322,19 @@ def _make_plotly_figure(self): return figure def _figure_with_one_or_two_y_axes(self): - has_secondary_y_axis = lambda series: isinstance(series, Series) and series.y2 if self._subplot_on: max_row = max(series.subplot for series in self) figure = subplots.make_subplots(rows=max_row, cols=1) figure.update_layout(showlegend=False) return figure - elif any(has_secondary_y_axis(series) for series in self): + elif self._any_have_secondary_y_axis(): return subplots.make_subplots(specs=[[{"secondary_y": True}]]) else: return go.Figure() + def _any_have_secondary_y_axis(self): + return any(isinstance(series, Series) and series.y2 for series in self) + def _set_legend(self, figure): default_colorbar_x = 1.02 colorbar_spacing = 0.18 @@ -351,14 +359,24 @@ def _set_legend(self, figure): cb.x = default_colorbar_x + i * colorbar_spacing max_colorbar_x = max(max_colorbar_x, cb.x) - # Step 3: Place legend safely to the right of all colorbars + # Step 3: Place legend safely to the right of all colorbars, or clear of the + # secondary y-axis title (which sits on the right) if there is one. The right + # margin auto-expands to fit the shifted legend. if max_colorbar_x > 1.0: legend_x = max_colorbar_x + colorbar_spacing figure.update_layout( legend=dict(x=legend_x, y=1.0, xanchor="left", yanchor="top") ) - - figure.update_layout(margin=dict(r=120)) + figure.update_layout(margin=dict(r=120)) + elif self._any_have_secondary_y_axis(): + # push the legend clear of the secondary y-axis title on the right and widen + # the margin so the (potentially long) labels are not clipped + figure.update_layout( + legend=dict(x=1.13, y=1.0, xanchor="left", yanchor="top") + ) + figure.update_layout(margin=dict(r=200)) + else: + figure.update_layout(margin=dict(r=120)) def _set_xaxis_options(self, figure): if self._subplot_on: @@ -374,6 +392,8 @@ def _set_xaxis_options(self, figure): figure.layout.xaxis.ticktext = self._xtick_labels() if self.xrange: figure.layout.xaxis.range = self.xrange + if self.xscale: + figure.layout.xaxis.type = self.xscale if self._all_are_contour(): figure.layout.xaxis.visible = False @@ -388,8 +408,12 @@ def _set_yaxis_options(self, figure): figure.update_yaxes(title_text=ylabel, row=row + 1, col=1) else: figure.layout.yaxis.title.text = self.ylabel + if self.yscale: + figure.layout.yaxis.type = self.yscale if self.y2label: figure.layout.yaxis2.title.text = self.y2label + if self.y2scale and self._any_have_secondary_y_axis(): + figure.layout.yaxis2.type = self.y2scale if self.yrange: figure.layout.yaxis.range = self.yrange if self._all_are_contour(): diff --git a/src/py4vasp/_third_party/graph/series.py b/src/py4vasp/_third_party/graph/series.py index 90d912fc..4c5f840b 100644 --- a/src/py4vasp/_third_party/graph/series.py +++ b/src/py4vasp/_third_party/graph/series.py @@ -1,7 +1,7 @@ # Copyright © VASP Software GmbH, # Licensed under the Apache License 2.0 (http://www.apache.org/licenses/LICENSE-2.0) from dataclasses import dataclass, fields -from typing import Generator, Optional, Tuple +from typing import Generator, Optional, Tuple, Union import numpy as np @@ -11,6 +11,37 @@ go = import_.optional("plotly.graph_objects") +# Map py4vasp's friendly marker names onto plotly symbols so that the plotly vocabulary +# does not leak into the public interface. Unknown names are passed through unchanged. +_SYMBOL_ALIASES = { + "o": "circle", + "*": "star", + "s": "square", + "^": "triangle-up", + "v": "triangle-down", + "x": "x", + "d": "diamond", + "+": "cross", +} + + +@dataclass +class Marker: + """Style of the markers drawn for a series. + + Parameters + ---------- + symbol : str + Which marker to draw. Use one of the friendly names "o" (circle), "*" (star), + "s" (square), "^" (triangle-up), "v" (triangle-down), "x" (cross), "d" (diamond), + or "+" (plus). + size : Optional[float] + The size of the marker. If not set, a default size is used. + """ + + symbol: str + size: Optional[float] = None + @dataclass class Series(trace.Trace): @@ -80,9 +111,12 @@ class Series(trace.Trace): color: Optional[str] = None """The color used for this series. Accepts any valid CSS color string, including rgba() for transparency in area plots.""" - marker: Optional[str] = None - """Marker style for scatter plots. If None, displays as a line plot. Common values - include 'circle', 'square', 'diamond', etc.""" + marker: Optional[Union[str, Marker]] = None + """Marker style for scatter plots. If None, displays as a line plot. May be a + :class:`Marker` or, as a shorthand, a symbol string (see :class:`Marker`).""" + show_legend: bool = True + """If False, this series is not given its own entry in the legend. Useful to let + several series share a single legend entry.""" _frozen = False def __post_init__(self): @@ -191,6 +225,13 @@ def _is_line(self): def _is_area(self): return (self.weight is not None) and (self.marker is None) + def _marker_symbol_and_size(self): + marker = ( + Marker(symbol=self.marker) if isinstance(self.marker, str) else self.marker + ) + symbol = _SYMBOL_ALIASES.get(marker.symbol, marker.symbol) + return symbol, marker.size + def _options_line(self, y): return { "x": self.x, @@ -211,19 +252,21 @@ def _options_area(self, y, weight): } def _options_scaled_points(self, y, weight): - return { - "x": self.x, - "y": y, - "mode": "markers", - "marker": {"size": weight, "sizemode": "area", "color": self.color}, - } + symbol, size = self._marker_symbol_and_size() + marker = {"symbol": symbol, "color": self.color} + if weight is not None: + marker.update(size=weight, sizemode="area") + elif size is not None: + marker["size"] = size + return {"x": self.x, "y": y, "mode": "markers", "marker": marker} def _options_colored_points(self, y, weight): + symbol, _ = self._marker_symbol_and_size() return { "x": self.x, "y": y, "mode": "markers", - "marker": {"color": weight, "coloraxis": "coloraxis"}, + "marker": {"symbol": symbol, "color": weight, "coloraxis": "coloraxis"}, } def _common_options(self, first_trace): @@ -231,7 +274,7 @@ def _common_options(self, first_trace): "name": self.label, "text": self._convert_annotations(), "legendgroup": self.label, - "showlegend": first_trace, + "showlegend": first_trace and self.show_legend, "yaxis": "y2" if self.y2 else "y", } diff --git a/tests/calculation/test_electronic_minimization.py b/tests/calculation/test_electronic_minimization.py index 2bbfece0..bef9d20a 100644 --- a/tests/calculation/test_electronic_minimization.py +++ b/tests/calculation/test_electronic_minimization.py @@ -7,12 +7,13 @@ import numpy as np import pytest -from py4vasp import exception +from py4vasp import exception, raw from py4vasp._calculation.electronic_minimization import ( ElectronicMinimization, ElectronicMinimizationHandler, ) from py4vasp._raw.data_db import ElectronicMinimization_DB +from py4vasp._third_party.graph import Marker @pytest.fixture @@ -27,7 +28,10 @@ def electronic_minimization(raw_data): electronic_minimization.ref.deps = convergence_data[:, 3] electronic_minimization.ref.ncg = convergence_data[:, 4] electronic_minimization.ref.rms = convergence_data[:, 5] - electronic_minimization.ref.rmsc = convergence_data[:, 6] + # read() nulls out the not-yet-computed (zero) rms(c) entries of the early steps + rms_c = np.array(convergence_data[:, 6], dtype=float) + rms_c[rms_c == 0.0] = np.nan + electronic_minimization.ref.rms_c = rms_c is_elmin_converged = [raw_elmin.is_elmin_converged == [0.0]] electronic_minimization.ref.is_elmin_converged = is_elmin_converged string_rep = "N\t\tE\t\tdE\t\tdeps\t\tncg\trms\t\trms(c)\n" @@ -54,19 +58,60 @@ def test_read(electronic_minimization, Assert): Assert.allclose(actual["deps"], expected.deps) Assert.allclose(actual["ncg"], expected.ncg) Assert.allclose(actual["rms"], expected.rms) - Assert.allclose(actual["rms(c)"], expected.rmsc) + Assert.allclose(actual["rms_c"], expected.rms_c) @pytest.mark.parametrize( - "quantity_name", ["N", "E", "dE", "deps", "ncg", "rms", "rms(c)"] + "quantity_name", ["N", "E", "dE", "deps", "ncg", "rms", "rms_c"] ) def test_read_selection(quantity_name, electronic_minimization, Assert): actual = electronic_minimization.read(quantity_name) - name_without_parenthesis = quantity_name.replace("(", "").replace(")", "") - expected = getattr(electronic_minimization.ref, name_without_parenthesis) + expected = getattr(electronic_minimization.ref, quantity_name) Assert.allclose(actual[quantity_name], expected) +def test_read_list(electronic_minimization, Assert): + actual = electronic_minimization.read("E, dE") + assert set(actual) == {"E", "dE"} + Assert.allclose(actual["E"], electronic_minimization.ref.E) + Assert.allclose(actual["dE"], electronic_minimization.ref.dE) + + +def test_read_addition(electronic_minimization, Assert): + actual = electronic_minimization.read("E + dE") + assert list(actual) == ["E + dE"] + expected = electronic_minimization.ref.E + electronic_minimization.ref.dE + Assert.allclose(actual["E + dE"], expected) + + +def test_read_nulls_early_rms_c(electronic_minimization, Assert): + # the first (zero) rms(c) entries should be reported as NaN, real ones kept + actual = electronic_minimization.read("rms_c")["rms_c"] + assert np.all(np.isnan(actual[:5])) + assert not np.any(np.isnan(actual[5:])) + Assert.allclose(actual, electronic_minimization.ref.rms_c) + + +def test_sanity_check_applies_to_any_column(Assert): + # values below the threshold are nulled regardless of which column they are in + convergence_data = np.array( + [ + [1, -8.0, 1e-30, 5e-1, 5, 3e-1, 0.0], + [2, -8.1, 1e-3, 1e-25, 6, 1e-1, 0.0], + ] + ) + raw_elmin = raw.ElectronicMinimization( + convergence_data=raw.VaspData(convergence_data), + label=raw.VaspData([b"N", b"E", b"dE", b"deps", b"ncg", b"rms", b"rms(c)"]), + is_elmin_converged=[0], + ) + data = ElectronicMinimizationHandler.from_data(raw_elmin).to_dict() + assert np.isnan(data["dE"][0]) and not np.isnan(data["dE"][1]) + assert not np.isnan(data["deps"][0]) and np.isnan(data["deps"][1]) + assert np.all(np.isnan(data["rms_c"])) + assert not np.any(np.isnan(data["E"])) + + def test_read_incorrect_selection(electronic_minimization): with pytest.raises(exception.RefinementError): electronic_minimization.read("forces") @@ -81,16 +126,129 @@ def test_slice(electronic_minimization, Assert): Assert.allclose(actual["deps"], expected.deps) Assert.allclose(actual["ncg"], expected.ncg) Assert.allclose(actual["rms"], expected.rms) - Assert.allclose(actual["rms(c)"], expected.rmsc) + Assert.allclose(actual["rms_c"], expected.rms_c) def test_plot(electronic_minimization, Assert): graph = electronic_minimization.plot() + ref = electronic_minimization.ref assert graph.xlabel == "Iteration number" - assert graph.ylabel == "E" - assert len(graph.series) == 1 - Assert.allclose(graph.series[0].x, electronic_minimization.ref.N) - Assert.allclose(graph.series[0].y, electronic_minimization.ref.E) + assert graph.ylabel == "Energy change (eV)" + assert graph.y2label == "Residual" + assert graph.yscale == "log" + assert graph.y2scale == "log" + assert len(graph.series) == 5 + by_label = {series.label: series for series in graph.series} + for series in graph.series: + Assert.allclose(series.x, ref.N) + # energy changes on the (log) left axis; the demo energy converges downward, so the + # distance E - E_final is positive and plotted (and labelled) as-is + energy_change = ref.E - ref.E[-1] + energy_change[-1] = np.nan # final point is zero and dropped on the log axis + Assert.allclose(by_label["E - E_final"].y, energy_change) + Assert.allclose(by_label["dE"].y, ref.dE) + Assert.allclose(by_label["d eps"].y, ref.deps) + assert not by_label["E - E_final"].y2 + assert not by_label["dE"].y2 + # residuals on the (log) secondary axis + Assert.allclose(by_label["rms"].y, ref.rms) + Assert.allclose(by_label["rms_c"].y, ref.rms_c) + assert by_label["rms"].y2 + assert by_label["rms_c"].y2 + + +def test_plot_with_ort_residual(Assert): + # ALGO=All labels the 7th column "ort" instead of "rms(c)"; the overview must derive + # the residual columns from the labels instead of hardcoding the name + convergence_data = np.array( + [ + [1, -8.0, 0.5, 0.2, 5, 0.3, -0.1], + [2, -8.2, 0.3, 0.4, 6, 0.2, -0.2], + ] + ) + raw_elmin = raw.ElectronicMinimization( + convergence_data=raw.VaspData(convergence_data), + label=raw.VaspData([b"N", b"E", b"dE", b"deps", b"ncg", b"rms", b"ort"]), + is_elmin_converged=[0], + ) + graph = ElectronicMinimizationHandler.from_data(raw_elmin).to_graph() + by_label = {series.label: series for series in graph.series} + # "ort" is negative, so it is shown as its magnitude on the (secondary) residual axis + assert "|ort|" in by_label + assert "rms" in by_label + assert by_label["|ort|"].y2 + Assert.allclose(by_label["|ort|"].y, np.abs(convergence_data[:, 6])) + Assert.allclose(by_label["rms"].y, convergence_data[:, 5]) + + +def test_plot_uses_absolute_value_for_signed_series(Assert): + # a directly plotted series with negative values (e.g. "ort", mislabelled as + # rms(c) for ALGO=All) is shown as its magnitude and labelled "|label|" + convergence_data = np.array( + [ + [1, -8.0, -0.5, 0.2, 5, 0.3, -0.1], + [2, -8.2, 0.3, -0.4, 6, 0.2, -0.2], + [3, -8.3, 0.1, 0.05, 7, 0.1, -0.05], + ] + ) + raw_elmin = raw.ElectronicMinimization( + convergence_data=raw.VaspData(convergence_data), + label=raw.VaspData([b"N", b"E", b"dE", b"deps", b"ncg", b"rms", b"rms(c)"]), + is_elmin_converged=[0], + ) + graph = ElectronicMinimizationHandler.from_data(raw_elmin).to_graph() + by_label = {series.label: series for series in graph.series} + # energy decreases here, so E - E_final stays positive and keeps its plain label + main_labels = {label for label in by_label if label != "unusual sign"} + assert main_labels == {"E - E_final", "|dE|", "|d eps|", "rms", "|rms_c|"} + Assert.allclose(by_label["|rms_c|"].y, np.abs(convergence_data[:, 6])) + Assert.allclose(by_label["|dE|"].y, np.abs(convergence_data[:, 2])) + Assert.allclose(by_label["rms"].y, convergence_data[:, 5]) + + +def test_plot_marks_unusual_sign_points(Assert): + # dE is usually negative with one positive blip (left axis); ort is usually negative + # with one positive blip (right axis). Only these atypical-sign points are flagged. + convergence_data = np.array( + [ + [1, -8.0, -0.5, 0.3, 5, 0.30, 0.1], + [2, -8.2, -0.3, 0.2, 6, 0.20, -0.2], + [3, -8.3, -0.1, 0.1, 7, 0.10, -0.3], + [4, -8.4, 0.2, 0.05, 8, 0.05, -0.15], + ] + ) + raw_elmin = raw.ElectronicMinimization( + convergence_data=raw.VaspData(convergence_data), + label=raw.VaspData([b"N", b"E", b"dE", b"deps", b"ncg", b"rms", b"ort"]), + is_elmin_converged=[0], + ) + graph = ElectronicMinimizationHandler.from_data(raw_elmin).to_graph() + main = [series for series in graph.series if series.label != "unusual sign"] + unusual = [series for series in graph.series if series.label == "unusual sign"] + assert len(main) == 5 + # every overlay is markers-only using the small "x" marker + assert all(series.marker == Marker(symbol="x", size=8) for series in unusual) + # a single shared legend entry even though the points span both axes + assert {series.y2 for series in unusual} == {False, True} + assert sum(series.show_legend for series in unusual) == 1 + left = next(series for series in unusual if not series.y2) + right = next(series for series in unusual if series.y2) + Assert.allclose(left.x, [4]) # dE > 0 only at step 4 (energy went up) + Assert.allclose(left.y, [0.2]) + Assert.allclose(right.x, [1]) # ort > 0 only at step 1 + Assert.allclose(right.y, [0.1]) + + +def test_plot_selection(electronic_minimization, Assert): + graph = electronic_minimization.plot("E, rms") + ref = electronic_minimization.ref + assert graph.yscale == "log" + # the total energy is negative, so on the log axis it is shown as its magnitude + assert {series.label for series in graph.series} == {"|E|", "rms"} + by_label = {series.label: series for series in graph.series} + Assert.allclose(by_label["|E|"].x, ref.N) + Assert.allclose(by_label["|E|"].y, np.abs(ref.E)) + Assert.allclose(by_label["rms"].y, ref.rms) def test_print(electronic_minimization, format_): diff --git a/tests/third_party/graph/test_graph.py b/tests/third_party/graph/test_graph.py index e526c22a..f31ec8ae 100644 --- a/tests/third_party/graph/test_graph.py +++ b/tests/third_party/graph/test_graph.py @@ -8,7 +8,7 @@ import pytest from py4vasp import _config, exception -from py4vasp._third_party.graph import Contour, Graph, Series +from py4vasp._third_party.graph import Contour, Graph, Marker, Series from py4vasp._util import import_, slicing px = import_.optional("plotly.express") @@ -292,6 +292,47 @@ def test_secondary_yaxis(parabola, sine, Assert): assert fig.layout.yaxis2.title.text == graph.y2label assert len(fig.data) == 2 assert fig.data[1].yaxis == "y2" + # the legend is pushed clear of the secondary y-axis title on the right, with a + # wider right margin so the labels are not clipped + assert fig.layout.legend.x > 1.05 + assert fig.layout.margin.r > 120 + + +def test_single_yaxis_uses_default_legend(parabola, sine): + pytest.importorskip("plotly") + fig = Graph([parabola, sine]).to_plotly() + # no secondary axis: keep plotly's default legend placement + assert fig.layout.legend.x is None + + +def test_default_axis_scale_is_linear(parabola): + pytest.importorskip("plotly") + fig = Graph(parabola).to_plotly() + assert fig.layout.xaxis.type in (None, "linear") + assert fig.layout.yaxis.type in (None, "linear") + + +def test_log_scale_axes(parabola): + pytest.importorskip("plotly") + graph = Graph(parabola, xscale="log", yscale="log") + fig = graph.to_plotly() + assert fig.layout.xaxis.type == "log" + assert fig.layout.yaxis.type == "log" + + +def test_secondary_yaxis_log_scale(parabola, sine): + pytest.importorskip("plotly") + sine.y2 = True + graph = Graph([parabola, sine], y2scale="log") + fig = graph.to_plotly() + assert fig.layout.yaxis2.type == "log" + + +def test_scale_preserved_on_merge(parabola, sine): + pytest.importorskip("plotly") + merged = Graph(parabola, yscale="log") + Graph(sine, yscale="log") + assert merged.yscale == "log" + assert merged.to_plotly().layout.yaxis.type == "log" def check_legend_group(converted, original, first_trace): @@ -366,6 +407,38 @@ def test_simple_with_marker(sine): assert fig.data[0].marker.size is None +def test_marker_symbol(parabola): + pytest.importorskip("plotly") + series = dataclasses.replace(parabola, marker=Marker(symbol="*")) + fig = Graph(series).to_plotly() + assert fig.data[0].mode == "markers" + assert fig.data[0].marker.symbol == "star" # friendly name mapped to plotly + + +def test_marker_symbol_and_size(parabola): + pytest.importorskip("plotly") + series = dataclasses.replace(parabola, marker=Marker(symbol="*", size=7)) + fig = Graph(series).to_plotly() + assert fig.data[0].marker.symbol == "star" + assert fig.data[0].marker.size == 7 + + +def test_marker_string_shorthand(parabola): + pytest.importorskip("plotly") + series = dataclasses.replace(parabola, marker="o") + fig = Graph(series).to_plotly() + assert fig.data[0].mode == "markers" + assert fig.data[0].marker.symbol == "circle" + + +def test_hide_legend_entry(parabola, sine): + pytest.importorskip("plotly") + hidden = dataclasses.replace(sine, show_legend=False) + fig = Graph([parabola, hidden]).to_plotly() + assert fig.data[0].showlegend is True + assert fig.data[1].showlegend is False + + def test_fatband_with_marker(fatband, Assert): pytest.importorskip("plotly") with_marker = dataclasses.replace(fatband, marker="o")