diff --git a/src/aiidalab_qe_vibroscopy/app/widgets/euphonicmodel.py b/src/aiidalab_qe_vibroscopy/app/widgets/euphonicmodel.py index 687b838..b2f5616 100644 --- a/src/aiidalab_qe_vibroscopy/app/widgets/euphonicmodel.py +++ b/src/aiidalab_qe_vibroscopy/app/widgets/euphonicmodel.py @@ -33,12 +33,14 @@ class EuphonicBaseResultsModel(Model): # 1. single crystal data: sc # 2. powder average: pa # 3. Q planes: qp - # TRY it on the detached app. - # AAA TOBE defined with respect to the type spectra = {} path = [] q_path = None + spectrum_type = "single_crystal" + x_label = None + y_label = "Energy (meV)" + detached_app = False # Settings for single crystal and powder average q_spacing = tl.Float(0.01) @@ -46,6 +48,27 @@ class EuphonicBaseResultsModel(Model): energy_bins = tl.Int(200) temperature = tl.Float(0) weighting = tl.Unicode("coherent") + + def set_model_state(self, parameters: dict): + for k, v in parameters.items(): + setattr(self, k, v) + + def _get_default(self, trait): + if trait in ["h_vec", "k_vec"]: + return [1, 1, 1, 100, 1] + elif trait == "Q0_vec": + return [0.0, 0.0, 0.0] + return self.traits()[trait].default_value + + def get_model_state(self): + return {trait: getattr(self, trait) for trait in self.traits()} + + def reset( + self, + ): + with self.hold_trait_notifications(): + for trait in self.traits(): + setattr(self, trait, self._get_default(trait)) def fetch_data(self): """Fetch the data from the database or from the uploaded files.""" @@ -65,18 +88,20 @@ def fetch_data(self): def _inject_single_crystal_settings( self, ): - self.parameters = copy.deepcopy( - par_dict - ) # need to be different if powder or q section. + # Case in which we want to inject the model into the single crystal widget + # we define specific parameters dictionary and callback function for the single crystal case + self.parameters = copy.deepcopy(parameters_single_crystal) self._callback_spectra_generation = produce_bands_weigthed_data + # Dynamically add a trait for single crystal settings self.add_traits(custom_kpath=tl.Unicode("")) def _inject_powder_settings( self, ): - self.parameters = copy.deepcopy(par_dict_powder) + self.parameters = copy.deepcopy(parameters_powder) self._callback_spectra_generation = produce_powder_data + # Dynamically add a trait for powder settings self.add_traits(q_min=tl.Float(0.0)) self.add_traits(q_max=tl.Float(1)) @@ -96,85 +121,77 @@ def _inject_qsection_settings( k_vec=tl.List(trait=tl.Float(), default_value=[1, 1, 1, 100, 1]) ) - def set_model_state(self, parameters: dict): - for k, v in parameters.items(): - setattr(self, k, v) - - def _get_default(self, trait): - if trait in ["h_vec", "k_vec"]: - return [1, 1, 1, 100, 1] - elif trait == "Q0_vec": - return [0.0, 0.0, 0.0] - return self.traits()[trait].default_value - - def get_model_state(self): - return {trait: getattr(self, trait) for trait in self.traits()} - - def reset( - self, - ): - with self.hold_trait_notifications(): - for trait in self.traits(): - setattr(self, trait, self._get_default(trait)) - - def _update_spectra( + def get_spectra( self, ): # This is used to update the spectra when the parameters are changed # and the if not hasattr(self, "parameters"): self._inject_single_crystal_settings() + + if self.spectrum_type == "q_planes": + self._get_qsection_spectra() + return self.parameters.update(self.get_model_state()) - parameters_ = AttrDict(self.parameters) # custom linear path custom_kpath = self.custom_kpath if hasattr(self, "custom_kpath") else "" if len(custom_kpath) > 1: - coordinates, labels = self.curate_path_and_labels() + coordinates, labels = self._curate_path_and_labels() qpath = { "coordinates": coordinates, "labels": labels, # ["$\Gamma$","X","X","(1,1,1)"], - "delta_q": parameters_["q_spacing"], + "delta_q": self.parameters["q_spacing"], } else: qpath = copy.deepcopy(self.q_path) if qpath: - qpath["delta_q"] = parameters_["q_spacing"] + qpath["delta_q"] = self.parameters["q_spacing"] spectra, parameters = self._callback_spectra_generation( - params=parameters_, + params=AttrDict(self.parameters), fc=self.fc, linear_path=qpath, plot=False, ) # curated spectra (labels and so on...) - if hasattr(self, "custom_kpath"): # single crystal case + if spectrum_type == "single_crystal": # single crystal case self.x, self.y = np.meshgrid( spectra[0].x_data.magnitude, spectra[0].y_data.magnitude ) ( - self.final_xspectra, - self.final_zspectra, + final_xspectra, + final_zspectra, self.ticks_positions, self.ticks_labels, ) = generated_curated_data(spectra) - else: + + self.z = final_zspectra.T + self.y = self.y[:,0] + self.x = None # we have the ticks positions and labels + + self.xlabel = "" + self.ylabel = "Energy (meV)" + + elif spectrum_type == "powder": # powder case # Spectrum2D as output of the powder data self.x, self.y = np.meshgrid( spectra.x_data.magnitude, spectra.y_data.magnitude ) - # we don't need to curate the powder data, - # we can directly use them: - self.final_xspectra = spectra.x_data.magnitude - self.final_zspectra = spectra.z_data.magnitude + # we don't need to curate the powder data, at variance with the single crystal case. + # We can directly use them: + self.x = spectra.x_data.magnitude[0] + self.y = self.y[:,0] + self.z = spectra.z_data.magnitude.T - def _update_qsection_spectra( + def _get_qsection_spectra( self, ): - parameters_ = AttrDict( + # This is used to update the spectra in the case we plot the Q planes (the third tab). + parameters_qplanes = AttrDict( { "h": np.array([i for i in self.h_vec[:-2]]), "k": np.array([i for i in self.k_vec[:-2]]), @@ -193,32 +210,32 @@ def _update_qsection_spectra( modes, q_array, h_array, k_array, labels, dw = produce_Q_section_modes( self.fc, - h=parameters_.h, - k=parameters_.k, - Q0=parameters_.Q0, - n_h=parameters_.n_h, - n_k=parameters_.n_k, - h_extension=parameters_.h_extension, - k_extension=parameters_.k_extension, - temperature=parameters_.temperature, + h=parameters_qplanes.h, + k=parameters_qplanes.k, + Q0=parameters_qplanes.Q0, + n_h=parameters_qplanes.n_h, + n_k=parameters_qplanes.n_k, + h_extension=parameters_qplanes.h_extension, + k_extension=parameters_qplanes.k_extension, + temperature=parameters_qplanes.temperature, ) - self.av_spec, self.q_array, self.h_array, self.k_array, self.labels = ( + self.av_spec, self.z, self.x, self.y, self.labels = ( produce_Q_section_spectrum( modes, q_array, h_array, k_array, - ecenter=parameters_.ecenter, - deltaE=parameters_.deltaE, - bins=parameters_.bins, - spectrum_type=parameters_.spectrum_type, + ecenter=parameters_qplanes.ecenter, + deltaE=parameters_qplanes.deltaE, + bins=parameters_qplanes.bins, + spectrum_type=parameters_qplanes.spectrum_type, dw=dw, labels=labels, ) ) - def curate_path_and_labels( + def _curate_path_and_labels( self, ): # This is used to curate the path and labels of the spectra if custom kpath is provided. @@ -240,9 +257,6 @@ def curate_path_and_labels( coordinates.append(scoords) return coordinates, labels - def _clone(self): - return copy.deepcopy(self) - def produce_phonopy_files(self): # This is used to produce the phonopy files from # PhonopyCalculation data. The files are phonopy.yaml and force_constants.hdf5 @@ -250,3 +264,12 @@ def produce_phonopy_files(self): self.node.phonon_bands.creator, mode="download" ) return phonopy_yaml, fc_hdf5 + + def prepare_data_for_download(self): + raise NotImplementedError("Need to implement the download of a CSV file") + + def _clone(self): + # in case we want to clone the model. + # This is the case when we have the same data and we inject in three + # different models: we don't need to fetch three times. + return copy.deepcopy(self) \ No newline at end of file diff --git a/src/aiidalab_qe_vibroscopy/utils/euphonic/data_manipulation/bands_pdos.py b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/bands_pdos.py similarity index 98% rename from src/aiidalab_qe_vibroscopy/utils/euphonic/data_manipulation/bands_pdos.py rename to src/aiidalab_qe_vibroscopy/utils/euphonic/data/bands_pdos.py index deaafce..58a0e80 100644 --- a/src/aiidalab_qe_vibroscopy/utils/euphonic/data_manipulation/bands_pdos.py +++ b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/bands_pdos.py @@ -21,7 +21,7 @@ (from euphonic, using the force constants instances as obtained from phonopy.yaml). These are then used in the widgets to plot the corresponding quantities. -NB: no more used in the QE-app. We use phonopy instead. +NB: not used in the QE-app. We use phonopy instead. """ ######################## diff --git a/src/aiidalab_qe_vibroscopy/utils/euphonic/data/export_vibronic_to_euphonic.py b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/export_vibronic_to_euphonic.py new file mode 100644 index 0000000..7254e93 --- /dev/null +++ b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/export_vibronic_to_euphonic.py @@ -0,0 +1,41 @@ +from aiida.orm import Dict +from aiida_phonopy.common.raw_parsers import get_force_constants_from_phonopy +from aiida_phonopy.workflows.phonon import generate_force_constant_instance + +def export_euphonic_data(output_vibronic, fermi_energy=None): + if "phonon_bands" not in output_vibronic: + return None + + output_set = output_vibronic.phonon_bands + + if any(not element for element in output_set.creator.caller.inputs.structure.pbc): + vibro_bands = output_set.creator.caller.inputs.phonopy_bands_dict.get_dict() + # Group the band and band_labels + band = vibro_bands["band"] + band_labels = vibro_bands["band_labels"] + + grouped_bands = [ + item + for sublist in [band_labels[i : i + 2] for i in range(len(band_labels) - 1)] + for item in sublist + ] + grouped_q = [ + [tuple(band[i : i + 3]), tuple(band[i + 3 : i + 6])] + for i in range(0, len(band) - 3, 3) + ] + q_path = { + "coordinates": grouped_q, + "labels": grouped_bands, + "delta_q": 0.01, # 1/A + } + else: + q_path = None + + phonopy_calc = output_set.creator + fc = generate_force_constant_instance(phonopy_calc) + # bands = compute_bands(fc) + # pdos = compute_pdos(fc) + return { + "fc": fc, + "q_path": q_path, + } # "bands": bands, "pdos": pdos, "thermal": None} diff --git a/src/aiidalab_qe_vibroscopy/utils/euphonic/data/parameters.py b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/parameters.py new file mode 100644 index 0000000..4db6364 --- /dev/null +++ b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/parameters.py @@ -0,0 +1,52 @@ +"""Set of parameters for given Euphonic calculation. + +We distinguish between parameters for a single crystal and for a powder calculation: the former requires a path in reciprocal space, +while the latter requires a range of q-points. +We have a set of common parameters that are shared between the two types of calculations. +""" +common_parameters = { + "weighting": "coherent", # Spectral weighting to plot: DOS, coherent inelastic neutron scattering (default: dos) + "grid": None, # FWHM of broadening on q axis in 1/LENGTH_UNIT (no broadening if unspecified). (default: None) + "grid_spacing": 0.1, # q-point spacing of Monkhorst-Pack grid. (default: 0.1) + "energy_unit": "THz", + "temperature": 0, # Temperature in K; enable Debye-Waller factor calculation. (Only applicable when --weighting=coherent). (default: None) + "shape": "gauss", # The broadening shape (default: gauss) + "length_unit": "angstrom", + "q_spacing": 0.01, # Target distance between q-point samples in 1/LENGTH_UNIT (default: 0.025) + "energy_broadening": 1, + "q_broadening": None, # FWHM of broadening on q axis in 1/LENGTH_UNIT (no broadening if unspecified). (default: None) + "ebins": 200, # Number of energy bins (default: 200) + "e_min": 0, + "e_max": None, + "title": None, + "ylabel": "THz", + "xlabel": "", + "save_json": False, + "no_base_style": False, + "style": False, + "vmin": None, + "vmax": None, + "save_to": None, + "asr": None, # Apply an acoustic-sum-rule (ASR) correction to the data: "realspace" applies the correction to the force constant matrix in real space. "reciprocal" applies the correction to the dynamical matrix at each q-point. (default: None) + "dipole_parameter": 1.0, # Set the cutoff in real/reciprocal space for the dipole Ewald sum; higher values use more reciprocal terms. If tuned correctly this can result in performance improvements. See euphonic-optimise-dipole-parameter program for help on choosing a good DIPOLE_PARAMETER. (default: 1.0) + "use_c": None, + "n_threads": None, +} + +parameters_single_crystal = { + **common_parameters, +} + +parameters_powder = { + **common_parameters, + "q_min": 0, + "q_max": 1, + "npts": 150, + "npts_density": None, + "pdos": None, + "e_i": None, + "sampling": "golden", + "jitter": True, + "e_f": None, + "disable_widgets": True, +} diff --git a/src/aiidalab_qe_vibroscopy/utils/euphonic/data/phonopy_interface.py b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/phonopy_interface.py new file mode 100644 index 0000000..048257f --- /dev/null +++ b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/phonopy_interface.py @@ -0,0 +1,118 @@ +import tempfile +import pathlib +import base64 +from typing import Optional +import euphonic +from euphonic.io.phonopy import write_force_constants_to_hdf5 + +def generate_force_constant_instance( + phonopy_calc=None, + path: str = None, + summary_name: str = None, + born_name: Optional[str] = None, + fc_name: str = "FORCE_CONSTANTS", + fc_format: Optional[str] = None, + mode="stream", # "download" to have the download of phonopy.yaml and fc.hdf5 . TOBE IMPLEMENTED. +): + """ + Basically allows to obtain the ForceConstants instance from phonopy, both via files (from the second + input parameters we have the same one of `euphonic.ForceConstants.from_phonopy`), or via a + PhonopyCalculation instance. Respectively, the two ways will support independent euphonic app and integration + of Euphonic into aiidalab. + """ + blockPrint() + + ####### This is done to support the detached app (from aiidalab) with the same code: + if path and summary_name: + fc = euphonic.ForceConstants.from_phonopy( + path=path, + summary_name=summary_name, + fc_name=fc_name, + ) + return fc + elif not phonopy_calc: + raise NotImplementedError( + "Please provide or the files or the phonopy calculation node." + ) + + ####### This is almost copied from PhonopyCalculation and is done to support functionalities in aiidalab env: + from phonopy.interface.phonopy_yaml import PhonopyYaml + + kwargs = {} + + if "settings" in phonopy_calc.inputs: + the_settings = phonopy_calc.inputs.settings.get_dict() + for key in ["symmetrize_nac", "factor_nac", "subtract_residual_forces"]: + if key in the_settings: + kwargs.update({key: the_settings[key]}) + + if "phonopy_data" in phonopy_calc.inputs: + ph = phonopy_calc.inputs.phonopy_data.get_phonopy_instance(**kwargs) + p2s_map = phonopy_calc.inputs.phonopy_data.get_cells_mappings()["primitive"][ + "p2s_map" + ] + ph.produce_force_constants() + elif "force_constants" in phonopy_calc.inputs: + ph = phonopy_calc.inputs.force_constants.get_phonopy_instance(**kwargs) + p2s_map = phonopy_calc.inputs.force_constants.get_cells_mappings()["primitive"][ + "p2s_map" + ] + ph.force_constants = phonopy_calc.inputs.force_constants.get_array( + "force_constants" + ) + + ####### + + # Create temporary directory + # + with tempfile.TemporaryDirectory() as dirpath: + # phonopy.yaml generation: + phpy_yaml = PhonopyYaml() + phpy_yaml.set_phonon_info(ph) + phpy_yaml_txt = str(phpy_yaml) + + with open( + pathlib.Path(dirpath) / "phonopy.yaml", "w", encoding="utf8" + ) as handle: + handle.write(phpy_yaml_txt) + + # Force constants hdf5 file generation: + # all this is needed to load the euphonic instance, in case no FC are written in phonopy.yaml + # which is the case + + write_force_constants_to_hdf5( + force_constants=ph.force_constants, + filename=pathlib.Path(dirpath) / "fc.hdf5", + p2s_map=p2s_map, + ) + + # Here below we trigger the download mode. Can be improved avoiding the repetitions of lines + if mode == "download": + with open( + pathlib.Path(dirpath) / "phonopy.yaml", "r", encoding="utf8" + ) as handle: + file_content = handle.read() + phonopy_yaml_bitstream = base64.b64encode(file_content.encode()).decode( + "utf-8" + ) + + with open( + pathlib.Path(dirpath) / "fc.hdf5", + "rb", + ) as handle: + file_content = handle.read() + fc_hdf5_bitstream = base64.b64encode(file_content).decode() + + return phonopy_yaml_bitstream, fc_hdf5_bitstream + + # Read force constants (fc.hdf5) and summary+NAC (phonopy.yaml) + + fc = euphonic.ForceConstants.from_phonopy( + path=dirpath, + summary_name="phonopy.yaml", + fc_name="fc.hdf5", + ) + # print(filename) + # print(dirpath) + enablePrint() + return fc \ No newline at end of file diff --git a/src/aiidalab_qe_vibroscopy/utils/euphonic/data_manipulation/intensity_maps.py b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/structure_factors.py similarity index 71% rename from src/aiidalab_qe_vibroscopy/utils/euphonic/data_manipulation/intensity_maps.py rename to src/aiidalab_qe_vibroscopy/utils/euphonic/data/structure_factors.py index 66095ca..edb65d9 100644 --- a/src/aiidalab_qe_vibroscopy/utils/euphonic/data_manipulation/intensity_maps.py +++ b/src/aiidalab_qe_vibroscopy/utils/euphonic/data/structure_factors.py @@ -85,6 +85,8 @@ def __init__(self, *args, **kwargs): and the powder maps(from euphonic, using the force constants instances as obtained from phonopy.yaml). These are then used in the widgets to plot the corresponding quantities. +NOTE: the two main functions here are produce_bands_weigthed_data and produce_powder_data. + PLEASE NOTE: the scattering lengths are tabulated (Euphonic/euphonic/data/sears-1992.json) and are from Sears (1992) Neutron News 3(3) pp26--37. """ @@ -174,38 +176,6 @@ def join_q_paths(coordinates: list, labels: list, delta_q=0.1, G=[0, 0, 0]): ################################ START INTENSITY PLOT GENERATOR ######################## -par_dict = { - "weighting": "coherent", # Spectral weighting to plot: DOS, coherent inelastic neutron scattering (default: dos) - "grid": None, # FWHM of broadening on q axis in 1/LENGTH_UNIT (no broadening if unspecified). (default: None) - "grid_spacing": 0.1, # q-point spacing of Monkhorst-Pack grid. (default: 0.1) - "energy_unit": "THz", - "temperature": 0, # Temperature in K; enable Debye-Waller factor calculation. (Only applicable when --weighting=coherent). (default: None) - #'btol':, - "shape": "gauss", # The broadening shape (default: gauss) - "length_unit": "angstrom", - "q_spacing": 0.01, # Target distance between q-point samples in 1/LENGTH_UNIT (default: 0.025) - "energy_broadening": 1, - "q_broadening": None, # FWHM of broadening on q axis in 1/LENGTH_UNIT (no broadening if unspecified). (default: None) - "ebins": 200, # Number of energy bins (default: 200) - "e_min": 0, - "e_max": None, - "title": None, - "ylabel": "THz", - "xlabel": "", - "save_json": False, - "no_base_style": False, - "style": False, - "vmin": None, - "vmax": None, - "save_to": None, - "asr": None, # Apply an acoustic-sum-rule (ASR) correction to the data: "realspace" applies the correction to the force constant matrix in real space. "reciprocal" applies the correction to the dynamical matrix at each q-point. (default: None) - "dipole_parameter": 1.0, # Set the cutoff in real/reciprocal space for the dipole Ewald sum; higher values use more reciprocal terms. If tuned correctly this can result in performance improvements. See euphonic-optimise-dipole-parameter program for help on choosing a good DIPOLE_PARAMETER. (default: 1.0) - "use_c": None, - "n_threads": None, -} - -parameters = par_dict - def produce_bands_weigthed_data( params: Optional[List[str]] = parameters, @@ -385,48 +355,8 @@ def produce_bands_weigthed_data( ################################ START POWDER ######################## -par_dict_powder = { - "weighting": "coherent", # Spectral weighting to plot: DOS, coherent inelastic neutron scattering (default: dos) - "grid": None, # FWHM of broadening on q axis in 1/LENGTH_UNIT (no broadening if unspecified). (default: None) - "grid_spacing": 0.1, # q-point spacing of Monkhorst-Pack grid. (default: 0.1) - "q_min": 0, - "q_max": 1, - "temperature": None, # Temperature in K; enable Debye-Waller factor calculation. (Only applicable when --weighting=coherent). (default: None) - "ebins": 200, # Number of energy bins (default: 200) - "q_spacing": 0.01, # Target distance between q-point samples in 1/LENGTH_UNIT (default: 0.025) - "energy_broadening": 1, - "npts": 150, - #'grid':, - "energy_unit": "THz", - #'btol':, - "shape": "gauss", # The broadening shape (default: gauss) - "length_unit": "angstrom", - "q_broadening": None, # FWHM of broadening on q axis in 1/LENGTH_UNIT (no broadening if unspecified). (default: None) - "e_min": 0, - "e_max": None, - "title": None, - "ylabel": "THz", - "xlabel": "", - "save_json": False, - "no_base_style": False, - "style": False, - "vmin": None, - "vmax": None, - "save_to": None, - "asr": None, - "dipole_parameter": 1.0, - "use_c": None, - "n_threads": None, - "npts_density": None, - "pdos": None, - "e_i": None, - "sampling": "golden", - "jitter": True, - "e_f": None, - "disable_widgets": True, -} - -parameters_powder = AttrDict(par_dict_powder) + +#parameters_powder = AttrDict(par_dict_powder) def produce_powder_data( @@ -436,6 +366,8 @@ def produce_powder_data( linear_path=None, ) -> None: blockPrint() + """Read the description of the produce_bands_weigthed_data function for more details. + """ # args = get_args(get_parser(), params) if not params: @@ -648,158 +580,6 @@ def update_max(max_val): matplotlib_save_or_show(save_filename=args.save_to) -def generate_force_constant_instance( - phonopy_calc=None, - path: str = None, - summary_name: str = None, - born_name: Optional[str] = None, - fc_name: str = "FORCE_CONSTANTS", - fc_format: Optional[str] = None, - mode="stream", # "download" to have the download of phonopy.yaml and fc.hdf5 . TOBE IMPLEMENTED. -): - """ - Basically allows to obtain the ForceConstants instance from phonopy, both via files (from the second - input parameters we have the same one of `euphonic.ForceConstants.from_phonopy`), or via a - PhonopyCalculation instance. Respectively, the two ways will support independent euphonic app and integration - of Euphonic into aiidalab. - """ - blockPrint() - - ####### This is done to support the detached app (from aiidalab) with the same code: - if path and summary_name: - fc = euphonic.ForceConstants.from_phonopy( - path=path, - summary_name=summary_name, - fc_name=fc_name, - ) - return fc - elif not phonopy_calc: - raise NotImplementedError( - "Please provide or the files or the phonopy calculation node." - ) - - ####### This is almost copied from PhonopyCalculation and is done to support functionalities in aiidalab env: - from phonopy.interface.phonopy_yaml import PhonopyYaml - - kwargs = {} - - if "settings" in phonopy_calc.inputs: - the_settings = phonopy_calc.inputs.settings.get_dict() - for key in ["symmetrize_nac", "factor_nac", "subtract_residual_forces"]: - if key in the_settings: - kwargs.update({key: the_settings[key]}) - - if "phonopy_data" in phonopy_calc.inputs: - ph = phonopy_calc.inputs.phonopy_data.get_phonopy_instance(**kwargs) - p2s_map = phonopy_calc.inputs.phonopy_data.get_cells_mappings()["primitive"][ - "p2s_map" - ] - ph.produce_force_constants() - elif "force_constants" in phonopy_calc.inputs: - ph = phonopy_calc.inputs.force_constants.get_phonopy_instance(**kwargs) - p2s_map = phonopy_calc.inputs.force_constants.get_cells_mappings()["primitive"][ - "p2s_map" - ] - ph.force_constants = phonopy_calc.inputs.force_constants.get_array( - "force_constants" - ) - - ####### - - # Create temporary directory - # - with tempfile.TemporaryDirectory() as dirpath: - # phonopy.yaml generation: - phpy_yaml = PhonopyYaml() - phpy_yaml.set_phonon_info(ph) - phpy_yaml_txt = str(phpy_yaml) - - with open( - pathlib.Path(dirpath) / "phonopy.yaml", "w", encoding="utf8" - ) as handle: - handle.write(phpy_yaml_txt) - - # Force constants hdf5 file generation: - # all this is needed to load the euphonic instance, in case no FC are written in phonopy.yaml - # which is the case - - write_force_constants_to_hdf5( - force_constants=ph.force_constants, - filename=pathlib.Path(dirpath) / "fc.hdf5", - p2s_map=p2s_map, - ) - - # Here below we trigger the download mode. Can be improved avoiding the repetitions of lines - if mode == "download": - with open( - pathlib.Path(dirpath) / "phonopy.yaml", "r", encoding="utf8" - ) as handle: - file_content = handle.read() - phonopy_yaml_bitstream = base64.b64encode(file_content.encode()).decode( - "utf-8" - ) - - with open( - pathlib.Path(dirpath) / "fc.hdf5", - "rb", - ) as handle: - file_content = handle.read() - fc_hdf5_bitstream = base64.b64encode(file_content).decode() - - return phonopy_yaml_bitstream, fc_hdf5_bitstream - - # Read force constants (fc.hdf5) and summary+NAC (phonopy.yaml) - - fc = euphonic.ForceConstants.from_phonopy( - path=dirpath, - summary_name="phonopy.yaml", - fc_name="fc.hdf5", - ) - # print(filename) - # print(dirpath) - enablePrint() - return fc - - -def export_euphonic_data(output_vibronic, fermi_energy=None): - if "phonon_bands" not in output_vibronic: - return None - - output_set = output_vibronic.phonon_bands - - if any(not element for element in output_set.creator.caller.inputs.structure.pbc): - vibro_bands = output_set.creator.caller.inputs.phonopy_bands_dict.get_dict() - # Group the band and band_labels - band = vibro_bands["band"] - band_labels = vibro_bands["band_labels"] - - grouped_bands = [ - item - for sublist in [band_labels[i : i + 2] for i in range(len(band_labels) - 1)] - for item in sublist - ] - grouped_q = [ - [tuple(band[i : i + 3]), tuple(band[i + 3 : i + 6])] - for i in range(0, len(band) - 3, 3) - ] - q_path = { - "coordinates": grouped_q, - "labels": grouped_bands, - "delta_q": 0.01, # 1/A - } - else: - q_path = None - - phonopy_calc = output_set.creator - fc = generate_force_constant_instance(phonopy_calc) - # bands = compute_bands(fc) - # pdos = compute_pdos(fc) - return { - "fc": fc, - "q_path": q_path, - } # "bands": bands, "pdos": pdos, "thermal": None} - - def generated_curated_data(spectra): # here we concatenate the bands groups and create the ticks and labels. diff --git a/src/aiidalab_qe_vibroscopy/utils/euphonic/structurefactorwidget.py b/src/aiidalab_qe_vibroscopy/utils/euphonic/structurefactorwidget.py new file mode 100644 index 0000000..3d9c9bf --- /dev/null +++ b/src/aiidalab_qe_vibroscopy/utils/euphonic/structurefactorwidget.py @@ -0,0 +1,458 @@ +class EuphonicStructureFactorWidget(ipw.VBox): + """Description. + + Collects all the button and widget used to define settings for Neutron dynamic structure factor, + in all the three cases: single crystal, powder, and q-section.... + """ + + def __init__(self, model, spectrum_type = "single_crystal", detached_app = False, **kwargs): + super().__init__() + node + self._model = model + self._model.spectrum_type = spectrum_type + self._model.detached_app = detached_app + self.rendered = False + + def render(self): + """Render the widget. + + This means render the plot button. + """ + if self.rendered: + return + + self.tab_widget = ipw.Tab() + self.tab_widget.layout.display = "none" + self.tab_widget.set_title(0, "Single crystal") + self.tab_widget.set_title(1, "Powder sample") + self.tab_widget.set_title(2, "Q-plane view") + self.tab_widget.children = () + + self.plot_button = ipw.Button( + description="Initialise INS data", + icon="pencil", + button_style="primary", + disabled=True, + layout=ipw.Layout(width="auto"), + ) + self.plot_button.on_click(self._render_for_real) + + self.loading_widget = LoadingWidget("Loading INS data") + self.loading_widget.layout.display = "none" + + if not self._model.detached_app: + self.plot_button.disabled = False + else: + self.upload_widget = UploadPhonopyWidget() + self.upload_widget.reset_uploads.on_click(self._on_reset_uploads_button_clicked) + self.upload_widget.children[0].observe(self._on_upload_yaml, "_counter") + self.upload_widget.children[1].observe(self._on_upload_hdf5, "_counter") + self.children += (self.upload_widget,) + + self.download_widget = DownloadYamlHdf5Widget(model=self._model) + self.download_widget.layout.display = "none" + + self.children += ( + self.plot_button, + self.loading_widget, + self.tab_widget, + self.download_widget, + ) + + # NOTE: we initialise here the figure widget, but we do not plot anything yet. + # this is useful to call the init_view method, which contains the update for the figure. + self.fig = go.FigureWidget() + + self.rendered = True + + def _render_for_real(self, change=None): + + self.plot_button.layout.display = "none" + self.loading_widget.layout.display = "block" + + self._init_view() + + slider_intensity = ipw.FloatRangeSlider( + value=[1, 100], # Default selected interval + min=1, + max=100, + step=1, + orientation="horizontal", + readout=True, + readout_format=".0f", + layout=ipw.Layout( + width="400px", + ), + ) + slider_intensity.observe(self._update_intensity_filter, "value") + specification_intensity = ipw.HTML( + "(Intensity is relative to the maximum intensity at T=0K)" + ) + + E_units_button = ipw.ToggleButtons( + options=[ + ("meV", "meV"), + ("THz", "THz"), + # ("cm-1", "cm-1"), + ], + value="meV", + description="Energy units:", + disabled=False, + layout=ipw.Layout( + width="auto", + ), + ) + E_units_button.observe(self._update_energy_units, "value") + # MAYBE WE LINK ALSO THIS TO THE MODEL? so we can download the data with the preferred units. + + q_spacing = ipw.FloatText( + value=self._model.q_spacing, + step=0.001, + description="q step (1/A)", + tooltip="q spacing in 1/A", + ) + ipw.link( + (self._model, "q_spacing"), + (self.q_spacing, "value"), + ) + q_spacing.observe(self._on_setting_change, names="value") + + energy_broadening = ipw.FloatText( + value=self._model.energy_broadening, + step=0.01, + description="ΔE (meV)", + tooltip="Energy broadening in meV", + ) + ipw.link( + (self._model, "energy_broadening"), + (energy_broadening, "value"), + ) + energy_broadening.observe(self._on_setting_change, names="value") + + energy_bins = ipw.IntText( + value=self._model.energy_bins, + description="#E bins", + tooltip="Number of energy bins", + ) + ipw.link( + (self._model, "energy_bins"), + (energy_bins, "value"), + ) + energy_bins.observe(self._on_setting_change, names="value") + + temperature = ipw.FloatText( + value=self._model.temperature, + step=0.01, + description="T (K)", + disabled=False, + ) + ipw.link( + (self._model, "temperature"), + (temperature, "value"), + ) + temperature.observe(self._on_setting_change, names="value") + + weight_button = ipw.ToggleButtons( + options=[ + ("Coherent", "coherent"), + ("DOS", "dos"), + ], + value=self._model.weighting, + description="weight:", + disabled=False, + style={"description_width": "initial"}, + ) + ipw.link( + (self._model, "weighting"), + (weight_button, "value"), + ) + weight_button.observe(self._on_weight_button_change, names="value") + + plot_button = ipw.Button( + description="Replot", + icon="pencil", + button_style="primary", + disabled=True, + layout=ipw.Layout(width="auto"), + ) + plot_button.observe(self._on_plot_button_change, names="disabled") + + reset_button = ipw.Button( + description="Reset", + icon="recycle", + button_style="primary", + disabled=False, + layout=ipw.Layout(width="auto"), + ) + reset_button.on_click(self._reset_settings) + + download_button = ipw.Button( + description="Download Data and Plot", + icon="download", + button_style="primary", + disabled=False, # Large files... + layout=ipw.Layout(width="auto"), + ) + download_button.on_click(self._download_data) + + if self._model.spectrum_type == "single_crystal": + self.custom_kpath_description = ipw.HTML( + """ +
+ Custom q-points path for the structure factor:
+ we can provide it via a specific format:
+ (1) each linear path should be divided by '|';
+ (2) each path is composed of 'qxi qyi qzi - qxf qyf qzf' where qxi and qxf are, respectively, + the start and end x-coordinate of the q direction, in reciprocal lattice units (rlu).
+ An example path is: '0 0 0 - 1 1 1 | 1 1 1 - 0.5 0.5 0.5'.
+ For now, we do not support fractions (i.e. we accept 0.5 but not 1/2). +
+ """ + ) + + self.custom_kpath_text = ipw.Text( + value="", + description="Custom path (rlu):", + style={"description_width": "initial"}, + ) + custom_style = '' + display(ipw.HTML(custom_style)) + self.custom_kpath_text.add_class("custom-font") + ipw.link( + (self._model, "custom_kpath"), + (self.custom_kpath_text, "value"), + ) + self.custom_kpath_text.observe(self._on_setting_changed, names="value") + # fi self._model.spectrum_type == "single_crystal" + elif self._model.spectrum_type == "powder": + self.qmin = ipw.FloatText( + value=0, + description="|q|min (1/A)", + ) + ipw.link( + (self._model, "q_min"), + (self.qmin, "value"), + ) + self.qmin.observe(self._on_setting_changed, names="value") + + self.qmax = ipw.FloatText( + step=0.01, + value=1, + description="|q|max (1/A)", + ) + ipw.link( + (self._model, "q_max"), + (self.qmax, "value"), + ) + self.qmax.observe(self._on_setting_changed, names="value") + + self.int_npts = ipw.IntText( + value=100, + description="npts", + tooltip="Number of points to be used in the average sphere.", + ) + ipw.link( + (self._model, "npts"), + (self.int_npts, "value"), + ) + self.int_npts.observe(self._on_setting_changed, names="value") + # fi self._model.spectrum_type == "powder" + elif self._model.spectrum_type == "q_planes": + self.ecenter = ipw.FloatText( + value=0, + description="E (meV)", + ) + ipw.link( + (self._model, "center_e"), + (self.ecenter, "value"), + ) + self.ecenter.observe(self._on_setting_changed, names="value") + + self.plane_description_widget = ipw.HTML( + """ +
+ Q-plane definition:
+ To define a plane in the reciprocal space,
+ you should define a point in the reciprocal space, Q0, + and two vectors h⃗ and k⃗. Then, each Q point is defined as: Q = Q0 + α*h⃗ + β*k⃗.
+ Then you can select the number of q points in both directions and the α and β parameters.
+ Coordinates are reciprocal lattice units (rlu). +
+ """ + ) + + self.Q0_vec = ipw.HBox( + [ipw.FloatText(value=0, layout={"width": "60px"}) for j in range(3)] + + [ + ipw.HTML( + "Nhq, Nkq ↓", + layout={"width": "60px"}, + ), + ipw.HTML(r"α, β ↓", layout={"width": "60px"}), + ] + ) + + self.h_vec = ipw.HBox( + [ + ipw.FloatText(value=1, layout={"width": "60px"}) # coordinates + for j in range(3) + ] + + [ + ipw.IntText(value=100, layout={"width": "60px"}), + ipw.IntText(value=1, layout={"width": "60px"}), + ] # number of points along this dir, i.e. n_h; and multiplicative factor alpha + ) + self.k_vec = ipw.HBox( + [ipw.FloatText(value=1, layout={"width": "60px"}) for j in range(3)] + + [ + ipw.IntText(value=100, layout={"width": "60px"}), + ipw.IntText(value=1, layout={"width": "60px"}), + ] + ) + + for vec in [self.Q0_vec, self.h_vec, self.k_vec]: + for child in vec.children: + child.observe(self._on_setting_changed, names="value") + child.observe(self._on_vector_changed, names="value") + + self.Q0_widget = ipw.HBox( + [ipw.HTML("Q0: ", layout={"width": "20px"}), self.Q0_vec] + ) + self.h_widget = ipw.HBox( + [ipw.HTML("h: ", layout={"width": "20px"}), self.h_vec] + ) + self.k_widget = ipw.HBox( + [ipw.HTML("k: ", layout={"width": "20px"}), self.k_vec] + ) + + self.energy_broadening = ipw.FloatText( + value=0.5, + description="ΔE (meV)", + tooltip="Energy window in eV", + ) + ipw.link( + (self._model, "energy_broadening"), + (self.energy_broadening, "value"), + ) + self.energy_broadening.observe(self._on_setting_changed, names="value") + + self.plot_button.disabled = False + self.plot_button.description = "Plot" + # self.reset_button.disabled = True + self.download_button.disabled = True + # fi self._model.spectrum_type == "q_planes" + + self.children += ( + ... + ) + + def _init_view(self, _=None): + self._model.fetch_data() + self._update_plot() + + def _on_plot_button_change(self, change): + self.download_button.disabled = not change["new"] + + def _on_weight_button_change(self, change): + self._model.temperature = 0 + self.temperature.disabled = True if change["new"] == "dos" else False + self.plot_button.disabled = False + + def _on_setting_change( + self, change + ): # think if we want to do something more evident... + self.plot_button.disabled = False + + def _update_plot(self): + # update the spectra, i.e. the data contained in the _model. + # TODO: we need to treat differently the update of intensity and units. + # they anyway need to modify the data, but no additional spectra re-generation is really needed. + # so the update_spectra need some more logic, or we call another method. + self._model.get_spectra() + + if not self.rendered: + # First time we render, we set several layout settings. + # Layout settings + self.fig["layout"]["xaxis"].update( + title=self._model.xlabel, + range=[min(self._model.x), max(self._model.x)], + ) + self.fig["layout"]["yaxis"].update( + title=self._model.ylabel, + range=[min(self._model.y), max(self._model.y)], + ) + + if self.fig.layout.images: + for image in self.fig.layout.images: + image["scl"] = 2 # Set the scale for each image + + self.fig.update_layout( + height=500, + width=700, + margin=dict(l=15, r=15, t=15, b=15), + ) + # Update x-axis and y-axis to enable autoscaling + self.fig.update_xaxes(autorange=True) + self.fig.update_yaxes(autorange=True) + + # Update the layout to enable autoscaling + self.fig.update_layout(autosize=True) + + + heatmap_trace = go.Heatmap( + z=self._model.z, + y=(self._model.y), + x=self._model.x, + colorbar=COLORBAR_DICT, + colorscale=COLORSCALE, # imported from euphonic_base_widgets + ) + + # change the path wants also a change in the labels + if "ticks_positions" in self._model and "ticks_labels" in self._model: + self.fig.update_layout( + xaxis=dict( + tickmode="array", + tickvals=self._model.ticks_positions, + ticktext=self._model.ticks_labels, + ) + ) + + # Add colorbar + colorbar = heatmap_trace.colorbar + colorbar.x = 1.05 # Move colorbar to the right + colorbar.y = 0.5 # Center colorbar vertically + + # Add heatmap trace to figure + self.fig.add_trace(heatmap_trace) + self.fig.data = [self.fig.data[1]] + + def _reset_settings(self, _): + self._model.reset() + + def _download_data(self, _=None): + data, filename = self._model.prepare_data_for_download() + self._download(data, filename) + + @staticmethod + def _download(payload, filename): + from IPython.display import Javascript + + javas = Javascript( + """ + var link = document.createElement('a'); + link.href = 'data:text/json;charset=utf-8;base64,{payload}' + link.download = "{filename}" + document.body.appendChild(link); + link.click(); + document.body.removeChild(link); + """.format(payload=payload, filename=filename) + ) + display(javas) + + def _on_vector_changed(self, change=None): + """ + Update the model. + """ + self._model.Q0_vec = [i.value for i in self.Q0_vec.children[:-2]] + self._model.h_vec = [i.value for i in self.h_vec.children] + self._model.k_vec = [i.value for i in self.k_vec.children]