recipes

__all__ module-attribute

__all__ = ['dft_ensemble_analyse', 'dft_ensemble_flow', 'get_final_autoSKZCAM_Hads', 'skzcam_analyse', 'skzcam_eint_flow', 'skzcam_initialise']

dft_ensemble_analyse

dft_ensemble_analyse(calc_dir: Path | str, xc_ensemble: list[str] | dict[str, str], geom_error_xc: str, vib_xc_ensemble: list[str], freeze_surface_vib: bool, temperature: float = 200.0) -> dict[str, list[float]]

Analyses the completed DFT ensemble calculations.

Parameters:

• calc_dir (Path or str) –

  The directory where the calculations were performed.

• xc_ensemble (dict[str, str]) –

  A dictionary containing the xc functionals to be used as keys and the corresponding settings as values.

• geom_error_xc (str) –

  The xc functional to be used for the geometry error calculation.

• vib_xc_ensemble (list[str]) –

  A list of xc functionals for which the vibrational calculations were performed.

• freeze_surface_vib (bool) –

  True if the vibrational calculations on the slab should be skipped.

• temperature (float, default: 200.0 ) –

  The temperature to get the vibrational contributions to.

Returns:

• dict[str, list[float]] –

  A dictionary containing the relaxation energy (and its geometry error) and DeltaH contributions from the DFT ensemble.

Source code in autoSKZCAM/recipes_dft.py

def dft_ensemble_analyse(
    calc_dir: Path | str,
    xc_ensemble: list[str] | dict[str, str],
    geom_error_xc: str,
    vib_xc_ensemble: list[str],
    freeze_surface_vib: bool,
    temperature: float = 200.0,
) -> dict[str, list[float]]:
    """
    Analyses the completed DFT ensemble calculations.

    Parameters
    ----------
    calc_dir : Path or str
        The directory where the calculations were performed.
    xc_ensemble : dict[str, str]
        A dictionary containing the xc functionals to be used as keys and the corresponding settings as values.
    geom_error_xc : str
        The xc functional to be used for the geometry error calculation.
    vib_xc_ensemble : list[str]
        A list of xc functionals for which the vibrational calculations were performed.
    freeze_surface_vib : bool
        True if the vibrational calculations on the slab should be skipped.
    temperature : float
        The temperature to get the vibrational contributions to.

    Returns
    -------
    dict[str, list[float]]
        A dictionary containing the relaxation energy (and its geometry error) and DeltaH contributions from the DFT ensemble.

    """

    # Ensure that all the functionals in vib_xc_ensemble are also in xc_ensemble
    for vib_xc in vib_xc_ensemble:
        if vib_xc not in xc_ensemble:
            raise ValueError(
                f"The functional {vib_xc} in vib_xc_ensemble is not in the xc_ensemble."
            )

    # Ensure that the geom_error_xc is in the xc_ensemble
    if geom_error_xc is not None and geom_error_xc not in xc_ensemble:
        raise ValueError(
            f"The functional {geom_error_xc} in geom_error_xc is not in the xc_ensemble."
        )

    dft_ensemble_results = read_completed_calculations(
        calc_dir, xc_ensemble, vib_xc_ensemble, freeze_surface_vib
    )

    # Confirm that all the calculations have been completed
    for job_type in dft_ensemble_results:
        for xc_func in xc_ensemble:
            if "vib" in job_type and (
                xc_func not in vib_xc_ensemble
                or (job_type == "07-slab_vib" and freeze_surface_vib)
            ):
                continue
            if dft_ensemble_results[job_type][xc_func] is None or (
                "results" in dft_ensemble_results[job_type]
                and "energy" not in dft_ensemble_results[job_type][xc_func]["results"]
            ):
                raise ValueError(
                    f"The {job_type} calculation for the functional {xc_func} has not been completed."
                )

    xc_eads_dict = dict.fromkeys(xc_ensemble, 0)
    xc_eint_dict = dict.fromkeys(xc_ensemble, 0)
    xc_vib_dict = dict.fromkeys(vib_xc_ensemble, 0)

    for xc_func in xc_ensemble:
        xc_eads_dict[xc_func] = (
            dft_ensemble_results["04-adsorbate_slab"][xc_func]["results"]["energy"]
            - dft_ensemble_results["03-slab"][xc_func]["results"]["energy"]
            - dft_ensemble_results["01-molecule"][xc_func]["results"]["energy"]
        )
        xc_eint_dict[xc_func] = (
            dft_ensemble_results["08-eint_adsorbate_slab"][xc_func]["results"]["energy"]
            - dft_ensemble_results["09-eint_adsorbate"][xc_func]["results"]["energy"]
            - dft_ensemble_results["10-eint_slab"][xc_func]["results"]["energy"]
        )

        if xc_func in vib_xc_ensemble:
            adsorbate_slab_dU, _, _, _ = get_quasi_rrho(
                dft_ensemble_results["05-adsorbate_slab_vib"][xc_func]["results"][
                    "real_vib_freqs"
                ],
                dft_ensemble_results["05-adsorbate_slab_vib"][xc_func]["results"][
                    "imag_vib_freqs"
                ],
                temperature,
            )
            adsorbate_dU, _, _, _ = get_quasi_rrho(
                dft_ensemble_results["06-molecule_vib"][xc_func]["results"][
                    "real_vib_freqs"
                ],
                dft_ensemble_results["06-molecule_vib"][xc_func]["results"][
                    "imag_vib_freqs"
                ],
                temperature,
            )

            if freeze_surface_vib is False:
                slab_dU, _, _, _ = get_quasi_rrho(
                    dft_ensemble_results["07-slab_vib"][xc_func]["results"][
                        "real_vib_freqs"
                    ],
                    dft_ensemble_results["07-slab_vib"][xc_func]["results"][
                        "imag_vib_freqs"
                    ],
                    temperature,
                )
            else:
                slab_dU = 0
            xc_vib_dict[xc_func] = adsorbate_slab_dU - adsorbate_dU - slab_dU

    erlx = xc_eads_dict[geom_error_xc] - xc_eint_dict[geom_error_xc]
    geom_error = 2 * np.sqrt(
        np.mean(
            [
                (xc_eads_dict[xc_func] - xc_eint_dict[xc_func] - erlx) ** 2
                for xc_func in xc_ensemble
                if xc_func != geom_error_xc
            ]
        )
    )
    return {
        "DFT Erlx": [erlx * 1000, geom_error * 1000],
        "DFT DeltaH": [
            np.mean([xc_vib_dict[xc_func] for xc_func in vib_xc_ensemble]),
            2 * np.std([xc_vib_dict[xc_func] for xc_func in vib_xc_ensemble]),
        ],
    }

