diff --git a/CHANGELOG.md b/CHANGELOG.md index da45ded..2627a1f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -12,6 +12,32 @@ The versioning scheme followed is `YYYY.MAJOR.MINOR`: - `MAJOR`: major feature (= change in API, no backward compatibility) - `MINOR`: minor feature (= bug fixes, minor changes not affecting the API, backward compatible) +## [2026.2.0 - not released yet] + +### `fdmnes` + + - Split `fdmnes.py` into a `fdmnes` package (`fdmnes/input.py`) with refactored input handling; `FdmnesXasInput` and `struct2fdmnes` remain importable from `larixite.fdmnes` as before + - Direct `Cif_file` structure type support: a CIF can now be passed straight to FDMNES instead of only an expanded `Crystal` description + - Added convolution parameters (`E_cut`, `Dec`, `Ecent`, `Elarg`, `Gamma_hole`, `Gamma_max`, `Gaussian`, `Estart`) and automatic `Conv_out` section in the generated input + - `green`/`scf` are now plain attributes kept in sync with `params["Green"]`/`params["SCF"]` (previously computed properties) + - Added `run_sbatch()` to submit jobs via `sbatch --parsable` and track status in `status.yaml` + - Added `dump_params()` / `from_yaml()` to store and restore a `FdmnesXasInput` configuration as YAML + - Extended `write_sbatch()` with `nnodes`, `mem_per_cpu`, `walltime`, `constraint`, `partition` options (default `ncpus` raised from 8 to 12) + - Removed the `get_vmax()`/`get_rself()` helper methods; `vmax`/`rself` are now rendered directly from their attribute values (also fixes `Vmax` previously always being written as `-6` regardless of the value set) + +### `struct` + + - Migration of coordination environment analysis and visualization from `struct2xas` into separate modules `analyze` and `visualize` (new `py3Dmol`-based 3D visualization, new `pyproject.toml` dependency) + +### Bug fixes + + - Fix `FdmnesXasInput(optimize=True)` not actually enabling `SCF` in the generated input despite logging "SCF enabled" (regression introduced during the `fdmnes` refactor, caught before release; affects the web app's FDMNES "optimize" option) + - Restore `from larixite.fdmnes import logger`, broken by the module split (used in `examples/fdmnes_xas_workflow.ipynb`) + +### Tests + + - Added tests for analysis, visualization, and fdmnes input modules + ## [2026.1.0 - 18-May-2026] ### web application diff --git a/larixite/fdmnes/__init__.py b/larixite/fdmnes/__init__.py new file mode 100644 index 0000000..e77f857 --- /dev/null +++ b/larixite/fdmnes/__init__.py @@ -0,0 +1 @@ +from .input import FdmnesXasInput, struct2fdmnes, logger diff --git a/larixite/fdmnes.py b/larixite/fdmnes/input.py similarity index 64% rename from larixite/fdmnes.py rename to larixite/fdmnes/input.py index f463382..c1272f6 100644 --- a/larixite/fdmnes.py +++ b/larixite/fdmnes/input.py @@ -9,23 +9,29 @@ spectroscopy (XAS, XES, RIXS) from the atomic structures """ + import time +import yaml from dataclasses import dataclass from pathlib import Path +import os +import shutil +import subprocess from pymatgen.core import __version__ as pymatgen_version, Element, Molecule from larixite.struct import get_structure, get_structure_from_text from larixite.struct.xas import XasStructure from larixite.utils import get_logger, strict_ascii, isotime, read_textfile from larixite.version import __version__ as larixite_version + logger = get_logger("larixite.fdmnes") -TEMPLATE_FOLDER = Path(Path(__file__).parent, "templates") +TEMPLATE_FOLDER = Path(Path(__file__).parent.parent, "templates") FDMNES_DEFAULT_PARAMS = { #: "FDMNES key name": True/False "Energpho": False, "Quadrupole": False, - "Polarize": False, #: -> TODO + "Polarize": False, #: -> TODO: implement giving polarization vectors "Density": False, "Density_all": False, "Green": True, @@ -34,6 +40,10 @@ "Spinorbit": False, "SCF": False, "SCF_exc": False, + "R_self": True, + "N_self": True, + "P_self": True, + "Vmax": False, "Screening": False, #: -> TODO "Dilatorb": False, #: -> TODO "TDDFT": False, @@ -42,9 +52,25 @@ "Atom_conf": False, #: preferred over `Atom` (permits to keep atomic number in the list of atoms) -> TODO "COOP": False, #: -> TODO "Convolution": True, + "E_cut": False, #: taken from E_fermi + "Dec": False, #: apply energy shift by E_cut + "Ecent": True, #: arctan center + "Elarg": True, #: arctan width + "Gamma_hole": False, #: use core-hole broadening + "Gamma_max": True, + "Gaussian": False, } +def _clean_env(): + """Return environment without SLURM_ variables.""" + env = dict(os.environ) + for key in list(env): + if key.startswith("SLURM_"): + del env[key] + return env + + @dataclass class FdmnesXasInput: """Input generator for a XAS calculation with FDMNES""" @@ -59,21 +85,26 @@ class FdmnesXasInput: 0 #: index of the frame inside the structure (e.g. for multi-frame XYZ files) ) radius: float = 7 #: radius of the calculation + green: bool = True #: use muffin-tin approximation + scf: bool = False #: enable SCF erange: str = "-10.0 0.25 60.0 1.0 100.0" #: energy range rself: float | None = None #: radius for the SCF calculation nself: int = 100 #: number of maximum iterations for the SCF calculation pself: float = 0.025 #: initial weight for the SCF calculation - vmax: float | None = None #: maximum potential value for molecules - fileout_prefix: str = ( - "job" #: prefix of the output filename for the FDMNES job (extension: .inp) - ) - tmplpath: str | Path | None = ( - None #: path to the templates directory (FDMNES inputs and sbatch) - ) + vmax: float | str | None = None #: maximum potential value for molecules + fileout_prefix: str | None = None #: prefix for the FDMNES job (deafault: "job" / extension: .inp) + tmplpath: str | Path | None = None #: path to the templates directory params: dict[str, bool] | None = None #: enable/disable parameters for FDMNES optimize: bool = False #: optimize the input parameters spacer: str = " " #: spacer string for the FDMNES input text outdir: str | Path | None = None #: path to the output directory of the FDMNES jobs + ecut: float | str | None = None #: energy cutoff for the convolution + ecent: float = 30. #: energy center for the convolution + elarg: float = 30. #: energy width for the convolution + gamma_hole: float | str | None = None #: start width of the energy convolution (None -> core-hole broadening) + gamma_max: float = 8 #: maximum energy for the convolution + gaussian: float | str | None = None #: gaussian broadening (= experimental braodening) + estart: float | str | None = None #: starting energy for the output spectrum def __post_init__(self): """Validate and optimize attributes""" @@ -90,11 +121,26 @@ def __post_init__(self): self.structpath = self._xs.filepath else: self._xs = get_structure(self.structpath, absorber=self.absorber) + if self.fileout_prefix is None: + self.fileout_prefix = "job" + #: parameters + if self.params is None: + self.params = FDMNES_DEFAULT_PARAMS + self.params["Green"] = self.green + self.params["SCF"] = self.scf #: radius self.set_radius(self.radius) #: R_self if self.rself is None: self.set_rself() + else: + self.params["R_self"] = True + #: atictivate Vmax if given by the user + if self.vmax is None: + self.params["Vmax"] = False + self.vmax = "!" + else: + self.params["Vmax"] = True #: structure type if self.struct_type is None: self.struct_type = self._xs.struct_type @@ -107,9 +153,6 @@ def __post_init__(self): self.tmplpath = Path(self.tmplpath) #: absorption edge self.validate_edge() - #: parameters - if self.params is None: - self.params = FDMNES_DEFAULT_PARAMS #: optimize params if self.optimize: self.params = self.optimize_params() @@ -118,30 +161,42 @@ def __post_init__(self): self.outdir = Path(self.outdir) if self.outdir is None: self.outdir = Path.home() / ".larixite" / "fdmnes" - + #: CONVOLUTION + if self.ecut is None: + self.params["E_cut"] = False + self.ecut = "!" + else: + self.params["E_cut"] = True + if self.gamma_hole is None: + self.params["Gamma_hole"] = False + self.gamma_hole = "!" + else: + self.params["Gamma_hole"] = True + if self.gaussian is None: + self.params["Gaussian"] = False + self.gaussian = "!" + else: + self.params["Gaussian"] = True + if self.estart is None: + self.params["Estart"] = False + self.estart = "!" + else: + self.params["Estart"] = True #: store a list of jobs self._jobs = [] + def __setattr__(self, name: str, value): + """Sync green/scf with params dict on assignment""" + super().__setattr__(name, value) + if name == "green" and hasattr(self, "params") and self.params is not None: + self.params["Green"] = value + elif name == "scf" and hasattr(self, "params") and self.params is not None: + self.params["SCF"] = value + @property def xs(self): return self._xs - @property - def green(self): - return self.params["Green"] - - @green.setter - def green(self, val): - self.params["Green"] = val - - @property - def scf(self): - return self.params["SCF"] - - @scf.setter - def scf(self, val): - self.params["SCF"] = val - def set_radius(self, value: float): self.radius = value self._xs.radius = value @@ -213,10 +268,12 @@ def optimize_params(self) -> dict: ) if "mol" in self.struct_type.lower(): + params["Vmax"] = True self.vmax = -6 - logger.info(f"Vmax set to {self.vmax} for molecule") + logger.info(f"Molecule: Vmax enablend and set to {self.vmax}") - #enable SCF + #: enable SCF + self.scf = True params["SCF"] = True logger.info("SCF enabled") @@ -237,12 +294,12 @@ def get_structure(self, struct_type: str | None = None) -> str: FDMNES supports various structure types: - Crystal -> Implemented (default) - Molecule -> TODO + - Cif_file -> Implemented - Film -> Not implemented yet - Surface -> Not implemented yet - Interface -> Not implemented yet - Pdb_file -> Not implemented yet - Film_Pdb_file -> Not implemented yet - - Cif_file -> Not implemented yet - Film_Cif_file -> Not implemented yet """ @@ -295,6 +352,11 @@ def get_structure(self, struct_type: str | None = None) -> str: ): sitestr = f"{elem.Z:>3d} {site.a:15.10f} {site.b:15.10f} {site.c:15.10f} {occupancy:>5.3f} !{site.label:>4s} {wyckoff:>4s} {elstr:>4s}" structout.append(sitestr) + elif "cif" in struct_type.lower(): + structout.append("Z_absorber") + structout.append(f"{self.spacer}{self.absorber.Z}") + structout.append("Cif_file") + structout.append(f"{self.spacer}{Path(self._xs.filepath).name}") elif "mol" in struct_type.lower(): #: build the cluster and map by distance from absorber at (0,0,0) mol = self._xs.cluster @@ -323,22 +385,6 @@ def get_structure(self, struct_type: str | None = None) -> str: raise AttributeError(errmsg) return "\n".join(structout) - def get_vmax(self) -> str: - """Get the Vmax section of the input""" - if self.vmax is not None: - vmax = ["Vmax"] - vmax.append(f"{self.spacer}-6") - return "\n".join(vmax) - else: - return "! Vmax" - - def get_rself(self): - """Get the R_self section of the input""" - if self.rself is not None: - return f"{self.spacer}{self.rself}" - else: - return "!" - def get_input( self, comment: str = "", @@ -370,12 +416,19 @@ def get_input( "edge": f"{self.spacer}{self.edge}", "radius": f"{self.spacer}{self.radius:.2f}", "erange": self.spacer + self.erange, - "vmax": self.get_vmax(), - "rself": self.get_rself(), + "rself": f"{self.spacer}{self.rself}", "nself": f"{self.spacer}{self.nself}", - "pself": f"{self.spacer}{self.pself:.3f}", + "pself": f"{self.spacer}{self.pself}", + "vmax": f"{self.spacer}{self.vmax}", "struct_type": self.struct_type, "structure": self.get_structure(struct_type=struct_type), + "ecut": f"{self.spacer}{self.ecut}", + "ecent": f"{self.spacer}{self.ecent}", + "elarg": f"{self.spacer}{self.elarg}", + "gamma_hole": f"{self.spacer}{self.gamma_hole}", + "gamma_max": f"{self.spacer}{self.gamma_max}", + "gaussian": f"{self.spacer}{self.gaussian}", + "estart": f"{self.spacer}{self.estart}", } for parkey, parval in params.items(): conf[parkey] = str(parkey) if parval is True else f"! {parkey}" @@ -402,15 +455,31 @@ def write_input( with open(jobdir / "fdmfile.txt", "w") as fp: fp.write(f"1\n{fileout_name}\n") logger.info(f"written {fp.name}") + if self.struct_type is not None and "cif" in self.struct_type.lower(): + src = Path(self._xs.filepath) + shutil.copy2(src, jobdir) + logger.info(f"copied {src.name} to {jobdir}") self._jobs.append(jobdir) + _ = self.dump_params(jobdir / f"{self.fileout_prefix}_params.yaml") return jobdir - def write_sbatch(self, jobdir: str | Path | None = None, template: str | Path | None = None, ncpus: int = 8): + def write_sbatch( + self, + jobdir: str | Path | None = None, + template: str | Path | None = None, + ncpus: int = 12, + nnodes: int = 1, + mem_per_cpu: str = "16GB", + walltime: str = "8:00:00", + constraint: str = "", + partition: str = "nice,nice-long", + **kwargs, + ) -> Path: """Generates a SBATCH file (SLURM workload manager) using a template Arguments --------- - jobdir: str, Path or None + jobdir: str, Path or None path to the job directory to write the SBATCH file if None: uses last job directory template: str, Path or None @@ -435,11 +504,131 @@ def write_sbatch(self, jobdir: str | Path | None = None, template: str | Path | if isinstance(jobdir, str): jobdir = Path(jobdir) sbatchout = jobdir / f"{self.fileout_prefix}.sbatch" - kwargs = {"jobname": self.fileout_prefix, - "ncpus": ncpus,} + kwargs = { + "jobname": self.fileout_prefix, + "nnodes": nnodes, + "ncpus": ncpus, + "mem_per_cpu": mem_per_cpu, + "walltime": walltime, + "constraint": constraint, + "partition": partition, + **kwargs, + } with open(sbatchout, "w") as fp, open(template) as tp: fp.write(tp.read().format(**kwargs)) logger.info(f"written {fp.name}") + return sbatchout + + def run_sbatch(self, jobdir: Path) -> str | None: + """Submit the sbatch script via `sbatch --parsable`. + + Writes `status.yaml` with status=``submitted`` and the SLURM job ID. + Returns the job ID on success, None on failure. + """ + sbatch_script = jobdir / f"{self.fileout_prefix}.sbatch" + if not sbatch_script.exists(): + logger.warning(f"No sbatch script found for job {self.fileout_prefix}") + return None + + try: + result = subprocess.run( + ["sbatch", "--parsable", str(sbatch_script)], + capture_output=True, + text=True, + check=True, + cwd=jobdir, + env=_clean_env(), + ) + slurm_job_id = result.stdout.strip() + status_file = jobdir / "status.yaml" + status_file.write_text( + yaml.dump( + { + "status": "submitted", + "slurm_job_id": slurm_job_id, + } + ) + ) + logger.info( + f"SUBMITTED job {self.fileout_prefix} (slurm job id: {slurm_job_id})" + ) + return slurm_job_id + except subprocess.CalledProcessError as e: + logger.error(f"FAILED to submit job {self.fileout_prefix}: {e}") + return None + + def dump_params(self, yamlpath: str | Path | None = None) -> Path: + """Dump input parameters to a YAML file + + Parameters + ---------- + yamlpath : str | Path | None + path to the output YAML file + if None: writes `{jobdir}/{fileout_prefix}_params.yaml` + if jobdir doesn't exist yet, writes next to the structural file + + Returns + ------- + Path + path to the written YAML file + """ + params_dict = { + "structpath": str(self.structpath), + "absorber": self.absorber.symbol, + "edge": self.edge, + "struct_type": self.struct_type, + "frame": self.frame, + "radius": self.radius, + "erange": self.erange, + "rself": self.rself, + "nself": self.nself, + "pself": self.pself, + "vmax": self.vmax, + "fileout_prefix": self.fileout_prefix, + "params": self.params, + "optimize": self.optimize, + "spacer": self.spacer, + "outdir": str(self.outdir) if self.outdir else None, + } + + if yamlpath is None: + if self._jobs: + yamlpath = self._jobs[-1] / f"{self.fileout_prefix}_params.yaml" + else: + yamlpath = ( + Path(self.structpath).parent / f"{self.fileout_prefix}_params.yaml" + ) + + if isinstance(yamlpath, str): + yamlpath = Path(yamlpath) + + with open(yamlpath, "w") as fp: + yaml.dump(params_dict, fp, default_flow_style=False) + logger.info(f"written {yamlpath.name}") + return yamlpath + + @classmethod + def from_yaml(cls, yamlpath: str | Path) -> "FdmnesXasInput": + """Restore a FdmnesXasInput object from a YAML parameters file + + Parameters + ---------- + yamlpath : str | Path + path to the YAML file previously written by `dump_params()` + + Returns + ------- + FdmnesXasInput + reconstructed input object + """ + if isinstance(yamlpath, str): + yamlpath = Path(yamlpath) + + with open(yamlpath) as fp: + params_dict = yaml.safe_load(fp) + + kwargs = {k: v for k, v in params_dict.items() if k != "structpath"} + return cls(params_dict["structpath"], **kwargs) def struct2fdmnes( diff --git a/larixite/struct/__init__.py b/larixite/struct/__init__.py index 2b7e24b..35750f0 100644 --- a/larixite/struct/__init__.py +++ b/larixite/struct/__init__.py @@ -22,6 +22,8 @@ from larixite.struct.xas_xyz import XasStructureXyz from larixite.utils import get_logger, read_textfile from larixite.amcsd_utils import PMG_CIF_OPTS +from larixite.struct.visualize import visualize +from larixite.struct.analyze import get_coord_env, show_coord_env logger = get_logger("larixite.struct") diff --git a/larixite/struct/analyze.py b/larixite/struct/analyze.py new file mode 100644 index 0000000..6a3e1a0 --- /dev/null +++ b/larixite/struct/analyze.py @@ -0,0 +1,323 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- + +""" +Coordination environment analysis for XAS structures +==================================================== +""" + +import numpy as np +from larixite.struct.xas import XasStructure +from larixite.utils import get_logger, pprint + +logger = get_logger("larixite.struct.analyze") + +# Mapping of chemenv coordination environment symbols to human-readable names +# Based on pymatgen chemenv ideal polyhedra definitions +CE_SYMBOL_DESCRIPTION = { + "L:2": "Linear", + "T:3": "Trigonal planar", + "B:3": "Bent", + "T:4": "Tetrahedral", + "S:4": "Square planar", + "F:4": "See-saw", + "T:5": "Trigonal bipyramidal", + "S:5": "Square pyramidal", + "B:5": "Bay", + "R:5": "Rocket", + "G:5": "Envelope", + "E:5": "Eclipsed pentagonal", + "S:6": "Hexagonal planar", + "A:6": "Hexagonal antiprismatic", + "T:6": "Trigonal prismatic", + "O:6": "Octahedral", + "Y:6": "Pentagonal pyramidal", + "M:6": "Pentagonal bisphenoidal", + "H:6": "Triangular cup", + "S:7": "Pentagonal bipyramidal", + "D:7": "Heptagonal planar", + "S:8": "Square antiprismatic", + "D:8": "Dodecahedral", + "A:8": "Hexagonal bipyramidal", + "A:9": "Monocapped square antiprismatic", + "A:10": "Bicapped square antiprismatic", + "S:12": "Icosahedral", +} + + +def _group_into_shells(dists, toler): + """Group sorted distances into shells. + + Distances within `toler` Angstrom of each other are considered + belonging to the same coordination shell. + + Parameters + ---------- + dists : list of float + Sorted distances from the absorber. + toler : float + Tolerance [Angstrom] to group distances into the same shell. + + Returns + ------- + list of list of float + Each inner list contains the distances belonging to one shell. + """ + if len(dists) == 0: + return [] + + shells = [[dists[0]]] + for d in dists[1:]: + if d - shells[-1][-1] <= toler: + shells[-1].append(d) + else: + shells.append([d]) + + return shells + + +def get_coord_env( + struct: XasStructure, shell: int = 1, toler: float = 0.1 +) -> list: + """Determine the coordination environment around the absorber atom. + + For crystal structures (CIF), uses pymatgen's chemenv toolkit + (LocalGeometryFinder, BVAnalyzer, MultiWeightsChemenvStrategy) to + identify coordination geometry symbols such as: + + - S:4 - Square Plane + - T:4 - Tetrahedral + - T:5 - Trigonal bipyramid + - S:5 - Square pyramidal + - O:6 - Octahedral + - T:6 - Trigonal prism + + For molecular structures (XYZ), uses pymatgen's CrystalNN to + identify neighboring atoms. + + Neighbors are grouped into coordination shells according to their + distances from the absorber: distances differing by less than + `toler` are considered the same shell. + + Parameters + ---------- + struct : XasStructure + The XAS structure object. + shell : int, optional + Number of coordination shells to return. Default is 1 (first + coordination shell only). + toler : float, optional + Distance tolerance [Angstrom] to group neighbors into the same + coordination shell. Default is 0.1. + + Returns + ------- + list + List of lists, each containing: + + [0]: Info string about which site is being analyzed. + [1]: Coordination environment dictionary with keys: + - ce_symbol: coordination geometry symbol + - ce_fraction: probability for given coordination env (0-1) + - csm: continuous symmetry measure (lower = more symmetric) + - permutation: atom arrangement order + [2]: List of shell dictionaries, one per requested shell, with keys: + - shell_idx: shell index (1-based) + - avg_dist: average distance of atoms in the shell [Angstrom] + - cn: coordination number (number of atoms in the shell) + - elems: dict mapping element symbol to count in this shell + - neighbors: list of neighbor site dicts with keys: + - site: pymatgen Site with element and coordinates + - site_index: structure index of the neighbor atom + - dist: distance from absorber [Angstrom] + - image: displacement vector from original site (CIF) + - weight: significance of the neighbor (CIF) + """ + idx_abs_site = struct.absorber_idx + abs_site = struct.struct[idx_abs_site] + coord_env_ce = {} + + # Determine geometry symbol for the 1st coordination shell + if struct.struct_type == "crystal": + try: + from pymatgen.analysis.bond_valence import BVAnalyzer + from pymatgen.analysis.chemenv.coordination_environments.coordination_geometry_finder import ( + LocalGeometryFinder, + ) + from pymatgen.analysis.chemenv.coordination_environments.structure_environments import ( + LightStructureEnvironments, + ) + from pymatgen.analysis.chemenv.coordination_environments.chemenv_strategies import ( + MultiWeightsChemenvStrategy, + ) + + lgf = LocalGeometryFinder() + lgf.setup_structure(struct.struct) + + bva = BVAnalyzer() + try: + valences = bva.get_valences(structure=struct.struct) + except ValueError: + valences = "undefined" + + se = lgf.compute_structure_environments( + max_cn=6, + valences=valences, + only_indices=[idx_abs_site], + only_symbols=["S:4", "T:4", "T:5", "S:5", "O:6", "T:6"], + ) + + logger.debug(se.voronoi.neighbors_distances[idx_abs_site][0]["max"]) + + strategy = MultiWeightsChemenvStrategy.stats_article_weights_parameters() + lse = LightStructureEnvironments.from_structure_environments( + strategy=strategy, structure_environments=se + ) + coord_env_ce = lse.coordination_environments[idx_abs_site] + except Exception: + logger.warning( + "chemenv analysis failed, falling back to basic neighbor search" + ) + coord_env_ce = {"ce_symbol": "chemenv analysis failed"} + + elif struct.struct_type == "molecule": + try: + from pymatgen.analysis.local_env import CrystalNN + + obj = CrystalNN() + coord_dict = obj.get_cn_dict(struct.struct, idx_abs_site) + coord_env_ce = {"ce_symbol": f"Elements Dict = {coord_dict}"} + except Exception: + logger.warning( + "CrystalNN analysis failed, falling back to basic neighbor search" + ) + coord_env_ce = {"ce_symbol": "CrystalNN analysis failed"} + + # Build full neighbor list using struct.get_neighbors() up to cluster_size + # This works for arbitrary distances, unlike chemenv/CrystalNN which only + # cover the 1st coordination shell + search_radius = struct.cluster_size + sphere = struct.struct.get_neighbors(abs_site, search_radius) + + neighbors = [] + for item in sphere: + # Structure: (PeriodicNeighbor, distance, site_index, jimage) + # Molecule: (PeriodicNeighbor, distance, site_index) + site = item[0] + dist = round(float(item[1]), 5) + si = item[2] if len(item) > 2 else None + if dist < 1e-6: + continue + neighbors.append({"site": site, "dist": dist, "site_index": si}) + + neighbors.sort(key=lambda x: x["dist"]) + + # Group into shells by distance clustering + dists = [nb["dist"] for nb in neighbors] + dist_shells = _group_into_shells(dists, toler) + + # Build shell dicts for requested number of shells + shell_dicts = [] + idx = 0 + for s in range(min(shell, len(dist_shells))): + shell_dists = dist_shells[s] + shell_len = len(shell_dists) + shell_nbs = neighbors[idx : idx + shell_len] + idx += shell_len + + avg_dist = round(float(np.mean(shell_dists)), 5) + elems = {} + for nb in shell_nbs: + sym = nb["site"].species_string + elems[sym] = elems.get(sym, 0) + 1 + + shell_dicts.append( + { + "shell_idx": s + 1, + "avg_dist": avg_dist, + "cn": shell_len, + "elems": elems, + "neighbors": shell_nbs, + } + ) + + coord_env_list = [ + [ + f"Coord. Env. for Site {struct.absorber_sites[0][0]}", + coord_env_ce, + shell_dicts, + ] + ] + return coord_env_list + + +def show_coord_env(struct: XasStructure, **kws): + """Pretty-print coordination environment information for the absorber atom. + + Calls get_coord_env() and formats the results as a table with shells, + distances, coordination numbers and element breakdown. + + Parameters + ---------- + struct : XasStructure + The XAS structure object. + **kws + Additional keyword arguments passed to get_coord_env(): + shell : int, default 1 + toler : float, default 0.1 + """ + coord_env = get_coord_env(struct, **kws)[0] + shells = coord_env[2] + + if len(shells) == 0: + logger.warning("Empty coordination environment") + return None + + # Extract coordination geometry symbol(s) + ce_info = coord_env[1] + if isinstance(ce_info, list): + # chemenv returns a list, extract symbol from the first entry + coord_sym = ce_info[0]["ce_symbol"] if ce_info else "not determined" + else: + # fallback: direct dict with ce_symbol key + coord_sym = ce_info.get("ce_symbol", "not determined") + + print( + f"Coord. Env. from absorber atom: {struct.absorber.symbol} " + f"at site {struct.absorber_sites[0][0]}" + ) + desc = CE_SYMBOL_DESCRIPTION.get(coord_sym, "") + if desc: + print(f"Geometry: {coord_sym} ({desc})") + else: + print(f"Geometry: {coord_sym}") + print() + + header = ["Shell", "Dist. (Ang)", "CN", "Elements"] + matrix = [header] + + for sh in shells: + elems_str = ", ".join( + f"{sym}: {cnt}" for sym, cnt in sorted(sh["elems"].items()) + ) + matrix.append( + [ + sh["shell_idx"], + sh["avg_dist"], + sh["cn"], + elems_str, + ] + ) + + pprint(matrix) + print() + + # Detailed neighbor table per shell + for sh in shells: + print(f"Shell {sh['shell_idx']} (avg dist: {sh['avg_dist']} A, CN: {sh['cn']})") + det_header = ["Element", "Distance (Ang)"] + det_matrix = [det_header] + for nb in sh["neighbors"]: + det_matrix.append([nb["site"].species_string, nb["dist"]]) + pprint(det_matrix) + print() diff --git a/larixite/struct/visualize.py b/larixite/struct/visualize.py new file mode 100644 index 0000000..8f41f83 --- /dev/null +++ b/larixite/struct/visualize.py @@ -0,0 +1,314 @@ +#!/usr/bin/env python +# -*- coding: utf-8 -*- + +""" +3D visualization of XasStructure local environment using py3Dmol +=============================================================== + +TODO (future improvements): + +- Add coordination polyhedra around the absorber atom. py3Dmol does not support + polyhedron rendering directly. Possible approaches: + + 1. **nglview** — the `nglview` Jupyter widget supports `add_polyhedron()` and + works directly with pymatgen structures. Could be offered as an alternative + backend (e.g. `visualize(..., backend="nglview")`). + + 2. **py3Dmol workaround with triangles** — compute the convex hull of the first + coordination shell with `scipy.spatial.ConvexHull`, then draw the facets as + semi-transparent plane surfaces via `py3Dmol`'s `addSurface()` or by adding + dummy atoms at facet centroids with connecting sticks to form a wireframe. + + 3. **pymatgen IPyWidgets** — `pymatgen.vis.ipyvolumetric` or `pymatgen.vis` can + render structures with coordination environments in interactive HTML. + + 4. **Export for VESTA** — write the cluster as a `.vesta` file (VESTA's native + format supports polyhedra definitions) so the user can open it in VESTA. + +- Expose the pymatgen Molecule cluster for downstream rendering. +- Support partial occupancy visualization (transparency, multi-snapshot). +- Allow custom color schemes beyond VESTA. +""" + +import sys + +from pymatgen.analysis.local_env import CovalentRadius + +from larixite.struct.xas import XasStructure +from larixite.utils import get_logger + +logger = get_logger("larixite.struct") + +try: + import py3Dmol + + HAS_PY3DMOL = True +except ImportError: + HAS_PY3DMOL = False + + +VESTA_COLORS = { + "Ac": ("#70abfa", "sky blue"), + "Ag": ("#c0c0c0", "silver"), + "Al": ("#81b2d6", "light grey blue"), + "Am": ("#545cf2", "warm blue"), + "Ar": ("#cffec4", "very light green"), + "As": ("#74d057", "fresh green"), + "At": ("#754f45", "brownish purple"), + "Au": ("#ffd123", "sun yellow"), + "B": ("#1fa20f", "green"), + "Ba": ("#00c900", "vibrant green"), + "Be": ("#5ed77b", "soft green"), + "Bh": ("#e00038", "cherry"), + "Bi": ("#9e4fb5", "purply"), + "Bk": ("#8a4fe3", "purple"), + "Br": ("#7e3102", "reddish brown"), + "C": ("#4c4c4c", "charcoal"), + "Ca": ("#5a96bd", "faded blue"), + "Cd": ("#ffd98f", "wheat"), + "Ce": ("#ffffc7", "ecru"), + "Cf": ("#a136d4", "dark orchid"), + "Cl": ("#31fc02", "vivid green"), + "Cm": ("#785ce3", "purple"), + "Co": ("#0000af", "royal blue"), + "Cr": ("#00009e", "royal blue"), + "Cs": ("#57178f", "indigo"), + "Cu": ("#2247dc", "true blue"), + "Db": ("#d1004f", "ruby"), + "Dy": ("#1fffc7", "greenish cyan"), + "Er": ("#00e675", "tealish green"), + "Es": ("#b31fd4", "magenta purple"), + "Eu": ("#61ffc7", "light greenish blue"), + "F": ("#b0b9e6", "light steel blue"), + "Fe": ("#b57100", "brown orange"), + "Fm": ("#b31fba", "magenta purple"), + "Fr": ("#420066", "royal purple"), + "Ga": ("#9ee373", "spring green"), + "Gd": ("#45ffc7", "greenish cyan"), + "Ge": ("#7e6ea6", "dark lavender"), + "H": ("#ffcccc", "light rose"), + "He": ("#fce8ce", "blanched almond"), + "Hf": ("#4dc2ff", "sky blue"), + "Hg": ("#b8b8d0", "cloudy blue"), + "Ho": ("#00ff9c", "medium spring green"), + "Hs": ("#e6002e", "cherry red"), + "I": ("#940094", "dark magenta"), + "In": ("#a67573", "reddish grey"), + "Ir": ("#175487", "light navy"), + "K": ("#a121f6", "electric purple"), + "Kr": ("#fac1f3", "light lavender"), + "La": ("#5ac449", "leafy green"), + "Li": ("#86df73", "light green"), + "Lr": ("#c70066", "deep pink"), + "Lu": ("#00ab24", "kelly green"), + "Md": ("#b30da6", "magenta"), + "Mg": ("#fb7b15", "pumpkin orange"), + "Mn": ("#a7089d", "barney purple"), + "Mo": ("#54b5b5", "medium aquamarine"), + "Mt": ("#eb0026", "cherry red"), + "N": ("#b0b9e6", "light steel blue"), + "Na": ("#f9dc3c", "off yellow"), + "Nb": ("#73c2c9", "seafoam blue"), + "Nd": ("#c7ffc7", "very pale green"), + "Ne": ("#fe37b5", "hot pink"), + "Ni": ("#b7bbbd", "silver"), + "No": ("#bd0d87", "medium violet red"), + "Np": ("#0080ff", "deep sky blue"), + "O": ("#fe0300", "red"), + "Os": ("#266696", "ugly blue"), + "P": ("#c09cc2", "pale purple"), + "Pa": ("#00a1ff", "azure"), + "Pb": ("#575961", "gunmetal"), + "Pd": ("#006985", "peacock blue"), + "Pm": ("#a3ffc7", "light seafoam"), + "Po": ("#ab5c00", "orange brown"), + "Pr": ("#d9ffc7", "very light green"), + "Pt": ("#d0d0e0", "light grey"), + "Pu": ("#006bff", "bright blue"), + "Ra": ("#007d00", "green"), + "Rb": ("#702eb0", "rebecca purple"), + "Re": ("#267dab", "bluish"), + "Rf": ("#cc0059", "deep pink"), + "Rh": ("#0a7d8c", "ocean"), + "Rn": ("#428296", "dirty blue"), + "Ru": ("#248f8f", "blue green"), + "S": ("#fffa00", "bright yellow"), + "Sb": ("#9e63b5", "amethyst"), + "Sc": ("#b563ab", "soft purple"), + "Se": ("#9aef0f", "acid green"), + "Sg": ("#d90045", "ruby"), + "Si": ("#1b3bfa", "vivid blue"), + "Sm": ("#8fffc7", "light seafoam"), + "Sn": ("#9a8eb9", "heather"), + "Sr": ("#00ff00", "lime"), + "Ta": ("#4da6ff", "cornflower blue"), + "Tb": ("#30ffc7", "greenish cyan"), + "Tc": ("#3b9e9e", "sea"), + "Te": ("#d47a00", "pumpkin"), + "Th": ("#00baff", "deep sky blue"), + "Ti": ("#78caff", "light blue"), + "Tl": ("#a6544d", "light maroon"), + "Tm": ("#00d452", "shamrock green"), + "U": ("#008fff", "azure"), + "V": ("#e51900", "tomato red"), + "W": ("#2194d6", "water blue"), + "Xe": ("#429eb0", "cool blue"), + "Y": ("#94ffff", "robin egg blue"), + "Yb": ("#00bf38", "shamrock green"), + "Zn": ("#8f8f81", "warm grey"), + "Zr": ("#00ff00", "lime"), +} + + +def _lighten_hex(hex_color: str, factor: float = 0.75) -> str: + """Return a darker hex color obtained by scaling each RGB channel by *factor*.""" + r = int(int(hex_color[1:3], 16) * factor) + g = int(int(hex_color[3:5], 16) * factor) + b = int(int(hex_color[5:7], 16) * factor) + return "#{:02x}{:02x}{:02x}".format(r, g, b) + + +def _print_color_legend(color_elems: dict[str, str], natoms: int) -> None: + """Print a human-readable color legend, adapting to the execution environment. + + In Jupyter: renders HTML color swatches. + In a terminal: uses ANSI 24-bit color escape codes. + Fallback: plain text with element names and color descriptions. + """ + try: + from IPython.display import display, HTML + + swatches = " ".join( + f'' + f'{e}' + for e, c in color_elems.items() + ) + display(HTML(f'
{natoms} atoms — {swatches}
')) + return + except ImportError: + pass + + if hasattr(sys.stdout, "isatty") and sys.stdout.isatty() and sys.stdout.encoding == "UTF-8": + for elem, hex_color in color_elems.items(): + r = int(hex_color[1:3], 16) + g = int(hex_color[3:5], 16) + b = int(hex_color[5:7], 16) + sys.stdout.write("\x1b[38;2;{};{};{};m{}\x1b[0m ".format(r, g, b, elem)) + sys.stdout.write(f" ({natoms} atoms)\n") + return + + legend = ", ".join(f"{e}: {VESTA_COLORS.get(e, ('#888888', 'grey'))[1]}" for e in color_elems) + logger.info(f"Colors ({natoms} atoms): {legend}") + + +def _adjust_font_color(hex_color: str, dark: str = "#222222", light: str = "#ffffff") -> str: + """Return dark or light text color for contrast against the given background hex color.""" + r = int(hex_color[1:3], 16) + g = int(hex_color[3:5], 16) + b = int(hex_color[5:7], 16) + luminance = (0.299 * r + 0.587 * g + 0.114 * b) / 255 + return dark if luminance > 0.55 else light + + +def _cluster_to_xyz(mol) -> tuple[str, list[str]]: + """Convert a pymatgen Molecule to an XYZ format string and the list of unique elements.""" + coords = [] + elements = [] + for site in mol: + elem = site.species_string + cart = site.coords + coords.append((elem, cart)) + if elem not in elements: + elements.append(elem) + + elements = sorted(elements) + + output_str = f"{len(coords)}\n\n" + for elem, coord in coords: + coord_str = " ".join(f"{c:.6f}" for c in coord) + output_str += f"{elem} {coord_str}\n" + + return output_str, elements + + +def visualize(struct: XasStructure, radius: float = 2.5, unitcell: bool = False): + """ + Display a 3D visualization of the local structure around the absorber atom. + + Args: + struct: XasStructure object with the absorber defined. + radius: radius in Angstroms for the visualization sphere around the central atom. + unitcell: if True, displays the structure unit cell (CIF structures only). + + Returns: + py3Dmol.view object, or None if py3Dmol is not installed. + """ + if not HAS_PY3DMOL: + logger.error("py3Dmol not installed! -> run `pip install py3Dmol`") + return + + mol = struct.build_cluster(radius=radius) + + xyz, elems = _cluster_to_xyz(mol) + + xyzview = py3Dmol.view(width=600, height=600) + xyzview.addModel(xyz, "xyz") + + if unitcell: + try: + lat = struct.struct.lattice + a, b, c = lat.a, lat.b, lat.c + alpha, beta, gamma = lat.alpha, lat.beta, lat.gamma + m = xyzview.getModel() + m.setCrystData(a, b, c, alpha, beta, gamma) + xyzview.addUnitCell() + except AttributeError: + logger.warning("Unit cell visualization not available for this structure type") + + color_elems = {} + for elem in elems: + color_elems[elem] = VESTA_COLORS.get(elem, ("#888888", "grey"))[0] + + covalent_r = CovalentRadius().radius + + absorber_sym = struct.absorber.symbol + absorber_color = color_elems.get(absorber_sym, "#888888") + absorber_r = max(0.2, min(covalent_r.get(absorber_sym, 0.5) / 2, 0.75)) + absorber_light = _lighten_hex(absorber_color, 0.6) + + for elem in elems: + color = color_elems[elem] + sphere_r = max(0.2, min(covalent_r.get(elem, 0.5) / 2, 0.75)) + xyzview.setStyle( + {"elem": elem}, + { + "stick": { + "radius": 0.1, + "opacity": 1, + "hidden": False, + "color": color, + }, + "sphere": {"color": color, "radius": sphere_r, "opacity": 1}, + }, + ) + + #: reset index 0 style to avoid duplication, then reapply with lighter color + xyzview.setStyle({"index": 0}, None) # clear + xyzview.setStyle( + {"index": 0}, + { + "stick": { + "radius": 0.1, + "opacity": 1, + "hidden": False, + "color": absorber_light, + }, + "sphere": {"color": absorber_light, "radius": absorber_r, "opacity": 1}, + }, + ) + + xyzview.zoomTo() + xyzview.show() + + _print_color_legend(color_elems, len(mol)) + return xyzview diff --git a/larixite/templates/fdmnes_sbatch_esrf.tmpl b/larixite/templates/fdmnes_sbatch_esrf.tmpl index 2cce9e3..a621e3b 100644 --- a/larixite/templates/fdmnes_sbatch_esrf.tmpl +++ b/larixite/templates/fdmnes_sbatch_esrf.tmpl @@ -1,13 +1,14 @@ #!/bin/bash -l -#SBATCH --nodes=1 #: Number of nodes / if more nodes, add `#SBATCH --constraint="hpc7"` (type of machines) +#SBATCH --nodes={nnodes} #: Number of nodes / if more nodes, add `#SBATCH --constraint="hpc7"` (type of machines) #SBATCH --ntasks-per-node={ncpus} #: Total number of tasks (=CPUs cores) per node #SBATCH --cpus-per-task=1 #: CPUs for each task / DON'T CHANGE THIS - required by FDMNES -#SBATCH --mem-per-cpu=16G #: Total memory for each CPU -#SBATCH --time=120:00:00 #: Wall-time limit -#SBATCH --partition=nice,nice-long #: ESRF cluster +#SBATCH --mem-per-cpu={mem_per_cpu} #: Total memory for each CPU +#SBATCH --time={walltime} #: Wall-time limit +#SBATCH --partition="{partition}" #: ESRF cluster +#SBATCH --constraint="{constraint}" #SBATCH --job-name="{jobname}" #: Name used to set the output files -#SBATCH --output=%x_output.txt -#SBATCH --error=%x_output.txt +#SBATCH --output="{jobname}_slurm.stdout" +#SBATCH --error="{jobname}_slurm.stderr" function compress {{ gzip -f {jobname}_bav.txt diff --git a/larixite/templates/fdmnes_xas.tmpl b/larixite/templates/fdmnes_xas.tmpl index 2e9d901..e0ed819 100755 --- a/larixite/templates/fdmnes_xas.tmpl +++ b/larixite/templates/fdmnes_xas.tmpl @@ -22,8 +22,8 @@ Radius ! {Quadrupole} ! -{Polarize} -! =>TODO +{Polarize} ! without giving a vector it calculates non equivalent directions +! ! =>TODO: specify polarization and wave vectors {TDDFT} {Memory_save} {Relativism} @@ -36,17 +36,37 @@ Radius {SpGr_atom} !`Full_atom` is default since March 2026 ! {SCF} -R_self +{R_self} {rself} -N_self +{N_self} {nself} -P_self +{P_self} {pself} {SCF_exc} +! +{Vmax} {vmax} +! {Atom_conf} ! {structure} ! +Conv_out +{filout}_conv {Convolution} +{E_cut} +{ecut} +{Dec} +{Gamma_hole} +{gamma_hole} +{Ecent} +{ecent} +{Elarg} +{elarg} +{Gamma_max} +{gamma_max} +{Gaussian} +{gaussian} +{Estart} +{estart} End diff --git a/larixite/utils.py b/larixite/utils.py index 1789339..44b5afe 100644 --- a/larixite/utils.py +++ b/larixite/utils.py @@ -184,6 +184,7 @@ def get_logger(name="larixite", level="INFO"): """Get a custom logger for larixite""" logger = logging.getLogger(name) logger.setLevel(log_levels[level]) + logger.propagate = False formatter = ColorFormatter() handler = logging.StreamHandler() handler.setFormatter(formatter) diff --git a/pyproject.toml b/pyproject.toml index a90dc92..f455754 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -27,6 +27,7 @@ dependencies = [ "mp-api>=0.45.8", "emmet-core>=0.84.9", "requests", + "py3Dmol", ] license = "MIT" diff --git a/tests/test_analyze.py b/tests/test_analyze.py new file mode 100644 index 0000000..6891552 --- /dev/null +++ b/tests/test_analyze.py @@ -0,0 +1,135 @@ +"""Tests for larixite.struct.analyze — get_coord_env and show_coord_env.""" + +from pathlib import Path +from io import StringIO +import sys +import warnings + +import pytest +import numpy as np + +from larixite.struct import get_structure +from larixite.struct.analyze import get_coord_env, show_coord_env + +STRUCTS = Path(__file__).parent / "structs" + + +# ── helpers ─────────────────────────────────────────────────────── + + +def cap_show(output, struct, **kws): + """Capture the stdout produced by show_coord_env and return it.""" + old_out = sys.stdout + sys.stdout = buf = StringIO() + with warnings.catch_warnings(): + warnings.simplefilter("ignore") + output = show_coord_env(struct, **kws) + sys.stdout = old_out + return buf.getvalue() + + +@pytest.fixture(scope="module") +def cif_zno(): + with warnings.catch_warnings(): + warnings.simplefilter("ignore", UserWarning) + return get_structure(STRUCTS / "ZnO_mp-2133.cif", "Zn") + + +@pytest.fixture(scope="module") +def xyz_cuo6(): + with warnings.catch_warnings(): + warnings.simplefilter("ignore", UserWarning) + return get_structure(STRUCTS / "CuO6_D4h.xyz", "Cu") + + +# ── CIF (crystal) tests ───────────────────────────────────────── + + +class TestCoordEnvCIF: + def test_shell1_returned(self, cif_zno): + ce = get_coord_env(cif_zno, shell=1) + shells = ce[0][2] + assert len(shells) == 1 + assert shells[0]["shell_idx"] == 1 + + def test_shell1_cn_positive(self, cif_zno): + ce = get_coord_env(cif_zno, shell=1) + shells = ce[0][2] + assert shells[0]["cn"] > 0 + + def test_shell1_has_elements(self, cif_zno): + ce = get_coord_env(cif_zno, shell=1) + shells = ce[0][2] + assert len(shells[0]["elems"]) > 0 + + def test_shell1_neighbors_sorted(self, cif_zno): + ce = get_coord_env(cif_zno, shell=1) + shells = ce[0][2] + dists = [nb["dist"] for nb in shells[0]["neighbors"]] + assert dists == sorted(dists) + + def test_multiple_shells_returned(self, cif_zno): + ce = get_coord_env(cif_zno, shell=3) + shells = ce[0][2] + assert len(shells) >= 2 + + def test_multiple_shells_increasing_distance(self, cif_zno): + ce = get_coord_env(cif_zno, shell=3) + shells = ce[0][2] + for i in range(1, len(shells)): + assert shells[i]["avg_dist"] > shells[i - 1]["avg_dist"] + + def test_cn_sum_matches_requested_shells(self, cif_zno): + shell = 3 + ce = get_coord_env(cif_zno, shell=shell) + shells = ce[0][2] + expected_cn = 0 + distances = sorted(nb["dist"] for s in shells for nb in s["neighbors"]) + assert len(distances) == sum(s["cn"] for s in shells) + + def test_show_coord_env_output(self, cif_zno): + output = cap_show(None, cif_zno, shell=3) + assert "Coord. Env. from absorber atom: Zn" in output + assert "Shell" in output + assert "CN" in output + + def test_toler_groups_closely_spaced_distances(self, cif_zno): + ce_strict = get_coord_env(cif_zno, shell=10, toler=0.01) + ce_loose = get_coord_env(cif_zno, shell=10, toler=1.0) + n_shells_strict = len(ce_strict[0][2]) + n_shells_loose = len(ce_loose[0][2]) + # Larger tolerance should group into fewer (or equal) shells + assert n_shells_loose <= n_shells_strict + + +# ── XYZ (molecule) tests ───────────────────────────────────────── + + +class TestCoordEnvXYZ: + def test_shell1_returned(self, xyz_cuo6): + ce = get_coord_env(xyz_cuo6, shell=1) + shells = ce[0][2] + assert len(shells) >= 1 + assert shells[0]["shell_idx"] == 1 + + def test_shell1_cn_positive(self, xyz_cuo6): + ce = get_coord_env(xyz_cuo6, shell=1) + shells = ce[0][2] + assert shells[0]["cn"] > 0 + + def test_shell1_has_elements(self, xyz_cuo6): + ce = get_coord_env(xyz_cuo6, shell=1) + shells = ce[0][2] + assert len(shells[0]["elems"]) > 0 + + def test_shell1_neighbors_sorted(self, xyz_cuo6): + ce = get_coord_env(xyz_cuo6, shell=1) + shells = ce[0][2] + dists = [nb["dist"] for nb in shells[0]["neighbors"]] + assert dists == sorted(dists) + + def test_show_coord_env_output(self, xyz_cuo6): + output = cap_show(None, xyz_cuo6, shell=2) + assert "Coord. Env. from absorber atom: Cu" in output + assert "Shell" in output + assert "CN" in output diff --git a/tests/test_fdmnes.py b/tests/test_fdmnes.py index a4d9f27..63c84ab 100644 --- a/tests/test_fdmnes.py +++ b/tests/test_fdmnes.py @@ -30,5 +30,24 @@ def test_fdmnes(): assert outdir.exists() +def test_green_scf_sync(): + db = get_amcsd() + cif = db.get_cif(4438) + outfile = cif.write_cif(verbose=True) + sg = get_structure(outfile, "Fe") + + f = FdmnesXasInput(sg, absorber="Fe", green=False, scf=True, optimize=False) + assert f.green is False + assert f.scf is True + assert f.params["Green"] is False + assert f.params["SCF"] is True + + f.green = True + f.scf = False + assert f.params["Green"] is True + assert f.params["SCF"] is False + + if __name__ == "__main__": test_fdmnes() + test_green_scf_sync() diff --git a/tests/test_visualize.py b/tests/test_visualize.py new file mode 100644 index 0000000..4e5201e --- /dev/null +++ b/tests/test_visualize.py @@ -0,0 +1,191 @@ +"""Tests for larixite.struct.visualize — _cluster_to_xyz, _lighten_hex, VESTA_COLORS, visualize().""" + +from pathlib import Path + +import pytest + +from larixite.struct import get_structure +from larixite.struct.visualize import ( + HAS_PY3DMOL, + VESTA_COLORS, + _cluster_to_xyz, + _lighten_hex, + visualize, +) + +STRUCTS = Path(__file__).parent / "structs" + +# ── fixtures ────────────────────────────────────────────────────── + + +@pytest.fixture(scope="module") +def cif_zno(): + return get_structure(STRUCTS / "ZnO_mp-2133.cif", "Zn") + + +@pytest.fixture(scope="module") +def xyz_cuo6(): + return get_structure(STRUCTS / "CuO6_D4h.xyz", "Cu") + + + + +# ── _cluster_to_xyz ───────────────────────────────────────────── + + +class TestClusterToXyz: + def test_returns_strings_and_elems(self, cif_zno): + mol = cif_zno.build_cluster(radius=2.5) + xyz_str, elems = _cluster_to_xyz(mol) + assert isinstance(xyz_str, str) + assert isinstance(elems, list) + assert all(isinstance(e, str) for e in elems) + + def test_xyz_header_matches_atom_count(self, cif_zno): + mol = cif_zno.build_cluster(radius=2.5) + xyz_str, elems = _cluster_to_xyz(mol) + header = int(xyz_str.strip().split("\n")[0]) + assert header == len(mol) + + def test_absorber_first(self, cif_zno): + mol = cif_zno.build_cluster(radius=2.5) + xyz_str, _ = _cluster_to_xyz(mol) + lines = xyz_str.strip().split("\n") + absorber = lines[2].split()[0] + assert absorber == cif_zno.absorber.symbol + + def test_coords_zero_origin(self, cif_zno): + mol = cif_zno.build_cluster(radius=2.5) + xyz_str, _ = _cluster_to_xyz(mol) + lines = xyz_str.strip().split("\n") + coords = lines[2].split()[1:] + for c in coords: + assert float(c) == pytest.approx(0.0, abs=1e-5) + + def test_elements_sorted(self, cif_zno): + mol = cif_zno.build_cluster(radius=2.5) + _, elems = _cluster_to_xyz(mol) + assert elems == sorted(elems) + + def test_elements_match_structure(self, cif_zno): + mol = cif_zno.build_cluster(radius=7.0) + _, elems = _cluster_to_xyz(mol) + assert "Zn" in elems + assert "O" in elems + + def test_xyz_single_element(self, xyz_cuo6): + mol = xyz_cuo6.build_cluster(radius=0.5) + xyz_str, elems = _cluster_to_xyz(mol) + assert len(elems) == 1 + assert elems == ["Cu"] + + +# ── _lighten_hex ──────────────────────────────────────────────── + + +class TestLightenHex: + def test_returns_hex(self): + assert _lighten_hex("#ff0000", 0.5).startswith("#") + + def test_darkens_correctly(self): + result = _lighten_hex("#ff0000", 0.5) + # int(255 * 0.5) == 127 (truncation), not 128 + assert result == "#7f0000" + + def test_identity(self): + assert _lighten_hex("#b57100", 1.0) == "#b57100" + + def test_factor_0_8(self): + result = _lighten_hex("#ff0000", 0.8) + r = int(int("ff", 16) * 0.8) + assert "#" in result and len(result) == 7 + + def test_fe_color(self): + fe_hex, fe_name = VESTA_COLORS["Fe"] + light = _lighten_hex(fe_hex, 0.8) + assert light != fe_hex + assert light.startswith("#") + assert fe_name == "brown orange" + + def test_all_vesta_colors_tuplified(self): + for sym, entry in VESTA_COLORS.items(): + assert isinstance(entry, tuple) and len(entry) == 2 + hex_color, name = entry + light = _lighten_hex(hex_color, 0.8) + assert light.startswith("#") and len(light) == 7 + assert isinstance(name, str) and name + + def test_invalid_hex_unchanged(self): + with pytest.raises(ValueError): + _lighten_hex("notahex", 0.5) + + +# ── VESTA_COLORS ──────────────────────────────────────────────── + + +class TestVestaColors: + COMMON = ["Fe", "O", "C", "N", "Si", "Al", "Ti", "Cu", "Zn", "S", "P", "Mg", "Ca"] + + def test_common_elements_present(self): + for sym in self.COMMON: + assert sym in VESTA_COLORS + + def test_colors_are_hex(self): + for hex_color, name in VESTA_COLORS.values(): + assert hex_color.startswith("#") and len(hex_color) == 7 + assert isinstance(name, str) and name + + def test_all_vesta_values_valid_int(self): + for hex_color, name in VESTA_COLORS.values(): + int(hex_color[1:3], 16) + int(hex_color[3:5], 16) + int(hex_color[5:7], 16) + + def test_no_black_in_palette(self): + assert "#000000" not in [h for h, n in VESTA_COLORS.values()] + + def test_color_count(self): + assert len(VESTA_COLORS) > 100 + + +# ── visualize ─────────────────────────────────────────────────── + + +class TestVisualize: + @pytest.mark.skipif(not HAS_PY3DMOL, reason="py3Dmol not installed") + def test_returns_view_cif(self, cif_zno): + result = visualize(cif_zno, radius=2.5) + assert result is not None + + @pytest.mark.skipif(not HAS_PY3DMOL, reason="py3Dmol not installed") + def test_returns_view_xyz(self, xyz_cuo6): + result = visualize(xyz_cuo6, radius=2.0) + assert result is not None + + @pytest.mark.skipif(not HAS_PY3DMOL, reason="py3Dmol not installed") + def test_default_radius(self, cif_zno): + result = visualize(cif_zno) + assert result is not None + + @pytest.mark.skipif(not HAS_PY3DMOL, reason="py3Dmol not installed") + def test_unitcell_flag(self, cif_zno): + result = visualize(cif_zno, radius=2.5, unitcell=True) + assert result is not None + + @pytest.mark.skipif(not HAS_PY3DMOL, reason="py3Dmol not installed") + def test_unitcell_xyz_no_error(self, xyz_cuo6): + result = visualize(xyz_cuo6, radius=2.0, unitcell=True) + assert result is not None + + @pytest.mark.skipif(not HAS_PY3DMOL, reason="py3Dmol not installed") + def test_radius_affects_cluster(self, cif_zno): + result = visualize(cif_zno, radius=1.5) + assert result is not None + + @pytest.mark.skipif(HAS_PY3DMOL, reason="py3Dmol is installed") + def test_no_py3dmol_returns_none(self, cif_zno, monkeypatch): + import larixite.struct.visualize as viz + monkeypatch.setattr(viz, "HAS_PY3DMOL", False) + monkeypatch.setattr(viz, "py3Dmol", None) + result = viz.visualize(cif_zno) + assert result is None