From a91372b199640ffbdd51fb6e642f86f272e3936b Mon Sep 17 00:00:00 2001 From: preetibayas Date: Tue, 27 Oct 2026 12:14:18 +0530 Subject: [PATCH 1/3] Harden Psi4 screening and charged-state workflows Add explicit dipole extraction, validated spin references, isolated calculation workers and artifacts, and bounded optimization recovery with focused regression tests and runtime smoke coverage. Co-authored-by: Copilot <223556219+Copilot@users.noreply.github.com> Copilot-Session: 87cf205e-4bdc-4491-9fdd-cbec1660935f --- agents/psi4/README.md | 160 ++++++++- agents/psi4/agent.yaml | 47 ++- agents/psi4/tools/psi4/Dockerfile | 6 +- agents/psi4/tools/psi4/psi4_state_job.py | 51 +++ agents/psi4/tools/psi4/psi4_utils.py | 240 +++++++++++-- agents/psi4/tools/psi4/run_psi4_job.py | 83 +++++ agents/psi4/tools/psi4/smoke_psi4_runtime.py | 44 +++ .../psi4/tools/psi4/test_psi4_reliability.py | 323 ++++++++++++++++++ agents/psi4/tools/psi4/tool.yaml | 2 +- 9 files changed, 912 insertions(+), 44 deletions(-) create mode 100644 agents/psi4/tools/psi4/psi4_state_job.py create mode 100644 agents/psi4/tools/psi4/run_psi4_job.py create mode 100644 agents/psi4/tools/psi4/smoke_psi4_runtime.py create mode 100644 agents/psi4/tools/psi4/test_psi4_reliability.py diff --git a/agents/psi4/README.md b/agents/psi4/README.md index 0002dabf..9727b61d 100644 --- a/agents/psi4/README.md +++ b/agents/psi4/README.md @@ -27,7 +27,7 @@ Psi4 features a Python API for scriptable workflows, density fitting for computa ```bash # Navigate to the psi4 directory -cd 6-solutions/tools-and-models/psi4/ +cd agents/psi4/tools/psi4/ # Build the Docker image docker build -t psi4:latest . @@ -50,7 +50,10 @@ psi4/ ├── psi4-tool-definition.yaml # Tool infrastructure definition ├── psi4-agent-definition.yaml # Agent instructions and capabilities ├── psi4_utils.py # Python utilities library -├── test_psi4_utils.py # Unit tests +├── run_psi4_job.py # Isolated script launcher +├── psi4_state_job.py # Private fresh-process four-point worker +├── test_psi4_reliability.py # Mocked-Psi4 unittest suite +├── smoke_psi4_runtime.py # Real-Psi4 runtime smoke ├── README.md # This file └── example-input-files/ ├── water.xyz # Water molecule @@ -89,6 +92,9 @@ psi4/ | Function | Description | |----------|-------------| | `compute_energy()` | Single-point energy calculation | +| `compute_electronic_properties()` | One wavefunction for energy, orbital gap, and explicit dipole | +| `compute_dipole_moment()` | Explicit one-electron dipole, atomic units and Debye | +| `compute_reorganization_energies()` | Gas-phase four-point hole/electron reorganization | | `compute_gradient()` | Energy gradient (forces) | | `optimize_geometry()` | Geometry optimization | | `compute_frequencies()` | Vibrational frequencies | @@ -151,14 +157,18 @@ psi4.set_options({ ### Local Testing (without Psi4) -```bash -cd 6-solutions/tools-and-models/psi4/ -python -m venv .venv-psi4 -source .venv-psi4/bin/activate # or .venv-psi4\Scripts\activate on Windows -pip install numpy scipy matplotlib pandas pytest -pytest test_psi4_utils.py -v +```console +cd agents/psi4/tools/psi4 +python -m unittest -v test_psi4_reliability +python -m py_compile psi4_utils.py run_psi4_job.py psi4_state_job.py smoke_psi4_runtime.py test_psi4_reliability.py ``` +The test runner is stdlib `unittest`; it mocks Psi4 and uses the library's existing +NumPy dependency. Tests include reference selection/validation, single-wavefunction +dipole conversion, nonfinite rejection, optimization retry/exhaustion, candidate-ID +and partial-failure preservation, four-point formulas, concurrent launcher isolation, +nested artifact copying, and child/copy failure exit codes. + ### Docker Testing ```bash @@ -256,7 +266,139 @@ Microsoft Discovery team contact: discovery-catalog@microsoft.com ## Known Limitations -No known limitations at this time. If you encounter issues, please report them via the support channel above. +Mocked tests do not validate numerical accuracy, SCF/optimizer behavior, basis/method +availability, or PSIO behavior in the deployed image. Python 3.14 host tests are not +a substitute for the container's Python 3.11/Psi4 runtime. The scratch/working-directory +changes are isolation hardening, **not a proven PSIO root-cause fix**. Reproduce the +original failing input in the approved runtime after the smoke below. + +## Gas-phase electronic properties and four-point reorganization + +Public helper signatures: + +```python +compute_dipole_moment(molecule, method='b3lyp', basis='def2-tzvp', **kwargs) +compute_electronic_properties(molecule, method='b3lyp', basis='def2-tzvp', **kwargs) +compute_reorganization_energies( + molecule, candidate_id, method='b3lyp-d3bj', + optimization_basis='def2-svp', energy_basis='def2-tzvp', + neutral_multiplicity=1, cation_multiplicity=2, anion_multiplicity=2, **kwargs +) +``` + +`compute_energy`, `compute_homo_lumo_gap`, and `optimize_geometry` retain their +signatures and select `rhf` for singlets or `uhf` for higher multiplicities on **each** +calculation, including DFT. Explicit `reference='rhf'|'uhf'|'rohf'` remains supported. +Closed-shell `rhf` with open-shell multiplicity is rejected before calculation. +Do not substitute `rks`/`uks` keywords. Unconverged SCF acceptance is disabled. +Open-shell orbital gaps retain the existing **alpha-spin** convention; they are not +spin-combined, fundamental, or optical gaps. + +The combined helper returns the existing `homo_energy_hartree/eV`, +`lumo_energy_hartree/eV`, `gap_hartree/eV`, and `total_energy_hartree/eV` keys, plus +`reference`, `dipole_vector_au`, `dipole_vector_debye`, and `dipole_moment_debye`. +It computes one energy/wavefunction and explicitly calls +`psi4.oeprop(wfn, 'DIPOLE', title='DISCOVERY')`, retrieving +`wfn.variable('DISCOVERY DIPOLE')` as a finite three-component vector. +The conversion is $1\,e a_0 = 2.5417464518950259\,\mathrm{D}$. +Charged-species dipoles depend on the coordinate origin. No polarizability or +optional-property artifact is required for these helpers to succeed. + +API/units references verified against upstream documentation and implementation: + +- [oeprop signature, wavefunction argument, and title prefix](https://psicode.org/psi4manual/master/oeprop.html) +- [Vector dipoles and atomic units](https://psicode.org/psi4manual/master/glossary_psivariables.html#psivar-SCF-DIPOLE) +- [Atomic-unit-to-Debye conversion](https://psicode.org/psi4manual/master/autodoc_physconst.html) +- [OEProp title-prefixed array storage on the wavefunction](https://github.com/psi4/psi4/blob/master/psi4/src/psi4/libmints/oeprop.cc) +- [OptimizationConvergenceError last wavefunction](https://psicode.org/psi4manual/master/_modules/psi4/driver/p4util/exceptions.html#OptimizationConvergenceError) + +Four-point calculations always reoptimize the neutral, cation, and anion; an +upstream partial neutral result is **never** assumed converged. Defaults represent +a singlet neutral and doublet singly charged ions; supply the actual multiplicities +for other spin states. Input must be neutral. Each state uses its own validated +reference. Every optimization has at most **one** retry, from the exception's last +wavefunction geometry, with doubled `geom_maxiter`. If no last wavefunction is +available it fails without restarting the original geometry. Neither convergence +tolerances nor method/basis are weakened. + +Let $Q_0$, $Q_+$, and $Q_-$ denote the optimized neutral, cation, and anion geometries. +At the **same final method/energy basis**, the two four-point sums are: + +$$ +\lambda_h = [E_+(Q_0)-E_+(Q_+)] + [E_0(Q_+)-E_0(Q_0)] +$$ +$$ +\lambda_e = [E_-(Q_0)-E_-(Q_-)] + [E_0(Q_-)-E_0(Q_0)]. +$$ + +Results preserve `candidate_id`, `status`, method/bases, `energy_units='hartree'`, +`states` (charges, multiplicities, references, convergence flags, geometries in +angstrom), and `components`. Component keys are `E_neutral_at_neutral`, +`E_cation_at_cation`, `E_anion_at_anion`, `E_cation_at_neutral`, +`E_neutral_at_cation`, `E_anion_at_neutral`, and `E_neutral_at_anion`. +Success adds `hole_reorganization_hartree/eV` and `electron_reorganization_hartree/eV`. +On failure `status='failed'`, error details and all completed components remain; +do not treat them as completed reorganization results. Negative values are not +clamped or hidden. The helper persists `reorganization_result.json` in its unique +`artifact_directory`, with per-state requests/results and Psi4/worker logs. + +### Invocation isolation and deployment use + +The Dockerfile **packages source only when an approved image is built**. Editing it +does not update a deployed image. Do not build, push, retag, replace an image, or +change live tool references without explicit approval. + +The execution command in the matching tool definition should use positional args: + +```console +/opt/conda/envs/psi4env/bin/python /usr/local/lib/run_psi4_job.py /path/to/script.py /path/to/output +``` + +The launcher executes `[sys.executable, script_path]` without shell reparsing. +It creates unique writable work/scratch directories under the output directory's +parent, preflights writes, sets `PSI_SCRATCH` **before imports**, streams combined +stdout/stderr to terminal and `logs.out`, and copies nested work artifacts even +after failure. Scratch and symlinks are excluded. A child nonzero exit is preserved; +copy failure makes an otherwise successful invocation fail. `launcher_failure.json` +records failures. Invocation directories are cleaned only after copying. + +`quick_setup(input_dir=..., output_dir=...)` retains the supplied data-handling +input/output paths; work resolves to invocation-specific `PSI4_WORK_DIR`, then an +explicit `work_dir`, otherwise a unique directory. `setup_psi4` binds the Psi4 +IOManager to the selected writable scratch path. Use these settings rather than +shared work/scratch names. Four-point state jobs launch seven fresh processes +(three optimizations plus final diagonal energies, four cross energies) with their +own scratch directories and filenames. Non-JSON wavefunctions remain process-local. + +```python +from psi4_utils import ( + quick_setup, setup_psi4, create_molecule, compute_electronic_properties, + compute_reorganization_energies, save_final_results, quick_finish, +) +# Supply actual paths from the task's dataHandlingContext. +quick_setup(input_dir=input_directory, output_dir=output_directory) +setup_psi4(memory='4 GB', nthreads=1) +molecule = create_molecule(geometry, charge=0, multiplicity=1) +properties = compute_electronic_properties(molecule, method='b3lyp', basis='def2-tzvp') +reorganization = compute_reorganization_energies(molecule, candidate_id='candidate-001') +save_final_results( + {'electronic_properties': properties, 'reorganization': reorganization}, + status='completed' if reorganization['status'] == 'success' else 'failed', +) +quick_finish() +``` + +In the approved Psi4 runtime, run the real smoke using the launcher: + +```console +python run_psi4_job.py smoke_psi4_runtime.py smoke-output +``` + +It checks a water dipole, charged doublet water SCF energies, and tiny LiH +four-point calculations with a diffuse basis. It writes `smoke_results.json` and +asserts convergence and finite values. **Not yet run successfully in a real Psi4 +runtime**: the local Windows Python has no Psi4 and Docker is unavailable. This is +the remaining verification gate before relying on L4b production numerical results. ## Contributing diff --git a/agents/psi4/agent.yaml b/agents/psi4/agent.yaml index 99364764..3928c8f2 100644 --- a/agents/psi4/agent.yaml +++ b/agents/psi4/agent.yaml @@ -37,6 +37,35 @@ instructions: | 1. **Plan**: Identify calculation type, method/basis, sequence 2. **Script**: ONE comprehensive Python script executed immediately + ## AUTOMATED SCREENING AND RECOVERY + + - For neutral gas-phase screening, use the packaged `compute_electronic_properties` + helper instead of inventing dipole variable names or parsing printed output. + Obtain orbital energies and dipole from the same converged property wavefunction. + Polarizability is optional unless explicitly required; its failure must not erase + successful mandatory properties. + - For four-point reorganization energies, use the packaged + `compute_reorganization_energies` helper. Preserve upstream candidate identifiers, + geometry references, all component energies, units, and separate neutral/cation/anion + convergence outcomes. A partially optimized neutral geometry must be converged + before reporting a successful adiabatic reorganization energy. + - Use the helper's bounded recovery within this execution. Reuse successful stages + rather than repeating every optimization for a missing property. Do not weaken + convergence tolerances, change the final method/basis, or substitute zero for a + missing value to make validation pass. + - The tool launcher provides invocation-isolated working and scratch directories. + Do not override them with shared paths, delete another job's scratch files, or + put scratch files in the published output directory. + - Before declaring completion, check that every required numeric value is finite, + each geometry file exists, and all requested component energies and units are + present. Distinguish a valid zero dipole from an absent dipole. + - Save `final_results.json` before `quick_finish()`. Always publish diagnostic logs + and explicitly labeled partial/failed results, even if every candidate fails. + Do not withhold failure artifacts from downstream validation. + - If bounded recovery is exhausted, report the exact failing stage and error. + Never claim completion, fabricate scientific results, or request disabling + validation merely to avoid human escalation. + # SYSTEM TYPE DETECTION | System Type | Method Recommendation | Basis Set | Notes | @@ -96,12 +125,14 @@ instructions: | # PSI4 SETUP setup_psi4, # setup_psi4(memory='4 GB', nthreads=8) - create_molecule, # create_molecule(geometry_str, charge=0, mult=1) + create_molecule, # create_molecule(geometry_str, charge=0, multiplicity=1) read_xyz_file, # read_xyz_file('molecule.xyz') -> geometry string # ENERGY CALCULATIONS compute_energy, # compute_energy(mol, method='b3lyp', basis='def2-tzvp') # -> {energy_hartree, energy_eV, energy_kcal} + compute_electronic_properties, # Orbital energies and explicit dipole from one wavefunction + compute_reorganization_energies, # Four-point workflow with charged-state isolation compute_homo_lumo_gap, # compute_homo_lumo_gap(mol, method='b3lyp', basis='def2-tzvp') # -> {homo_energy_eV, lumo_energy_eV, gap_eV, caveat, recommendations} # Returns KS orbital gap (NOT the fundamental or optical gap) @@ -202,8 +233,8 @@ instructions: | results['error'] = str(e) finally: - quick_finish() save_final_results(results, status=results.get('status', 'completed')) + quick_finish() ``` # EXAMPLE WORKFLOWS @@ -254,8 +285,8 @@ instructions: | results['status'] = 'failed' results['error'] = str(e) finally: - quick_finish() save_final_results(results, status=results.get('status', 'completed')) + quick_finish() ``` ## Example 2: SAPT Interaction Energy @@ -311,8 +342,8 @@ instructions: | results['status'] = 'failed' results['error'] = str(e) finally: - quick_finish() save_final_results(results, status=results.get('status', 'completed')) + quick_finish() ``` ## Example 3: Excited States @@ -360,15 +391,15 @@ instructions: | results['status'] = 'failed' results['error'] = str(e) finally: - quick_finish() save_final_results(results, status=results.get('status', 'completed')) + quick_finish() ``` # COMMON ERRORS AND SOLUTIONS | Error | Cause | Solution | |-------|-------|----------| - | `SCF did not converge` | Difficult electronic structure | Use `scf_type='pk'`, reduce `d_convergence` | + | `SCF did not converge` | Difficult electronic structure | Use bounded recovery without loosening the requested convergence criteria; retain diagnostics | | `Linear dependency` | Large basis with small molecule | Use smaller basis or `S_orthogonalize` | | `Memory exceeded` | Insufficient memory | Increase memory or use smaller basis | | `Open-shell singlet` | Wrong reference | Set `reference='uhf'` or `'rohf'` | @@ -404,10 +435,10 @@ instructions: | # CRITICAL REMINDERS - 1. Generate ONE script - NEVER multiple scripts + 1. Generate ONE orchestration script. The packaged workflow may use isolated worker processes and bounded internal retries; do not restart the entire campaign for a failed property. 2. COPY inputs using quick_setup(input_dir=..., output_dir=...) with paths from dataHandlingContext 3. Always call setup_psi4() before calculations - 4. For open-shell systems, set correct multiplicity + 4. Set charge, multiplicity, and reference for every calculation. Neutral singlets use a restricted reference; doublet cations/anions use an unrestricted reference. Do not rely on options left over from the previous molecule. 5. Optimize geometry BEFORE frequencies 6. Pass actual status to save_final_results: `save_final_results(results, status=results.get('status', 'completed'))` 7. DFT frequency calculations: compute_frequencies() defaults to dertype='gradient' (numerical Hessian from analytical gradients) which works for ALL methods. Do NOT pass dertype='hessian' with DFT methods (B3LYP, PBE, etc.) - analytical Hessians are only available for HF. diff --git a/agents/psi4/tools/psi4/Dockerfile b/agents/psi4/tools/psi4/Dockerfile index bd03adfc..e8411e09 100644 --- a/agents/psi4/tools/psi4/Dockerfile +++ b/agents/psi4/tools/psi4/Dockerfile @@ -94,7 +94,6 @@ RUN ${MAMBA_ROOT_PREFIX}/envs/psi4env/bin/pip install --no-cache-dir \ || true # Set Psi4 environment variables -ENV PSI_SCRATCH=/tmp/psi4_scratch ENV OMP_NUM_THREADS=1 ENV MKL_NUM_THREADS=1 # Disable OpenMP CPU affinity - improves compatibility in containerized/virtualized environments @@ -102,12 +101,15 @@ ENV KMP_AFFINITY=disabled # Copy utilities library COPY psi4_utils.py /usr/local/lib/psi4_utils.py +COPY run_psi4_job.py /usr/local/lib/run_psi4_job.py +COPY psi4_state_job.py /usr/local/lib/psi4_state_job.py # Install to Python site-packages RUN cp /usr/local/lib/psi4_utils.py ${MAMBA_ROOT_PREFIX}/envs/psi4env/lib/python3.11/site-packages/ +RUN cp /usr/local/lib/psi4_state_job.py ${MAMBA_ROOT_PREFIX}/envs/psi4env/lib/python3.11/site-packages/ # Create workspace directories -RUN mkdir -p /app/workdir /input /output /tmp/psi4_scratch +RUN mkdir -p /app/workdir /input /output WORKDIR /app/workdir diff --git a/agents/psi4/tools/psi4/psi4_state_job.py b/agents/psi4/tools/psi4/psi4_state_job.py new file mode 100644 index 00000000..2d8150f9 --- /dev/null +++ b/agents/psi4/tools/psi4/psi4_state_job.py @@ -0,0 +1,51 @@ +#!/usr/bin/env python3 +"""Private JSON worker: one fresh Psi4 process per optimization/electronic state.""" +import json +import os +from pathlib import Path +import sys +import traceback + + +def main(): + request_path = Path(sys.argv[1]).resolve() + request = json.loads(request_path.read_text(encoding='utf-8')) + scratch = request_path.parent / 'scratch' + scratch.mkdir(exist_ok=True) + os.environ['PSI_SCRATCH'] = str(scratch) + from psi4_utils import create_molecule, setup_psi4, optimize_geometry, compute_energy + setup_psi4(memory=os.environ.get('PSI4_MEMORY', '4 GB'), + nthreads=int(os.environ.get('OMP_NUM_THREADS', '1')), + output_file=str(request_path.parent / 'psi4.out')) + molecule = create_molecule(request['geometry'], charge=request['charge'], + multiplicity=request['multiplicity'], symmetry='c1') + options = dict(request['options']) + if request['optimize']: + result = optimize_geometry(molecule, request['method'], request['optimization_basis'], + **options) + molecule = result.pop('molecule') + if not result['converged']: + result['status'] = 'failed' + return result + else: + result = {} + options.pop('geom_maxiter', None) + options.pop('full_hess_every', None) + result.update(compute_energy(molecule, request['method'], request['energy_basis'], **options)) + result['status'] = 'success' + return result + + +if __name__ == '__main__': + result_path = Path(sys.argv[1]).resolve().with_name('result.json') + try: + result = main() + result_path.write_text(json.dumps(result, indent=2, allow_nan=False), encoding='utf-8') + sys.exit(0 if result['status'] == 'success' else 1) + except Exception as error: + # Process boundary: persist failure evidence, but never turn an exception into success. + traceback.print_exc() + result_path.write_text(json.dumps({'status': 'failed', 'error': str(error), + 'error_type': type(error).__name__}, indent=2), + encoding='utf-8') + sys.exit(1) diff --git a/agents/psi4/tools/psi4/psi4_utils.py b/agents/psi4/tools/psi4/psi4_utils.py index dba5b0ff..bfb0dc66 100644 --- a/agents/psi4/tools/psi4/psi4_utils.py +++ b/agents/psi4/tools/psi4/psi4_utils.py @@ -33,6 +33,7 @@ import logging import subprocess import multiprocessing +import tempfile from pathlib import Path from typing import Dict, List, Optional, Tuple, Any, Union @@ -51,12 +52,13 @@ _DIRS_CONFIGURED = False # Psi4 specific directories -PSI4_SCRATCH = os.environ.get('PSI_SCRATCH', '/tmp/psi4_scratch') +PSI4_SCRATCH = os.environ.get('PSI_SCRATCH', '') # Physical constants HARTREE_TO_EV = 27.211386245988 HARTREE_TO_KCAL = 627.5094740631 HARTREE_TO_KJ = 2625.4996394799 +DIPOLE_AU_TO_DEBYE = 2.5417464518950259 BOHR_TO_ANGSTROM = 0.529177210903 ANGSTROM_TO_BOHR = 1.8897259886 @@ -105,13 +107,18 @@ def setup_directories(input_dir, output_dir, work_dir=None, copy_input=True): Args: input_dir: Path to the input directory (required). output_dir: Path to the output directory (required). - work_dir: Path to the working directory. If None, works directly in output_dir. - copy_input: If True and work_dir is set, copy input files to work_dir. + work_dir: Working directory; PSI4_WORK_DIR takes precedence. Otherwise unique. + copy_input: If True, copy matching input files to the working directory. """ - global INPUT_DIR, WORK_DIR, OUTPUT_DIR, _DIRS_CONFIGURED - INPUT_DIR = input_dir - OUTPUT_DIR = output_dir - WORK_DIR = work_dir if work_dir is not None else '/app/workdir' + global INPUT_DIR, WORK_DIR, OUTPUT_DIR, PSI4_SCRATCH, _DIRS_CONFIGURED + INPUT_DIR = os.path.abspath(input_dir) + OUTPUT_DIR = os.path.abspath(output_dir) + WORK_DIR = os.path.abspath(os.environ.get('PSI4_WORK_DIR') or work_dir or + tempfile.mkdtemp(prefix='psi4-work-', dir=OUTPUT_DIR + if os.path.isdir(OUTPUT_DIR) else os.getcwd())) + PSI4_SCRATCH = os.path.abspath( + os.environ.get('PSI_SCRATCH') or os.path.join(WORK_DIR, 'scratch')) + os.environ['PSI_SCRATCH'] = PSI4_SCRATCH _DIRS_CONFIGURED = True os.makedirs(WORK_DIR, exist_ok=True) os.makedirs(OUTPUT_DIR, exist_ok=True) @@ -122,7 +129,7 @@ def setup_directories(input_dir, output_dir, work_dir=None, copy_input=True): logging.info(f"Output directory: {OUTPUT_DIR}") logging.info(f"Detected {NUM_CORES} CPU cores available") logging.info(f"Psi4 scratch directory: {PSI4_SCRATCH}") - if copy_input and work_dir is not None: + if copy_input and Path(INPUT_DIR) != Path(WORK_DIR): copy_input_files() @@ -296,8 +303,15 @@ def setup_psi4(memory: str = '4 GB', nthreads: int = None, output_file: str = No Returns: psi4 module reference """ + scratch = (os.environ.get('PSI_SCRATCH') or PSI4_SCRATCH or + tempfile.mkdtemp(prefix='psi4-scratch-', dir=os.getcwd())) + os.makedirs(scratch, exist_ok=True) + with tempfile.TemporaryFile(dir=scratch) as probe: + probe.write(b'psi4 scratch preflight') + os.environ['PSI_SCRATCH'] = os.path.abspath(scratch) import psi4 + psi4.core.IOManager.shared_object().set_default_path(os.path.abspath(scratch)) psi4.set_memory(memory) if nthreads is None: @@ -309,7 +323,6 @@ def setup_psi4(memory: str = '4 GB', nthreads: int = None, output_file: str = No psi4.set_options({ 'scf_type': 'df', # Density fitting for efficiency - 'reference': 'rhf', # Closed-shell reference }) logging.info(f"Psi4 initialized: {memory} memory, {nthreads} threads") @@ -491,6 +504,22 @@ def write_xyz_file(molecule, filename: str, comment: str = ""): # ============= ENERGY CALCULATIONS ============= +def _calculation_options(molecule, basis, options): + """Select a reference per calculation, never inheriting another state's reference.""" + options = {str(key).lower(): value for key, value in options.items()} + multiplicity = molecule.multiplicity() + if not isinstance(multiplicity, int) or isinstance(multiplicity, bool) or multiplicity < 1: + raise ValueError("Multiplicity must be positive.") + reference = str(options.pop('reference', 'rhf' if multiplicity == 1 else 'uhf')).lower() + if reference not in {'rhf', 'uhf', 'rohf'}: + raise ValueError("reference must be rhf, uhf, or rohf (also for DFT).") + if multiplicity != 1 and reference == 'rhf': + raise ValueError("Restricted closed-shell RHF is incompatible with open-shell multiplicity.") + if options.get('fail_on_maxiter') is False: + raise ValueError("fail_on_maxiter=False would accept unconverged SCF states.") + return {'basis': basis, **options, 'reference': reference, 'fail_on_maxiter': True} + + def compute_energy(molecule, method: str = 'hf', basis: str = 'cc-pvdz', return_wfn: bool = False, **kwargs) -> Dict: """Compute single-point energy. @@ -517,7 +546,8 @@ def compute_energy(molecule, method: str = 'hf', basis: str = 'cc-pvdz', """ import psi4 - psi4.set_options({'basis': basis, **kwargs}) + options = _calculation_options(molecule, basis, kwargs) + psi4.set_options(options) logging.info(f"Computing {method.upper()}/{basis} energy...") @@ -526,6 +556,8 @@ def compute_energy(molecule, method: str = 'hf', basis: str = 'cc-pvdz', else: energy = psi4.energy(method, molecule=molecule) wfn = None + if not np.isfinite(energy): + raise ValueError("Psi4 returned a nonfinite electronic energy.") result = { 'energy_hartree': energy, @@ -533,6 +565,7 @@ def compute_energy(molecule, method: str = 'hf', basis: str = 'cc-pvdz', 'energy_kcal': energy * HARTREE_TO_KCAL, 'method': method.upper(), 'basis': basis, + 'reference': options['reference'], } if return_wfn: @@ -624,16 +657,19 @@ def compute_homo_lumo_gap(molecule, method: str = 'b3lyp', basis: str = 'def2-tz """ import psi4 - psi4.set_options({'basis': basis, **kwargs}) + result = compute_energy(molecule, method, basis, return_wfn=True, **kwargs) + return _orbital_gap_result(result['energy_hartree'], result['wfn'], method, basis) - logging.info(f"Computing {method.upper()}/{basis} HOMO-LUMO gap...") - energy, wfn = psi4.energy(method, molecule=molecule, return_wfn=True) +def _orbital_gap_result(energy, wfn, method, basis): + """Format existing alpha-spin orbital-gap contract from a computed wavefunction.""" # Extract orbital energies (handles multi-irrep wavefunctions like D2h benzene) occ_energies, virt_energies = _extract_orbital_energies(wfn) nocc = len(occ_energies) nvirt = len(virt_energies) + if not nocc or not nvirt: + raise ValueError("Orbital gap requires occupied and virtual alpha orbitals.") homo_energy = float(occ_energies[-1]) lumo_energy = float(virt_energies[0]) @@ -690,6 +726,146 @@ def compute_homo_lumo_gap(molecule, method: str = 'b3lyp', basis: str = 'def2-tz return result +def _dipole_result(wfn): + import psi4 + + # A title scopes the documented vector variable to this wavefunction/property call. + psi4.oeprop(wfn, 'DIPOLE', title='DISCOVERY') + vector = np.asarray(wfn.variable('DISCOVERY DIPOLE'), dtype=float).reshape(-1) + if vector.size != 3 or not np.all(np.isfinite(vector)): + raise ValueError("Psi4 returned an invalid dipole vector; expected three finite atomic units.") + debye = vector * DIPOLE_AU_TO_DEBYE + magnitude = float(np.linalg.norm(debye)) + if not np.all(np.isfinite(debye)) or not np.isfinite(magnitude): + raise ValueError("Dipole conversion returned a nonfinite vector or magnitude.") + return { + 'dipole_vector_au': vector.tolist(), + 'dipole_vector_debye': debye.tolist(), + 'dipole_moment_debye': magnitude, + 'dipole_origin_note': 'Psi4 molecular coordinate origin; charged-species dipoles depend on origin.', + } + + +def compute_dipole_moment(molecule, method='b3lyp', basis='def2-tzvp', **kwargs): + """Compute an explicit one-electron dipole from the returned wavefunction.""" + energy = compute_energy(molecule, method, basis, return_wfn=True, **kwargs) + return {**{key: value for key, value in energy.items() if key != 'wfn'}, + **_dipole_result(energy['wfn'])} + + +def compute_electronic_properties(molecule, method='b3lyp', basis='def2-tzvp', **kwargs): + """Compute energy, alpha-spin orbital gap and dipole using a single wavefunction.""" + energy = compute_energy(molecule, method, basis, return_wfn=True, **kwargs) + return {**_orbital_gap_result(energy['energy_hartree'], energy['wfn'], method, basis), + 'reference': energy['reference'], **_dipole_result(energy['wfn'])} + + +def _geometry_text(molecule): + return '\n'.join( + f'{molecule.symbol(i)} {molecule.x(i) * BOHR_TO_ANGSTROM:.12f} ' + f'{molecule.y(i) * BOHR_TO_ANGSTROM:.12f} {molecule.z(i) * BOHR_TO_ANGSTROM:.12f}' + for i in range(molecule.natom())) + + +def _run_reorganization_state(directory, name, request): + """Use argv (never shell interpolation), isolated scratch and persistent state evidence.""" + state_dir = directory / name + state_dir.mkdir() + request_path = state_dir / 'request.json' + request_path.write_text(json.dumps(request, indent=2, allow_nan=False), encoding='utf-8') + with (state_dir / 'worker.out').open('w', encoding='utf-8') as log: + process = subprocess.run([sys.executable, str(Path(__file__).with_name('psi4_state_job.py')), + str(request_path)], cwd=state_dir, stdout=log, + stderr=subprocess.STDOUT, check=False) + result_path = state_dir / 'result.json' + result = json.loads(result_path.read_text(encoding='utf-8')) if result_path.exists() else {} + if process.returncode or result.get('status') != 'success': + raise RuntimeError(f"{name} failed (exit {process.returncode}): " + f"{result.get('error', 'not converged; inspect state artifacts')}") + energy = result.get('energy_hartree') + if energy is None or not np.isfinite(energy): + raise ValueError(f'{name} returned a missing/nonfinite energy.') + return result + + +def compute_reorganization_energies(molecule, candidate_id, method='b3lyp-d3bj', + optimization_basis='def2-svp', energy_basis='def2-tzvp', + neutral_multiplicity=1, cation_multiplicity=2, + anion_multiplicity=2, **kwargs): + """Gas-phase four-point hole/electron reorganization; reoptimize every state. + + Each optimization/state energy runs in a fresh process with its own scratch. + Return JSON-safe candidate-tagged status, all available components, and artifact path. + A failed state stops this candidate without claiming a lambda or converged result. + """ + directory = Path(tempfile.mkdtemp(prefix='reorganization-', dir=WORK_DIR or os.getcwd())) + artifact_parent = OUTPUT_DIR or os.environ.get('PSI4_OUTPUT_DIR') or str(directory.parent) + result = { + 'candidate_id': candidate_id, 'status': 'failed', 'method': method, + 'optimization_basis': optimization_basis, 'energy_basis': energy_basis, + 'energy_units': 'hartree', 'components': {}, 'states': {}, + 'artifact_directory': str(Path(artifact_parent) / directory.name), + } + try: + if molecule.molecular_charge() != 0: + raise ValueError('Four-point input must be a neutral molecule.') + geometry = _geometry_text(molecule) + states = { + 'neutral': (0, neutral_multiplicity), + 'cation': (1, cation_multiplicity), + 'anion': (-1, anion_multiplicity), + } + # Validate all spin/reference choices before spending computation on any state. + for name, (charge, multiplicity) in states.items(): + state = molecule.clone() + state.set_molecular_charge(charge) + state.set_multiplicity(multiplicity) + options = _calculation_options(state, optimization_basis, kwargs) + result['states'][name] = {'charge': charge, 'multiplicity': multiplicity, + 'reference': options['reference'], 'converged': False} + common = {'method': method, 'optimization_basis': optimization_basis, + 'energy_basis': energy_basis, + 'options': {str(key).lower(): value for key, value in kwargs.items()}} + for name, (charge, multiplicity) in states.items(): + state_result = _run_reorganization_state(directory, f'opt_{name}', { + **common, 'geometry': geometry, 'charge': charge, + 'multiplicity': multiplicity, 'optimize': True, + }) + if not state_result.get('converged'): + raise RuntimeError(f'{name} optimization did not converge.') + result['states'][name].update({ + 'converged': True, 'geometry': state_result['optimized_geometry'], + 'geometry_units': 'angstrom', + 'optimization_attempts': state_result['optimization_attempts'], + 'optimization_energy_hartree': state_result['optimized_energy_hartree'], + }) + result['components'][f'E_{name}_at_{name}'] = state_result['energy_hartree'] + for state, at in [('cation', 'neutral'), ('neutral', 'cation'), + ('anion', 'neutral'), ('neutral', 'anion')]: + charge, multiplicity = states[state] + state_result = _run_reorganization_state(directory, f'{state}_at_{at}', { + **common, 'geometry': result['states'][at]['geometry'], + 'charge': charge, 'multiplicity': multiplicity, 'optimize': False, + }) + result['components'][f'E_{state}_at_{at}'] = state_result['energy_hartree'] + energies = result['components'] + for carrier, ion in [('hole', 'cation'), ('electron', 'anion')]: + ion_term = energies[f'E_{ion}_at_neutral'] - energies[f'E_{ion}_at_{ion}'] + neutral_term = energies[f'E_neutral_at_{ion}'] - energies['E_neutral_at_neutral'] + value = ion_term + neutral_term + if not np.isfinite(value): + raise ValueError(f'{carrier} reorganization energy is nonfinite.') + result[f'{carrier}_reorganization_hartree'] = value + result[f'{carrier}_reorganization_eV'] = value * HARTREE_TO_EV + result['status'] = 'success' + except (OSError, ValueError, TypeError, KeyError, RuntimeError, subprocess.SubprocessError) as error: + result['error'] = str(error) + result['error_type'] = type(error).__name__ + (directory / 'reorganization_result.json').write_text( + json.dumps(result, indent=2, allow_nan=False), encoding='utf-8') + return result + + def compute_fundamental_gap(molecule, method: str = 'b3lyp', basis: str = 'def2-tzvp', **kwargs) -> Dict: """Compute the fundamental gap (IP - EA) using the ΔSCF method. @@ -879,11 +1055,14 @@ def optimize_geometry(molecule, method: str = 'hf', basis: str = 'cc-pvdz', """ import psi4 + if isinstance(geom_maxiter, bool) or not isinstance(geom_maxiter, int) or geom_maxiter < 1: + raise ValueError("geom_maxiter must be a positive integer.") + options = _calculation_options(molecule, basis, kwargs) psi4.set_options({ 'basis': basis, 'geom_maxiter': geom_maxiter, 'full_hess_every': full_hess_every, - **kwargs + **options }) logging.info(f"Optimizing geometry with {method.upper()}/{basis}...") @@ -893,14 +1072,25 @@ def optimize_geometry(molecule, method: str = 'hf', basis: str = 'cc-pvdz', initial_energy = psi4.energy(method, molecule=molecule) # Run optimization - try: - opt_energy, opt_wfn = psi4.optimize(method, molecule=molecule, return_wfn=True) - converged = True - except psi4.OptimizationConvergenceError as e: - logging.warning(f"Optimization did not converge: {e}") - opt_energy = psi4.core.variable('CURRENT ENERGY') - converged = False - opt_wfn = None + converged = False + opt_energy = float('nan') + for attempt in range(2): + psi4.set_options({'geom_maxiter': geom_maxiter * (attempt + 1)}) + try: + opt_energy, opt_wfn = psi4.optimize(method, molecule=molecule, return_wfn=True) + molecule = opt_wfn.molecule() + converged = True + break + except psi4.OptimizationConvergenceError as error: + logging.warning("Optimization attempt %s did not converge: %s", attempt + 1, error) + # Psi4 exposes the last iterate on the exception wavefunction. + last_wfn = getattr(error, 'wfn', None) + if last_wfn is None: + break + molecule = last_wfn.molecule() + opt_energy = last_wfn.energy() + if converged and not np.isfinite(opt_energy): + raise ValueError("Optimization returned a nonfinite energy.") # Get optimized geometry opt_mol = molecule # Psi4 modifies molecule in place @@ -918,6 +1108,8 @@ def optimize_geometry(molecule, method: str = 'hf', basis: str = 'cc-pvdz', 'initial_energy_hartree': initial_energy, 'energy_change_kcal': (opt_energy - initial_energy) * HARTREE_TO_KCAL, 'converged': converged, + 'optimization_attempts': attempt + 1, + 'reference': options['reference'], 'optimized_geometry': geometry_str.strip(), 'molecule': opt_mol, } @@ -1582,8 +1774,8 @@ def quick_setup(input_dir, output_dir, work_dir=None, copy_input=True) -> List[s Args: input_dir: Path to the input directory (required). output_dir: Path to the output directory (required). - work_dir: Path to the working directory. If None, works directly in output_dir. - copy_input: If True and work_dir is set, copy input files to work_dir. + work_dir: Working directory; PSI4_WORK_DIR takes precedence. Otherwise unique. + copy_input: If True, copy matching input files to the working directory. Returns: List of copied input files (empty if copy_input is False or work_dir is None) diff --git a/agents/psi4/tools/psi4/run_psi4_job.py b/agents/psi4/tools/psi4/run_psi4_job.py new file mode 100644 index 00000000..73585bce --- /dev/null +++ b/agents/psi4/tools/psi4/run_psi4_job.py @@ -0,0 +1,83 @@ +#!/usr/bin/env python3 +"""Run one script in an isolated workspace; preserve logs/artifacts even on failure.""" +import argparse +import json +import os +from pathlib import Path +import shutil +import subprocess +import sys +import tempfile + + +def copy_artifacts(work, output): + """Copy nested regular files, excluding scratch and symlink traversal.""" + for source in work.rglob('*'): + relative = source.relative_to(work) + if 'scratch' in relative.parts or any(parent.is_symlink() for parent in + [source, *source.parents]): + continue + if source.is_file(): + # Keep the launcher transcript intact if the script also writes logs.out. + target = output / ('script-logs.out' if relative == Path('logs.out') else relative) + target.parent.mkdir(parents=True, exist_ok=True) + shutil.copy2(source, target) + + +def main(argv=None): + parser = argparse.ArgumentParser(description=__doc__) + parser.add_argument('script_path', type=Path) + parser.add_argument('output_dir', type=Path) + args = parser.parse_args(argv) + script = args.script_path.resolve(strict=True) + output = args.output_dir.resolve() + output.mkdir(parents=True, exist_ok=True) + with tempfile.TemporaryDirectory(prefix='psi4-job-', dir=output.parent) as directory: + root = Path(directory) + work, scratch = root / 'work', root / 'scratch' + work.mkdir() + scratch.mkdir() + # Preflight write permissions before importing Psi4 in the child process. + for folder in (work, scratch): + probe = folder / '.write-probe' + probe.write_bytes(b'psi4') + probe.unlink() + env = dict(os.environ, PSI4_WORK_DIR=str(work), PSI4_OUTPUT_DIR=str(output), + PSI_SCRATCH=str(scratch), + TMPDIR=str(scratch), TEMP=str(scratch), TMP=str(scratch), + PYTHONUNBUFFERED='1') + exit_code = 1 + failure = None + with (output / 'logs.out').open('w', encoding='utf-8') as log: + try: + process = subprocess.Popen( + [sys.executable, str(script)], cwd=work, env=env, + stdout=subprocess.PIPE, stderr=subprocess.STDOUT, + text=True, encoding='utf-8', errors='replace', bufsize=1) + for line in process.stdout: + log.write(line) + log.flush() + sys.stdout.write(line) + sys.stdout.flush() + process.stdout.close() + exit_code = process.wait() + except OSError as error: + failure = str(error) + log.write(f'Launcher failed: {error}\n') + finally: + try: + copy_artifacts(work, output) + except (OSError, shutil.Error) as error: + failure = f'Artifact copy failed: {error}' + print(failure, file=sys.stderr) + if exit_code == 0: + exit_code = 1 + if exit_code != 0 or failure: + (output / 'launcher_failure.json').write_text(json.dumps({ + 'status': 'failed', 'exit_code': exit_code, 'error': failure, + }, indent=2), encoding='utf-8') + return exit_code + + +if __name__ == '__main__': + sys.exit(main()) diff --git a/agents/psi4/tools/psi4/smoke_psi4_runtime.py b/agents/psi4/tools/psi4/smoke_psi4_runtime.py new file mode 100644 index 00000000..e42532e1 --- /dev/null +++ b/agents/psi4/tools/psi4/smoke_psi4_runtime.py @@ -0,0 +1,44 @@ +"""Real-Psi4 smoke; run through run_psi4_job.py in the approved built runtime.""" +import json +import math +import os +from pathlib import Path + +from psi4_utils import ( + setup_psi4, create_molecule, compute_energy, compute_electronic_properties, + compute_reorganization_energies, +) + + +def main(): + setup_psi4(memory='1 GB', nthreads=1, output_file='smoke.out') + water = create_molecule('O 0 0 0\nH 0 0.757 0.587\nH 0 -0.757 0.587', symmetry='c1') + properties = compute_electronic_properties(water, method='hf', basis='sto-3g') + assert 0.1 < properties['dipole_moment_debye'] < 5 + assert len(properties['dipole_vector_debye']) == 3 + assert all(math.isfinite(v) for v in properties['dipole_vector_debye']) + # Water doublet ions exercise unrestricted charged SCF with a diffuse basis. + for charge in (1, -1): + ion = water.clone() + ion.set_molecular_charge(charge) + ion.set_multiplicity(2) + energy = compute_energy(ion, method='hf', basis='6-31+g*') + assert energy['reference'] == 'uhf' + assert math.isfinite(energy['energy_hartree']) + # LiH has bound neutral/cation/anion states and a tiny nuclear geometry. + molecule = create_molecule('Li 0 0 0\nH 0 0 1.6', symmetry='c1') + reorganization = compute_reorganization_energies( + molecule, 'smoke-LiH', method='hf', optimization_basis='6-31+g*', + energy_basis='6-31+g*', geom_maxiter=50) + Path('smoke_results.json').write_text(json.dumps( + {'properties': properties, 'reorganization': reorganization}, indent=2, allow_nan=False)) + assert reorganization['status'] == 'success', reorganization + assert all(state['converged'] for state in reorganization['states'].values()) + assert len(reorganization['components']) == 7 + assert all(math.isfinite(value) for value in reorganization['components'].values()) + assert all(math.isfinite(reorganization[f'{carrier}_reorganization_eV']) + for carrier in ('hole', 'electron')) + + +if __name__ == '__main__': + main() diff --git a/agents/psi4/tools/psi4/test_psi4_reliability.py b/agents/psi4/tools/psi4/test_psi4_reliability.py new file mode 100644 index 00000000..5a429522 --- /dev/null +++ b/agents/psi4/tools/psi4/test_psi4_reliability.py @@ -0,0 +1,323 @@ +"""Targeted stdlib unittest suite; Psi4 is mocked, NumPy is the existing dependency.""" +import copy +import json +import os +from pathlib import Path +import subprocess +import sys +import tempfile +import unittest +from unittest.mock import MagicMock, patch + +import numpy as np +import psi4_utils as utils +import run_psi4_job as launcher + + +class Molecule: + def __init__(self, charge=0, mult=1, x=0.0): + self.charge, self.mult, self.position = charge, mult, x + + def multiplicity(self): + return self.mult + + def molecular_charge(self): + return self.charge + + def set_molecular_charge(self, value): + self.charge = value + + def set_multiplicity(self, value): + self.mult = value + + def clone(self): + return copy.copy(self) + + def natom(self): + return 1 + + def symbol(self, index): + return 'He' + + def x(self, index): + return self.position + + def y(self, index): + return 0.0 + + def z(self, index): + return 0.0 + + +class OptimizationError(Exception): + def __init__(self, wfn=None): + self.wfn = wfn + + +class UtilityTests(unittest.TestCase): + def setUp(self): + self.psi4 = MagicMock() + self.psi4.OptimizationConvergenceError = OptimizationError + self.patch = patch.dict(sys.modules, {'psi4': self.psi4}) + self.patch.start() + self.addCleanup(self.patch.stop) + output_patch = patch.object(utils, 'OUTPUT_DIR', '') + output_patch.start() + self.addCleanup(output_patch.stop) + self.molecule = Molecule() + self.wfn = MagicMock() + self.wfn.nirrep.return_value = 1 + self.wfn.epsilon_a.return_value = [-0.5, 0.1] + self.wfn.nalpha.return_value = 1 + self.wfn.variable.return_value = np.array([[0.0, 0.0, 1.0]]) + self.wfn.molecule.return_value = self.molecule + self.wfn.energy.return_value = -1.0 + self.psi4.energy.return_value = (-1.0, self.wfn) + + def test_reference_selected_for_each_energy_including_dft(self): + for mult, reference in [(1, 'rhf'), (2, 'uhf'), (1, 'rhf')]: + utils.compute_energy(Molecule(mult=mult), 'b3lyp', return_wfn=True) + self.assertEqual(self.psi4.set_options.call_args.args[0]['reference'], reference) + + def test_explicit_reference_validated_before_energy(self): + for reference in ['rhf', 'rks', 'invalid']: + with self.assertRaises(ValueError): + utils.compute_energy(Molecule(mult=2), reference=reference) + self.psi4.energy.assert_not_called() + utils.compute_energy(Molecule(mult=2), reference='ROHF', return_wfn=True) + self.assertEqual(self.psi4.set_options.call_args.args[0]['reference'], 'rohf') + + def test_combined_properties_one_wavefunction_and_debye_conversion(self): + result = utils.compute_electronic_properties(self.molecule) + self.psi4.energy.assert_called_once() + self.psi4.oeprop.assert_called_once_with(self.wfn, 'DIPOLE', title='DISCOVERY') + self.wfn.variable.assert_called_once_with('DISCOVERY DIPOLE') + self.assertAlmostEqual(result['dipole_moment_debye'], utils.DIPOLE_AU_TO_DEBYE) + self.assertEqual(result['dipole_vector_debye'], [0.0, 0.0, utils.DIPOLE_AU_TO_DEBYE]) + self.assertAlmostEqual(result['gap_hartree'], 0.6) + + def test_gap_automatically_selects_uhf(self): + result = utils.compute_homo_lumo_gap(Molecule(mult=2)) + self.assertAlmostEqual(result['gap_hartree'], 0.6) + self.assertEqual(self.psi4.set_options.call_args.args[0]['reference'], 'uhf') + + def test_dipole_invalid_vectors_rejected(self): + for vector in [[1, 2], [0, float('nan'), 0], [float('inf'), 0, 0]]: + self.wfn.variable.return_value = vector + with self.assertRaises(ValueError): + utils.compute_dipole_moment(self.molecule) + + def test_empty_orbitals_rejected(self): + self.wfn.nalpha.return_value = 0 + with self.assertRaises(ValueError): + utils.compute_electronic_properties(self.molecule) + + def test_optimization_retries_once_from_last_geometry_without_weakening(self): + self.psi4.energy.return_value = -0.9 + last = Molecule(mult=2, x=1.0) + self.wfn.molecule.return_value = last + self.psi4.optimize.side_effect = [OptimizationError(self.wfn), (-1.1, self.wfn)] + result = utils.optimize_geometry(Molecule(mult=2), 'b3lyp', geom_maxiter=5, + g_convergence='gau_tight') + self.assertTrue(result['converged']) + self.assertEqual(result['optimization_attempts'], 2) + self.assertIs(self.psi4.optimize.call_args.kwargs['molecule'], last) + settings = [call.args[0] for call in self.psi4.set_options.call_args_list] + self.assertEqual(settings[-2:], [{'geom_maxiter': 5}, {'geom_maxiter': 10}]) + self.assertEqual(settings[0]['g_convergence'], 'gau_tight') + self.assertEqual(settings[0]['reference'], 'uhf') + + def test_optimization_exhaustion_is_not_success(self): + self.psi4.energy.return_value = -0.9 + self.psi4.optimize.side_effect = OptimizationError(self.wfn) + result = utils.optimize_geometry(self.molecule) + self.assertFalse(result['converged']) + self.assertEqual(self.psi4.optimize.call_count, 2) + + def test_optimization_without_last_wavefunction_does_not_retry_initial_geometry(self): + self.psi4.energy.return_value = -0.9 + self.psi4.optimize.side_effect = OptimizationError() + self.assertFalse(utils.optimize_geometry(self.molecule)['converged']) + self.assertEqual(self.psi4.optimize.call_count, 1) + + def test_scf_failure_not_swallowed(self): + self.psi4.energy.side_effect = RuntimeError('SCF failed') + with self.assertRaisesRegex(RuntimeError, 'SCF failed'): + utils.compute_electronic_properties(self.molecule) + + def test_nonfinite_energy_and_unconverged_acceptance_rejected(self): + self.psi4.energy.return_value = float('nan') + with self.assertRaises(ValueError): + utils.compute_energy(self.molecule) + with self.assertRaises(ValueError): + utils.compute_energy(self.molecule, fail_on_maxiter=False) + + def test_uppercase_restricted_override_cannot_bypass_validation(self): + with self.assertRaises(ValueError): + utils.compute_energy(Molecule(mult=2), REFERENCE='RHF') + + def test_setup_binds_iomanager_scratch(self): + with tempfile.TemporaryDirectory(dir=Path.cwd()) as directory: + with patch.dict(os.environ, {'PSI_SCRATCH': directory}): + utils.setup_psi4() + self.psi4.core.IOManager.shared_object().set_default_path.assert_called_with( + os.path.abspath(directory)) + + +class ReorganizationTests(unittest.TestCase): + def setUp(self): + self.directory = tempfile.TemporaryDirectory(dir=Path.cwd()) + self.addCleanup(self.directory.cleanup) + self.patch = patch.object(utils, 'WORK_DIR', self.directory.name) + self.patch.start() + self.addCleanup(self.patch.stop) + output_patch = patch.object(utils, 'OUTPUT_DIR', '') + output_patch.start() + self.addCleanup(output_patch.stop) + + @staticmethod + def state_result(directory, name, request): + energies = {'opt_neutral': -10.0, 'opt_cation': -9.0, 'opt_anion': -10.2, + 'cation_at_neutral': -8.9, 'neutral_at_cation': -9.95, + 'anion_at_neutral': -10.1, 'neutral_at_anion': -9.97} + return {'status': 'success', 'converged': True, 'energy_hartree': energies[name], + 'optimized_geometry': 'He 0 0 0', 'optimization_attempts': 1, + 'optimized_energy_hartree': energies[name]} + + def test_four_point_formulas_ids_units_states_and_seven_fresh_jobs(self): + with patch.object(utils, '_run_reorganization_state', side_effect=self.state_result) as run: + result = utils.compute_reorganization_energies(Molecule(), 'candidate-A') + self.assertEqual(result['status'], 'success') + self.assertEqual(result['candidate_id'], 'candidate-A') + self.assertAlmostEqual(result['hole_reorganization_hartree'], 0.15) + self.assertAlmostEqual(result['electron_reorganization_hartree'], 0.13) + self.assertAlmostEqual(result['electron_reorganization_eV'], 0.13 * utils.HARTREE_TO_EV) + self.assertEqual(result['energy_units'], 'hartree') + self.assertEqual(run.call_count, 7) + self.assertEqual([call.args[2]['charge'] for call in run.call_args_list], + [0, 1, -1, 1, 0, -1, 0]) + self.assertEqual([call.args[2]['multiplicity'] for call in run.call_args_list], + [1, 2, 2, 2, 1, 2, 1]) + self.assertEqual([state['reference'] for state in result['states'].values()], + ['rhf', 'uhf', 'uhf']) + self.assertEqual(len(result['components']), 7) + self.assertTrue((Path(result['artifact_directory']) / 'reorganization_result.json').exists()) + + def test_failure_preserves_id_partial_components_and_json(self): + def run(directory, name, request): + if name == 'opt_cation': + raise RuntimeError('unconverged') + return self.state_result(directory, name, request) + with patch.object(utils, '_run_reorganization_state', side_effect=run): + result = utils.compute_reorganization_energies(Molecule(), 'candidate-B') + self.assertEqual(result['status'], 'failed') + self.assertEqual(result['components'], {'E_neutral_at_neutral': -10.0}) + self.assertNotIn('hole_reorganization_eV', result) + self.assertIn('unconverged', result['error']) + saved = json.loads((Path(result['artifact_directory']) / 'reorganization_result.json').read_text()) + self.assertEqual(saved['candidate_id'], 'candidate-B') + + def test_upstream_partial_neutral_never_assumed_converged(self): + with patch.object(utils, '_run_reorganization_state', + return_value={'converged': False}) as run: + result = utils.compute_reorganization_energies(Molecule(), 'candidate-C') + self.assertEqual(run.call_args.args[1], 'opt_neutral') + self.assertEqual(result['status'], 'failed') + self.assertFalse(result['states']['neutral']['converged']) + + def test_restricted_open_shell_rejected_before_jobs(self): + with patch.object(utils, '_run_reorganization_state') as run: + result = utils.compute_reorganization_energies(Molecule(), 'candidate-D', reference='rhf') + self.assertEqual(result['status'], 'failed') + run.assert_not_called() + + def test_worker_nonzero_exit_with_result_still_fails(self): + directory = Path(self.directory.name) + with patch.object(utils.subprocess, 'run', return_value=subprocess.CompletedProcess([], 7)): + with self.assertRaisesRegex(RuntimeError, 'exit 7'): + utils._run_reorganization_state(directory, 'bad_state', {'options': {}}) + self.assertTrue((directory / 'bad_state' / 'worker.out').exists()) + + +class LauncherTests(unittest.TestCase): + def setUp(self): + self.directory = tempfile.TemporaryDirectory(dir=Path.cwd(), prefix='launcher-tests-') + self.root = Path(self.directory.name) + self.addCleanup(self.directory.cleanup) + + def start(self, script, output): + return subprocess.Popen([sys.executable, str(Path(launcher.__file__)), str(script), str(output)], + stdout=subprocess.PIPE, stderr=subprocess.STDOUT, text=True) + + def test_concurrent_invocations_isolate_scratch_work_and_copy_nested_artifacts(self): + script = self.root / 'script ; & spaces.py' + script.write_text( + "import json,os,time\nfrom pathlib import Path\n" + "Path('nested').mkdir()\n" + "Path('nested/scratch').mkdir()\n" + "Path('nested/scratch/secret.dat').write_text('nested scratch')\n" + "Path('logs.out').write_text('script own transcript')\n" + "Path('nested/result.json').write_text(json.dumps(dict(work=os.getcwd()," + "scratch=os.environ['PSI_SCRATCH'])))\n" + "Path(os.environ['PSI_SCRATCH'],'secret.dat').write_text('scratch')\n" + "print('stdout evidence', flush=True)\n" + "import sys; print('stderr evidence',file=sys.stderr,flush=True)\n" + "time.sleep(0.1)\n", encoding='utf-8') + outputs = [self.root / 'output one', self.root / 'output two'] + processes = [self.start(script, output) for output in outputs] + for process in processes: + stdout, _ = process.communicate(timeout=30) + self.assertEqual(process.returncode, 0, stdout) + results = [json.loads((output / 'nested' / 'result.json').read_text()) for output in outputs] + self.assertNotEqual(results[0]['work'], results[1]['work']) + self.assertNotEqual(results[0]['scratch'], results[1]['scratch']) + for output in outputs: + self.assertIn('stderr evidence', (output / 'logs.out').read_text()) + self.assertEqual((output / 'script-logs.out').read_text(), 'script own transcript') + self.assertFalse(list(output.rglob('secret.dat'))) + + def test_child_exit_preserved_and_partial_artifacts_retained(self): + script = self.root / 'fail.py' + script.write_text("from pathlib import Path\nPath('partial.json').write_text('{}')\n" + "print('failing evidence')\nraise SystemExit(7)\n") + output = self.root / 'output' + process = self.start(script, output) + stdout, _ = process.communicate(timeout=30) + self.assertEqual(process.returncode, 7, stdout) + self.assertTrue((output / 'partial.json').exists()) + self.assertEqual(json.loads((output / 'launcher_failure.json').read_text())['exit_code'], 7) + + def test_copy_failure_converts_success_to_failure_and_is_recorded(self): + script = self.root / 'ok.py' + script.write_text("print('ok')\n") + output = self.root / 'output' + with patch.object(launcher, 'copy_artifacts', side_effect=OSError('copy denied')): + self.assertEqual(launcher.main([str(script), str(output)]), 1) + self.assertIn('copy denied', (output / 'launcher_failure.json').read_text()) + + def test_copy_failure_does_not_hide_child_exit(self): + script = self.root / 'bad.py' + script.write_text('raise SystemExit(7)\n') + output = self.root / 'output' + with patch.object(launcher, 'copy_artifacts', side_effect=OSError('copy denied')): + self.assertEqual(launcher.main([str(script), str(output)]), 7) + + def test_setup_honors_work_env_without_replacing_input_output_paths(self): + original = Path.cwd() + input_dir, output_dir, work_dir = [self.root / name for name in ('input', 'out', 'work')] + input_dir.mkdir() + (input_dir / 'molecule.xyz').write_text('input') + try: + with patch.dict(os.environ, {'PSI4_WORK_DIR': str(work_dir), + 'PSI_SCRATCH': str(self.root / 'scratch')}): + utils.setup_directories(str(input_dir), str(output_dir), work_dir='ignored') + self.assertEqual(utils.INPUT_DIR, str(input_dir)) + self.assertEqual(utils.OUTPUT_DIR, str(output_dir)) + self.assertEqual(utils.WORK_DIR, str(work_dir)) + self.assertEqual((work_dir / 'molecule.xyz').read_text(), 'input') + finally: + os.chdir(original) + + +if __name__ == '__main__': + unittest.main() diff --git a/agents/psi4/tools/psi4/tool.yaml b/agents/psi4/tools/psi4/tool.yaml index 4b15f365..1593eec2 100644 --- a/agents/psi4/tools/psi4/tool.yaml +++ b/agents/psi4/tools/psi4/tool.yaml @@ -37,7 +37,7 @@ infra: code_environments: - language: python3 - command: "python3 -u \"/{{scriptName}}\" > {{output_directory}}/logs.out 2>&1; cp -r /app/workdir/* {{output_directory}}/ 2>/dev/null || true; tail -n 100 {{output_directory}}/logs.out" + command: "/opt/conda/envs/psi4env/bin/python -u /usr/local/lib/run_psi4_job.py \"/{{scriptName}}\" \"{{output_directory}}\"" description: | Python code environment with Psi4 for ab initio quantum chemistry calculations. Use this for complete automated workflows including single-point energies, From f48063a4c1ff936b8ad1467462f15737584ad1fa Mon Sep 17 00:00:00 2001 From: preetibayas Date: Wed, 28 Oct 2026 18:26:41 +0530 Subject: [PATCH 2/3] Add genuine UMAP support to molToolkit workflows Install umap-learn in molenv, require a real ten-candidate precomputed-distance embedding during image builds, and document algorithm fidelity, reproducibility, and coordinated image/agent updates. Co-authored-by: Copilot <223556219+Copilot@users.noreply.github.com> Copilot-Session: 87cf205e-4bdc-4491-9fdd-cbec1660935f --- agents/mol-toolkit/README.md | 35 +++++++++++++++++++ agents/mol-toolkit/agent.yaml | 4 ++- .../mol-toolkit/tools/mol-toolkit/Dockerfile | 3 ++ .../tools/mol-toolkit/app/smoke_umap.py | 33 +++++++++++++++++ 4 files changed, 74 insertions(+), 1 deletion(-) create mode 100644 agents/mol-toolkit/tools/mol-toolkit/app/smoke_umap.py diff --git a/agents/mol-toolkit/README.md b/agents/mol-toolkit/README.md index 3d122d8d..d145ec30 100644 --- a/agents/mol-toolkit/README.md +++ b/agents/mol-toolkit/README.md @@ -20,6 +20,41 @@ Before starting the deployment, ensure you have: 3. Docker installed locally for image building 4. Azure CLI or PowerShell for resource management +## UMAP dependency and verification + +The `molenv` environment includes `umap-learn` (imported as `from umap import +UMAP`). The Dockerfile runs `app/smoke_umap.py` during the build: it computes real +UMAP coordinates from a ten-row fingerprint Tanimoto-distance matrix and rejects +non-finite coordinates or an output shape other than `(10, 2)`. This is an +embedding execution check, not just an import check or an MDS/PCA substitute. + +In the matching tool environment, the check can also be run without rebuilding: + +```bash +python /app/smoke_umap.py +``` + +Success prints `UMAP_SMOKE_OK`, the installed package version, and `(10, 2)`. +The check requires the image's Python dependencies; compiling the script alone +does not validate UMAP execution. It does not validate scientific cluster quality. + +| Prompt behavior before | Prompt behavior after | +|---|---| +| UMAP was not in the agent's available-package list | The agent explicitly uses `umap-learn` | +| Missing UMAP could lead to an MDS fallback | MDS, PCA, and t-SNE cannot be reported as UMAP | +| Precomputed similarities could be confused with distances | Tanimoto similarities are validated and converted to distances before UMAP | + +For example, "visualize ten candidate fingerprints using UMAP" now requires a +real two-dimensional UMAP embedding, a plot, saved coordinates, package version, +and parameters. A failed import or embedding must be reported as a failure with +diagnostics, not as a completed substitute visualization. + +These are source changes, not an update to any existing image or deployed agent. +Obtain approval before building, publishing, retagging, or replacing an image. +After approval, build and scan the matching image first, verify the smoke check, +then update the deployed tool and agent together before rerunning affected tasks. +Changing agent instructions alone cannot add the missing dependency to an image. + ## Build Docker Image ### Step 1: Build and Publish Docker Image diff --git a/agents/mol-toolkit/agent.yaml b/agents/mol-toolkit/agent.yaml index ab25cdf7..d792554f 100644 --- a/agents/mol-toolkit/agent.yaml +++ b/agents/mol-toolkit/agent.yaml @@ -23,11 +23,13 @@ instructions: | The container image has pre-built molecular chemistry functions exposed as Python modules in `/app`. You should **ALWAYS prioritize using these pre-built functions** by importing them directly from their modules to achieve user goals. Only generate custom Python scripts when absolutely necessary - when no combination of pre-built functions can accomplish the task. ***Open source Python packages available in the MolToolkit container:*** - ase, biopython, matplotlib, MDAnalysis, numpy, openbabel, openpyxl, pandas, pyarrow, pymatgen, pymol-open-source, rdkit, scanpy, scikit-bio, scikit-learn, scipy, seaborn, sympy + ase, biopython, matplotlib, MDAnalysis, numpy, openbabel, openpyxl, pandas, pyarrow, pymatgen, pymol-open-source, rdkit, scanpy, scikit-bio, scikit-learn, scipy, seaborn, sympy, umap-learn You cannot use any other open-source packages that are not listed here as they are not available in your environment. **IMPORTANT scikit-learn API note:** In scikit-learn >= 1.4, `AgglomerativeClustering` renamed the `affinity` parameter to `metric`. Always use `metric=` (e.g., `AgglomerativeClustering(metric='precomputed', linkage='average')`). Using `affinity=` will raise a TypeError. + **UMAP requirements:** For a requested UMAP embedding, use `from umap import UMAP` from the installed `umap-learn` package. Never substitute MDS, PCA, or t-SNE and label it UMAP. For a Tanimoto similarity matrix, validate symmetry, unit diagonal, finite values, and values in [0, 1], then use `distances = 1.0 - similarities` with a zero diagonal and `metric="precomputed"`; UMAP expects distances, not similarities. For at least three candidates, set `n_components=2`, `n_neighbors=min(15, n_candidates - 1)`, `init="random"`, `random_state=42`, and `n_jobs=1`. For fewer than three candidates, report insufficient data instead of inventing coordinates. Validate that the embedding contains one finite 2D coordinate per candidate, preserve candidate ordering, and save the coordinates, requested plot, package version, and parameters in the workflow outputs. If the import or fit fails, preserve diagnostics and report failure; do not silently switch algorithms or declare UMAP complete. + **CODE EXECUTION METHOD**: Import functions directly from `/app` modules (e.g., `from molecular_utils import function_name`) for maximum flexibility, seamless function composition, better error handling, and custom workflow creation. Unless the user explicitly specifies working with a single file or specific files, scripts should iterate through ALL relevant files in the input directory (path from dataHandlingContext) using functions from `app/io_utils.py` (like `find_smiles_files()` and `read_smiles_from_file()`). # MANDATORY INSTRUCTIONS: diff --git a/agents/mol-toolkit/tools/mol-toolkit/Dockerfile b/agents/mol-toolkit/tools/mol-toolkit/Dockerfile index d4f8aa0f..7e6791f9 100644 --- a/agents/mol-toolkit/tools/mol-toolkit/Dockerfile +++ b/agents/mol-toolkit/tools/mol-toolkit/Dockerfile @@ -48,6 +48,7 @@ RUN mamba create -y -n molenv -c conda-forge \ scipy \ seaborn \ sympy \ + umap-learn \ && mamba clean -a -y # ---- Patch base conda env (Trivy scans all site-packages dirs) ---- @@ -83,6 +84,8 @@ ENV PATH=/opt/conda/envs/molenv/bin:/opt/conda/bin:$PATH ENV PYTHONPATH=/app ENV PYTHONUNBUFFERED=1 +RUN mamba run -n molenv python /app/smoke_umap.py + # Create common mount points / working dirs (single layer) RUN mkdir -p /workdir diff --git a/agents/mol-toolkit/tools/mol-toolkit/app/smoke_umap.py b/agents/mol-toolkit/tools/mol-toolkit/app/smoke_umap.py new file mode 100644 index 00000000..63e7f2ff --- /dev/null +++ b/agents/mol-toolkit/tools/mol-toolkit/app/smoke_umap.py @@ -0,0 +1,33 @@ +"""Verify real UMAP execution on the ten-candidate Tanimoto-distance workflow.""" + +from importlib.metadata import version + +import numpy as np +from umap import UMAP + + +def main() -> None: + fingerprints = np.random.RandomState(42).randint(0, 2, size=(10, 64)) + intersections = fingerprints @ fingerprints.T + counts = np.diag(intersections) + unions = counts[:, None] + counts[None, :] - intersections + similarities = intersections / unions + distances = 1.0 - similarities + np.fill_diagonal(distances, 0.0) + + embedding = UMAP( + n_components=2, + metric="precomputed", + n_neighbors=min(15, distances.shape[0] - 1), + init="random", + random_state=42, + n_jobs=1, + ).fit_transform(distances) + + if embedding.shape != (10, 2) or not np.isfinite(embedding).all(): + raise RuntimeError("UMAP must produce ten finite two-dimensional coordinates") + print(f"UMAP_SMOKE_OK: umap-learn={version('umap-learn')}, shape={embedding.shape}") + + +if __name__ == "__main__": + main() From 2e9a3e1a921bf197f98f28cfac5028506efbfd82 Mon Sep 17 00:00:00 2001 From: preetibayas Date: Wed, 28 Oct 2026 18:28:32 +0530 Subject: [PATCH 3/3] Normalize new Psi4 Python scripts to LF line endings Keep Linux launcher shebangs and the full branch whitespace check valid without changing script behavior. Co-authored-by: Copilot <223556219+Copilot@users.noreply.github.com> Copilot-Session: 87cf205e-4bdc-4491-9fdd-cbec1660935f --- agents/psi4/tools/psi4/psi4_state_job.py | 102 +-- agents/psi4/tools/psi4/run_psi4_job.py | 166 ++--- agents/psi4/tools/psi4/smoke_psi4_runtime.py | 88 +-- .../psi4/tools/psi4/test_psi4_reliability.py | 646 +++++++++--------- 4 files changed, 501 insertions(+), 501 deletions(-) diff --git a/agents/psi4/tools/psi4/psi4_state_job.py b/agents/psi4/tools/psi4/psi4_state_job.py index 2d8150f9..cd3bf70e 100644 --- a/agents/psi4/tools/psi4/psi4_state_job.py +++ b/agents/psi4/tools/psi4/psi4_state_job.py @@ -1,51 +1,51 @@ -#!/usr/bin/env python3 -"""Private JSON worker: one fresh Psi4 process per optimization/electronic state.""" -import json -import os -from pathlib import Path -import sys -import traceback - - -def main(): - request_path = Path(sys.argv[1]).resolve() - request = json.loads(request_path.read_text(encoding='utf-8')) - scratch = request_path.parent / 'scratch' - scratch.mkdir(exist_ok=True) - os.environ['PSI_SCRATCH'] = str(scratch) - from psi4_utils import create_molecule, setup_psi4, optimize_geometry, compute_energy - setup_psi4(memory=os.environ.get('PSI4_MEMORY', '4 GB'), - nthreads=int(os.environ.get('OMP_NUM_THREADS', '1')), - output_file=str(request_path.parent / 'psi4.out')) - molecule = create_molecule(request['geometry'], charge=request['charge'], - multiplicity=request['multiplicity'], symmetry='c1') - options = dict(request['options']) - if request['optimize']: - result = optimize_geometry(molecule, request['method'], request['optimization_basis'], - **options) - molecule = result.pop('molecule') - if not result['converged']: - result['status'] = 'failed' - return result - else: - result = {} - options.pop('geom_maxiter', None) - options.pop('full_hess_every', None) - result.update(compute_energy(molecule, request['method'], request['energy_basis'], **options)) - result['status'] = 'success' - return result - - -if __name__ == '__main__': - result_path = Path(sys.argv[1]).resolve().with_name('result.json') - try: - result = main() - result_path.write_text(json.dumps(result, indent=2, allow_nan=False), encoding='utf-8') - sys.exit(0 if result['status'] == 'success' else 1) - except Exception as error: - # Process boundary: persist failure evidence, but never turn an exception into success. - traceback.print_exc() - result_path.write_text(json.dumps({'status': 'failed', 'error': str(error), - 'error_type': type(error).__name__}, indent=2), - encoding='utf-8') - sys.exit(1) +#!/usr/bin/env python3 +"""Private JSON worker: one fresh Psi4 process per optimization/electronic state.""" +import json +import os +from pathlib import Path +import sys +import traceback + + +def main(): + request_path = Path(sys.argv[1]).resolve() + request = json.loads(request_path.read_text(encoding='utf-8')) + scratch = request_path.parent / 'scratch' + scratch.mkdir(exist_ok=True) + os.environ['PSI_SCRATCH'] = str(scratch) + from psi4_utils import create_molecule, setup_psi4, optimize_geometry, compute_energy + setup_psi4(memory=os.environ.get('PSI4_MEMORY', '4 GB'), + nthreads=int(os.environ.get('OMP_NUM_THREADS', '1')), + output_file=str(request_path.parent / 'psi4.out')) + molecule = create_molecule(request['geometry'], charge=request['charge'], + multiplicity=request['multiplicity'], symmetry='c1') + options = dict(request['options']) + if request['optimize']: + result = optimize_geometry(molecule, request['method'], request['optimization_basis'], + **options) + molecule = result.pop('molecule') + if not result['converged']: + result['status'] = 'failed' + return result + else: + result = {} + options.pop('geom_maxiter', None) + options.pop('full_hess_every', None) + result.update(compute_energy(molecule, request['method'], request['energy_basis'], **options)) + result['status'] = 'success' + return result + + +if __name__ == '__main__': + result_path = Path(sys.argv[1]).resolve().with_name('result.json') + try: + result = main() + result_path.write_text(json.dumps(result, indent=2, allow_nan=False), encoding='utf-8') + sys.exit(0 if result['status'] == 'success' else 1) + except Exception as error: + # Process boundary: persist failure evidence, but never turn an exception into success. + traceback.print_exc() + result_path.write_text(json.dumps({'status': 'failed', 'error': str(error), + 'error_type': type(error).__name__}, indent=2), + encoding='utf-8') + sys.exit(1) diff --git a/agents/psi4/tools/psi4/run_psi4_job.py b/agents/psi4/tools/psi4/run_psi4_job.py index 73585bce..22b752dd 100644 --- a/agents/psi4/tools/psi4/run_psi4_job.py +++ b/agents/psi4/tools/psi4/run_psi4_job.py @@ -1,83 +1,83 @@ -#!/usr/bin/env python3 -"""Run one script in an isolated workspace; preserve logs/artifacts even on failure.""" -import argparse -import json -import os -from pathlib import Path -import shutil -import subprocess -import sys -import tempfile - - -def copy_artifacts(work, output): - """Copy nested regular files, excluding scratch and symlink traversal.""" - for source in work.rglob('*'): - relative = source.relative_to(work) - if 'scratch' in relative.parts or any(parent.is_symlink() for parent in - [source, *source.parents]): - continue - if source.is_file(): - # Keep the launcher transcript intact if the script also writes logs.out. - target = output / ('script-logs.out' if relative == Path('logs.out') else relative) - target.parent.mkdir(parents=True, exist_ok=True) - shutil.copy2(source, target) - - -def main(argv=None): - parser = argparse.ArgumentParser(description=__doc__) - parser.add_argument('script_path', type=Path) - parser.add_argument('output_dir', type=Path) - args = parser.parse_args(argv) - script = args.script_path.resolve(strict=True) - output = args.output_dir.resolve() - output.mkdir(parents=True, exist_ok=True) - with tempfile.TemporaryDirectory(prefix='psi4-job-', dir=output.parent) as directory: - root = Path(directory) - work, scratch = root / 'work', root / 'scratch' - work.mkdir() - scratch.mkdir() - # Preflight write permissions before importing Psi4 in the child process. - for folder in (work, scratch): - probe = folder / '.write-probe' - probe.write_bytes(b'psi4') - probe.unlink() - env = dict(os.environ, PSI4_WORK_DIR=str(work), PSI4_OUTPUT_DIR=str(output), - PSI_SCRATCH=str(scratch), - TMPDIR=str(scratch), TEMP=str(scratch), TMP=str(scratch), - PYTHONUNBUFFERED='1') - exit_code = 1 - failure = None - with (output / 'logs.out').open('w', encoding='utf-8') as log: - try: - process = subprocess.Popen( - [sys.executable, str(script)], cwd=work, env=env, - stdout=subprocess.PIPE, stderr=subprocess.STDOUT, - text=True, encoding='utf-8', errors='replace', bufsize=1) - for line in process.stdout: - log.write(line) - log.flush() - sys.stdout.write(line) - sys.stdout.flush() - process.stdout.close() - exit_code = process.wait() - except OSError as error: - failure = str(error) - log.write(f'Launcher failed: {error}\n') - finally: - try: - copy_artifacts(work, output) - except (OSError, shutil.Error) as error: - failure = f'Artifact copy failed: {error}' - print(failure, file=sys.stderr) - if exit_code == 0: - exit_code = 1 - if exit_code != 0 or failure: - (output / 'launcher_failure.json').write_text(json.dumps({ - 'status': 'failed', 'exit_code': exit_code, 'error': failure, - }, indent=2), encoding='utf-8') - return exit_code - - -if __name__ == '__main__': - sys.exit(main()) +#!/usr/bin/env python3 +"""Run one script in an isolated workspace; preserve logs/artifacts even on failure.""" +import argparse +import json +import os +from pathlib import Path +import shutil +import subprocess +import sys +import tempfile + + +def copy_artifacts(work, output): + """Copy nested regular files, excluding scratch and symlink traversal.""" + for source in work.rglob('*'): + relative = source.relative_to(work) + if 'scratch' in relative.parts or any(parent.is_symlink() for parent in + [source, *source.parents]): + continue + if source.is_file(): + # Keep the launcher transcript intact if the script also writes logs.out. + target = output / ('script-logs.out' if relative == Path('logs.out') else relative) + target.parent.mkdir(parents=True, exist_ok=True) + shutil.copy2(source, target) + + +def main(argv=None): + parser = argparse.ArgumentParser(description=__doc__) + parser.add_argument('script_path', type=Path) + parser.add_argument('output_dir', type=Path) + args = parser.parse_args(argv) + script = args.script_path.resolve(strict=True) + output = args.output_dir.resolve() + output.mkdir(parents=True, exist_ok=True) + with tempfile.TemporaryDirectory(prefix='psi4-job-', dir=output.parent) as directory: + root = Path(directory) + work, scratch = root / 'work', root / 'scratch' + work.mkdir() + scratch.mkdir() + # Preflight write permissions before importing Psi4 in the child process. + for folder in (work, scratch): + probe = folder / '.write-probe' + probe.write_bytes(b'psi4') + probe.unlink() + env = dict(os.environ, PSI4_WORK_DIR=str(work), PSI4_OUTPUT_DIR=str(output), + PSI_SCRATCH=str(scratch), + TMPDIR=str(scratch), TEMP=str(scratch), TMP=str(scratch), + PYTHONUNBUFFERED='1') + exit_code = 1 + failure = None + with (output / 'logs.out').open('w', encoding='utf-8') as log: + try: + process = subprocess.Popen( + [sys.executable, str(script)], cwd=work, env=env, + stdout=subprocess.PIPE, stderr=subprocess.STDOUT, + text=True, encoding='utf-8', errors='replace', bufsize=1) + for line in process.stdout: + log.write(line) + log.flush() + sys.stdout.write(line) + sys.stdout.flush() + process.stdout.close() + exit_code = process.wait() + except OSError as error: + failure = str(error) + log.write(f'Launcher failed: {error}\n') + finally: + try: + copy_artifacts(work, output) + except (OSError, shutil.Error) as error: + failure = f'Artifact copy failed: {error}' + print(failure, file=sys.stderr) + if exit_code == 0: + exit_code = 1 + if exit_code != 0 or failure: + (output / 'launcher_failure.json').write_text(json.dumps({ + 'status': 'failed', 'exit_code': exit_code, 'error': failure, + }, indent=2), encoding='utf-8') + return exit_code + + +if __name__ == '__main__': + sys.exit(main()) diff --git a/agents/psi4/tools/psi4/smoke_psi4_runtime.py b/agents/psi4/tools/psi4/smoke_psi4_runtime.py index e42532e1..e9095d1f 100644 --- a/agents/psi4/tools/psi4/smoke_psi4_runtime.py +++ b/agents/psi4/tools/psi4/smoke_psi4_runtime.py @@ -1,44 +1,44 @@ -"""Real-Psi4 smoke; run through run_psi4_job.py in the approved built runtime.""" -import json -import math -import os -from pathlib import Path - -from psi4_utils import ( - setup_psi4, create_molecule, compute_energy, compute_electronic_properties, - compute_reorganization_energies, -) - - -def main(): - setup_psi4(memory='1 GB', nthreads=1, output_file='smoke.out') - water = create_molecule('O 0 0 0\nH 0 0.757 0.587\nH 0 -0.757 0.587', symmetry='c1') - properties = compute_electronic_properties(water, method='hf', basis='sto-3g') - assert 0.1 < properties['dipole_moment_debye'] < 5 - assert len(properties['dipole_vector_debye']) == 3 - assert all(math.isfinite(v) for v in properties['dipole_vector_debye']) - # Water doublet ions exercise unrestricted charged SCF with a diffuse basis. - for charge in (1, -1): - ion = water.clone() - ion.set_molecular_charge(charge) - ion.set_multiplicity(2) - energy = compute_energy(ion, method='hf', basis='6-31+g*') - assert energy['reference'] == 'uhf' - assert math.isfinite(energy['energy_hartree']) - # LiH has bound neutral/cation/anion states and a tiny nuclear geometry. - molecule = create_molecule('Li 0 0 0\nH 0 0 1.6', symmetry='c1') - reorganization = compute_reorganization_energies( - molecule, 'smoke-LiH', method='hf', optimization_basis='6-31+g*', - energy_basis='6-31+g*', geom_maxiter=50) - Path('smoke_results.json').write_text(json.dumps( - {'properties': properties, 'reorganization': reorganization}, indent=2, allow_nan=False)) - assert reorganization['status'] == 'success', reorganization - assert all(state['converged'] for state in reorganization['states'].values()) - assert len(reorganization['components']) == 7 - assert all(math.isfinite(value) for value in reorganization['components'].values()) - assert all(math.isfinite(reorganization[f'{carrier}_reorganization_eV']) - for carrier in ('hole', 'electron')) - - -if __name__ == '__main__': - main() +"""Real-Psi4 smoke; run through run_psi4_job.py in the approved built runtime.""" +import json +import math +import os +from pathlib import Path + +from psi4_utils import ( + setup_psi4, create_molecule, compute_energy, compute_electronic_properties, + compute_reorganization_energies, +) + + +def main(): + setup_psi4(memory='1 GB', nthreads=1, output_file='smoke.out') + water = create_molecule('O 0 0 0\nH 0 0.757 0.587\nH 0 -0.757 0.587', symmetry='c1') + properties = compute_electronic_properties(water, method='hf', basis='sto-3g') + assert 0.1 < properties['dipole_moment_debye'] < 5 + assert len(properties['dipole_vector_debye']) == 3 + assert all(math.isfinite(v) for v in properties['dipole_vector_debye']) + # Water doublet ions exercise unrestricted charged SCF with a diffuse basis. + for charge in (1, -1): + ion = water.clone() + ion.set_molecular_charge(charge) + ion.set_multiplicity(2) + energy = compute_energy(ion, method='hf', basis='6-31+g*') + assert energy['reference'] == 'uhf' + assert math.isfinite(energy['energy_hartree']) + # LiH has bound neutral/cation/anion states and a tiny nuclear geometry. + molecule = create_molecule('Li 0 0 0\nH 0 0 1.6', symmetry='c1') + reorganization = compute_reorganization_energies( + molecule, 'smoke-LiH', method='hf', optimization_basis='6-31+g*', + energy_basis='6-31+g*', geom_maxiter=50) + Path('smoke_results.json').write_text(json.dumps( + {'properties': properties, 'reorganization': reorganization}, indent=2, allow_nan=False)) + assert reorganization['status'] == 'success', reorganization + assert all(state['converged'] for state in reorganization['states'].values()) + assert len(reorganization['components']) == 7 + assert all(math.isfinite(value) for value in reorganization['components'].values()) + assert all(math.isfinite(reorganization[f'{carrier}_reorganization_eV']) + for carrier in ('hole', 'electron')) + + +if __name__ == '__main__': + main() diff --git a/agents/psi4/tools/psi4/test_psi4_reliability.py b/agents/psi4/tools/psi4/test_psi4_reliability.py index 5a429522..edd7057c 100644 --- a/agents/psi4/tools/psi4/test_psi4_reliability.py +++ b/agents/psi4/tools/psi4/test_psi4_reliability.py @@ -1,323 +1,323 @@ -"""Targeted stdlib unittest suite; Psi4 is mocked, NumPy is the existing dependency.""" -import copy -import json -import os -from pathlib import Path -import subprocess -import sys -import tempfile -import unittest -from unittest.mock import MagicMock, patch - -import numpy as np -import psi4_utils as utils -import run_psi4_job as launcher - - -class Molecule: - def __init__(self, charge=0, mult=1, x=0.0): - self.charge, self.mult, self.position = charge, mult, x - - def multiplicity(self): - return self.mult - - def molecular_charge(self): - return self.charge - - def set_molecular_charge(self, value): - self.charge = value - - def set_multiplicity(self, value): - self.mult = value - - def clone(self): - return copy.copy(self) - - def natom(self): - return 1 - - def symbol(self, index): - return 'He' - - def x(self, index): - return self.position - - def y(self, index): - return 0.0 - - def z(self, index): - return 0.0 - - -class OptimizationError(Exception): - def __init__(self, wfn=None): - self.wfn = wfn - - -class UtilityTests(unittest.TestCase): - def setUp(self): - self.psi4 = MagicMock() - self.psi4.OptimizationConvergenceError = OptimizationError - self.patch = patch.dict(sys.modules, {'psi4': self.psi4}) - self.patch.start() - self.addCleanup(self.patch.stop) - output_patch = patch.object(utils, 'OUTPUT_DIR', '') - output_patch.start() - self.addCleanup(output_patch.stop) - self.molecule = Molecule() - self.wfn = MagicMock() - self.wfn.nirrep.return_value = 1 - self.wfn.epsilon_a.return_value = [-0.5, 0.1] - self.wfn.nalpha.return_value = 1 - self.wfn.variable.return_value = np.array([[0.0, 0.0, 1.0]]) - self.wfn.molecule.return_value = self.molecule - self.wfn.energy.return_value = -1.0 - self.psi4.energy.return_value = (-1.0, self.wfn) - - def test_reference_selected_for_each_energy_including_dft(self): - for mult, reference in [(1, 'rhf'), (2, 'uhf'), (1, 'rhf')]: - utils.compute_energy(Molecule(mult=mult), 'b3lyp', return_wfn=True) - self.assertEqual(self.psi4.set_options.call_args.args[0]['reference'], reference) - - def test_explicit_reference_validated_before_energy(self): - for reference in ['rhf', 'rks', 'invalid']: - with self.assertRaises(ValueError): - utils.compute_energy(Molecule(mult=2), reference=reference) - self.psi4.energy.assert_not_called() - utils.compute_energy(Molecule(mult=2), reference='ROHF', return_wfn=True) - self.assertEqual(self.psi4.set_options.call_args.args[0]['reference'], 'rohf') - - def test_combined_properties_one_wavefunction_and_debye_conversion(self): - result = utils.compute_electronic_properties(self.molecule) - self.psi4.energy.assert_called_once() - self.psi4.oeprop.assert_called_once_with(self.wfn, 'DIPOLE', title='DISCOVERY') - self.wfn.variable.assert_called_once_with('DISCOVERY DIPOLE') - self.assertAlmostEqual(result['dipole_moment_debye'], utils.DIPOLE_AU_TO_DEBYE) - self.assertEqual(result['dipole_vector_debye'], [0.0, 0.0, utils.DIPOLE_AU_TO_DEBYE]) - self.assertAlmostEqual(result['gap_hartree'], 0.6) - - def test_gap_automatically_selects_uhf(self): - result = utils.compute_homo_lumo_gap(Molecule(mult=2)) - self.assertAlmostEqual(result['gap_hartree'], 0.6) - self.assertEqual(self.psi4.set_options.call_args.args[0]['reference'], 'uhf') - - def test_dipole_invalid_vectors_rejected(self): - for vector in [[1, 2], [0, float('nan'), 0], [float('inf'), 0, 0]]: - self.wfn.variable.return_value = vector - with self.assertRaises(ValueError): - utils.compute_dipole_moment(self.molecule) - - def test_empty_orbitals_rejected(self): - self.wfn.nalpha.return_value = 0 - with self.assertRaises(ValueError): - utils.compute_electronic_properties(self.molecule) - - def test_optimization_retries_once_from_last_geometry_without_weakening(self): - self.psi4.energy.return_value = -0.9 - last = Molecule(mult=2, x=1.0) - self.wfn.molecule.return_value = last - self.psi4.optimize.side_effect = [OptimizationError(self.wfn), (-1.1, self.wfn)] - result = utils.optimize_geometry(Molecule(mult=2), 'b3lyp', geom_maxiter=5, - g_convergence='gau_tight') - self.assertTrue(result['converged']) - self.assertEqual(result['optimization_attempts'], 2) - self.assertIs(self.psi4.optimize.call_args.kwargs['molecule'], last) - settings = [call.args[0] for call in self.psi4.set_options.call_args_list] - self.assertEqual(settings[-2:], [{'geom_maxiter': 5}, {'geom_maxiter': 10}]) - self.assertEqual(settings[0]['g_convergence'], 'gau_tight') - self.assertEqual(settings[0]['reference'], 'uhf') - - def test_optimization_exhaustion_is_not_success(self): - self.psi4.energy.return_value = -0.9 - self.psi4.optimize.side_effect = OptimizationError(self.wfn) - result = utils.optimize_geometry(self.molecule) - self.assertFalse(result['converged']) - self.assertEqual(self.psi4.optimize.call_count, 2) - - def test_optimization_without_last_wavefunction_does_not_retry_initial_geometry(self): - self.psi4.energy.return_value = -0.9 - self.psi4.optimize.side_effect = OptimizationError() - self.assertFalse(utils.optimize_geometry(self.molecule)['converged']) - self.assertEqual(self.psi4.optimize.call_count, 1) - - def test_scf_failure_not_swallowed(self): - self.psi4.energy.side_effect = RuntimeError('SCF failed') - with self.assertRaisesRegex(RuntimeError, 'SCF failed'): - utils.compute_electronic_properties(self.molecule) - - def test_nonfinite_energy_and_unconverged_acceptance_rejected(self): - self.psi4.energy.return_value = float('nan') - with self.assertRaises(ValueError): - utils.compute_energy(self.molecule) - with self.assertRaises(ValueError): - utils.compute_energy(self.molecule, fail_on_maxiter=False) - - def test_uppercase_restricted_override_cannot_bypass_validation(self): - with self.assertRaises(ValueError): - utils.compute_energy(Molecule(mult=2), REFERENCE='RHF') - - def test_setup_binds_iomanager_scratch(self): - with tempfile.TemporaryDirectory(dir=Path.cwd()) as directory: - with patch.dict(os.environ, {'PSI_SCRATCH': directory}): - utils.setup_psi4() - self.psi4.core.IOManager.shared_object().set_default_path.assert_called_with( - os.path.abspath(directory)) - - -class ReorganizationTests(unittest.TestCase): - def setUp(self): - self.directory = tempfile.TemporaryDirectory(dir=Path.cwd()) - self.addCleanup(self.directory.cleanup) - self.patch = patch.object(utils, 'WORK_DIR', self.directory.name) - self.patch.start() - self.addCleanup(self.patch.stop) - output_patch = patch.object(utils, 'OUTPUT_DIR', '') - output_patch.start() - self.addCleanup(output_patch.stop) - - @staticmethod - def state_result(directory, name, request): - energies = {'opt_neutral': -10.0, 'opt_cation': -9.0, 'opt_anion': -10.2, - 'cation_at_neutral': -8.9, 'neutral_at_cation': -9.95, - 'anion_at_neutral': -10.1, 'neutral_at_anion': -9.97} - return {'status': 'success', 'converged': True, 'energy_hartree': energies[name], - 'optimized_geometry': 'He 0 0 0', 'optimization_attempts': 1, - 'optimized_energy_hartree': energies[name]} - - def test_four_point_formulas_ids_units_states_and_seven_fresh_jobs(self): - with patch.object(utils, '_run_reorganization_state', side_effect=self.state_result) as run: - result = utils.compute_reorganization_energies(Molecule(), 'candidate-A') - self.assertEqual(result['status'], 'success') - self.assertEqual(result['candidate_id'], 'candidate-A') - self.assertAlmostEqual(result['hole_reorganization_hartree'], 0.15) - self.assertAlmostEqual(result['electron_reorganization_hartree'], 0.13) - self.assertAlmostEqual(result['electron_reorganization_eV'], 0.13 * utils.HARTREE_TO_EV) - self.assertEqual(result['energy_units'], 'hartree') - self.assertEqual(run.call_count, 7) - self.assertEqual([call.args[2]['charge'] for call in run.call_args_list], - [0, 1, -1, 1, 0, -1, 0]) - self.assertEqual([call.args[2]['multiplicity'] for call in run.call_args_list], - [1, 2, 2, 2, 1, 2, 1]) - self.assertEqual([state['reference'] for state in result['states'].values()], - ['rhf', 'uhf', 'uhf']) - self.assertEqual(len(result['components']), 7) - self.assertTrue((Path(result['artifact_directory']) / 'reorganization_result.json').exists()) - - def test_failure_preserves_id_partial_components_and_json(self): - def run(directory, name, request): - if name == 'opt_cation': - raise RuntimeError('unconverged') - return self.state_result(directory, name, request) - with patch.object(utils, '_run_reorganization_state', side_effect=run): - result = utils.compute_reorganization_energies(Molecule(), 'candidate-B') - self.assertEqual(result['status'], 'failed') - self.assertEqual(result['components'], {'E_neutral_at_neutral': -10.0}) - self.assertNotIn('hole_reorganization_eV', result) - self.assertIn('unconverged', result['error']) - saved = json.loads((Path(result['artifact_directory']) / 'reorganization_result.json').read_text()) - self.assertEqual(saved['candidate_id'], 'candidate-B') - - def test_upstream_partial_neutral_never_assumed_converged(self): - with patch.object(utils, '_run_reorganization_state', - return_value={'converged': False}) as run: - result = utils.compute_reorganization_energies(Molecule(), 'candidate-C') - self.assertEqual(run.call_args.args[1], 'opt_neutral') - self.assertEqual(result['status'], 'failed') - self.assertFalse(result['states']['neutral']['converged']) - - def test_restricted_open_shell_rejected_before_jobs(self): - with patch.object(utils, '_run_reorganization_state') as run: - result = utils.compute_reorganization_energies(Molecule(), 'candidate-D', reference='rhf') - self.assertEqual(result['status'], 'failed') - run.assert_not_called() - - def test_worker_nonzero_exit_with_result_still_fails(self): - directory = Path(self.directory.name) - with patch.object(utils.subprocess, 'run', return_value=subprocess.CompletedProcess([], 7)): - with self.assertRaisesRegex(RuntimeError, 'exit 7'): - utils._run_reorganization_state(directory, 'bad_state', {'options': {}}) - self.assertTrue((directory / 'bad_state' / 'worker.out').exists()) - - -class LauncherTests(unittest.TestCase): - def setUp(self): - self.directory = tempfile.TemporaryDirectory(dir=Path.cwd(), prefix='launcher-tests-') - self.root = Path(self.directory.name) - self.addCleanup(self.directory.cleanup) - - def start(self, script, output): - return subprocess.Popen([sys.executable, str(Path(launcher.__file__)), str(script), str(output)], - stdout=subprocess.PIPE, stderr=subprocess.STDOUT, text=True) - - def test_concurrent_invocations_isolate_scratch_work_and_copy_nested_artifacts(self): - script = self.root / 'script ; & spaces.py' - script.write_text( - "import json,os,time\nfrom pathlib import Path\n" - "Path('nested').mkdir()\n" - "Path('nested/scratch').mkdir()\n" - "Path('nested/scratch/secret.dat').write_text('nested scratch')\n" - "Path('logs.out').write_text('script own transcript')\n" - "Path('nested/result.json').write_text(json.dumps(dict(work=os.getcwd()," - "scratch=os.environ['PSI_SCRATCH'])))\n" - "Path(os.environ['PSI_SCRATCH'],'secret.dat').write_text('scratch')\n" - "print('stdout evidence', flush=True)\n" - "import sys; print('stderr evidence',file=sys.stderr,flush=True)\n" - "time.sleep(0.1)\n", encoding='utf-8') - outputs = [self.root / 'output one', self.root / 'output two'] - processes = [self.start(script, output) for output in outputs] - for process in processes: - stdout, _ = process.communicate(timeout=30) - self.assertEqual(process.returncode, 0, stdout) - results = [json.loads((output / 'nested' / 'result.json').read_text()) for output in outputs] - self.assertNotEqual(results[0]['work'], results[1]['work']) - self.assertNotEqual(results[0]['scratch'], results[1]['scratch']) - for output in outputs: - self.assertIn('stderr evidence', (output / 'logs.out').read_text()) - self.assertEqual((output / 'script-logs.out').read_text(), 'script own transcript') - self.assertFalse(list(output.rglob('secret.dat'))) - - def test_child_exit_preserved_and_partial_artifacts_retained(self): - script = self.root / 'fail.py' - script.write_text("from pathlib import Path\nPath('partial.json').write_text('{}')\n" - "print('failing evidence')\nraise SystemExit(7)\n") - output = self.root / 'output' - process = self.start(script, output) - stdout, _ = process.communicate(timeout=30) - self.assertEqual(process.returncode, 7, stdout) - self.assertTrue((output / 'partial.json').exists()) - self.assertEqual(json.loads((output / 'launcher_failure.json').read_text())['exit_code'], 7) - - def test_copy_failure_converts_success_to_failure_and_is_recorded(self): - script = self.root / 'ok.py' - script.write_text("print('ok')\n") - output = self.root / 'output' - with patch.object(launcher, 'copy_artifacts', side_effect=OSError('copy denied')): - self.assertEqual(launcher.main([str(script), str(output)]), 1) - self.assertIn('copy denied', (output / 'launcher_failure.json').read_text()) - - def test_copy_failure_does_not_hide_child_exit(self): - script = self.root / 'bad.py' - script.write_text('raise SystemExit(7)\n') - output = self.root / 'output' - with patch.object(launcher, 'copy_artifacts', side_effect=OSError('copy denied')): - self.assertEqual(launcher.main([str(script), str(output)]), 7) - - def test_setup_honors_work_env_without_replacing_input_output_paths(self): - original = Path.cwd() - input_dir, output_dir, work_dir = [self.root / name for name in ('input', 'out', 'work')] - input_dir.mkdir() - (input_dir / 'molecule.xyz').write_text('input') - try: - with patch.dict(os.environ, {'PSI4_WORK_DIR': str(work_dir), - 'PSI_SCRATCH': str(self.root / 'scratch')}): - utils.setup_directories(str(input_dir), str(output_dir), work_dir='ignored') - self.assertEqual(utils.INPUT_DIR, str(input_dir)) - self.assertEqual(utils.OUTPUT_DIR, str(output_dir)) - self.assertEqual(utils.WORK_DIR, str(work_dir)) - self.assertEqual((work_dir / 'molecule.xyz').read_text(), 'input') - finally: - os.chdir(original) - - -if __name__ == '__main__': - unittest.main() +"""Targeted stdlib unittest suite; Psi4 is mocked, NumPy is the existing dependency.""" +import copy +import json +import os +from pathlib import Path +import subprocess +import sys +import tempfile +import unittest +from unittest.mock import MagicMock, patch + +import numpy as np +import psi4_utils as utils +import run_psi4_job as launcher + + +class Molecule: + def __init__(self, charge=0, mult=1, x=0.0): + self.charge, self.mult, self.position = charge, mult, x + + def multiplicity(self): + return self.mult + + def molecular_charge(self): + return self.charge + + def set_molecular_charge(self, value): + self.charge = value + + def set_multiplicity(self, value): + self.mult = value + + def clone(self): + return copy.copy(self) + + def natom(self): + return 1 + + def symbol(self, index): + return 'He' + + def x(self, index): + return self.position + + def y(self, index): + return 0.0 + + def z(self, index): + return 0.0 + + +class OptimizationError(Exception): + def __init__(self, wfn=None): + self.wfn = wfn + + +class UtilityTests(unittest.TestCase): + def setUp(self): + self.psi4 = MagicMock() + self.psi4.OptimizationConvergenceError = OptimizationError + self.patch = patch.dict(sys.modules, {'psi4': self.psi4}) + self.patch.start() + self.addCleanup(self.patch.stop) + output_patch = patch.object(utils, 'OUTPUT_DIR', '') + output_patch.start() + self.addCleanup(output_patch.stop) + self.molecule = Molecule() + self.wfn = MagicMock() + self.wfn.nirrep.return_value = 1 + self.wfn.epsilon_a.return_value = [-0.5, 0.1] + self.wfn.nalpha.return_value = 1 + self.wfn.variable.return_value = np.array([[0.0, 0.0, 1.0]]) + self.wfn.molecule.return_value = self.molecule + self.wfn.energy.return_value = -1.0 + self.psi4.energy.return_value = (-1.0, self.wfn) + + def test_reference_selected_for_each_energy_including_dft(self): + for mult, reference in [(1, 'rhf'), (2, 'uhf'), (1, 'rhf')]: + utils.compute_energy(Molecule(mult=mult), 'b3lyp', return_wfn=True) + self.assertEqual(self.psi4.set_options.call_args.args[0]['reference'], reference) + + def test_explicit_reference_validated_before_energy(self): + for reference in ['rhf', 'rks', 'invalid']: + with self.assertRaises(ValueError): + utils.compute_energy(Molecule(mult=2), reference=reference) + self.psi4.energy.assert_not_called() + utils.compute_energy(Molecule(mult=2), reference='ROHF', return_wfn=True) + self.assertEqual(self.psi4.set_options.call_args.args[0]['reference'], 'rohf') + + def test_combined_properties_one_wavefunction_and_debye_conversion(self): + result = utils.compute_electronic_properties(self.molecule) + self.psi4.energy.assert_called_once() + self.psi4.oeprop.assert_called_once_with(self.wfn, 'DIPOLE', title='DISCOVERY') + self.wfn.variable.assert_called_once_with('DISCOVERY DIPOLE') + self.assertAlmostEqual(result['dipole_moment_debye'], utils.DIPOLE_AU_TO_DEBYE) + self.assertEqual(result['dipole_vector_debye'], [0.0, 0.0, utils.DIPOLE_AU_TO_DEBYE]) + self.assertAlmostEqual(result['gap_hartree'], 0.6) + + def test_gap_automatically_selects_uhf(self): + result = utils.compute_homo_lumo_gap(Molecule(mult=2)) + self.assertAlmostEqual(result['gap_hartree'], 0.6) + self.assertEqual(self.psi4.set_options.call_args.args[0]['reference'], 'uhf') + + def test_dipole_invalid_vectors_rejected(self): + for vector in [[1, 2], [0, float('nan'), 0], [float('inf'), 0, 0]]: + self.wfn.variable.return_value = vector + with self.assertRaises(ValueError): + utils.compute_dipole_moment(self.molecule) + + def test_empty_orbitals_rejected(self): + self.wfn.nalpha.return_value = 0 + with self.assertRaises(ValueError): + utils.compute_electronic_properties(self.molecule) + + def test_optimization_retries_once_from_last_geometry_without_weakening(self): + self.psi4.energy.return_value = -0.9 + last = Molecule(mult=2, x=1.0) + self.wfn.molecule.return_value = last + self.psi4.optimize.side_effect = [OptimizationError(self.wfn), (-1.1, self.wfn)] + result = utils.optimize_geometry(Molecule(mult=2), 'b3lyp', geom_maxiter=5, + g_convergence='gau_tight') + self.assertTrue(result['converged']) + self.assertEqual(result['optimization_attempts'], 2) + self.assertIs(self.psi4.optimize.call_args.kwargs['molecule'], last) + settings = [call.args[0] for call in self.psi4.set_options.call_args_list] + self.assertEqual(settings[-2:], [{'geom_maxiter': 5}, {'geom_maxiter': 10}]) + self.assertEqual(settings[0]['g_convergence'], 'gau_tight') + self.assertEqual(settings[0]['reference'], 'uhf') + + def test_optimization_exhaustion_is_not_success(self): + self.psi4.energy.return_value = -0.9 + self.psi4.optimize.side_effect = OptimizationError(self.wfn) + result = utils.optimize_geometry(self.molecule) + self.assertFalse(result['converged']) + self.assertEqual(self.psi4.optimize.call_count, 2) + + def test_optimization_without_last_wavefunction_does_not_retry_initial_geometry(self): + self.psi4.energy.return_value = -0.9 + self.psi4.optimize.side_effect = OptimizationError() + self.assertFalse(utils.optimize_geometry(self.molecule)['converged']) + self.assertEqual(self.psi4.optimize.call_count, 1) + + def test_scf_failure_not_swallowed(self): + self.psi4.energy.side_effect = RuntimeError('SCF failed') + with self.assertRaisesRegex(RuntimeError, 'SCF failed'): + utils.compute_electronic_properties(self.molecule) + + def test_nonfinite_energy_and_unconverged_acceptance_rejected(self): + self.psi4.energy.return_value = float('nan') + with self.assertRaises(ValueError): + utils.compute_energy(self.molecule) + with self.assertRaises(ValueError): + utils.compute_energy(self.molecule, fail_on_maxiter=False) + + def test_uppercase_restricted_override_cannot_bypass_validation(self): + with self.assertRaises(ValueError): + utils.compute_energy(Molecule(mult=2), REFERENCE='RHF') + + def test_setup_binds_iomanager_scratch(self): + with tempfile.TemporaryDirectory(dir=Path.cwd()) as directory: + with patch.dict(os.environ, {'PSI_SCRATCH': directory}): + utils.setup_psi4() + self.psi4.core.IOManager.shared_object().set_default_path.assert_called_with( + os.path.abspath(directory)) + + +class ReorganizationTests(unittest.TestCase): + def setUp(self): + self.directory = tempfile.TemporaryDirectory(dir=Path.cwd()) + self.addCleanup(self.directory.cleanup) + self.patch = patch.object(utils, 'WORK_DIR', self.directory.name) + self.patch.start() + self.addCleanup(self.patch.stop) + output_patch = patch.object(utils, 'OUTPUT_DIR', '') + output_patch.start() + self.addCleanup(output_patch.stop) + + @staticmethod + def state_result(directory, name, request): + energies = {'opt_neutral': -10.0, 'opt_cation': -9.0, 'opt_anion': -10.2, + 'cation_at_neutral': -8.9, 'neutral_at_cation': -9.95, + 'anion_at_neutral': -10.1, 'neutral_at_anion': -9.97} + return {'status': 'success', 'converged': True, 'energy_hartree': energies[name], + 'optimized_geometry': 'He 0 0 0', 'optimization_attempts': 1, + 'optimized_energy_hartree': energies[name]} + + def test_four_point_formulas_ids_units_states_and_seven_fresh_jobs(self): + with patch.object(utils, '_run_reorganization_state', side_effect=self.state_result) as run: + result = utils.compute_reorganization_energies(Molecule(), 'candidate-A') + self.assertEqual(result['status'], 'success') + self.assertEqual(result['candidate_id'], 'candidate-A') + self.assertAlmostEqual(result['hole_reorganization_hartree'], 0.15) + self.assertAlmostEqual(result['electron_reorganization_hartree'], 0.13) + self.assertAlmostEqual(result['electron_reorganization_eV'], 0.13 * utils.HARTREE_TO_EV) + self.assertEqual(result['energy_units'], 'hartree') + self.assertEqual(run.call_count, 7) + self.assertEqual([call.args[2]['charge'] for call in run.call_args_list], + [0, 1, -1, 1, 0, -1, 0]) + self.assertEqual([call.args[2]['multiplicity'] for call in run.call_args_list], + [1, 2, 2, 2, 1, 2, 1]) + self.assertEqual([state['reference'] for state in result['states'].values()], + ['rhf', 'uhf', 'uhf']) + self.assertEqual(len(result['components']), 7) + self.assertTrue((Path(result['artifact_directory']) / 'reorganization_result.json').exists()) + + def test_failure_preserves_id_partial_components_and_json(self): + def run(directory, name, request): + if name == 'opt_cation': + raise RuntimeError('unconverged') + return self.state_result(directory, name, request) + with patch.object(utils, '_run_reorganization_state', side_effect=run): + result = utils.compute_reorganization_energies(Molecule(), 'candidate-B') + self.assertEqual(result['status'], 'failed') + self.assertEqual(result['components'], {'E_neutral_at_neutral': -10.0}) + self.assertNotIn('hole_reorganization_eV', result) + self.assertIn('unconverged', result['error']) + saved = json.loads((Path(result['artifact_directory']) / 'reorganization_result.json').read_text()) + self.assertEqual(saved['candidate_id'], 'candidate-B') + + def test_upstream_partial_neutral_never_assumed_converged(self): + with patch.object(utils, '_run_reorganization_state', + return_value={'converged': False}) as run: + result = utils.compute_reorganization_energies(Molecule(), 'candidate-C') + self.assertEqual(run.call_args.args[1], 'opt_neutral') + self.assertEqual(result['status'], 'failed') + self.assertFalse(result['states']['neutral']['converged']) + + def test_restricted_open_shell_rejected_before_jobs(self): + with patch.object(utils, '_run_reorganization_state') as run: + result = utils.compute_reorganization_energies(Molecule(), 'candidate-D', reference='rhf') + self.assertEqual(result['status'], 'failed') + run.assert_not_called() + + def test_worker_nonzero_exit_with_result_still_fails(self): + directory = Path(self.directory.name) + with patch.object(utils.subprocess, 'run', return_value=subprocess.CompletedProcess([], 7)): + with self.assertRaisesRegex(RuntimeError, 'exit 7'): + utils._run_reorganization_state(directory, 'bad_state', {'options': {}}) + self.assertTrue((directory / 'bad_state' / 'worker.out').exists()) + + +class LauncherTests(unittest.TestCase): + def setUp(self): + self.directory = tempfile.TemporaryDirectory(dir=Path.cwd(), prefix='launcher-tests-') + self.root = Path(self.directory.name) + self.addCleanup(self.directory.cleanup) + + def start(self, script, output): + return subprocess.Popen([sys.executable, str(Path(launcher.__file__)), str(script), str(output)], + stdout=subprocess.PIPE, stderr=subprocess.STDOUT, text=True) + + def test_concurrent_invocations_isolate_scratch_work_and_copy_nested_artifacts(self): + script = self.root / 'script ; & spaces.py' + script.write_text( + "import json,os,time\nfrom pathlib import Path\n" + "Path('nested').mkdir()\n" + "Path('nested/scratch').mkdir()\n" + "Path('nested/scratch/secret.dat').write_text('nested scratch')\n" + "Path('logs.out').write_text('script own transcript')\n" + "Path('nested/result.json').write_text(json.dumps(dict(work=os.getcwd()," + "scratch=os.environ['PSI_SCRATCH'])))\n" + "Path(os.environ['PSI_SCRATCH'],'secret.dat').write_text('scratch')\n" + "print('stdout evidence', flush=True)\n" + "import sys; print('stderr evidence',file=sys.stderr,flush=True)\n" + "time.sleep(0.1)\n", encoding='utf-8') + outputs = [self.root / 'output one', self.root / 'output two'] + processes = [self.start(script, output) for output in outputs] + for process in processes: + stdout, _ = process.communicate(timeout=30) + self.assertEqual(process.returncode, 0, stdout) + results = [json.loads((output / 'nested' / 'result.json').read_text()) for output in outputs] + self.assertNotEqual(results[0]['work'], results[1]['work']) + self.assertNotEqual(results[0]['scratch'], results[1]['scratch']) + for output in outputs: + self.assertIn('stderr evidence', (output / 'logs.out').read_text()) + self.assertEqual((output / 'script-logs.out').read_text(), 'script own transcript') + self.assertFalse(list(output.rglob('secret.dat'))) + + def test_child_exit_preserved_and_partial_artifacts_retained(self): + script = self.root / 'fail.py' + script.write_text("from pathlib import Path\nPath('partial.json').write_text('{}')\n" + "print('failing evidence')\nraise SystemExit(7)\n") + output = self.root / 'output' + process = self.start(script, output) + stdout, _ = process.communicate(timeout=30) + self.assertEqual(process.returncode, 7, stdout) + self.assertTrue((output / 'partial.json').exists()) + self.assertEqual(json.loads((output / 'launcher_failure.json').read_text())['exit_code'], 7) + + def test_copy_failure_converts_success_to_failure_and_is_recorded(self): + script = self.root / 'ok.py' + script.write_text("print('ok')\n") + output = self.root / 'output' + with patch.object(launcher, 'copy_artifacts', side_effect=OSError('copy denied')): + self.assertEqual(launcher.main([str(script), str(output)]), 1) + self.assertIn('copy denied', (output / 'launcher_failure.json').read_text()) + + def test_copy_failure_does_not_hide_child_exit(self): + script = self.root / 'bad.py' + script.write_text('raise SystemExit(7)\n') + output = self.root / 'output' + with patch.object(launcher, 'copy_artifacts', side_effect=OSError('copy denied')): + self.assertEqual(launcher.main([str(script), str(output)]), 7) + + def test_setup_honors_work_env_without_replacing_input_output_paths(self): + original = Path.cwd() + input_dir, output_dir, work_dir = [self.root / name for name in ('input', 'out', 'work')] + input_dir.mkdir() + (input_dir / 'molecule.xyz').write_text('input') + try: + with patch.dict(os.environ, {'PSI4_WORK_DIR': str(work_dir), + 'PSI_SCRATCH': str(self.root / 'scratch')}): + utils.setup_directories(str(input_dir), str(output_dir), work_dir='ignored') + self.assertEqual(utils.INPUT_DIR, str(input_dir)) + self.assertEqual(utils.OUTPUT_DIR, str(output_dir)) + self.assertEqual(utils.WORK_DIR, str(work_dir)) + self.assertEqual((work_dir / 'molecule.xyz').read_text(), 'input') + finally: + os.chdir(original) + + +if __name__ == '__main__': + unittest.main()