dft_ensemble_flow

dft_ensemble_flow(xc_ensemble: dict[str, dict[str, Any]], vib_xc_ensemble: list[str] | None = None, geom_error_xc: str | None = None, freeze_surface_vib: bool = True, job_params: dict[str, dict[str, Any]] | None = None, adsorbate: Atoms | None = None, unit_cell: Atoms | None = None, calc_dir: str | Path = './dft_calc_dir', slab_gen_func: Callable[[Atoms], Atoms] | None = None, adsorbate_slab_gen_func: Callable[[Atoms], Atoms] | None = None)

Workflow to perform the DFT ensemble calculations to obtain the geometry error and get the DeltaH contribution. The workflow consists of the following steps:

1. Relax the gas-phase molecule for each functional in the ensemble.

2. Relax the unit cell for each functional in the ensemble.

3. Generate and relax the slab from the relaxed solid for each functional in the ensemble.

4. Generate and relax the adsorbate-slab complex from the relaxed adsorbate and slab for each functional in the ensemble.

5. Perform the vibrational calculation for each functional in the ensemble.

6. Perform the eint calculation on the chosen functional for each functional in the ensemble.

Parameters:

• xc_ensemble (dict[str, dict[str, Any]]) –

  A dictionary containing the xc functionals to be used as keys and the corresponding settings as values.

• job_params (dict[str, dict[str, Any]], default: None ) –

  A dictionary containing the job parameters to be used for each functional in the ensemble. If not provided, the default job parameters will be used.

• adsorbate (Atoms, default: None ) –

  The adsorbate molecule to be used for the calculations. If not provided, will attempt to read in the adsorbate from the calc_dir.

