Repository: Quantum-Accelerators/quacc
Pull request #2042 (2024-04-22 23:23:12)

Add perturbation-based quasi-IRC workflows for ORCA and Q-Chem
## Summary of Changes

Title says most of it. I also added a frequency workflow for ORCA; felt strange that there wasn't one already.

### Checklist

- [X] I have read the ["Guidelines" section](https://quantum-accelerators.github.io/quacc/dev/contributing.html#guidelines) of the contributing guide. Don't lie! 😉
- [X] My PR is on a custom branch and is _not_ named `main`.
- [X] I have added relevant, comprehensive [unit tests](https://quantum-accelerators.github.io/quacc/dev/contributing.html#unit-tests).

### Notes

- Your PR will likely not be merged without proper and thorough tests.
- If you are an external contributor, you will see a comment from [@buildbot-princeton](https://github.com/buildbot-princeton). This is solely for the maintainers.
- When your code is ready for review, ping one of the [active maintainers](https://quantum-accelerators.github.io/quacc/about/contributors.html#active-maintainers).


buildbot-princeton: Can one of the admins verify this patch?
espottesmith: Pinging @Andrew-S-Rosen
codecov[bot]: ## [Codecov](https://app.codecov.io/gh/Quantum-Accelerators/quacc/pull/2042?dropdown=coverage&src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=Quantum-Accelerators) Report
Attention: Patch coverage is `97.36842%` with `1 lines` in your changes are missing coverage. Please review.
> Project coverage is 99.09%. Comparing base [(`487376b`)](https://app.codecov.io/gh/Quantum-Accelerators/quacc/commit/487376b5450eed2c720fb6ab76bcdbb11b7921d7?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=Quantum-Accelerators) to head [(`f2bbce7`)](https://app.codecov.io/gh/Quantum-Accelerators/quacc/pull/2042?dropdown=coverage&src=pr&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=Quantum-Accelerators).
> Report is 22 commits behind head on main.

| [Files](https://app.codecov.io/gh/Quantum-Accelerators/quacc/pull/2042?dropdown=coverage&src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=Quantum-Accelerators) | Patch % | Lines |
|---|---|---|
| [src/quacc/recipes/orca/core.py](https://app.codecov.io/gh/Quantum-Accelerators/quacc/pull/2042?src=pr&el=tree&filepath=src%2Fquacc%2Frecipes%2Forca%2Fcore.py&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=Quantum-Accelerators#diff-c3JjL3F1YWNjL3JlY2lwZXMvb3JjYS9jb3JlLnB5) | 94.11% | [1 Missing :warning: ](https://app.codecov.io/gh/Quantum-Accelerators/quacc/pull/2042?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=Quantum-Accelerators) |

<details><summary>Additional details and impacted files</summary>


```diff
@@            Coverage Diff             @@
##             main    #2042      +/-   ##
==========================================
- Coverage   99.13%   99.09%   -0.04%     
==========================================
  Files          81       82       +1     
  Lines        3336     3314      -22     
==========================================
- Hits         3307     3284      -23     
- Misses         29       30       +1     
```



</details>

[:umbrella: View full report in Codecov by Sentry](https://app.codecov.io/gh/Quantum-Accelerators/quacc/pull/2042?dropdown=coverage&src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=Quantum-Accelerators).   
:loudspeaker: Have feedback on the report? [Share it here](https://about.codecov.io/codecov-pr-comment-feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=Quantum-Accelerators).

espottesmith: Alright, think we're looking good!
Andrew-S-Rosen: Thank you very much, @espottesmith, for the contribution and addressing the comments! I can confirm this is good to go.

At some point, we will probably want to add a specific test for `atoms.core.perturb`, but I won't let that hold up this PR.

Patch:
diff --git a/docs/about/contributors.md b/docs/about/contributors.md
index 654a7887..1321e9c5 100644
--- a/docs/about/contributors.md
+++ b/docs/about/contributors.md
@@ -26,6 +26,7 @@ Additional contributions were made by the individuals listed [here](https://gith
 - [@ViktoriiaBaib](https://github.com/ViktoriiaBaib): Initial testing of quacc.
 - [@zulissimeta](https://github.com/zulissimeta): Dask support
 - [@yw-fang](https://github.com/yw-fang): VASP Non-SCF recipe
+- [@espottesmith](http://github.com/espottesmith): Some ORCA and Q-Chem workflows, mostly re: IRC
 
 ## Inspiration
 
diff --git a/docs/user/recipes/recipes_list.md b/docs/user/recipes/recipes_list.md
index 24fc4bcb..9085e094 100644
--- a/docs/user/recipes/recipes_list.md
+++ b/docs/user/recipes/recipes_list.md
@@ -140,11 +140,14 @@ The list of available quacc recipes is shown below. The "Req'd Extras" column sp
 
 <center>
 
-| Name           | Decorator       | Documentation                             | Req'd Extras |
-| -------------- | --------------- | ----------------------------------------- | ------------ |
-| ORCA Static    | `#!Python @job` | [quacc.recipes.orca.core.static_job][]    |              |
-| ORCA Relax     | `#!Python @job` | [quacc.recipes.orca.core.relax_job][]     |              |
-| ORCA ASE Relax | `#!Python @job` | [quacc.recipes.orca.core.ase_relax_job][] |              |
+| Name                       | Decorator       | Documentation                                         | Req'd Extras |
+| -------------------------- | --------------- | ----------------------------------------------------- | ------------ |
+| ORCA Static                | `#!Python @job` | [quacc.recipes.orca.core.static_job][]                |              |
+| ORCA Relax                 | `#!Python @job` | [quacc.recipes.orca.core.relax_job][]                 |              |
+| ORCA Freq                  | `#!Python @job` | [quacc.recipes.orca.core.freq_job][]                  |              |
+| ORCA ASE Relax             | `#!Python @job` | [quacc.recipes.orca.core.ase_relax_job][]             |              |
+| ORCA ASE Quasi-IRC Perturb | `#!Python @job` | [quacc.recipes.orca.core.ase_quasi_irc_perturb_job][] |              |
+
 
 </center>
 
@@ -170,14 +173,15 @@ The list of available quacc recipes is shown below. The "Req'd Extras" column sp
 
 <center>
 
-| Name             | Decorator       | Documentation                            | Req'd Extras   |
-| ---------------- | --------------- | ---------------------------------------- | -------------- |
-| Q-Chem Static    | `#!Python @job` | [quacc.recipes.qchem.core.static_job][]  |                |
-| Q-Chem Relax     | `#!Python @job` | [quacc.recipes.qchem.core.relax_job][]   |                |
-| Q-Chem Frequency | `#!Python @job` | [quacc.recipes.qchem.core.freq_job][]    |                |
-| Q-Chem TS        | `#!Python @job` | [quacc.recipes.qchem.ts.ts_job][]        | `quacc[sella]` |
-| Q-Chem IRC       | `#!Python @job` | [quacc.recipes.qchem.ts.irc_job][]       | `quacc[sella]` |
-| Q-Chem Quasi IRC | `#!Python @job` | [quacc.recipes.qchem.ts.quasi_irc_job][] | `quacc[sella]` |
+| Name                     | Decorator       | Documentation                                    | Req'd Extras   |
+| ------------------------ | --------------- | ------------------------------------------------ | -------------- |
+| Q-Chem Static            | `#!Python @job` | [quacc.recipes.qchem.core.static_job][]          |                |
+| Q-Chem Relax             | `#!Python @job` | [quacc.recipes.qchem.core.relax_job][]           |                |
+| Q-Chem Frequency         | `#!Python @job` | [quacc.recipes.qchem.core.freq_job][]            |                |
+| Q-Chem TS                | `#!Python @job` | [quacc.recipes.qchem.ts.ts_job][]                | `quacc[sella]` |
+| Q-Chem IRC               | `#!Python @job` | [quacc.recipes.qchem.ts.irc_job][]               | `quacc[sella]` |
+| Q-Chem Quasi IRC         | `#!Python @job` | [quacc.recipes.qchem.ts.quasi_irc_job][]         | `quacc[sella]` |
+| Q-Chem Quasi IRC Perturb | `#!Python @job` | [quacc.recipes.qchem.ts.quasi_irc_perturb_job][] | `quacc[sella]` |
 
 </center>
 
diff --git a/src/quacc/atoms/core.py b/src/quacc/atoms/core.py
index 95ec212b..f57b0c7c 100644
--- a/src/quacc/atoms/core.py
+++ b/src/quacc/atoms/core.py
@@ -15,6 +15,7 @@ from pymatgen.io.ase import AseAtomsAdaptor
 if TYPE_CHECKING:
     from ase.atoms import Atoms
     from ase.optimize.optimize import Dynamics
+    from numpy.typing import NDArray
 
 logger = logging.getLogger(__name__)
 
@@ -254,3 +255,34 @@ def get_final_atoms_from_dynamics(dynamics: Dynamics) -> Atoms:
     return (
         dynamics.atoms.atoms if isinstance(dynamics.atoms, Filter) else dynamics.atoms
     )
+
+
+def perturb(mol: Atoms, matrix: list[list[float]] | NDArray, scale: float) -> Atoms:
+    """
+    Perturb each atom in a molecule by a (scaled) 1x3 vector, reflecting e.g. a vibrational normal mode.
+
+    Parameters
+    ----------
+    mol
+        ASE Atoms object representing a molecule
+    matrix
+        Nx3 matrix, where N is the numebr of atoms. This means that there is potentially a different translation
+        vector for each atom in the molecule.
+    scale
+        Scaling factor for perturbation
+
+    Returns
+    -------
+    Atoms
+        The input molecule after perturbation
+    """
+
+    mol_copy = copy_atoms(mol)
+    mode = np.asarray(matrix)
+
+    orig_pos = mol_copy.get_positions()
+
+    pos = orig_pos + mode * scale
+    mol_copy.set_positions(pos)
+
+    return mol_copy
diff --git a/src/quacc/recipes/orca/core.py b/src/quacc/recipes/orca/core.py
index 53a4e71b..60d44489 100644
--- a/src/quacc/recipes/orca/core.py
+++ b/src/quacc/recipes/orca/core.py
@@ -7,12 +7,14 @@ from typing import TYPE_CHECKING
 import psutil
 
 from quacc import job
+from quacc.atoms.core import perturb
 from quacc.recipes.orca._base import run_and_summarize, run_and_summarize_opt
 
 if TYPE_CHECKING:
     from typing import Any, Literal
 
     from ase.atoms import Atoms
+    from numpy.typing import NDArray
 
     from quacc.schemas._aliases.cclib import cclibSchema
     from quacc.utils.files import Filenames, SourceDirectory
@@ -81,6 +83,74 @@ def static_job(
     )
 
 
+@job
+def freq_job(
+    atoms: Atoms,
+    charge: int = 0,
+    spin_multiplicity: int = 1,
+    xc: str = "wb97x-d3bj",
+    basis: str = "def2-tzvp",
+    numerical: bool = False,
+    orcasimpleinput: list[str] | None = None,
+    orcablocks: list[str] | None = None,
+    nprocs: int | Literal["max"] = "max",
+    copy_files: SourceDirectory | dict[SourceDirectory, Filenames] | None = None,
+) -> cclibSchema:
+    """
+    Carry out a vibrational frequency analysis calculation.
+
+    Parameters
+    ----------
+    atoms
+        Atoms object
+    charge
+        Charge of the system.
+    spin_multiplicity
+        Multiplicity of the system.
+    xc
+        Exchange-correlation functional
+    basis
+        Basis set
+    numerical
+        If True (default False), a numeric frequency calculation will be requested
+    orcasimpleinput
+        List of `orcasimpleinput` swaps for the calculator. To remove entries
+        from the defaults, put a `#` in front of the name. Refer to the
+        [ase.calculators.orca.ORCA][] calculator for details on `orcasimpleinput`.
+    orcablocks
+        List of `orcablocks` swaps for the calculator. To remove entries
+        from the defaults, put a `#` in front of the name. Refer to the
+        [ase.calculators.orca.ORCA][] calculator for details on `orcablocks`.
+    nprocs
+        Number of processors to use. Defaults to the number of physical cores.
+    copy_files
+        Files to copy (and decompress) from source to the runtime directory.
+
+    Returns
+    -------
+    cclibSchema
+        Dictionary of results from [quacc.schemas.cclib.cclib_summarize_run][].
+        See the type-hint for the data structure.
+    """
+    nprocs = psutil.cpu_count(logical=False) if nprocs == "max" else nprocs
+
+    default_inputs = [xc, basis, "normalprint", "numfreq" if numerical else "freq"]
+
+    default_blocks = [f"%pal nprocs {nprocs} end"]
+
+    return run_and_summarize(
+        atoms,
+        charge,
+        spin_multiplicity,
+        default_inputs=default_inputs,
+        default_blocks=default_blocks,
+        input_swaps=orcasimpleinput,
+        block_swaps=orcablocks,
+        additional_fields={"name": "ORCA vibrational frequency analysis"},
+        copy_files=copy_files,
+    )
+
+
 @job
 def relax_job(
     atoms: Atoms,
@@ -215,3 +285,84 @@ def ase_relax_job(
         additional_fields={"name": "ORCA ASE Relax"},
         copy_files=copy_files,
     )
+
+
+@job
+def ase_quasi_irc_perturb_job(
+    atoms: Atoms,
+    mode: list[list[float]] | NDArray,
+    charge: int = 0,
+    spin_multiplicity: int = 1,
+    perturb_magnitude: float = 0.6,
+    direction: Literal["forward", "reverse"] = "forward",
+    xc: str = "wb97x-d3bj",
+    basis: str = "def2-tzvp",
+    orcasimpleinput: list[str] | None = None,
+    orcablocks: list[str] | None = None,
+    opt_params: dict[str, Any] | None = None,
+    nprocs: int | Literal["max"] = "max",
+    copy_files: SourceDirectory | dict[SourceDirectory, Filenames] | None = None,
+) -> cclibSchema:
+    """
+    Quasi-IRC to optimize a reaction endpoint from a transition-state with known vibrational frequency modes.
+    Perturbs the structure of `atoms` by a finite amount (0.6 * the normalized mode magnitude) along the specified
+    vibrational frequency mode (assumed to be the transition mode), and then performs a `relax_job` on the perturbed
+    structure.
+
+    Parameters
+    ----------
+    atoms
+        Atoms object
+    mode
+        Transition mode. This should be an Nx3 matrix, where N is the number of atoms in `atoms`.
+    charge
+        Charge of the system.
+    spin_multiplicity
+        Multiplicity of the system.
+    perturb_magnitude
+        Factor to multiply the transition mode. Default is 0.6. In some cases, it may be advisable to increase this
+        factor, perhaps to 1.0 or 1.1. Lowering it is not generally found to be helpful.
+    direction
+        Direction of the (Quasi)IRC. Should be "forward" or "reverse".
+    xc
+        Exchange-correlation functional
+    basis
+        Basis set.
+    orcasimpleinput
+        List of `orcasimpleinput` swaps for the calculator. To remove entries
+        from the defaults, put a `#` in front of the name. Refer to the
+        [ase.calculators.orca.ORCA][] calculator for details on `orcasimpleinput`.
+    orcablocks
+        List of `orcablocks` swaps for the calculator. To remove entries
+        from the defaults, put a `#` in front of the name. Refer to the
+        [ase.calculators.orca.ORCA][] calculator for details on `orcablocks`.
+    nprocs
+        Number of processors to use. Defaults to the number of physical cores.
+    copy_files
+        Files to copy (and decompress) from source to the runtime directory.
+
+    Returns
+    -------
+    cclibASEOptSchema
+        Dictionary of results from [quacc.schemas.cclib.cclib_summarize_run][] merged with
+        the results from [quacc.schemas.ase.summarize_opt_run][].
+        See the type-hint for the data structure.
+    """
+    nprocs = psutil.cpu_count(logical=False) if nprocs == "max" else nprocs
+    default_inputs = [xc, basis, "engrad", "normalprint"]
+    default_blocks = [f"%pal nprocs {nprocs} end"]
+
+    scale = perturb_magnitude if direction == "forward" else perturb_magnitude * -1
+
+    return run_and_summarize_opt(
+        perturb(atoms, mode, scale),
+        charge=charge,
+        spin_multiplicity=spin_multiplicity,
+        default_inputs=default_inputs,
+        default_blocks=default_blocks,
+        input_swaps=orcasimpleinput,
+        block_swaps=orcablocks,
+        opt_params=opt_params,
+        additional_fields={"name": "ORCA ASE Quasi-IRC perturbed optimization"},
+        copy_files=copy_files,
+    )
diff --git a/src/quacc/recipes/qchem/ts.py b/src/quacc/recipes/qchem/ts.py
index fd5ed019..9f5a6cb6 100644
--- a/src/quacc/recipes/qchem/ts.py
+++ b/src/quacc/recipes/qchem/ts.py
@@ -7,6 +7,7 @@ from typing import TYPE_CHECKING
 from monty.dev import requires
 
 from quacc import SETTINGS, job, strip_decorator
+from quacc.atoms.core import perturb
 from quacc.recipes.qchem._base import run_and_summarize_opt
 from quacc.recipes.qchem.core import _BASE_SET, relax_job
 from quacc.utils.dicts import recursive_dict_merge
@@ -14,20 +15,23 @@ from quacc.utils.dicts import recursive_dict_merge
 try:
     from sella import IRC, Sella
 
+    has_sella = True
+
 except ImportError:
-    Sella = False
+    has_sella = False
 
 if TYPE_CHECKING:
     from typing import Any, Literal
 
     from ase.atoms import Atoms
+    from numpy.typing import NDArray
 
     from quacc.schemas._aliases.ase import OptSchema
     from quacc.utils.files import Filenames, SourceDirectory
 
 
 @job
-@requires(Sella, "Sella must be installed. Refer to the quacc documentation.")
+@requires(has_sella, "Sella must be installed. Refer to the quacc documentation.")
 def ts_job(
     atoms: Atoms,
     charge: int,
@@ -92,7 +96,7 @@ def ts_job(
 
 
 @job
-@requires(Sella, "Sella must be installed. Refer to the quacc documentation.")
+@requires(has_sella, "Sella must be installed. Refer to the quacc documentation.")
 def irc_job(
     atoms: Atoms,
     charge: int,
@@ -162,7 +166,7 @@ def irc_job(
 
 
 @job
-@requires(Sella, "Sella must be installed. Refer to the quacc documentation.")
+@requires(has_sella, "Sella must be installed. Refer to the quacc documentation.")
 def quasi_irc_job(
     atoms: Atoms,
     charge: int,
@@ -229,3 +233,81 @@ def quasi_irc_job(
     relax_summary["initial_irc"] = irc_summary
 
     return relax_summary
+
+
+@job
+@requires(has_sella, "Sella must be installed. Refer to the quacc documentation.")
+def quasi_irc_perturb_job(
+    atoms: Atoms,
+    charge: int,
+    spin_multiplicity: int,
+    mode: list[list[float]] | NDArray,
+    perturb_magnitude: float = 0.6,
+    direction: Literal["forward", "reverse"] = "forward",
+    method: str = "wb97mv",
+    basis: str = "def2-svpd",
+    opt_params: dict[str, Any] | None = None,
+    copy_files: SourceDirectory | dict[SourceDirectory, Filenames] | None = None,
+    **calc_kwargs,
+) -> OptSchema:
+    """
+    Quasi-IRC to optimize a reaction endpoint from a transition-state with known vibrational frequency modes.
+    Perturbs the structure of `atoms` by a finite amount (0.6 * the normalized mode magnitude) along the specified
+    vibrational frequency mode (assumed to be the transition mode), and then performs a `relax_job` on the perturbed
+    structure.
+
+    Parameters
+    ----------
+    atoms
+        Atoms object.
+    charge
+        Charge of the system.
+    spin_multiplicity
+        Multiplicity of the system.
+    mode
+        Transition mode. This should be an Nx3 matrix, where N is the number of atoms in `atoms`.
+    perturb_magnitude
+        Factor to multiply the transition mode. Default is 0.6. In some cases, it may be advisable to increase this
+        factor, perhaps to 1.0 or 1.1. Lowering it is not generally found to be helpful.
+    direction
+        Direction of the (Quasi)IRC. Should be "forward" or "reverse".
+    method
+        DFT exchange-correlation functional or other electronic structure
+        method.
+    basis
+        Basis set.
+    opt_params
+        Dictionary of custom kwargs for the optimization process. For a list
+        of available keys, refer to [quacc.runners.ase.run_opt][].
+    copy_files
+        Files to copy (and decompress) from source to the runtime directory.
+    **calc_kwargs
+        Custom kwargs for the calculator. Set a value to `quacc.Remove` to remove
+        a pre-existing key entirely. See [quacc.calculators.qchem.qchem.QChem][] for more
+        details.
+
+
+    Returns
+    -------
+    OptSchema
+        Dictionary of results from [quacc.schemas.ase.summarize_opt_run][]
+    """
+
+    calc_defaults = recursive_dict_merge(
+        _BASE_SET, {"rem": {"job_type": "force", "method": method, "basis": basis}}
+    )
+    opt_defaults = {"optimizer": Sella} if has_sella else {}
+
+    scale = perturb_magnitude if direction == "forward" else perturb_magnitude * -1
+
+    return run_and_summarize_opt(
+        perturb(atoms, mode, scale),
+        charge=charge,
+        spin_multiplicity=spin_multiplicity,
+        calc_defaults=calc_defaults,
+        calc_swaps=calc_kwargs,
+        opt_defaults=opt_defaults,
+        opt_params=opt_params,
+        additional_fields={"name": "Q-Chem Quasi-IRC perturbed optimization"},
+        copy_files=copy_files,
+    )


Test patch:
diff --git a/tests/core/recipes/orca_recipes/test_orca_recipes.py b/tests/core/recipes/orca_recipes/test_orca_recipes.py
index b919e0a7..555bbb5d 100644
--- a/tests/core/recipes/orca_recipes/test_orca_recipes.py
+++ b/tests/core/recipes/orca_recipes/test_orca_recipes.py
@@ -5,8 +5,22 @@ from pathlib import Path
 
 import pytest
 from ase.build import molecule
+from ase.io import read
 
-from quacc.recipes.orca.core import ase_relax_job, relax_job, static_job
+from quacc.recipes.orca.core import (
+    ase_quasi_irc_perturb_job,
+    ase_relax_job,
+    freq_job,
+    relax_job,
+    static_job,
+)
+
+FILE_DIR = Path(__file__).parent
+
+
+@pytest.fixture()
+def test_atoms():
+    return read(FILE_DIR / "xyz" / "ts_test.xyz")
 
 
 def test_static_job(tmp_path, monkeypatch):
@@ -141,3 +155,66 @@ def test_ase_relax_job_store(tmp_path, monkeypatch):
         assert "orca.xyz.gz" in os.listdir(Path(output["dir_name"], f"step{i}"))
     assert len(output["steps"]) == nsteps
     assert "attributes" in output["steps"][0]
+
+
+@pytest.mark.skipif(os.name == "nt", reason="mpirun not available on Windows")
+def test_freq_job(tmp_path, monkeypatch):
+    monkeypatch.chdir(tmp_path)
+
+    atoms = molecule("H2")
+
+    output = freq_job(
+        atoms,
+        xc="hf",
+        basis="def2-svp",
+        charge=0,
+        spin_multiplicity=1,
+        orcasimpleinput=["#normalprint"],
+    )
+    assert output["natoms"] == len(atoms)
+    assert output["parameters"]["charge"] == 0
+    assert output["parameters"]["mult"] == 1
+    assert output["parameters"]["orcasimpleinput"] == "def2-svp freq hf xyzfile"
+    assert output.get("attributes")
+
+
+@pytest.mark.skipif(os.name == "nt", reason="mpirun not available on Windows")
+def test_ase_quasi_irc_perturb_job(test_atoms, tmp_path, monkeypatch):
+    monkeypatch.chdir(tmp_path)
+
+    mode = [
+        [-0.164243, 0.289973, 0.026016],
+        [0.111755, -0.021282, -0.002864],
+        [0.011866, -0.071662, -0.040428],
+        [-0.087153, 0.038657, -0.038008],
+        [-0.017777, 0.013758, 0.001151],
+        [0.030409, -0.18767, 0.028847],
+        [0.750577, -0.378458, 0.179569],
+        [0.042058, 0.035122, 0.024727],
+        [-0.008014, -0.000956, -0.009153],
+        [-0.052972, -0.176662, -0.077928],
+        [0.037381, 0.036482, 0.028861],
+        [0.043279, 0.040925, 0.022537],
+        [0.035434, 0.032613, 0.019516],
+        [-0.002674, -0.03398, 0.011123],
+        [-0.006118, -0.009193, -0.122432],
+        [0.014124, -0.035613, 0.097518],
+    ]
+
+    output = ase_quasi_irc_perturb_job(
+        atoms=test_atoms,
+        mode=mode,
+        charge=0,
+        spin_multiplicity=1,
+        perturb_magnitude=1.0,
+        direction="reverse",
+        xc="hf",
+        basis="def2-svp",
+        orcasimpleinput=["#normalprint"],
+    )
+    assert output["natoms"] == len(test_atoms)
+    assert output["atoms"] != test_atoms
+    assert output["parameters"]["charge"] == 0
+    assert output["parameters"]["mult"] == 1
+    assert output["parameters"]["orcasimpleinput"] == "def2-svp engrad hf xyzfile"
+    assert output.get("attributes")
diff --git a/tests/core/recipes/orca_recipes/xyz/ts_test.xyz b/tests/core/recipes/orca_recipes/xyz/ts_test.xyz
new file mode 100644
index 00000000..6f68b038
--- /dev/null
+++ b/tests/core/recipes/orca_recipes/xyz/ts_test.xyz
@@ -0,0 +1,18 @@
+16
+
+C     -2.1846949292      0.1898391258      0.1964826227
+O     -0.4848515165     -0.2917846139     -0.1275732483
+C      0.5394541212      0.3658912018      0.2429641752
+O      0.4798059249      1.5225609297      0.7410807853
+O     -2.8764467052     -0.1509181676     -0.6925854737
+O     -1.7935252197      1.8615485866      0.1068277631
+H     -0.7532722964      1.8974038694      0.4912008424
+C      1.8832273355     -0.2733121087      0.0407321421
+C     -2.2825562627     -0.1155835245      1.6624047756
+H     -2.4183125325      2.3939897026      0.6146768439
+H      1.7851969133     -1.3556014160     -0.0789960465
+H      2.5408551922     -0.0219920239      0.8792458113
+H      2.3250034882      0.1492594553     -0.8712799313
+H     -1.4650505694      0.3352661303      2.2318827567
+H     -2.2601855390     -1.2048429495      1.7751129189
+H     -3.2478474045      0.2577558028      2.0283932628
diff --git a/tests/core/recipes/qchem_recipes/mocked/qchem_examples/mol.qin.qirc_forward b/tests/core/recipes/qchem_recipes/mocked/qchem_examples/mol.qin.qirc_forward
new file mode 100644
index 00000000..4deffc31
--- /dev/null
+++ b/tests/core/recipes/qchem_recipes/mocked/qchem_examples/mol.qin.qirc_forward
@@ -0,0 +1,34 @@
+$molecule
+ 0 1
+ C     -2.2327219082      0.3791393294      0.2143702282
+ O     -0.4001040911     -0.3209141098     -0.1334138299
+ C      0.5423733258      0.2977790640      0.2206586158
+ O      0.3924968984      1.4585347680      0.7033089993
+ O     -2.8376665474     -0.0973158207     -0.6674694144
+ O     -1.7015955091      1.7258053313      0.1289163517
+ H     -0.5810710172      1.5729449388      0.5275661907
+ C      1.8860908742     -0.2457908005      0.0738796537
+ C     -2.2383734832     -0.0934887378      1.6277789473
+ H     -2.4617721877      2.3579851511      0.6488511867
+ H      1.8010691866     -1.3457375140     -0.0557562494
+ H      2.5419651710      0.0138395515      0.9320914408
+ H      2.3352038754      0.1921439154     -0.8423775733
+ H     -1.3918401900      0.3251951235      2.2131306207
+ H     -2.2407267711     -1.2034273556      1.6755719861
+ H     -3.1861276262      0.2835871654      2.0670628463
+$end
+
+$rem
+   basis = def2-svpd
+   gen_scfman = true
+   job_type = force
+   max_scf_cycles = 100
+   method = wb97mv
+   resp_charges = true
+   s2thresh = 16
+   scf_algorithm = diis
+   sym_ignore = true
+   symmetry = false
+   thresh = 14
+   xc_grid = 3
+$end
diff --git a/tests/core/recipes/qchem_recipes/mocked/qchem_examples/mol.qin.qirc_reverse b/tests/core/recipes/qchem_recipes/mocked/qchem_examples/mol.qin.qirc_reverse
new file mode 100644
index 00000000..48875c1f
--- /dev/null
+++ b/tests/core/recipes/qchem_recipes/mocked/qchem_examples/mol.qin.qirc_reverse
@@ -0,0 +1,34 @@
+$molecule
+ -1 2
+ C     -2.0137704788     -0.1270734858      0.1642010188
+ O     -0.5995953068     -0.2921342976     -0.1231428019
+ C      0.5061642245      0.4275172940      0.2730711793
+ O      0.5334437588      1.4727239578      0.7474905023
+ O     -2.8412706626     -0.2099655589     -0.6873946438
+ O     -1.8288928518      2.0165855272      0.0417843943
+ H     -1.3429919048      2.3740186237      0.3845744402
+ C      1.8193658177     -0.3037461746      0.0198838218
+ C     -2.2672351085     -0.1210357673      1.6552929612
+ H     -2.4429538836      2.6113929006      0.7209957866
+ H      1.7385359612     -1.3967031272     -0.0862781690
+ H      2.4755055076     -0.0345738984      0.8626997493
+ H      2.2676650535      0.1131019024     -0.8976730188
+ H     -1.4518547674      0.3782355629      2.2010561136
+ H     -2.2382992481     -1.2010725621      1.8946348321
+ H     -3.2610161107      0.2842091034      1.9133738342
+$end
+
+$rem
+   basis = def2-tzvpd
+   gen_scfman = true
+   job_type = force
+   max_scf_cycles = 200
+   method = wb97mv
+   resp_charges = true
+   s2thresh = 16
+   scf_algorithm = gdm
+   sym_ignore = true
+   symmetry = false
+   thresh = 14
+   xc_grid = 3
+$end
diff --git a/tests/core/recipes/qchem_recipes/mocked/test_qchem_recipes.py b/tests/core/recipes/qchem_recipes/mocked/test_qchem_recipes.py
index a3aee3a6..89919f4b 100644
--- a/tests/core/recipes/qchem_recipes/mocked/test_qchem_recipes.py
+++ b/tests/core/recipes/qchem_recipes/mocked/test_qchem_recipes.py
@@ -14,7 +14,7 @@ from quacc import SETTINGS
 from quacc.atoms.core import check_charge_and_spin
 from quacc.calculators.qchem import QChem
 from quacc.recipes.qchem.core import freq_job, relax_job, static_job
-from quacc.recipes.qchem.ts import irc_job, quasi_irc_job, ts_job
+from quacc.recipes.qchem.ts import irc_job, quasi_irc_job, quasi_irc_perturb_job, ts_job
 
 try:
     import sella
@@ -32,6 +32,11 @@ def test_atoms():
     return read(FILE_DIR / "xyz" / "test.xyz")
 
 
+@pytest.fixture()
+def test_qirc_atoms():
+    return read(FILE_DIR / "xyz" / "ts_test.xyz")
+
+
 @pytest.fixture()
 def os_atoms():
     return read(FILE_DIR / "xyz" / "OS_test.xyz")
@@ -589,3 +594,79 @@ def test_quasi_irc_job(monkeypatch, tmp_path, test_atoms):
     qcin = QCInput.from_file(str(Path(output["dir_name"], "mol.qin.gz")))
     ref_qcin = QCInput.from_file(str(QCHEM_DIR / "mol.qin.quasi_irc_reverse"))
     qcinput_nearly_equal(qcin, ref_qcin)
+
+
+@pytest.mark.skipif(sella is None, reason="Does not have Sella")
+def test_quasi_irc_perturb_job(monkeypatch, tmp_path, test_qirc_atoms):
+    monkeypatch.chdir(tmp_path)
+
+    monkeypatch.setattr(QChem, "read_results", mock_read)
+    monkeypatch.setattr(QChem, "execute", mock_execute4)
+
+    # Transition mode for this transition-state
+    mode = [
+        [-0.164, 0.289, 0.027],
+        [0.112, -0.02, -0.004],
+        [0.012, -0.072, -0.042],
+        [-0.087, 0.039, -0.038],
+        [-0.017, 0.013, 0.001],
+        [0.028, -0.186, 0.028],
+        [0.751, -0.378, 0.186],
+        [0.042, 0.034, 0.025],
+        [-0.007, -0.001, -0.009],
+        [-0.056, -0.179, -0.076],
+        [0.036, 0.035, 0.027],
+        [0.043, 0.037, 0.023],
+        [0.036, 0.032, 0.021],
+        [-0.003, -0.032, 0.011],
+        [-0.006, -0.009, -0.118],
+        [0.014, -0.034, 0.094],
+    ]
+
+    charge, spin_multiplicity = check_charge_and_spin(test_qirc_atoms)
+    output = quasi_irc_perturb_job(
+        atoms=test_qirc_atoms,
+        charge=charge,
+        spin_multiplicity=spin_multiplicity,
+        mode=mode,
+        direction="forward",
+        method="wb97mv",
+        opt_params={"max_steps": 5},
+        basis="def2-svpd",
+    )
+
+    assert output["atoms"] != test_qirc_atoms
+    assert output["charge"] == 0
+    assert output["spin_multiplicity"] == 1
+    assert output["formula_alphabetical"] == "C4 H8 O4"
+    assert output["nelectrons"] == 64
+    assert output["parameters"]["charge"] == 0
+    assert output["parameters"]["spin_multiplicity"] == 1
+
+    qcin = QCInput.from_file(str(Path(output["dir_name"], "mol.qin.gz")))
+    ref_qcin = QCInput.from_file(str(QCHEM_DIR / "mol.qin.qirc_forward"))
+    qcinput_nearly_equal(qcin, ref_qcin)
+
+    output = quasi_irc_perturb_job(
+        test_qirc_atoms,
+        charge=-1,
+        spin_multiplicity=2,
+        mode=mode,
+        perturb_magnitude=1.0,
+        direction="reverse",
+        basis="def2-tzvpd",
+        opt_params={"max_steps": 6},
+        rem={"scf_algorithm": "gdm"},
+    )
+
+    assert output["atoms"] != test_qirc_atoms
+    assert output["charge"] == -1
+    assert output["spin_multiplicity"] == 2
+    assert output["formula_alphabetical"] == "C4 H8 O4"
+    assert output["nelectrons"] == 65
+    assert output["parameters"]["charge"] == -1
+    assert output["parameters"]["spin_multiplicity"] == 2
+
+    qcin = QCInput.from_file(str(Path(output["dir_name"], "mol.qin.gz")))
+    ref_qcin = QCInput.from_file(str(QCHEM_DIR / "mol.qin.qirc_reverse"))
+    qcinput_nearly_equal(qcin, ref_qcin)
diff --git a/tests/core/recipes/qchem_recipes/mocked/xyz/ts_test.xyz b/tests/core/recipes/qchem_recipes/mocked/xyz/ts_test.xyz
new file mode 100644
index 00000000..6f68b038
--- /dev/null
+++ b/tests/core/recipes/qchem_recipes/mocked/xyz/ts_test.xyz
@@ -0,0 +1,18 @@
+16
+
+C     -2.1846949292      0.1898391258      0.1964826227
+O     -0.4848515165     -0.2917846139     -0.1275732483
+C      0.5394541212      0.3658912018      0.2429641752
+O      0.4798059249      1.5225609297      0.7410807853
+O     -2.8764467052     -0.1509181676     -0.6925854737
+O     -1.7935252197      1.8615485866      0.1068277631
+H     -0.7532722964      1.8974038694      0.4912008424
+C      1.8832273355     -0.2733121087      0.0407321421
+C     -2.2825562627     -0.1155835245      1.6624047756
+H     -2.4183125325      2.3939897026      0.6146768439
+H      1.7851969133     -1.3556014160     -0.0789960465
+H      2.5408551922     -0.0219920239      0.8792458113
+H      2.3250034882      0.1492594553     -0.8712799313
+H     -1.4650505694      0.3352661303      2.2318827567
+H     -2.2601855390     -1.2048429495      1.7751129189
+H     -3.2478474045      0.2577558028      2.0283932628