diff --git a/doc/velph-subcommands.md b/doc/velph-subcommands.md index 3ed9433..8798a1b 100644 --- a/doc/velph-subcommands.md +++ b/doc/velph-subcommands.md @@ -216,7 +216,8 @@ VASP input files were generated in "transport". #### `velph transport plot transport` The following command opens a graphical window that contains plots of transport -properties in `vaspout.h5`. +properties. By default, it uses the `transport/vaspout.h5` file unless a +different file is specified as the first argument. ```bash % velph transport plot transport @@ -224,10 +225,12 @@ properties in `vaspout.h5`. #### `velph transport plot eigenvalues` -The following command opens a graphical window that displays eigenvalues in the -Brillouin zone with non-zero and non-one occupations. These eigenvalues are -obtained from `vaspout.h5`, and their occupations are computed based on the -temperature and Fermi level specified through the command options. +The following command opens a graphical window to display eigenvalues in the +Brillouin zone with occupations that are neither zero nor one. By default, these +eigenvalues are extracted from the `transport/vaspout.h5` file unless a +different file is specified as the first argument. Their occupations are +determined based on the temperature and Fermi level set through the command-line +options. ```bash % velph transport plot eigenvalues @@ -246,7 +249,9 @@ efficiently. #### `velph transport plot selfenergy` -The following command opens a graphical window that displays Fan self-energies. +The following command opens a graphical window to display Fan self-energies. By +default, it uses the `transport/vaspout.h5` file unless a different file is +specified as the first argument. ```bash % velph transport plot selfenergy diff --git a/src/phelel/velph/cli/transport/plot/__init__.py b/src/phelel/velph/cli/transport/plot/__init__.py index 9720018..92bc14b 100644 --- a/src/phelel/velph/cli/transport/plot/__init__.py +++ b/src/phelel/velph/cli/transport/plot/__init__.py @@ -126,14 +126,25 @@ def cmd_plot_eigenvalues( args = _get_f_h5py_and_plot_filename( "transport", vaspout_filename=pathlib.Path(vaspout_filename) ) - if args[0] is not None: - plot_eigenvalues( - args[0], - tid=tid, - temperature=temperature, - cutoff_occupancy=cutoff_occupancy, - mu=mu, - ) + if args[0] is None: + return + + retvals = plot_eigenvalues( + args[0], + tid=tid, + temperature=temperature, + cutoff_occupancy=cutoff_occupancy, + mu=mu, + ) + + if retvals is not None: + with open("transport/bz.dat", "w") as w: + for i, (e, wt, rk) in enumerate(zip(*retvals)): + print( + f"{i + 1} {e:.6f} {wt:.6f} [{rk[0]:.6f} {rk[1]:.6f} {rk[2]:.6f}]", + file=w, + ) + click.echo('"transport/bz.dat" file was created.') def _get_f_h5py_and_plot_filename( @@ -146,11 +157,13 @@ def _get_f_h5py_and_plot_filename( click.echo("Please specify vaspout.h5 file path.") return None, None + dir_name = vaspout_filename.parent + if plot_filename is None: - _plot_filename = vaspout_filename.parent / f"{property_name}.pdf" + _plot_filename = dir_name / f"{property_name}.pdf" else: _plot_filename = plot_filename f_h5py = h5py.File(vaspout_filename) - return f_h5py, _plot_filename + return f_h5py, _plot_filename, dir_name diff --git a/src/phelel/velph/cli/transport/plot/plot_eigenvalues.py b/src/phelel/velph/cli/transport/plot/plot_eigenvalues.py index 036445d..cbd3fac 100644 --- a/src/phelel/velph/cli/transport/plot/plot_eigenvalues.py +++ b/src/phelel/velph/cli/transport/plot/plot_eigenvalues.py @@ -36,7 +36,7 @@ def plot_eigenvalues( cutoff_occupancy: float = 1e-2, mu: Optional[float] = None, time_reversal: bool = True, -): +) -> Optional[tuple[np.ndarray, np.ndarray, np.ndarray]]: """Show eigenvalues, occupation, k-points and Fermi-Dirac distribution. Parameters @@ -114,13 +114,6 @@ def plot_eigenvalues( all_weights = np.array(all_weights) all_eigenvals = np.array(all_eigenvals) - with open("bz.dat", "w") as w: - for i, (e, wt, rk) in enumerate(zip(all_eigenvals, all_weights, all_kpoints)): - print( - f"{i + 1} {e:.6f} {wt:.6f} [{rk[0]:.6f} {rk[1]:.6f} {rk[2]:.6f}]", - file=w, - ) - _plot_eigenvalues_in_BZ( all_kpoints, all_weights, @@ -128,6 +121,8 @@ def plot_eigenvalues( title=f"mu={_mu:.6f} eV, temperature={_temperature:.1f} K", ) + return all_eigenvals, all_weights, all_kpoints + def _plot_eigenvalues_in_BZ( data: np.ndarray, diff --git a/src/phelel/velph/cli/transport/plot/plot_selfenergy.py b/src/phelel/velph/cli/transport/plot/plot_selfenergy.py index 9ffa473..03e723c 100644 --- a/src/phelel/velph/cli/transport/plot/plot_selfenergy.py +++ b/src/phelel/velph/cli/transport/plot/plot_selfenergy.py @@ -2,11 +2,18 @@ from __future__ import annotations +import pathlib + import click import h5py -def plot_selfenergy(f_h5py: h5py.File, plot_filename: str, save_plot: bool = False): +def plot_selfenergy( + f_h5py: h5py.File, + plot_filename: str, + dir_name: pathlib.Path, + save_plot: bool = False, +): """Plot imaginary part of self-energies. Number of "self_energy_*" is @@ -37,11 +44,18 @@ def plot_selfenergy(f_h5py: h5py.File, plot_filename: str, save_plot: bool = Fal nrows = len(selfens) // 2 fig, axs = plt.subplots(nrows, 2, figsize=(8, 4 * nrows), squeeze=True) + lines = [] for i in range(len(selfens)): selfen = selfens[i + 1] - _show(selfen, i + 1) + lines += _get_text_lines(selfen, i + 1) _plot(axs[i], selfen) + with open(dir_name / "selfenergy.txt", "w") as w: + print("\n".join(lines), file=w) + click.echo( + f'Summary of data structure was saved in "{dir_name / "selfenergy.txt"}".' + ) + plt.tight_layout() if save_plot: plt.rcParams["pdf.fonttype"] = 42 @@ -62,7 +76,7 @@ def _plot(ax, selfen): ) -def _show(selfen: h5py.Group, index: int): +def _get_text_lines(selfen: h5py.Group, index: int) -> list[str]: """Show self-energy properties. ['band_start', 'band_stop', 'bks_idx', 'carrier_per_cell', @@ -71,27 +85,31 @@ def _show(selfen: h5py.Group, index: int): 'selfen_dw', 'selfen_fan', 'static', 'tetrahedron'] """ - print(f"- parameters: # {index}") - print( - " scattering_approximation:", - selfen["scattering_approximation"][()].decode("utf-8"), - ) - print(f" static_approximation: {bool(selfen['static'][()])}") - print(f" use_tetrahedron_method: {bool(selfen['tetrahedron'][()])}") + lines = [ + f"- parameters: # {index}", + " scattering_approximation: {}".format( + selfen["scattering_approximation"][()].decode("utf-8") + ), + f" static_approximation: {bool(selfen['static'][()])}", + f" use_tetrahedron_method: {bool(selfen['tetrahedron'][()])}", + ] if not selfen["tetrahedron"][()]: - print(f" smearing_width: {selfen['delta'][()]}") - print( - f" band_start_stop: [{selfen['band_start'][()]}, {selfen['band_stop'][()]}]" - ) - print(f" nbands: {selfen['nbands'][()]}") - print(f" nbands_sum: {selfen['nbands_sum'][()]}") - print(f" nw: {selfen['nw'][()]}") - print(" temperatures:") + lines.append(f" smearing_width: {selfen['delta'][()]}") + lines += [ + f" band_start_stop: [{selfen['band_start'][()]}, {selfen['band_stop'][()]}]", + f" nbands: {selfen['nbands'][()]}", + f" nbands_sum: {selfen['nbands_sum'][()]}", + f" nw: {selfen['nw'][()]}", + " temperatures:", + ] for i, t in enumerate(selfen["temps"]): - print(f" - {t} # {i + 1}") - - print(" data_array_shapes:") - print(f" carrier_per_cell: {list(selfen['carrier_per_cell'].shape)}") - print(f" Fan_self_energy: {list(selfen['selfen_fan'].shape)}") - print(f" sampling_energy_points: {list(selfen['energies'].shape)}") - print(f" Fermi_energies: {list(selfen['efermi'].shape)}") + lines.append(f" - {t} # {i + 1}") + + lines += [ + " data_array_shapes:", + f" carrier_per_cell: {list(selfen['carrier_per_cell'].shape)}", + f" Fan_self_energy: {list(selfen['selfen_fan'].shape)}", + f" sampling_energy_points: {list(selfen['energies'].shape)}", + f" Fermi_energies: {list(selfen['efermi'].shape)}", + ] + return lines diff --git a/src/phelel/velph/cli/transport/plot/plot_transport.py b/src/phelel/velph/cli/transport/plot/plot_transport.py index 24ccb0d..a45071e 100644 --- a/src/phelel/velph/cli/transport/plot/plot_transport.py +++ b/src/phelel/velph/cli/transport/plot/plot_transport.py @@ -2,12 +2,19 @@ from __future__ import annotations +import pathlib + import click import h5py import numpy as np -def plot_transport(f_h5py: h5py.File, plot_filename: str, save_plot: bool = False): +def plot_transport( + f_h5py: h5py.File, + plot_filename: str, + dir_name: pathlib.Path, + save_plot: bool = False, +): """Plot transport properties. Number of "transport_*" is @@ -46,8 +53,14 @@ def plot_transport(f_h5py: h5py.File, plot_filename: str, save_plot: bool = Fals for n in range(n_transport) ] + lines = [] for i, transport in enumerate(transports): - _show(transport, i + 1) + lines += _show(transport, i + 1) + with open(dir_name / "transport.txt", "w") as w: + print("\n".join(lines), file=w) + click.echo( + f'Summary of data structure was saved in "{dir_name / "transport.txt"}".' + ) n_props = len(property_names) _, axs = plt.subplots( @@ -65,6 +78,7 @@ def plot_transport(f_h5py: h5py.File, plot_filename: str, save_plot: bool = Fals i, transport["id_idx"][:], last_transport_idx, + transport_dir_name=dir_name, ) plt.tight_layout() @@ -84,6 +98,7 @@ def _plot( index: int, transport_idx: np.ndarray, last_transport_idx: np.ndarray, + transport_dir_name: pathlib.Path = pathlib.Path("transport"), ): """Plot transport properties. @@ -93,7 +108,7 @@ def _plot( import matplotlib.ticker as ticker temps = transport["temps"][:] - dat_filename = f"transport-{index + 1}.dat" + dat_filename = transport_dir_name / f"transport-{index + 1}.dat" properties = [] label_property_names = [] @@ -109,20 +124,23 @@ def _plot( else: properties.append([np.trace(tensor) / 3 for tensor in transport[key][:]]) - with open(dat_filename, "w") as w: - print("# temperature", *label_property_names, file=w) - for temp, props in zip(temps, np.transpose(properties)): - print(temp, *props, file=w) - click.echo(f'Transport data {index + 1} was saved in "{dat_filename}".') - names = [name.decode("utf-8") for name in transport["id_name"][:]] labels = [] - click.echo([(name, int(i)) for name, i in zip(names, transport_idx)]) + click.echo( + f"{index + 1}. " + + " / ".join([f"{name} {int(i)}" for name, i in zip(names, transport_idx)]) + ) for i, name in enumerate(names): if last_transport_idx[i] > 1: if name == "selfen_approx": labels.append(transport["scattering_approximation"][()].decode("utf-8")) + with open(dat_filename, "w") as w: + print("# temperature", *label_property_names, file=w) + for temp, props in zip(temps, np.transpose(properties)): + print(temp, *props, file=w) + click.echo(f'Transport data {index + 1} was saved in "{dat_filename}".') + for i, label_property in enumerate(label_property_names): if label_property == "e_conductivity": axs[i].semilogy(temps, properties[i], ".-") @@ -143,7 +161,7 @@ def _plot( ax.set_yticklabels([]) -def _show(transport: h5py.Group, index: int): +def _show(transport: h5py.Group, index: int) -> list[str]: """Show transport properties. ['cell_volume', 'dfermi_tol', 'e_conductivity', 'e_t_conductivity', 'emax', @@ -152,13 +170,15 @@ def _show(transport: h5py.Group, index: int): 'tau_average', 'temperature', 'transport_function'] """ - print(f"- parameters: # {index}") + lines = [f"- parameters: # {index}"] for i, temp in enumerate(transport["temps"][:]): - print(f" - temperature: {temp}") - print( - " scattering_approximation:", - transport["scattering_approximation"][()].decode("utf-8"), - ) + lines += [ + f" - temperature: {temp}", + " scattering_approximation: {}".format( + transport["scattering_approximation"][()].decode("utf-8") + ), + ] + for key in ( "cell_volume", "dfermi_tol", @@ -166,12 +186,12 @@ def _show(transport: h5py.Group, index: int): "nedos", ): if isinstance(key, tuple): - print(f" {key[1]}: {transport[key[0]][()]}") + lines.append(f" {key[1]}: {transport[key[0]][()]}") else: - print(f" {key}: {transport[key][()]}") - print(f" static_approximation: {bool(transport['static'][()])}") + lines.append(f" {key}: {transport[key][()]}") + lines.append(f" static_approximation: {bool(transport['static'][()])}") - print(" data_scalar:") + lines.append(" data_scalar:") for key in ( ("emax", "emax_for_transport_function"), ("emin", "emin_for_transport_function"), @@ -182,11 +202,11 @@ def _show(transport: h5py.Group, index: int): "tau_average", ): if isinstance(key, tuple): - print(f" {key[1]}: {transport[key[0]][i]}") + lines.append(f" {key[1]}: {transport[key[0]][i]}") else: - print(f" {key}: {transport[key][i]}") + lines.append(f" {key}: {transport[key][i]}") - print(" data_array_diagonal:") + lines.append(" data_array_diagonal:") for key in ( "e_conductivity", "e_t_conductivity", @@ -195,18 +215,20 @@ def _show(transport: h5py.Group, index: int): "seebeck", ): v = transport[key][:].ravel() - print( + lines.append( f" {key}: [{v[0]:.3e}, {v[4]:.3e}, {v[8]:.3e}, " f"{v[5]:.3e}, {v[2]:.3e}, {v[1]:.3e}]" ) - print(" data_array_shapes:") + lines.append(" data_array_shapes:") for key in ( "energy", ("lab", "Onsager_coefficients"), "transport_function", ): if isinstance(key, tuple): - print(f" {key[1]}: {list(transport[key[0]][i].shape)}") + lines.append(f" {key[1]}: {list(transport[key[0]][i].shape)}") else: - print(f" {key}: {list(transport[key][i].shape)}") + lines.append(f" {key}: {list(transport[key][i].shape)}") + + return lines