• unit_cell (Atoms, default: None ) –

  The unit cell of the solid to be used for the calculations. If not provided, will attempt to read in the unit cell from the calc_dir.

• calc_dir (str or Path, default: './dft_calc_dir' ) –

  The directory where the calculations will be performed. Defaults to './calc_dir'.

• slab_gen_func (Callable[[Atoms], Atoms], default: None ) –

  The function to generate the slab from the unit cell.

• adsorbate_slab_gen_func (Callable[[Atoms], Atoms], default: None ) –

  The function to generate the adsorbate molecule. It is important that the indices of the adsorbates are always the first indices in the Atoms object, followed by the slab Atoms object.

Returns:

• dict[str, dic[str, VaspSchema]] –

  A dictionary containing the results of the DFT ensemble calculations for each functional in the ensemble.

Source code in autoSKZCAM/recipes_dft.py

@flow
def dft_ensemble_flow(
    xc_ensemble: dict[str, dict[str, Any]],
    vib_xc_ensemble: list[str] | None = None,
    geom_error_xc: str | None = None,
    freeze_surface_vib: bool = True,
    job_params: dict[str, dict[str, Any]] | None = None,
    adsorbate: Atoms | None = None,
    unit_cell: Atoms | None = None,
    calc_dir: str | Path = "./dft_calc_dir",
    slab_gen_func: Callable[[Atoms], Atoms] | None = None,
    adsorbate_slab_gen_func: Callable[[Atoms], Atoms] | None = None,
):
    """
    Workflow to perform the DFT ensemble calculations to obtain the geometry error and get the DeltaH contribution.
    The workflow consists of the following steps:

    1. Relax the gas-phase molecule for each functional in the ensemble.

    2. Relax the unit cell for each functional in the ensemble.

    3. Generate and relax the slab from the relaxed solid for each functional in the ensemble.

    4. Generate and relax the adsorbate-slab complex from the relaxed adsorbate and slab for each functional in the ensemble.

    5. Perform the vibrational calculation for each functional in the ensemble.

    6. Perform the eint calculation on the chosen functional for each functional in the ensemble.

    Parameters
    ----------
    xc_ensemble : dict[str, dict[str, Any]]
        A dictionary containing the xc functionals to be used as keys and the corresponding settings as values.
    job_params : dict[str, dict[str, Any]], optional
        A dictionary containing the job parameters to be used for each functional in the ensemble. If not provided, the default job parameters will be used.
    adsorbate : Atoms, optional
        The adsorbate molecule to be used for the calculations. If not provided, will attempt to read in the adsorbate from the calc_dir.
    unit_cell : Atoms, optional
        The unit cell of the solid to be used for the calculations. If not provided, will attempt to read in the unit cell from the calc_dir.
    calc_dir : str or Path, optional
        The directory where the calculations will be performed. Defaults to './calc_dir'.
    slab_gen_func : Callable[[Atoms], Atoms]
        The function to generate the slab from the unit cell.
    adsorbate_slab_gen_func : Callable[[Atoms], Atoms]
        The function to generate the adsorbate molecule. It is important that the indices of the adsorbates are always the first indices in the Atoms object, followed by the slab Atoms object.

    Returns
    -------
    dict[str, dic[str,VaspSchema]]
        A dictionary containing the results of the DFT ensemble calculations for each functional in the ensemble.

    """

    if vib_xc_ensemble is None:
        vib_xc_ensemble = []

    job_list = [
        "01-molecule",
        "02-unit_cell",
        "03-slab",
        "04-adsorbate_slab",
        "05-adsorbate_slab_vib",
        "06-molecule_vib",
        "07-slab_vib",
        "08-eint_adsorbate_slab",
        "09-eint_adsorbate",
        "10-eint_slab",
    ]

    # Ensure that all the functionals in vib_xc_ensemble are also in xc_ensemble
    for vib_xc in vib_xc_ensemble:
        if vib_xc not in xc_ensemble:
            raise ValueError(
                f"The functional {vib_xc} in vib_xc_ensemble is not in the xc_ensemble."
            )

    # Ensure that the geom_error_xc is in the xc_ensemble
    if geom_error_xc is not None and geom_error_xc not in xc_ensemble:
        raise ValueError(
            f"The functional {geom_error_xc} in geom_error_xc is not in the xc_ensemble."
        )

    if job_params is not None:
        for job_type in job_params:
            if job_type not in job_list:
                raise ValueError(
                    f"The {job_type} key in job_params is not valid. Please choose from the following: '01-molecule', '02-unit_cell', '03-slab', '04-adsorbate_slab', '05-adsorbate_slab_vib', '06-molecule_vib', '07-slab_vib', '08-eint_adsorbate_slab', '09-eint_adsorbate', '10-eint_slab'."
                )
    else:
        job_params = {}

    # Try to read in completed calculations from the calc_dir
    dft_ensemble_results = read_completed_calculations(
        calc_dir, xc_ensemble, vib_xc_ensemble, freeze_surface_vib
    )

    for xc_func, xc_func_kwargs in xc_ensemble.items():
        # relax molecule
        if dft_ensemble_results["01-molecule"][xc_func] is None:
            calc_kwargs = {**job_params.get("01-molecule", {}), **xc_func_kwargs}
            dft_ensemble_results["01-molecule"][xc_func] = relax_job(
                adsorbate,
                additional_fields={
                    "calc_results_dir": Path(calc_dir, "01-molecule", xc_func)
                },
                pmg_kpts={"kppvol": 1},
                **calc_kwargs,
            )

        # relax solid
        if dft_ensemble_results["02-unit_cell"][xc_func] is None:
            calc_kwargs = {**job_params.get("02-unit_cell", {}), **xc_func_kwargs}
            relax_job1 = relax_job(
                unit_cell,
                additional_fields={
                    "calc_results_dir": Path(calc_dir, "02-unit_cell", xc_func)
                },
                relax_cell=True,
                unique_dir=True,
                **calc_kwargs,
            )
            dft_ensemble_results["02-unit_cell"][xc_func] = relax_job(
                relax_job1["atoms"],
                relax_cell=True,
                additional_fields={
                    "calc_results_dir": Path(calc_dir, "02-unit_cell", xc_func),
                    "relax1": relax_job1,
                },
                **calc_kwargs,
            )

        # bulk to slab
        if dft_ensemble_results["03-slab"][xc_func] is None:
            calc_kwargs = {**job_params.get("03-slab", {}), **xc_func_kwargs}
            initial_slab = slab_gen_func(
                dft_ensemble_results["02-unit_cell"][xc_func]["atoms"]
            )
            dft_ensemble_results["03-slab"][xc_func] = relax_job(
                initial_slab,
                additional_fields={
                    "calc_results_dir": Path(calc_dir, "03-slab", xc_func)
                },
                pmg_kpts={"length_densities": [50, 50, 1]},
                **calc_kwargs,
            )

        # slab to ads_slab
        if dft_ensemble_results["04-adsorbate_slab"][xc_func] is None:
            calc_kwargs = {**job_params.get("04-adsorbate_slab", {}), **xc_func_kwargs}
            initial_adsorbate_slab = adsorbate_slab_gen_func(
                dft_ensemble_results["01-molecule"][xc_func]["atoms"],
                dft_ensemble_results["03-slab"][xc_func]["atoms"],
            )

            dft_ensemble_results["04-adsorbate_slab"][xc_func] = relax_job(
                initial_adsorbate_slab,
                additional_fields={
                    "calc_results_dir": Path(calc_dir, "04-adsorbate_slab", xc_func)
                },
                pmg_kpts={"length_densities": [50, 50, 1]},
                **calc_kwargs,
            )

        if xc_func in vib_xc_ensemble:
            # vibrational calculation
            if dft_ensemble_results["05-adsorbate_slab_vib"][xc_func] is None:
                calc_kwargs = {
                    **job_params.get("05-adsorbate_slab_vib", {}),
                    **xc_func_kwargs,
                }
                adsorbate_slab_vib_atoms = dft_ensemble_results["04-adsorbate_slab"][
                    xc_func
                ]["atoms"]

                if freeze_surface_vib:
                    adsorbate_len = len(
                        dft_ensemble_results["01-molecule"][xc_func]["atoms"]
                    )
                    slab_indices = list(range(len(adsorbate_slab_vib_atoms)))[
                        adsorbate_len:
                    ]
                    adsorbate_slab_vib_atoms.set_constraint(
                        FixAtoms(indices=slab_indices)
                    )

                dft_ensemble_results["05-adsorbate_slab_vib"][xc_func] = freq_job(
                    adsorbate_slab_vib_atoms,
                    additional_fields={
                        "calc_results_dir": Path(
                            calc_dir, "05-adsorbate_slab_vib", xc_func
                        )
                    },
                    pmg_kpts={"length_densities": [50, 50, 1]},
                    **calc_kwargs,
                )

            if dft_ensemble_results["06-molecule_vib"][xc_func] is None:
                calc_kwargs = {
                    **job_params.get("06-molecule_vib", {}),
                    **xc_func_kwargs,
                }
                dft_ensemble_results["06-molecule_vib"][xc_func] = freq_job(
                    dft_ensemble_results["01-molecule"][xc_func]["atoms"],
                    additional_fields={
                        "calc_results_dir": Path(calc_dir, "06-molecule_vib", xc_func)
                    },
                    pmg_kpts={"kppvol": 1},
                    **calc_kwargs,
                )

            if (
                dft_ensemble_results["07-slab_vib"][xc_func] is None
                and freeze_surface_vib is False
            ):
                calc_kwargs = {**job_params.get("07-slab_vib", {}), **xc_func_kwargs}
                dft_ensemble_results["07-slab_vib"][xc_func] = freq_job(
                    dft_ensemble_results["03-slab"][xc_func]["atoms"],
                    additional_fields={
                        "calc_results_dir": Path(calc_dir, "07-slab_vib", xc_func)
                    },
                    pmg_kpts={"length_densities": [50, 50, 1]},
                    **calc_kwargs,
                )

    if geom_error_xc is not None:
        adsorbate_len = len(dft_ensemble_results["01-molecule"][geom_error_xc]["atoms"])
        adsorbate_slab_atoms = dft_ensemble_results["04-adsorbate_slab"][geom_error_xc][
            "atoms"
        ]
        fixed_adsorbate_atoms = dft_ensemble_results["04-adsorbate_slab"][
            geom_error_xc
        ]["atoms"][:adsorbate_len]
        fixed_slab_atoms = dft_ensemble_results["04-adsorbate_slab"][geom_error_xc][
            "atoms"
        ][adsorbate_len:]

        for xc_func in xc_ensemble:
            # eint on adsorbate_slab
            if dft_ensemble_results["08-eint_adsorbate_slab"][xc_func] is None:
                calc_kwargs = {
                    **job_params.get("08-eint_adsorbate_slab", {}),
                    **xc_func_kwargs,
                }
                dft_ensemble_results["08-eint_adsorbate_slab"][xc_func] = static_job(
                    adsorbate_slab_atoms,
                    additional_fields={
                        "calc_results_dir": Path(
                            calc_dir, "08-eint_adsorbate_slab", xc_func
                        )
                    },
                    pmg_kpts={"length_densities": [50, 50, 1]},
                    **calc_kwargs,
                )

            # eint on adsorbate
            if dft_ensemble_results["09-eint_adsorbate"][xc_func] is None:
                calc_kwargs = {
                    **job_params.get("09-eint_adsorbate", {}),
                    **xc_func_kwargs,
                }
                dft_ensemble_results["09-eint_adsorbate"][xc_func] = static_job(
                    fixed_adsorbate_atoms,
                    additional_fields={
                        "calc_results_dir": Path(calc_dir, "09-eint_adsorbate", xc_func)
                    },
                    pmg_kpts={"length_densities": [50, 50, 1]},
                    **calc_kwargs,
                )

            # eint on slab
            if dft_ensemble_results["10-eint_slab"][xc_func] is None:
                calc_kwargs = {**job_params.get("10-eint_slab", {}), **xc_func_kwargs}
                dft_ensemble_results["10-eint_slab"][xc_func] = static_job(
                    fixed_slab_atoms,
                    additional_fields={
                        "calc_results_dir": Path(calc_dir, "10-eint_slab", xc_func)
                    },
                    pmg_kpts={"length_densities": [50, 50, 1]},
                    **calc_kwargs,
                )

    return dft_ensemble_results

