Skip to content

Commit

Permalink
dos --> phdos
Browse files Browse the repository at this point in the history
  • Loading branch information
gmatteo committed Dec 21, 2024
1 parent a29d307 commit cf13770
Show file tree
Hide file tree
Showing 18 changed files with 464 additions and 400 deletions.
15 changes: 12 additions & 3 deletions abipy/abio/inputs.py
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,7 @@
from monty.collections import dict2namedtuple
from monty.string import is_string, list_strings
from monty.json import MontyDecoder, MSONable
from monty.termcolor import cprint
from pymatgen.core.units import Energy
from pymatgen.symmetry.bandstructure import HighSymmKpath
from abipy.tools.numtools import is_diagonal
Expand Down Expand Up @@ -3916,7 +3917,11 @@ def phbands_and_dos(cls, structure, ngqpt, nqsmall, qppa=None, ndivsm=20, line_d
prtdos = 1
i = dos_method.find(":")
if i != -1:
value, eunit = dos_method[i+1:].split()
try:
value, eunit = dos_method[i+1:].split()
except Exception as exc:
raise ValueError(f"Invalid {dos_method=}") from exc

dossmear = Energy(float(value), eunit).to("Ha")
else:
raise NotImplementedError("Wrong value for dos_method: %s" % str(dos_method))
Expand Down Expand Up @@ -4234,8 +4239,12 @@ def to_string(self, sortmode=None, mode="text", verbose=0, files_file=False) ->
#if mode == "html": vname = root + "#%s" % vname
if mode == "html": vname = var_database[vname].html_link(label=vname)
value = format_string_abivars(vname, value, "anaddb")
#print("vname:", vname, "value:", value)
app(str(InputVariable(vname, value)))

try:
app(str(InputVariable(vname, value)))
except Exception as exc:
cprint(f"{vname=}, {value=}", color="red")
raise exc

return "\n".join(lines) if mode == "text" else "\n".join(lines).replace("\n", "<br>")

Expand Down
45 changes: 15 additions & 30 deletions abipy/dfpt/ddb.py
Original file line number Diff line number Diff line change
Expand Up @@ -1194,7 +1194,10 @@ def anaget_phbst_and_phdos_files(self, nqsmall=10, qppa=None, ndivsm=20, line_de
|PhbstFile| with the phonon band structure.
|PhdosFile| with the the phonon DOS.
"""
if ngqpt is None: ngqpt = self.guessed_ngqpt
if ngqpt is None:
ngqpt = self.guessed_ngqpt
if ngqpt is None:
raise RuntimeError(f"Not able to autodetect q-mesh associated to DDB file {self.filepath=}, {self.guessed_ngqpt=}")

if lo_to_splitting == "automatic":
lo_to_splitting = self.has_lo_to_data() and dipdip != 0
Expand Down Expand Up @@ -1240,13 +1243,11 @@ def anaget_phbst_and_phdos_files(self, nqsmall=10, qppa=None, ndivsm=20, line_de

def get_coarse(self, ngqpt_coarse, filepath=None) -> DdbFile:
"""
Get a version of this file on a coarse mesh
Return a |DdbFile| on a coarse q-mesh
Args:
ngqpt_coarse: list of ngqpt indexes that must be a sub-mesh of the original ngqpt
filepath: Filename for coarse DDB. If None, temporary filename is used.
Return: |DdbFile| on coarse mesh.
"""
# Check if ngqpt is a sub-mesh of ngqpt
ngqpt_fine = self.guessed_ngqpt
Expand Down Expand Up @@ -1339,6 +1340,7 @@ def anacompare_dipdip(self, chneut_list=(1,), asr=2, dipquad=1, quadquad=1,
Invoke anaddb to compute the phonon band structure and the phonon DOS with different
values of the ``dipdip`` input variable (dipole-dipole treatment).
Build and return |PhononDosPlotter| object.
Client code can use ``plotter.combiplot()`` or ``plotter.gridplot()`` to visualize the results.
Args:
chneut_list: List of ``chneut`` values to test (used for dipdip == 1).
Expand All @@ -1358,11 +1360,6 @@ def anacompare_dipdip(self, chneut_list=(1,), asr=2, dipquad=1, quadquad=1,
verbose: Verbosity level.
mpi_procs: Number of MPI processes used by anaddb.
pre_label: String to prepen to the default label used by the Plotter.
Return:
|PhononDosPlotter| object.
Client code can use ``plotter.combiplot()`` or ``plotter.gridplot()`` to visualize the results.
"""
phbands_plotter = PhononBandsPlotter()

Expand Down Expand Up @@ -1484,6 +1481,7 @@ def anacompare_rifcsph(self, rifcsph_list, asr=2, chneut=1, dipdip=1, dipquad=1,
Invoke anaddb to compute the phonon band structure and the phonon DOS with different
values of the ``asr`` input variable (acoustic sum rule treatment).
Build and return |PhononBandsPlotter| object.
Client code can use ``plotter.combiplot()`` or ``plotter.gridplot()`` to visualize the results.
Args:
rifcsph_list: List of rifcsph to analyze.
Expand All @@ -1498,11 +1496,6 @@ def anacompare_rifcsph(self, rifcsph_list, asr=2, chneut=1, dipdip=1, dipquad=1,
ngqpt: Number of divisions for the ab-initio q-mesh in the DDB file. Auto-detected if None (default)
verbose: Verbosity level.
mpi_procs: Number of MPI processes used by anaddb.
Return:
|PhononBandsPlotter| object.
Client code can use ``plotter.combiplot()`` or ``plotter.gridplot()`` to visualize the results.
"""
phbands_plotter = PhononBandsPlotter()

Expand All @@ -1526,6 +1519,7 @@ def anacompare_quad(self, asr=2, chneut=1, dipdip=1, lo_to_splitting="automatic"
Invoke anaddb to compute the phonon band structure and the phonon DOS by including
dipole-quadrupole and quadrupole-quadrupole terms in the dynamical matrix
Build and return |PhononBandsPlotter| object.
Client code can use ``plotter.combiplot()`` or ``plotter.gridplot()`` to visualize the results.
Args:
asr: Acoustic sum rule input variable.
Expand All @@ -1544,11 +1538,6 @@ def anacompare_quad(self, asr=2, chneut=1, dipdip=1, lo_to_splitting="automatic"
ngqpt: Number of divisions for the ab-initio q-mesh in the DDB file. Auto-detected if None (default)
verbose: Verbosity level.
mpi_procs: Number of MPI processes used by anaddb.
Return:
|PhononBandsPlotter| object.
Client code can use ``plotter.combiplot()`` or ``plotter.gridplot()`` to visualize the results.
"""
phbands_plotter = PhononBandsPlotter()

Expand Down Expand Up @@ -1929,8 +1918,6 @@ def anaget_raman(self, asr=2, chneut=1, ramansr=1, alphon=1, workdir=None, mpi_p
directions: list of 3D directions along which the non analytical contribution will be calculated.
If None the three cartesian direction will be used.
anaddb_kwargs: additional kwargs for anaddb.
Return: Raman object.
"""
#if not self.has_raman_terms():
# raise ValueError('The DDB file does not contain Raman terms.')
Expand Down Expand Up @@ -1997,7 +1984,7 @@ def write(self, filepath: str, filter_blocks=None) -> None:
with open(filepath, "wt") as f:
f.write("\n".join(lines))

def get_block_for_qpoint(self, qpt):
def get_block_for_qpoint(self, qpt) -> list[str]:
"""
Extracts the block data for the selected qpoint.
Returns a list of lines containing the block information
Expand All @@ -2008,7 +1995,7 @@ def get_block_for_qpoint(self, qpt):
if b['qpt'] is not None and np.allclose(b['qpt'], qpt):
return b["data"]

def replace_block_for_qpoint(self, qpt, data):
def replace_block_for_qpoint(self, qpt, data) -> bool:
"""
Change the block data for the selected qpoint. Object is modified in-place.
Data should be a list of strings representing the whole data block.
Expand All @@ -2026,7 +2013,7 @@ def replace_block_for_qpoint(self, qpt, data):

return False

def insert_block(self, data, replace=True):
def insert_block(self, data, replace=True) -> bool:
"""
Inserts a block in the list. Can replace a block if already present.
Expand All @@ -2038,8 +2025,7 @@ def insert_block(self, data, replace=True):
replace: if True and an equivalent block is already present it will be replaced,
otherwise the block will not be inserted.
Returns:
bool: True if the block was inserted.
Returns: True if the block was inserted.
"""
dord = data["dord"]
for i, b in enumerate(self.blocks):
Expand All @@ -2056,7 +2042,7 @@ def insert_block(self, data, replace=True):
self.blocks.append(data)
return True

def remove_block(self, dord, qpt=None, qpt3=None):
def remove_block(self, dord, qpt=None, qpt3=None) -> bool:
"""
Removes one block from the list of blocks in the ddb
Expand All @@ -2067,8 +2053,7 @@ def remove_block(self, dord, qpt=None, qpt3=None):
qpt3: a 3x3 matrix with the coordinates for the third order perturbations.
Should be present in dord=3.
Returns:
bool: True if a matching block was found and removed.
Returns: True if a matching block was found and removed.
"""
if dord == 2 and qpt is None:
raise ValueError("if dord==2 the qpt should be set")
Expand Down Expand Up @@ -2640,7 +2625,7 @@ def plotly(self, w_min=0, w_max=None, gamma_ev=1e-4, num=500, component='diag',

for reimf, reims in reimfs:
if isinstance(component, (list, tuple)):
label = reims % r'ε%s%s' % (SUBSCRIPT_UNICODE[str(component[0])],SUBSCRIPT_UNICODE[str(component[1])])
label = reims % r'ε%s%s' % (SUBSCRIPT_UNICODE[str(component[0])], SUBSCRIPT_UNICODE[str(component[1])])
fig.add_scatter(x=wmesh, y=reimf(t[:,component[0], component[1]]), mode='lines', showlegend=True,
name=label, row=ply_row, col=ply_col, **kwargs)
elif component == 'diag':
Expand Down
9 changes: 4 additions & 5 deletions abipy/dfpt/deformation_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,13 +2,11 @@
from __future__ import annotations

import numpy as np
import re

from pymatgen.core import Lattice, Element
from pymatgen.core import Lattice
#from abipy.core.structure import Structure
from abipy.core.symmetries import AbinitSpaceGroup
from pymatgen.symmetry.analyzer import SpacegroupAnalyzer
from abipy.abio.inputs import AbinitInput
#from pymatgen.symmetry.analyzer import SpacegroupAnalyzer


def generate_deformations_volumic(structure, eps_V=0.02, scales=None):
Expand All @@ -31,6 +29,7 @@ def generate_deformations_volumic(structure, eps_V=0.02, scales=None):

return structures_new


def generate_deformations(structure, eps=0.005):
spgrp = AbinitSpaceGroup.from_structure(structure )
print (spgrp)
Expand All @@ -48,7 +47,7 @@ def generate_deformations(structure, eps=0.005):
[0,1,0,1,1,1], [0,1,1,0,1,1], [0,1,1,1,0,1], [0,1,1,1,1,0], [1,0,1,0,1,1], [1,0,1,1,0,1],
[1,0,1,1,1,0], [1,1,0,1,0,1], [1,1,0,1,1,0], [1,1,1,0,1,0] , [0 ,0,0,0,0,0]]
if abs(rprim[1, 0]) > 1e-9 or abs(rprim[2, 0]) > 1e-9 or abs(rprim[2, 1]) > 1e-9:
print("Warning: The lattice is oriented such that xz =xy =yz =0 .")
print("Warning: The lattice is oriented such that xz = xy = yz = 0.")
rprim0 = np.copy(rprim)
a=rprim[0, :]
b=rprim[1, :]
Expand Down
Loading

0 comments on commit cf13770

Please sign in to comment.