From 52c038a0f0d06e2c139d5d4235d8aba44271c6da Mon Sep 17 00:00:00 2001 From: jiyuang Date: Wed, 27 Jul 2022 19:18:10 +0800 Subject: [PATCH 01/13] add projected band of different contributions --- tools/plot-tools/abacus_plot/band.py | 183 ++++++++++++++++++++++++-- tools/plot-tools/abacus_plot/utils.py | 4 +- 2 files changed, 173 insertions(+), 14 deletions(-) diff --git a/tools/plot-tools/abacus_plot/band.py b/tools/plot-tools/abacus_plot/band.py index 1915d297b5..b62601468a 100644 --- a/tools/plot-tools/abacus_plot/band.py +++ b/tools/plot-tools/abacus_plot/band.py @@ -10,7 +10,7 @@ from os import PathLike from typing import Sequence, Tuple, Union, Any, List, Dict from matplotlib.figure import Figure -from matplotlib.colors import Normalize +from matplotlib.colors import Normalize, ListedColormap from matplotlib.collections import LineCollection from matplotlib import axes import matplotlib.pyplot as plt @@ -22,15 +22,16 @@ class Band: """Parse Bands data""" - def __init__(self, bandfile: Union[PathLike, Sequence[PathLike]] = None, kptfile: PathLike = '') -> None: + def __init__(self, bandfile: Union[PathLike, Sequence[PathLike]] = None, kptfile: PathLike = '', old_ver=False) -> None: self.bandfile = bandfile + self.old_ver = old_ver if isinstance(bandfile, list) or isinstance(bandfile, tuple): self.energy = [] for file in self.bandfile: - self.k_index, e = self.read(file) + self.k_index, e, self.k_lines = self.read(file, self.old_ver) self.energy.append(e) else: - self.k_index, self.energy = self.read(self.bandfile) + self.k_index, self.energy, self.k_lines = self.read(self.bandfile, self.old_ver) self.kptfile = kptfile self.kpt = None if self.kptfile: @@ -42,17 +43,24 @@ def __init__(self, bandfile: Union[PathLike, Sequence[PathLike]] = None, kptfile self._kzip = self.k_index @classmethod - def read(cls, filename: PathLike): + def read(cls, filename: PathLike, old_ver=True): """Read band data file and return k-points and energy :params filename: string of band data file """ - data = np.loadtxt(filename, dtype=float) - X, y = np.split(data, (1, ), axis=1) - x = X.flatten() + z = None + if old_ver: + data = np.loadtxt(filename, dtype=float) + X, y = np.split(data, (1, ), axis=1) + x = X.flatten() + else: + data = np.loadtxt(filename, dtype=float) + X, Z, y = np.split(data, [1, 2, data.shape[-1]], axis=1) + x = X.flatten() + z = Z.flatten() - return x, y + return x, y, z @classmethod def direct_bandgap(cls, vb: namedtuple, cb: namedtuple, klength: int): @@ -713,9 +721,9 @@ def plot(self, Returns: BandPlot object: for manually plotting picture with bandplot.ax """ - nums = len(self.bandfile) if isinstance(self.bandfile, list): + nums = len(self.bandfile) if not efermi: efermi = [0.0 for i in range(nums)] _linestyle = kwargs.pop( @@ -751,6 +759,155 @@ def plot(self, return bandplot + def _plot_contributions(self, + fig: Figure, + ax: axes.Axes, + energy: np.ndarray, + species: Union[Sequence[Any], Dict[Any, List[int]], + Dict[Any, Dict[int, List[int]]]] = [], + efermi: float = 0, + energy_range: Sequence[float] = [], + shift: bool = False, + keyname: str = '', + colors: list = [], + scale_width_factor: int = 10, + **kwargs): + """Plot parsed projected bands data of different contributions + + Args: + fig (Figure): object of matplotlib.figure.Figure + ax (Union[axes.Axes, Sequence[axes.Axes]]): object of matplotlib.axes.Axes or a list of this objects + species (Union[Sequence[Any], Dict[Any, List[int]], Dict[Any, Dict[int, List[int]]]], optional): list of atomic species(index or atom index) or dict of atomic species(index or atom index) and its angular momentum list. Defaults to []. + efermi (float, optional): fermi level in unit eV. Defaults to 0. + energy_range (Sequence[float], optional): energy range in unit eV for plotting. Defaults to []. + shift (bool, optional): if shift energy by fermi level and set the VBM to zero, or not. Defaults to False. + keyname (str, optional): the keyword that extracts the PBANDS. Defaults to ''. + + Returns: + BandPlot object: for manually plotting picture with bandplot.ax + """ + + wei, totnum = parse_projected_data(self.orbitals, species, keyname) + energy = self._shift_energy(energy, efermi, shift) + + if not species: + bandplot = BandPlot(fig, ax, **kwargs) + bandplot = super().plot(fig, ax, efermi, energy_range, shift, **kwargs) + bandplot._set_figure(self._kzip, energy_range) + + return bandplot + + whole_data_parsed = [] + whole_label_parsed = [] + bandplot = BandPlot(fig, ax, **kwargs) + if isinstance(species, (list, tuple)): + for i, elem in enumerate(wei.keys()): + whole_label_parsed.append(elem) + whole_data_parsed.append(wei[elem]) + + elif isinstance(species, dict): + for i, elem in enumerate(wei.keys()): + for ang in wei[elem].keys(): + l_index = int(ang) + if isinstance(wei[elem][ang], dict): + for mag in wei[elem][ang].keys(): + m_index = int(mag) + whole_label_parsed.append(f"{elem}-{get_angular_momentum_name(l_index, m_index)}") + whole_data_parsed.append(wei[elem][ang][mag]) + + else: + whole_label_parsed.append(f"{elem}-{get_angular_momentum_label(l_index)}") + whole_data_parsed.append(wei[elem][ang]) + + argmax_index = np.argmax(whole_data_parsed, axis=0) + max_array = np.max(whole_data_parsed, axis=0) + + if len(colors) == 0: + cmap = plt.cm.get_cmap("tab10") + colors = [cmap(c) for c in np.linspace(0, 1, len(whole_label_parsed))] + + norm = Normalize(vmin=0, vmax=1) + cmaps = ListedColormap(colors) + for ib in range(self.nbands): + bandplot.ax.scatter(self.k_index, energy[:, ib], c=argmax_index[:,ib], s=max_array[:, ib]*scale_width_factor, norm=norm, cmap=cmaps) + + clb = plt.colorbar( + plt.cm.ScalarMappable(norm=norm, cmap=cmaps), ax=bandplot.ax + ) + clb.set_ticks(np.linspace(0, 1, len(whole_label_parsed))) + clb.set_ticklabels(whole_label_parsed) + bandplot._set_figure(self._kzip, energy_range) + + return bandplot + + def plot_contributions(self, + fig: Figure, + ax: Union[axes.Axes, Sequence[axes.Axes]], + index: Union[Sequence[int], Dict[int, List[int]], + Dict[int, Dict[int, List[int]]]] = [], + atom_index: Union[Sequence[int], Dict[int, List[int]], + Dict[int, Dict[int, List[int]]]] = [], + species: Union[Sequence[str], Dict[str, List[int]], + Dict[str, Dict[int, List[int]]]] = [], + efermi: Union[float, Sequence[float]] = [], + energy_range: Sequence[float] = [], + shift: bool = False, + colors: list = [], + **kwargs): + """Plot parsed projected band data of different contributions + + Args: + fig (Figure): object of matplotlib.figure.Figure + ax (Union[axes.Axes, Sequence[axes.Axes]]): object of matplotlib.axes.Axes or a list of this objects + index (Union[Sequence[int], Dict[int, List[int]], Dict[int, Dict[int, List[int]]]], optional): extract PBAND of each atom. Defaults to []. + atom_index (Union[Sequence[int], Dict[int, List[int]], Dict[int, Dict[int, List[int]]]], optional): extract PBAND of each atom with same atom_index. Defaults to []. + species (Union[Sequence[str], Dict[str, List[int]], Dict[str, Dict[int, List[int]]]], optional): extract PBAND of each atom with same species. Defaults to []. + efermi (float, optional): fermi level in unit eV. Defaults to 0. + energy_range (Sequence[float], optional): energy range in unit eV for plotting. Defaults to []. + shift (bool, optional): if shift energy by fermi level and set the VBM to zero, or not. Defaults to False. + colors (list, optional): Default:[] + + Returns: + BandPlot object: for manually plotting picture with bandplot.ax + """ + + if isinstance(self.bandfile, list): + nums = len(self.bandfile) + if not efermi: + efermi = [0.0 for i in range(nums)] + _linestyle = kwargs.pop( + 'linestyle', ['solid' for i in range(nums)]) + + for i, band in enumerate(self.energy): + if not index and not atom_index and not species: + bandplot = self._plot_contributions(fig=fig, ax=ax, energy=band, species=[ + ], efermi=efermi[i], energy_range=energy_range, shift=shift, keyname='', linestyle=_linestyle[i], **kwargs) + if index: + bandplot = self._plot_contributions(fig=fig, ax=ax, energy=band, species=index, efermi=efermi[i], + energy_range=energy_range, shift=shift, keyname='index', linestyle=_linestyle[i], colors=colors, **kwargs) + if atom_index: + bandplot = self._plot_contributions(fig=fig, ax=ax, energy=band, species=atom_index, efermi=efermi[i], + energy_range=energy_range, shift=shift, keyname='atom_index', linestyle=_linestyle[i], colors=colors, **kwargs) + if species: + bandplot = self._plot_contributions(fig=fig, ax=ax, energy=band, species=species, efermi=efermi[i], + energy_range=energy_range, shift=shift, keyname='species', linestyle=_linestyle[i], colors=colors, **kwargs) + + else: + if not index and not atom_index and not species: + bandplot = self._plot_contributions(fig=fig, ax=ax, energy=self.energy, species=[ + ], efermi=efermi, energy_range=energy_range, shift=shift, keyname='', colors=colors, **kwargs) + if index: + bandplot = self._plot_contributions(fig=fig, ax=ax, energy=self.energy, species=index, efermi=efermi, + energy_range=energy_range, shift=shift, keyname='index', colors=colors, **kwargs) + if atom_index: + bandplot = self._plot_contributions(fig=fig, ax=ax, energy=self.energy, species=atom_index, efermi=efermi, + energy_range=energy_range, shift=shift, keyname='atom_index', colors=colors, **kwargs) + if species: + bandplot = self._plot_contributions(fig=fig, ax=ax, energy=self.energy, species=species, efermi=efermi, + energy_range=energy_range, shift=shift, keyname='species', colors=colors, **kwargs) + + return bandplot + if __name__ == "__main__": import matplotlib.pyplot as plt @@ -768,6 +925,8 @@ def plot(self, # if you want to specify `species` or `index`, you need to # set `species=species` or `index=index` in the following two functions - pband.plot(fig, ax, atom_index=atom_index, efermi=efermi, + pband.plot_contributions(fig, ax, atom_index=atom_index, efermi=efermi, energy_range=energy_range, shift=shift) - pband.write(atom_index=atom_index) + #pband.write(atom_index=atom_index) + #pband.plot(fig, ax, atom_index=atom_index, efermi=efermi, + # energy_range=energy_range, shift=shift) diff --git a/tools/plot-tools/abacus_plot/utils.py b/tools/plot-tools/abacus_plot/utils.py index 09d82f4361..7fd3679f1c 100644 --- a/tools/plot-tools/abacus_plot/utils.py +++ b/tools/plot-tools/abacus_plot/utils.py @@ -97,8 +97,8 @@ class Kpt: def __init__(self, mode: str, numbers: list = [], special_k: list = [], offset: list = [0.0, 0.0, 0.0], klabel: list = []) -> None: """Initialize Kpt object - :params mode: ‘Gamma’, ‘MP’ or 'Line' - :params numbers: for ‘Gamma’ and ‘MP’ mode, list of three integers, for 'Line' mode, list of number of k points between two adjacent special k points. + :params mode: 'sGamma', 'MP' or 'Line' + :params numbers: for 'Gamma' and 'MP' mode, list of three integers, for 'Line' mode, list of number of k points between two adjacent special k points. :params special_k: list of special k-points :params offset: offset of the k grid. Default: [0.0, 0.0, 0.0] """ From e06be7af86f3b5025aae5e0db48997598d692101 Mon Sep 17 00:00:00 2001 From: jiyuang Date: Wed, 27 Jul 2022 19:20:19 +0800 Subject: [PATCH 02/13] add python package in .gitignore --- .gitignore | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/.gitignore b/.gitignore index 7299d770ff..677b17cc63 100644 --- a/.gitignore +++ b/.gitignore @@ -12,3 +12,7 @@ result.out html *.log STRU_READIN_ADJUST.cif +*.egg +*.egg-info +build +dist \ No newline at end of file From 90febaa20ea6189847963a9a3934bb44ff06cdda Mon Sep 17 00:00:00 2001 From: jiyuang Date: Wed, 27 Jul 2022 20:31:57 +0800 Subject: [PATCH 03/13] add weights check --- tools/plot-tools/abacus_plot/band.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/tools/plot-tools/abacus_plot/band.py b/tools/plot-tools/abacus_plot/band.py index b62601468a..acf853b1de 100644 --- a/tools/plot-tools/abacus_plot/band.py +++ b/tools/plot-tools/abacus_plot/band.py @@ -43,7 +43,7 @@ def __init__(self, bandfile: Union[PathLike, Sequence[PathLike]] = None, kptfile self._kzip = self.k_index @classmethod - def read(cls, filename: PathLike, old_ver=True): + def read(cls, filename: PathLike, old_ver=False): """Read band data file and return k-points and energy :params filename: string of band data file @@ -367,6 +367,8 @@ def __init__(self, bandfile: Union[PathLike, Sequence[PathLike]] = None, kptfile else: self._kzip = self.k_index + self._check_weights(self.weights) + def _check_energy(self, energy): assert energy.shape[0] == self.nkpoints, "The dimension of band structure dismatches with the number of k-points." assert energy.shape[1] == self.nbands, "The dimension of band structure dismatches with the number of bands." @@ -384,7 +386,6 @@ def weights(self): data = np.empty((self.norbitals, self.nkpoints, self.nbands)) for i, orb in enumerate(self.orbitals): data[i] = orb['data'] - self._check_weights(data) return data @classmethod From 6e04f0875462cb429651d8476ffa376bf08270db Mon Sep 17 00:00:00 2001 From: jiyuang Date: Wed, 27 Jul 2022 21:11:55 +0800 Subject: [PATCH 04/13] modify example --- tools/plot-tools/abacus_plot/band.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/tools/plot-tools/abacus_plot/band.py b/tools/plot-tools/abacus_plot/band.py index acf853b1de..6d1119df3b 100644 --- a/tools/plot-tools/abacus_plot/band.py +++ b/tools/plot-tools/abacus_plot/band.py @@ -916,13 +916,14 @@ def plot_contributions(self, parent = Path(r"../examples/Si") name = "PBANDS_1" path = parent/name + kptfile = parent/'KLINES' fig, ax = plt.subplots(figsize=(12, 6)) energy_range = [-5, 7] efermi = 6.585653952007503 shift = False #species = {"Ag": [2], "Cl": [1], "In": [0]} atom_index = {1: {1: [0, 1]}} - pband = PBand(str(path)) + pband = PBand(str(path), kptfile) # if you want to specify `species` or `index`, you need to # set `species=species` or `index=index` in the following two functions From 567a6773be4765d8693197fbaa2a22aced56bcc5 Mon Sep 17 00:00:00 2001 From: jiyuang Date: Fri, 29 Jul 2022 11:49:13 +0800 Subject: [PATCH 05/13] plot different contributions --- tools/plot-tools/abacus_plot/band.py | 148 ++++++++++++++++----------- 1 file changed, 90 insertions(+), 58 deletions(-) diff --git a/tools/plot-tools/abacus_plot/band.py b/tools/plot-tools/abacus_plot/band.py index 6d1119df3b..f66d2857e1 100644 --- a/tools/plot-tools/abacus_plot/band.py +++ b/tools/plot-tools/abacus_plot/band.py @@ -1,7 +1,7 @@ ''' Date: 2021-12-29 10:27:01 LastEditors: jiyuyang -LastEditTime: 2022-04-26 00:00:41 +LastEditTime: 2022-07-29 11:48:18 Mail: jiyuyang@mail.ustc.edu.cn, 1041176461@qq.com ''' @@ -10,7 +10,8 @@ from os import PathLike from typing import Sequence, Tuple, Union, Any, List, Dict from matplotlib.figure import Figure -from matplotlib.colors import Normalize, ListedColormap +from matplotlib.lines import Line2D +from matplotlib.colors import Normalize, ListedColormap, LinearSegmentedColormap from matplotlib.collections import LineCollection from matplotlib import axes import matplotlib.pyplot as plt @@ -31,7 +32,8 @@ def __init__(self, bandfile: Union[PathLike, Sequence[PathLike]] = None, kptfile self.k_index, e, self.k_lines = self.read(file, self.old_ver) self.energy.append(e) else: - self.k_index, self.energy, self.k_lines = self.read(self.bandfile, self.old_ver) + self.k_index, self.energy, self.k_lines = self.read( + self.bandfile, self.old_ver) self.kptfile = kptfile self.kpt = None if self.kptfile: @@ -250,9 +252,10 @@ def plot(self, class BandPlot: """Plot band structure""" - def __init__(self, fig: Figure, ax: axes.Axes, **kwargs) -> None: + def __init__(self, fig: Figure, ax: axes.Axes, point_to_line=False, **kwargs) -> None: self.fig = fig self.ax = ax + self._point_to_line = point_to_line self._lw = kwargs.pop('lw', 2) self._bwidth = kwargs.pop('bwdith', 3) self._label = kwargs.pop('label', None) @@ -270,12 +273,12 @@ def _set_figure(self, kzip, range: Sequence): keys = [] values = [] for t in kzip: - if isinstance(t, tuple): + if isinstance(t, tuple) or isinstance(t, list): keys.append(t[0]) values.append(t[1]) - elif isinstance(t, int): + else: keys.append('') - values.append(t) + values.append(int(t)) # x-axis self.ax.set_xticks(values) @@ -339,6 +342,23 @@ def _set_figure(self, kzip, range: Sequence): self.ax.legend(by_label.values(), by_label.keys(), prop={'size': 15}) + if self._point_to_line: + handles = [] + for c, l in zip(self._color, self._label): + handles.append(Line2D([0], [0], color=c, label=l)) + if "legend_prop" in self.plot_params.keys(): + self.ax.legend(handles=handles, + prop=self.plot_params["legend_prop"]) + else: + self.ax.legend(handles=handles, prop={'size': 15}) + + def _color_to_alpha_cmap(self, color): + cmap = LinearSegmentedColormap.from_list("chaos", ["white", color]) + my_cmap = cmap(np.arange(cmap.N)) + my_cmap[:, -1] = np.linspace(0, 1, cmap.N) # this adds alpha + my_cmap = ListedColormap(my_cmap) + return my_cmap + class PBand(Band): def __init__(self, bandfile: Union[PathLike, Sequence[PathLike]] = None, kptfile: str = '') -> None: @@ -761,18 +781,18 @@ def plot(self, return bandplot def _plot_contributions(self, - fig: Figure, - ax: axes.Axes, - energy: np.ndarray, - species: Union[Sequence[Any], Dict[Any, List[int]], - Dict[Any, Dict[int, List[int]]]] = [], - efermi: float = 0, - energy_range: Sequence[float] = [], - shift: bool = False, - keyname: str = '', - colors: list = [], - scale_width_factor: int = 10, - **kwargs): + fig: Figure, + ax: axes.Axes, + energy: np.ndarray, + species: Union[Sequence[Any], Dict[Any, List[int]], + Dict[Any, Dict[int, List[int]]]] = [], + efermi: float = 0, + energy_range: Sequence[float] = [], + shift: bool = False, + keyname: str = '', + colors: list = [], + scale_width_factor: int = 5, + **kwargs): """Plot parsed projected bands data of different contributions Args: @@ -800,7 +820,6 @@ def _plot_contributions(self, whole_data_parsed = [] whole_label_parsed = [] - bandplot = BandPlot(fig, ax, **kwargs) if isinstance(species, (list, tuple)): for i, elem in enumerate(wei.keys()): whole_label_parsed.append(elem) @@ -813,48 +832,54 @@ def _plot_contributions(self, if isinstance(wei[elem][ang], dict): for mag in wei[elem][ang].keys(): m_index = int(mag) - whole_label_parsed.append(f"{elem}-{get_angular_momentum_name(l_index, m_index)}") + whole_label_parsed.append( + f"{elem}-{get_angular_momentum_name(l_index, m_index)}") whole_data_parsed.append(wei[elem][ang][mag]) else: - whole_label_parsed.append(f"{elem}-{get_angular_momentum_label(l_index)}") + whole_label_parsed.append( + f"{elem}-{get_angular_momentum_label(l_index)}") whole_data_parsed.append(wei[elem][ang]) - argmax_index = np.argmax(whole_data_parsed, axis=0) - max_array = np.max(whole_data_parsed, axis=0) - if len(colors) == 0: cmap = plt.cm.get_cmap("tab10") - colors = [cmap(c) for c in np.linspace(0, 1, len(whole_label_parsed))] + colors = [cmap(c) + for c in np.linspace(0, 1, len(whole_label_parsed))] - norm = Normalize(vmin=0, vmax=1) - cmaps = ListedColormap(colors) - for ib in range(self.nbands): - bandplot.ax.scatter(self.k_index, energy[:, ib], c=argmax_index[:,ib], s=max_array[:, ib]*scale_width_factor, norm=norm, cmap=cmaps) + bandplot = BandPlot( + fig, ax, color=colors, label=whole_label_parsed, point_to_line=True, **kwargs) - clb = plt.colorbar( - plt.cm.ScalarMappable(norm=norm, cmap=cmaps), ax=bandplot.ax - ) - clb.set_ticks(np.linspace(0, 1, len(whole_label_parsed))) - clb.set_ticklabels(whole_label_parsed) + norm = Normalize(vmin=0, vmax=1) + cmaps = [bandplot._color_to_alpha_cmap(c) for c in colors] + swidth = np.array(whole_data_parsed)*scale_width_factor + for ib in range(self.nbands): + for i, con in enumerate(swidth): + bandplot.ax.scatter( + self.k_index, energy[:, ib], c=con[:, ib], s=con[:, ib], norm=norm, cmap=cmaps[i], label=whole_label_parsed[i]) + + # clb = plt.colorbar( + # plt.cm.ScalarMappable(norm=norm, cmap=cmaps), ax=bandplot.ax + # ) + # clb.set_ticks(np.linspace(0, 1, len(whole_label_parsed))) + # clb.set_ticklabels(whole_label_parsed) bandplot._set_figure(self._kzip, energy_range) return bandplot def plot_contributions(self, - fig: Figure, - ax: Union[axes.Axes, Sequence[axes.Axes]], - index: Union[Sequence[int], Dict[int, List[int]], - Dict[int, Dict[int, List[int]]]] = [], - atom_index: Union[Sequence[int], Dict[int, List[int]], - Dict[int, Dict[int, List[int]]]] = [], - species: Union[Sequence[str], Dict[str, List[int]], - Dict[str, Dict[int, List[int]]]] = [], - efermi: Union[float, Sequence[float]] = [], - energy_range: Sequence[float] = [], - shift: bool = False, - colors: list = [], - **kwargs): + fig: Figure, + ax: Union[axes.Axes, Sequence[axes.Axes]], + index: Union[Sequence[int], Dict[int, List[int]], + Dict[int, Dict[int, List[int]]]] = [], + atom_index: Union[Sequence[int], Dict[int, List[int]], + Dict[int, Dict[int, List[int]]]] = [], + species: Union[Sequence[str], Dict[str, List[int]], + Dict[str, Dict[int, List[int]]]] = [], + efermi: Union[float, Sequence[float]] = [], + energy_range: Sequence[float] = [], + shift: bool = False, + colors: list = [], + **kwargs): """Plot parsed projected band data of different contributions Args: @@ -885,13 +910,13 @@ def plot_contributions(self, ], efermi=efermi[i], energy_range=energy_range, shift=shift, keyname='', linestyle=_linestyle[i], **kwargs) if index: bandplot = self._plot_contributions(fig=fig, ax=ax, energy=band, species=index, efermi=efermi[i], - energy_range=energy_range, shift=shift, keyname='index', linestyle=_linestyle[i], colors=colors, **kwargs) + energy_range=energy_range, shift=shift, keyname='index', linestyle=_linestyle[i], colors=colors, **kwargs) if atom_index: bandplot = self._plot_contributions(fig=fig, ax=ax, energy=band, species=atom_index, efermi=efermi[i], - energy_range=energy_range, shift=shift, keyname='atom_index', linestyle=_linestyle[i], colors=colors, **kwargs) + energy_range=energy_range, shift=shift, keyname='atom_index', linestyle=_linestyle[i], colors=colors, **kwargs) if species: bandplot = self._plot_contributions(fig=fig, ax=ax, energy=band, species=species, efermi=efermi[i], - energy_range=energy_range, shift=shift, keyname='species', linestyle=_linestyle[i], colors=colors, **kwargs) + energy_range=energy_range, shift=shift, keyname='species', linestyle=_linestyle[i], colors=colors, **kwargs) else: if not index and not atom_index and not species: @@ -899,13 +924,13 @@ def plot_contributions(self, ], efermi=efermi, energy_range=energy_range, shift=shift, keyname='', colors=colors, **kwargs) if index: bandplot = self._plot_contributions(fig=fig, ax=ax, energy=self.energy, species=index, efermi=efermi, - energy_range=energy_range, shift=shift, keyname='index', colors=colors, **kwargs) + energy_range=energy_range, shift=shift, keyname='index', colors=colors, **kwargs) if atom_index: bandplot = self._plot_contributions(fig=fig, ax=ax, energy=self.energy, species=atom_index, efermi=efermi, - energy_range=energy_range, shift=shift, keyname='atom_index', colors=colors, **kwargs) + energy_range=energy_range, shift=shift, keyname='atom_index', colors=colors, **kwargs) if species: bandplot = self._plot_contributions(fig=fig, ax=ax, energy=self.energy, species=species, efermi=efermi, - energy_range=energy_range, shift=shift, keyname='species', colors=colors, **kwargs) + energy_range=energy_range, shift=shift, keyname='species', colors=colors, **kwargs) return bandplot @@ -927,8 +952,15 @@ def plot_contributions(self, # if you want to specify `species` or `index`, you need to # set `species=species` or `index=index` in the following two functions + + # 1. write different contributions to files + # pband.write(atom_index=atom_index) + + # 2. plot different contributions in single picture pband.plot_contributions(fig, ax, atom_index=atom_index, efermi=efermi, - energy_range=energy_range, shift=shift) - #pband.write(atom_index=atom_index) - #pband.plot(fig, ax, atom_index=atom_index, efermi=efermi, - # energy_range=energy_range, shift=shift) + energy_range=energy_range, shift=shift) + plt.show() + + # 3. plot different contributions to different pictures with colobar denoting weightes + # pband.plot(fig, ax, atom_index=atom_index, efermi=efermi, + # energy_range=energy_range, shift=shift) energy_range=energy_range, shift=shift) From bfb2da32d109443d6de7d0c0962324949ca8531b Mon Sep 17 00:00:00 2001 From: jiyuang Date: Tue, 2 Aug 2022 11:36:49 +0800 Subject: [PATCH 06/13] fix bug in plot new BANDS_*.dat --- tools/plot-tools/abacus_plot/band.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tools/plot-tools/abacus_plot/band.py b/tools/plot-tools/abacus_plot/band.py index f66d2857e1..4a4783ddf0 100644 --- a/tools/plot-tools/abacus_plot/band.py +++ b/tools/plot-tools/abacus_plot/band.py @@ -1,7 +1,7 @@ ''' Date: 2021-12-29 10:27:01 LastEditors: jiyuyang -LastEditTime: 2022-07-29 11:48:18 +LastEditTime: 2022-08-02 11:36:22 Mail: jiyuyang@mail.ustc.edu.cn, 1041176461@qq.com ''' @@ -58,7 +58,7 @@ def read(cls, filename: PathLike, old_ver=False): x = X.flatten() else: data = np.loadtxt(filename, dtype=float) - X, Z, y = np.split(data, [1, 2, data.shape[-1]], axis=1) + X, Z, y, _ = np.split(data, [1, 2, data.shape[-1]], axis=1) x = X.flatten() z = Z.flatten() From 8b4573657fef5a9307d96a3528ae411daa694091 Mon Sep 17 00:00:00 2001 From: jiyuyang <1041176461@qq.com> Date: Wed, 3 Aug 2022 17:31:05 +0800 Subject: [PATCH 07/13] Update band.py --- tools/plot-tools/abacus_plot/band.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/tools/plot-tools/abacus_plot/band.py b/tools/plot-tools/abacus_plot/band.py index 4a4783ddf0..3b9487e359 100644 --- a/tools/plot-tools/abacus_plot/band.py +++ b/tools/plot-tools/abacus_plot/band.py @@ -58,9 +58,12 @@ def read(cls, filename: PathLike, old_ver=False): x = X.flatten() else: data = np.loadtxt(filename, dtype=float) - X, Z, y, _ = np.split(data, [1, 2, data.shape[-1]], axis=1) - x = X.flatten() - z = Z.flatten() + #X, Z, y, _ = np.split(data, [1, 2, data.shape[-1]], axis=1) + #x = X.flatten() + #z = Z.flatten() + x = data[:, 0] # k-index + z = data[:, 1] # k-points + y = data[:, 2:] # band return x, y, z From 32652fb6cebd59b2c2dc548ec42477c5f044fb8d Mon Sep 17 00:00:00 2001 From: jiyuyang <1041176461@qq.com> Date: Wed, 3 Aug 2022 18:12:44 +0800 Subject: [PATCH 08/13] Update input-main.md --- docs/input-main.md | 19 ++++++++++++++++--- 1 file changed, 16 insertions(+), 3 deletions(-) diff --git a/docs/input-main.md b/docs/input-main.md index 8f91e11acc..0fd518ce35 100644 --- a/docs/input-main.md +++ b/docs/input-main.md @@ -1014,13 +1014,26 @@ This part of variables are relevant when using hybrid functionals #### exx_hse_omega -- **Type**: +- **Type**: Real - **Description**: range-separation parameter in HSE functional, such that $1/r=erfc(\omega r)/r+erf(\omega r)/r$. - **Default**: 0.11 -adial integration for pseudopotentials, in Bohr. -@@ -214,6 +279,13 @@ This part of variables are used to control general system para +#### exx_separate_loop + +- **Type**: Boolean +- **Description**: There are two types of iterative approach provided by ABACUS to evaluate Fock exchange. If this parameter is set to 0, it will start with a GGA-Loop, and then Hybrid-Loop, in which EXX Hamiltonian ($H_{exx}$) is updated with electronic iterations. If this parameter is set to 1, a two-step method is employed, i.e. in the inner iterations, density matrix is updated, while in the outer iterations, $H_{exx}$ is calculated based on density matrix that converges in the inner iteration. +- **Default**: 1 + +#### exx_hybrid_step +- **Type**: Integer +- **Description**: This variable indicates the maximal electronic iteration number in the evaluation of Fock exchange. +- **Default**: 100 + +#### exx_lambda + +- **Type**: Real +- **Description**: - **Default**: 0.3 #### exx_pca_threshold From c48b4e734b16b5ec9ba5ec7dd4bf993035fdfe41 Mon Sep 17 00:00:00 2001 From: jiyuyang <1041176461@qq.com> Date: Wed, 3 Aug 2022 18:30:19 +0800 Subject: [PATCH 09/13] Update input-main.md --- docs/input-main.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/docs/input-main.md b/docs/input-main.md index 0fd518ce35..0e5ecbdb8f 100644 --- a/docs/input-main.md +++ b/docs/input-main.md @@ -1021,7 +1021,7 @@ This part of variables are relevant when using hybrid functionals #### exx_separate_loop - **Type**: Boolean -- **Description**: There are two types of iterative approach provided by ABACUS to evaluate Fock exchange. If this parameter is set to 0, it will start with a GGA-Loop, and then Hybrid-Loop, in which EXX Hamiltonian ($H_{exx}$) is updated with electronic iterations. If this parameter is set to 1, a two-step method is employed, i.e. in the inner iterations, density matrix is updated, while in the outer iterations, $H_{exx}$ is calculated based on density matrix that converges in the inner iteration. +- **Description**: There are two types of iterative approach provided by ABACUS to evaluate Fock exchange. If this parameter is set to 0, it will start with a GGA-Loop, and then Hybrid-Loop, in which EXX Hamiltonian $H_{exx}$ is updated with electronic iterations. If this parameter is set to 1, a two-step method is employed, i.e. in the inner iterations, density matrix is updated, while in the outer iterations, $H_{exx}$ is calculated based on density matrix that converges in the inner iteration. - **Default**: 1 #### exx_hybrid_step @@ -1033,7 +1033,7 @@ This part of variables are relevant when using hybrid functionals #### exx_lambda - **Type**: Real -- **Description**: +- **Description**: It is used to compensate for divergence points at G=0 in the evaluation of Fock exchange using *lcao_in_pw* method. - **Default**: 0.3 #### exx_pca_threshold From 38fa9c464883f00e324e955fc8cc10e5b8be74da Mon Sep 17 00:00:00 2001 From: linpz Date: Wed, 3 Aug 2022 18:46:21 +0800 Subject: [PATCH 10/13] 1. update BlasConnector and LapackConnector --- source/module_base/blas_connector.h | 5 ++-- source/module_base/lapack_connector.h | 26 ++++++++++++++----- .../src_ri/exx_abfs-inverse_matrix_double.cpp | 4 +-- 3 files changed, 25 insertions(+), 10 deletions(-) diff --git a/source/module_base/blas_connector.h b/source/module_base/blas_connector.h index e0dfd77a50..ab5bac84be 100644 --- a/source/module_base/blas_connector.h +++ b/source/module_base/blas_connector.h @@ -36,8 +36,9 @@ extern "C" double dznrm2_( const int *n, const std::complex *X, const int *incX ); // level 2: matrix-std::vector operations, O(n^2) data and O(n^2) work. - void dgemv_(const char *transa, const int *m, const int *n, const double *alpha, const double *a, - const int *lda, const double *x, const int *incx, const double *beta, double *y, const int *incy); + void dgemv_(const char*const transa, const int*const m, const int*const n, + const double*const alpha, const double*const a, const int*const lda, const double*const x, const int*const incx, + const double*const beta, double*const y, const int*const incy); void zgemv_(const char *trans, const int *m, const int *n, const std::complex *alpha, const std::complex *a, const int *lda, const std::complex *x, const int *incx, diff --git a/source/module_base/lapack_connector.h b/source/module_base/lapack_connector.h index 90efd5b87e..b812f1d508 100644 --- a/source/module_base/lapack_connector.h +++ b/source/module_base/lapack_connector.h @@ -60,8 +60,8 @@ extern "C" // Peize Lin add dsptrf and dsptri 2016-06-21, to compute inverse real symmetry indefinit matrix. // dpotrf computes the Cholesky factorization of a real symmetric positive definite matrix // while dpotri taks its output to perform matrix inversion - void dpotrf_(const char* uplo,const int* n, double* A, const int* lda, int *info); - void dpotri_(const char* uplo,const int* n, double* A, const int* lda, int *info); + void dpotrf_(const char*const uplo, const int*const n, double*const A, const int*const lda, int*const info); + void dpotri_(const char*const uplo, const int*const n, double*const A, const int*const lda, int*const info); // zgetrf computes the LU factorization of a general matrix // while zgetri takes its output to perform matrix inversion @@ -317,18 +317,32 @@ class LapackConnector // Peize Lin add 2016-07-09 static inline - void dpotrf( char uplo, const int n, ModuleBase::matrix &a, const int lda, int *info ) + void dpotrf( const char &uplo, const int &n, double*const A, const int &lda, int &info ) { const char uplo_changed = change_uplo(uplo); - dpotrf_( &uplo_changed, &n, a.c, &lda, info ); + dpotrf_( &uplo_changed, &n, A, &lda, &info ); } // Peize Lin add 2016-07-09 static inline - void dpotri( char uplo, const int n, ModuleBase::matrix &a, const int lda, int *info ) + void dpotri( const char &uplo, const int &n, double*const A, const int &lda, int &info ) { const char uplo_changed = change_uplo(uplo); - dpotri_( &uplo_changed, &n, a.c, &lda, info); + dpotri_( &uplo_changed, &n, A, &lda, &info); + } + + // Peize Lin add 2016-07-09 + static inline + void dpotrf( const char &uplo, const int &n, ModuleBase::matrix &A, const int &lda, int &info ) + { + dpotrf( uplo, n, A.c, lda, info ); + } + + // Peize Lin add 2016-07-09 + static inline + void dpotri( const char &uplo, const int &n, ModuleBase::matrix &A, const int &lda, int &info ) + { + dpotri( uplo, n, A.c, lda, info); } // Peize Lin add 2019-04-14 diff --git a/source/src_ri/exx_abfs-inverse_matrix_double.cpp b/source/src_ri/exx_abfs-inverse_matrix_double.cpp index 51900dd49b..1ef39df37f 100644 --- a/source/src_ri/exx_abfs-inverse_matrix_double.cpp +++ b/source/src_ri/exx_abfs-inverse_matrix_double.cpp @@ -35,7 +35,7 @@ void Exx_Abfs::Inverse_Matrix_Double::cal_inverse( const Method &method ) void Exx_Abfs::Inverse_Matrix_Double::using_dpotrf() { - LapackConnector::dpotrf('U',dim,A,dim,&info); + LapackConnector::dpotrf('U',dim,A,dim,info); if(info!=0) { @@ -45,7 +45,7 @@ void Exx_Abfs::Inverse_Matrix_Double::using_dpotrf() ModuleBase::QUIT(); } - LapackConnector::dpotri('U',dim,A,dim,&info); + LapackConnector::dpotri('U',dim,A,dim,info); if(info!=0) { From e91d1c69b9181dd0325971fb80503bcd76e76522 Mon Sep 17 00:00:00 2001 From: linpz Date: Wed, 3 Aug 2022 18:53:23 +0800 Subject: [PATCH 11/13] 1. change row_set and col_set in LocalMatrix from int* to vector --- source/module_orbital/ORB_control.cpp | 6 ++---- source/module_orbital/parallel_orbitals.cpp | 10 ++-------- source/src_pdiag/pdiag_common.h | 5 +++-- 3 files changed, 7 insertions(+), 14 deletions(-) diff --git a/source/module_orbital/ORB_control.cpp b/source/module_orbital/ORB_control.cpp index d12e6b5cd6..3332b146db 100644 --- a/source/module_orbital/ORB_control.cpp +++ b/source/module_orbital/ORB_control.cpp @@ -416,8 +416,7 @@ int ORB_control::mat_2d(MPI_Comm vu, // (5) row_set, it's a global index : // save explicitly : every row in this processor // belongs to which row in the global matrix. - delete[] LM.row_set; - LM.row_set = new int[LM.row_num]; + LM.row_set.resize(LM.row_num); j = 0; for (i = 0; i < LM.row_b; i++) { @@ -475,8 +474,7 @@ int ORB_control::mat_2d(MPI_Comm vu, if (pv->testpb)ModuleBase::GlobalFunc::OUT(ofs_running, "Local columns (including nb)", LM.row_num); - delete[] LM.col_set; - LM.col_set = new int[LM.col_num]; + LM.col_set.resize(LM.col_num); j = 0; for (i = 0; i < LM.col_b; i++) diff --git a/source/module_orbital/parallel_orbitals.cpp b/source/module_orbital/parallel_orbitals.cpp index 1013caa9d9..97e5367f9d 100644 --- a/source/module_orbital/parallel_orbitals.cpp +++ b/source/module_orbital/parallel_orbitals.cpp @@ -19,8 +19,6 @@ Parallel_Orbitals::Parallel_Orbitals() // default value of nb is 1, // but can change to larger value from input. nb = 1; - MatrixInfo.row_set = nullptr; - MatrixInfo.col_set = nullptr; //in multi-k, 2D-block-division variables for FT (R<->k) nnr = 1; @@ -43,8 +41,6 @@ Parallel_Orbitals::~Parallel_Orbitals() } delete[] Z_LOC; } - delete[] MatrixInfo.row_set; - delete[] MatrixInfo.col_set; delete[] nlocdim; delete[] nlocstart; @@ -270,8 +266,7 @@ void ORB_control::divide_HS_2d this->set_parameters(ofs_running, ofs_warning); pv->MatrixInfo.row_b = 1; pv->MatrixInfo.row_num = nlocal; - delete[] pv->MatrixInfo.row_set; - pv->MatrixInfo.row_set = new int[nlocal]; + pv->MatrixInfo.row_set.resize(nlocal); for(int i=0; iMatrixInfo.row_set[i]=i; @@ -280,8 +275,7 @@ void ORB_control::divide_HS_2d pv->MatrixInfo.col_b = 1; pv->MatrixInfo.col_num = nlocal; - delete[] pv->MatrixInfo.col_set; - pv->MatrixInfo.col_set = new int[nlocal]; + pv->MatrixInfo.col_set.resize(nlocal); for(int i=0; iMatrixInfo.col_set[i]=i; diff --git a/source/src_pdiag/pdiag_common.h b/source/src_pdiag/pdiag_common.h index d6d2f306f1..b8908ec410 100644 --- a/source/src_pdiag/pdiag_common.h +++ b/source/src_pdiag/pdiag_common.h @@ -6,11 +6,12 @@ #endif #include "../module_base/blas_connector.h" #include "../module_base/lapack_connector.h" // Peize Lin add 2016-08-04 +#include struct LocalMatrix { - int *row_set; - int *col_set; + std::vector row_set; // Peize Lin change int* to vector 2022.08.03 + std::vector col_set; int col_num; int row_num; From ba2f1df9bcff2f8a9d762f4d31a50355158cc6ac Mon Sep 17 00:00:00 2001 From: linpz Date: Sat, 13 Aug 2022 19:02:54 +0800 Subject: [PATCH 12/13] 1. add openmp in ORB_table_phi::cal_ST_Phi12_R() --- source/module_orbital/ORB_table_phi.cpp | 31 +++++++++++++------------ 1 file changed, 16 insertions(+), 15 deletions(-) diff --git a/source/module_orbital/ORB_table_phi.cpp b/source/module_orbital/ORB_table_phi.cpp index 0d01a7c64f..b2c7008c9b 100644 --- a/source/module_orbital/ORB_table_phi.cpp +++ b/source/module_orbital/ORB_table_phi.cpp @@ -5,6 +5,10 @@ #include "../module_base/constants.h" #include "../module_base/timer.h" +#ifdef _OPENMP +#include +#endif + double ORB_table_phi::dr = -1.0; ORB_table_phi::ORB_table_phi() @@ -128,9 +132,7 @@ void ORB_table_phi::cal_ST_Phi12_R { ModuleBase::timer::tick("ORB_table_phi", "cal_ST_Phi12_R"); - double* k1_dot_k2 = new double[kmesh]; - double* k1_dot_k2_dot_kpoint = new double[kmesh]; - + std::vector k1_dot_k2(kmesh); // Peize Lin change 2017-12-12 switch(job) { @@ -166,6 +168,7 @@ void ORB_table_phi::cal_ST_Phi12_R break; } + std::vector k1_dot_k2_dot_kpoint(kmesh); for (int ik = 0; ik < kmesh; ik++) { k1_dot_k2_dot_kpoint[ik] = k1_dot_k2[ik] * this->kpoint[ik]; @@ -176,14 +179,18 @@ void ORB_table_phi::cal_ST_Phi12_R //previous version - double* integrated_func = new double[kmesh]; + //double* integrated_func = new double[kmesh]; const std::vector> &jlm1 = pSB->get_jlx()[l-1]; const std::vector> &jl = pSB->get_jlx()[l]; const std::vector> &jlp1 = pSB->get_jlx()[l+1]; +#ifdef _OPENMP + #pragma omp parallel for schedule(static) +#endif for (int ir = 0; ir < rmesh; ir++) { + std::vector integrated_func(kmesh); const std::vector &jl_r = jl[ir]; for (int ik=0; ik 0) { - ModuleBase::GlobalFunc::ZEROS(integrated_func,kmesh); + std::vector integrated_func(kmesh); double temp = 0.0; for (int ik = 0; ik < kmesh; ik++) { - integrated_func[ik] = k1_dot_k2[ik] * pow (kpoint[ik], l); + integrated_func[ik] = k1_dot_k2[ik] * std::pow (kpoint[ik], l); } - ModuleBase::Integral::Simpson_Integral(kmesh,integrated_func,kab,temp); + ModuleBase::Integral::Simpson_Integral(kmesh, integrated_func.data(), kab, temp); rs[0] = ModuleBase::FOUR_PI / ModuleBase::Mathzone_Add1::dualfac (2*l+1) * temp; } - delete [] integrated_func; - delete [] k1_dot_k2; - delete [] k1_dot_k2_dot_kpoint; - ModuleBase::timer::tick("ORB_table_phi", "cal_ST_Phi12_R"); - - return; } #include "../module_base/constants.h" From 36950ccd2f98c2e14ce54cdfefcdfffac69422b6 Mon Sep 17 00:00:00 2001 From: "H.Lu" <43399551+PoloTier@users.noreply.github.com> Date: Sun, 21 Aug 2022 15:10:28 +0800 Subject: [PATCH 13/13] Update install.md --- docs/install.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/install.md b/docs/install.md index 1123f6e150..85003e2193 100644 --- a/docs/install.md +++ b/docs/install.md @@ -98,7 +98,7 @@ If environment variable `MKLROOT` exists, `cmake` will take MKL as a preference, You can also choose to build with which components. ```bash -cmake -B build -DUSE_LIBXC=1 -DUSE_CUDA=1 +cmake -B build -DENABLE_LIBXC=1 -DUSE_CUDA=1 ``` If Libxc is not installed in standard path (i.e. installed with a custom prefix path), you may add the installation prefix of `FindLibxc.cmake` to `CMAKE_MODULE_PATH` environment variable, or set `Libxc_DIR` to the directory containing the file.