get_final_autoSKZCAM_Hads

get_final_autoSKZCAM_Hads(skzcam_eint_analysis: dict[str, list[float]], dft_ensemble_analysis: dict[str, list[float]]) -> dict[str, list[float]]

Gets the final Hads from the autoSKZCAM workflow after dft_ensemble and skzcam analysis.

Parameters:

• skzcam_eint_analysis (dict[str, list[float]]) –

  The dictionary of the SKZCAM Eint analysis.

• dft_ensemble_analysis (dict[str, list[float]]) –

  The dictionary of the DFT ensemble analysis.

Returns:

• dict[str, list[float]] –

  The final Hads dictionary containing the Hads contributions from the DFT ensemble and the SKZCAM calculations.

Source code in autoSKZCAM/recipes.py

def get_final_autoSKZCAM_Hads(
    skzcam_eint_analysis: dict[str, list[float]],
    dft_ensemble_analysis: dict[str, list[float]],
) -> dict[str, list[float]]:
    """
    Gets the final Hads from the autoSKZCAM workflow after dft_ensemble and skzcam analysis.

    Parameters
    ----------
    skzcam_eint_analysis
        The dictionary of the SKZCAM Eint analysis.
    dft_ensemble_analysis
        The dictionary of the DFT ensemble analysis.

    Returns
    -------
    dict[str, list[float]]
        The final Hads dictionary containing the Hads contributions from the DFT ensemble and the SKZCAM calculations.
    """
    final_Hads = skzcam_eint_analysis.copy()

    final = skzcam_eint_analysis["Overall Eint"].copy()
    for key, contribution in dft_ensemble_analysis.items():
        final_Hads[key] = contribution
        final[0] += contribution[0]
        final[1] = np.sqrt(final[1] ** 2 + contribution[1] ** 2)

    final_Hads["Final Hads"] = final
    return final_Hads

skzcam_analyse

skzcam_analyse(calc_dir: str | Path, embedded_cluster_npy_path: Path | str | None = None, OniomInfo: dict[str, OniomLayerInfo] | None = None, EmbeddedCluster: CreateEmbeddedCluster | None = None) -> dict[str, tuple[float, float]]

Analyze the completed SKZCAM calculations and compute the final ONIOM contributions.

Parameters:

• calc_dir (str | Path) –

  The directory containing the calculations.

• embedded_cluster_npy_path (Path | str | None, default: None ) –

  The path to the embedded cluster .npy object.

• EmbeddedCluster (CreateEmbeddedCluster | None, default: None ) –

  The CreateEmbeddedCluster object containing the embedded cluster.

• OniomInfo (dict[str, OniomLayerInfo] | None, default: None ) –

  A dictionary containing the information about the ONIOM layers.

Returns:

• dict[str, tuple[float, float]] –

  A dictionary containing the ONIOM layer as key and a tuple containing the contribution to the final interaction energy as well as its error.

Source code in autoSKZCAM/recipes_skzcam.py

def skzcam_analyse(
    calc_dir: str | Path,
    embedded_cluster_npy_path: Path | str | None = None,
    OniomInfo: dict[str, OniomLayerInfo] | None = None,
    EmbeddedCluster: CreateEmbeddedCluster | None = None,
) -> dict[str, tuple[float, float]]:
    """
    Analyze the completed SKZCAM calculations and compute the final ONIOM contributions.

    Parameters
    ----------
    calc_dir
        The directory containing the calculations.
    embedded_cluster_npy_path
        The path to the embedded cluster .npy object.
    EmbeddedCluster
        The CreateEmbeddedCluster object containing the embedded cluster.
    OniomInfo
        A dictionary containing the information about the ONIOM layers.

    Returns
    -------
    dict[str, tuple[float, float]]
        A dictionary containing the ONIOM layer as key and a tuple containing the contribution to the final interaction energy as well as its error.
    """

    # Initialize the EmbeddedCluster object if it is not provided
    if EmbeddedCluster is None and embedded_cluster_npy_path is None:
        # Check if the embedded_cluster.npy file exists in the calc_dir
        if not Path(calc_dir, "embedded_cluster.npy").exists():
            raise ValueError(
                "Either the EmbeddedCluster object must be provided or embedded_cluster_npy_path is set or embedded_cluster.npy is provided in calc_dir."
            )
        EmbeddedCluster = np.load(
            Path(calc_dir, "embedded_cluster.npy"), allow_pickle=True
        ).item()

    elif EmbeddedCluster is None and embedded_cluster_npy_path is not None:
        EmbeddedCluster = np.load(embedded_cluster_npy_path, allow_pickle=True).item()

    # Load the OniomInfo dictionary if it is not provided
    if EmbeddedCluster.OniomInfo is None:
        if OniomInfo is None:
            raise ValueError(
                "The OniomInfo dictionary must be provided in EmbeddedCluster or as an argument."
            )
        EmbeddedCluster.OniomInfo = OniomInfo

    if EmbeddedCluster.OniomInfo is not None and OniomInfo is None:
        OniomInfo = EmbeddedCluster.OniomInfo

    skzcam_calcs_analysis = analyze_calculations(
        calc_dir=calc_dir,
        embedded_cluster_path=embedded_cluster_npy_path,
        EmbeddedCluster=EmbeddedCluster,
    )

    return compute_skzcam_int_ene(
        skzcam_calcs_analysis=skzcam_calcs_analysis, OniomInfo=OniomInfo
    )

skzcam_eint_flow

skzcam_eint_flow(EmbeddedCluster: CreateEmbeddedCluster, OniomInfo: dict[str, OniomLayerInfo], capped_ecp: dict[Literal['mrcc', 'orca'], str] | None = None, **kwargs) -> None

The complete SKZCAM protocol to generate the embedded clusters, perform the calculations, and analyze the results.

Parameters:

• EmbeddedCluster (CreateEmbeddedCluster) –

  The CreateEmbeddedCluster object containing the embedded cluster. This is initialised using the skzcam_initialise() function.

• OniomInfo (dict[str, OniomLayerInfo]) –

  A dictionary containing the information about the ONIOM layers.

• **kwargs –

  Additional keyword arguments to pass to the skzcam_generate_job() and skzcam_calculate_job() functions.

Returns:

• None –

Source code in autoSKZCAM/recipes_skzcam.py

def skzcam_eint_flow(
    EmbeddedCluster: CreateEmbeddedCluster,
    OniomInfo: dict[str, OniomLayerInfo],
    capped_ecp: dict[Literal["mrcc", "orca"], str] | None = None,
    **kwargs,
) -> None:
    """
    The complete SKZCAM protocol to generate the embedded clusters, perform the calculations, and analyze the results.

    Parameters
    ----------
    EmbeddedCluster
        The CreateEmbeddedCluster object containing the embedded cluster. This is initialised using the skzcam_initialise() function.
    OniomInfo
        A dictionary containing the information about the ONIOM layers.
    **kwargs
        Additional keyword arguments to pass to the skzcam_generate_job() and skzcam_calculate_job() functions.

    Returns
    -------
    None
    """

    # Generate the skzcam embedded clusters
    skzcam_generate_job(EmbeddedCluster, **kwargs)

    # Perform the calculations on the embedded clusters
    skzcam_calculate_job(EmbeddedCluster, OniomInfo, capped_ecp=capped_ecp, **kwargs)

skzcam_initialise

skzcam_initialise(adsorbate_indices: list[int], slab_center_indices: list[int], atom_oxi_states: dict[ElementStr, int], adsorbate_slab_file: str | Path, pun_filepath: str | Path = './ChemShell_EmbeddedCluster.pun', run_chemshell: bool = False, chemsh_radius_active: float = 40.0, chemsh_radius_cluster: float = 60.0, chemsh_bq_layer: float = 6.0, write_xyz_file: bool = False, **kwargs) -> CreateEmbeddedCluster

Parameters to initialise the SKZCAM protocol to generate the embedded clusters.

Parameters:

• adsorbate_indices (list[int]) –

  The indices of the atoms that make up the adsorbate molecule.

• slab_center_indices (list[int]) –

  The indices of the atoms that make up the 'center' of the slab right beneath the adsorbate.

• atom_oxi_states (dict[ElementStr, int]) –

  A dictionary with the element symbol as the key and its oxidation state as the value.

• adsorbate_slab_file (str | Path) –

  The path to the file containing the adsorbate molecule on the surface slab. It can be in any format that ASE can read.

• pun_filepath (str | Path, default: './ChemShell_EmbeddedCluster.pun' ) –

  The path to the .pun file containing the atomic coordinates and charges of the adsorbate-slab complex if it has already been generated by ChemShell. If it is None, then ChemShell will need to be used to create this file.

• run_chemshell (bool, default: False ) –

  If True, the ChemShell calculations will be run to generate the .pun file.

• chemsh_radius_active (float, default: 40.0 ) –

  The radius of the active region in Angstroms. This 'active' region is simply region where the charge fitting is performed to ensure correct Madelung potential; it can be a relatively large value.

• chemsh_radius_cluster (float, default: 60.0 ) –

  The radius of the total embedded cluster in Angstroms.

• chemsh_bq_layer (float, default: 6.0 ) –

  The height above the surface to place some additional fitting point charges in Angstroms; simply for better reproduction of the electrostatic potential close to the adsorbate.

• write_xyz_file (bool, default: False ) –

  If True, the .xyz file will be written containing the atomic coordinates of the adsorbate-slab complex.

Returns:

• CreateEmbeddedCluster –

  The CreateEmbeddedCluster object containing the embedded cluster.

Source code in autoSKZCAM/recipes_skzcam.py

def skzcam_initialise(
    adsorbate_indices: list[int],
    slab_center_indices: list[int],
    atom_oxi_states: dict[ElementStr, int],
    adsorbate_slab_file: str | Path,
    pun_filepath: str | Path = "./ChemShell_EmbeddedCluster.pun",
    run_chemshell: bool = False,
    chemsh_radius_active: float = 40.0,
    chemsh_radius_cluster: float = 60.0,
    chemsh_bq_layer: float = 6.0,
    write_xyz_file: bool = False,
    **kwargs,  # noqa ARG001
) -> CreateEmbeddedCluster:
    """
    Parameters to initialise the SKZCAM protocol to generate the embedded clusters.

    Parameters
    ----------
    adsorbate_indices
        The indices of the atoms that make up the adsorbate molecule.
    slab_center_indices
        The indices of the atoms that make up the 'center' of the slab right beneath the adsorbate.
    atom_oxi_states
        A dictionary with the element symbol as the key and its oxidation state as the value.
    adsorbate_slab_file
        The path to the file containing the adsorbate molecule on the surface slab. It can be in any format that ASE can read.
    pun_filepath
        The path to the .pun file containing the atomic coordinates and charges of the adsorbate-slab complex if it has already been generated by ChemShell. If it is None, then ChemShell will need to be used to create this file.
    run_chemshell
        If True, the ChemShell calculations will be run to generate the .pun file.
    chemsh_radius_active
        The radius of the active region in Angstroms. This 'active' region is simply region where the charge fitting is performed to ensure correct Madelung potential; it can be a relatively large value.
    chemsh_radius_cluster
        The radius of the total embedded cluster in Angstroms.
    chemsh_bq_layer
        The height above the surface to place some additional fitting point charges in Angstroms; simply for better reproduction of the electrostatic potential close to the adsorbate.
    write_xyz_file
        If True, the .xyz file will be written containing the atomic coordinates of the adsorbate-slab complex.

    Returns
    -------
    CreateEmbeddedCluster
        The CreateEmbeddedCluster object containing the embedded cluster.

    """

    EmbeddedCluster = CreateEmbeddedCluster(
        adsorbate_indices=adsorbate_indices,
        slab_center_indices=slab_center_indices,
        atom_oxi_states=atom_oxi_states,
        adsorbate_slab_file=adsorbate_slab_file,
        pun_filepath=pun_filepath,
    )

    # Check that pun_filepath exists if run_chemshell is False
    if run_chemshell is False and not Path(pun_filepath).exists():
        raise ValueError(
            "The path to the .pun file from ChemShell must be provided in EmbeddedCluster if run_chemshell is False."
        )

    if run_chemshell:
        # Create the ChemShell input file
        EmbeddedCluster.run_chemshell(
            filepath=pun_filepath,
            chemsh_radius_active=chemsh_radius_active,
            chemsh_radius_cluster=chemsh_radius_cluster,
            chemsh_bq_layer=chemsh_bq_layer,
            write_xyz_file=write_xyz_file,
        )

    return EmbeddedCluster