Cage builders¶

CageBuilder is the base class every topology-specific builder inherits from; it holds the full build, optimise, enumerate, functionalise, save, and statistics API. Topology-specific classes (Tri4Di6CageBuilder, M6L12CubeCageBuilder, ...) only override the topology graph and the point-group transformation library - see Topologies for the full list.

CageBuilder ¶

CageBuilder(linker: MoleculeInput = None, node: MoleculeInput = None, linker_fg=None, node_fg=None, node_fg_count: int = None, smarts_pattern: str = None, optimizer_step: float = 0.01, distance_threshold: float = 5.0, expansion_factor: float = 1.0, scale_multiplier: float = 1.0, linker_smiles: str = None, node_smiles: str = None)

Bases: ABC

Abstract base class for building cage structures with isomer enumeration.

Each subclass hard-codes one STK topology graph and its symmetry transformation library. The shared workflow (build → functionalise → enumerate → statistics → SMILES) lives here.

Typical usage

cage = Tet2Di4CageBuilder(linker_smiles='Nc1ccc(N)cc1') cage.build() cage.functionalise() isomers = cage.enumerate_isomers(output_path='Tet2Di4_isomers') stats = cage.get_statistics() smi = cage.to_smiles()

For cages without an STK builder (e.g. Tet6Di12):

cage = Tet6Di12CageBuilder.load_from_file('base_structure/Tet6Di12.pdb') cage.functionalise()

Initialize a cage builder.

The linker and node may be provided as a SMILES string, an ASE Atoms object, or a molecular structure file readable by ASE (e.g. .xyz, .mol, .mol2, .sdf, .pdb or .cif).

Functional groups are assigned automatically using the following priority:

1. **Single atom** (nodes only)
If the node represents a single atom (e.g. ``"[Zn]"`` or an
ASE object containing one atom), an ``stk.SingleAtom``
building block is created automatically.

2. **Dummy atoms**
If dummy atoms (``[*]`` or ``[#0]``) are present, the SMARTS
factory

    ``[#0][*]``

is used automatically.

3. **User-specified factory**
If ``linker_fg`` or ``node_fg`` is supplied, that factory is
used.

4. **Automatic factory detection**
If no factory is supplied, all factories in
``FG_FACTORY_MAP`` are searched until a compatible functional
group is found.

Parameters:

Name	Type	Description	Default
`linker`	`MoleculeInput`	Linker building block. May be one of a SMILES string, an ASE `Atoms` object, or a molecular structure file readable by ASE. Defaults to `default_linker`.	`None`
`node`	`MoleculeInput`	Cage node. May be one of a SMILES string, an ASE `Atoms` object, or a molecular structure file readable by ASE. Defaults to `default_node`.	`None`
`linker_fg`	`str or FunctionalGroupFactory`	Functional group factory assigned to the linker. May be a key in `FG_FACTORY_MAP`, an `stk.FunctionalGroupFactory` instance, or `None` to enable automatic detection. Defaults to `default_linker_fg`.	`None`
`node_fg`	`str or FunctionalGroupFactory`	Functional group factory assigned to the node. Use `"single_atom"` for isolated metal nodes or `None` for automatic detection. Defaults to `default_node_fg` only when the default node is used. User-supplied nodes default to automatic detection.	`None`
`node_fg_count`	`int`	Number of functional groups assigned to a single-atom node. Ignored for multi-atom nodes. Defaults to `default_node_fg_count`.	`None`
`smarts_pattern`	`str`	SMARTS pattern used when `linker_fg` or `node_fg` is `"smarts"`.	`None`
`optimizer_step`	`float`	Step size used by the STK `Collapser` optimizer.	`0.01`
`distance_threshold`	`float`	Distance threshold (Å) supplied to the STK `Collapser` optimizer.	`5.0`
`scale_multiplier`	`float`	Forwarded to the STK topology graph's own `scale_multiplier`. STK sizes the entire topology graph - node-to-node distances included - by multiplying its unit vertex positions by `scale_multiplier * max(bb.get_maximum_diameter() for bb in (node_bb, linker_bb))`. For a bulky multi-atom node next to a much shorter linker (e.g. a UiO-66 Zr6 SBU with a short dicarboxylate), the node's diameter alone sets the scale for every edge, including ones spanned by the short linker - so the initial placement comes out stretched to several times the linker's natural length, forcing :meth:`optimise` into a drastic global contraction instead of a local bond-length correction. Lowering `scale_multiplier` shrinks the initial placement before construction, at the source, rather than compensating for it after the fact via `expansion_factor`. Rule of thumb: build once at the default 1.0, note the resulting node-to-node distance `d0` (e.g. via each node subunit's centroid), and the distance `d_target` you actually want (a reasonable estimate is the linker's own end-to-end length plus twice a typical M-L bond length); then set `scale_multiplier = d_target / d0`.	`1.0`
`expansion_factor`	`float`	After STK assembly, each building-block subunit is translated rigidly so that its distance from the cage centroid is multiplied by this factor. Internal geometries are preserved exactly. Use values slightly above 1.0 when bulky multi-atom nodes (e.g. Zr6 clusters with CH3-capped peripheral connections) physically overlap their neighbours in the initial STK placement. A value too large will leave a visible gap between the node and linker bond sites; use the smallest value that eliminates the clash. Rule of thumb: if the node's outer radius is `r_max` and the initial vertex-to-vertex distance is `d`, the minimum factor is `2r_max / d`. For `uio66_tri_node.xyz` with CH3 caps at 5.97 Å and d ≈ 11.66 Å, a factor of 1.15* works well.	`1.0`
`linker_smiles`	`str`	Deprecated alias for `linker`. Included for backwards compatibility with earlier releases.	`None`
`node_smiles`	`str`	Deprecated alias for `node`. Included for backwards compatibility with earlier releases.	`None`

Source code in cage_isomer_builder/cage/base.py

def __init__(
    self,
    linker: MoleculeInput = None,
    node: MoleculeInput = None,
    linker_fg=None,
    node_fg=None,
    node_fg_count: int = None,
    smarts_pattern: str = None,
    optimizer_step: float = 0.01,
    distance_threshold: float = 5.0,
    expansion_factor: float = 1.0,
    scale_multiplier: float = 1.0,
    linker_smiles: str = None,
    node_smiles: str = None,
):
    """
    Initialize a cage builder.

    The linker and node may be provided as a SMILES string, an ASE
    ``Atoms`` object, or a molecular structure file readable by ASE
    (e.g. ``.xyz``, ``.mol``, ``.mol2``, ``.sdf``, ``.pdb`` or ``.cif``).

    Functional groups are assigned automatically using the following
    priority:

        1. **Single atom** (nodes only)
        If the node represents a single atom (e.g. ``"[Zn]"`` or an
        ASE object containing one atom), an ``stk.SingleAtom``
        building block is created automatically.

        2. **Dummy atoms**
        If dummy atoms (``[*]`` or ``[#0]``) are present, the SMARTS
        factory

            ``[#0][*]``

        is used automatically.

        3. **User-specified factory**
        If ``linker_fg`` or ``node_fg`` is supplied, that factory is
        used.

        4. **Automatic factory detection**
        If no factory is supplied, all factories in
        ``FG_FACTORY_MAP`` are searched until a compatible functional
        group is found.

    Parameters
    ----------
    linker : MoleculeInput, optional
        Linker building block.

        May be one of

        * a SMILES string,
        * an ASE ``Atoms`` object, or
        * a molecular structure file readable by ASE.

        Defaults to ``default_linker``.

    node : MoleculeInput, optional
        Cage node.

        May be one of

        * a SMILES string,
        * an ASE ``Atoms`` object, or
        * a molecular structure file readable by ASE.

        Defaults to ``default_node``.

    linker_fg : str or stk.FunctionalGroupFactory, optional
        Functional group factory assigned to the linker.

        May be

        * a key in ``FG_FACTORY_MAP``,
        * an ``stk.FunctionalGroupFactory`` instance, or
        * ``None`` to enable automatic detection.

        Defaults to ``default_linker_fg``.

    node_fg : str or stk.FunctionalGroupFactory, optional
        Functional group factory assigned to the node.

        Use ``"single_atom"`` for isolated metal nodes or ``None`` for
        automatic detection.

        Defaults to ``default_node_fg`` only when the default node is
        used. User-supplied nodes default to automatic detection.

    node_fg_count : int, optional
        Number of functional groups assigned to a single-atom node.

        Ignored for multi-atom nodes.

        Defaults to ``default_node_fg_count``.

    smarts_pattern : str, optional
        SMARTS pattern used when ``linker_fg`` or ``node_fg`` is
        ``"smarts"``.

    optimizer_step : float, default=0.01
        Step size used by the STK ``Collapser`` optimizer.

    distance_threshold : float, default=5.0
        Distance threshold (Å) supplied to the STK ``Collapser``
        optimizer.

    scale_multiplier : float, default=1.0
        Forwarded to the STK topology graph's own ``scale_multiplier``.
        STK sizes the *entire* topology graph - node-to-node distances
        included - by multiplying its unit vertex positions by
        ``scale_multiplier * max(bb.get_maximum_diameter() for bb in
        (node_bb, linker_bb))``. For a bulky multi-atom node next to a
        much shorter linker (e.g. a UiO-66 Zr6 SBU with a short
        dicarboxylate), the *node's* diameter alone sets the scale for
        every edge, including ones spanned by the short linker - so the
        initial placement comes out stretched to several times the
        linker's natural length, forcing :meth:`optimise` into a drastic
        global contraction instead of a local bond-length correction.
        Lowering ``scale_multiplier`` shrinks the initial placement
        before construction, at the source, rather than compensating for
        it after the fact via ``expansion_factor``.

        Rule of thumb: build once at the default 1.0, note the resulting
        node-to-node distance ``d0`` (e.g. via each node subunit's
        centroid), and the distance ``d_target`` you actually want (a
        reasonable estimate is the linker's own end-to-end length plus
        twice a typical M-L bond length); then set
        ``scale_multiplier = d_target / d0``.

    expansion_factor : float, default=1.0
        After STK assembly, each building-block subunit is translated
        rigidly so that its distance from the cage centroid is multiplied
        by this factor. Internal geometries are preserved exactly.

        Use values slightly above 1.0 when bulky multi-atom nodes (e.g.
        Zr6 clusters with CH3-capped peripheral connections) physically
        overlap their neighbours in the initial STK placement.  A value
        too large will leave a visible gap between the node and linker
        bond sites; use the smallest value that eliminates the clash.

        Rule of thumb: if the node's outer radius is ``r_max`` and the
        initial vertex-to-vertex distance is ``d``, the minimum factor is
        ``2*r_max / d``.  For ``uio66_tri_node.xyz`` with CH3 caps at
        5.97 Å and d ≈ 11.66 Å, a factor of **1.15** works well.

    linker_smiles : str, optional
        Deprecated alias for ``linker``.

        Included for backwards compatibility with earlier releases.

    node_smiles : str, optional
        Deprecated alias for ``node``.

        Included for backwards compatibility with earlier releases.
    """

    # ------------------------------------------------------------------
    # Input normalization
    # ------------------------------------------------------------------
    linker_input = linker if linker is not None else linker_smiles
    node_input = node if node is not None else node_smiles

    self.linker = normalize_molecule_input(
        mol=linker_input,
        default=self.default_linker,
    )

    self.node = normalize_molecule_input(
        mol=node_input,
        default=self.default_node,
    )

    self.linker_fg = (
        linker_fg
        if linker_fg is not None
        else self.default_linker_fg
    )

    self.node_fg = (
        node_fg
        if node_fg is not None
        else self.default_node_fg
        if node_input is None
        else None
    )

    self.node_fg_count = (
        node_fg_count
        if node_fg_count is not None
        else self.default_node_fg_count
    )

    self.smarts_pattern = smarts_pattern
    self.optimizer_step = optimizer_step
    self.distance_threshold = distance_threshold
    self.expansion_factor = expansion_factor
    self.scale_multiplier = scale_multiplier

    # Internal state
    self._stk_cage = None
    self._ase_cage = None
    self._fg_anchors = None
    self._fg_anchor_indices = None
    self._fg_anchor_h_positions = None
    self._fg_per_linker = None
    self._linker_atom_indices = None
    self._node_centers = None

build ¶

build() -> stk.ConstructedMolecule

Build the cage using the configured node and linker.

Returns:

Type	Description
`ConstructedMolecule`	Constructed STK cage.

Notes

Building blocks are generated automatically from the supplied node and linker inputs before constructing the topology graph. If expansion_factor > 1, each subunit is translated rigidly outward after assembly (see :meth:_expand_cage_rigid).

Source code in cage_isomer_builder/cage/base.py

def build(self) -> stk.ConstructedMolecule:
    """
    Build the cage using the configured node and linker.

    Returns
    -------
    stk.ConstructedMolecule
        Constructed STK cage.

    Notes
    -----
    Building blocks are generated automatically from the supplied node
    and linker inputs before constructing the topology graph.
    If ``expansion_factor > 1``, each subunit is translated rigidly
    outward after assembly (see :meth:`_expand_cage_rigid`).
    """
    node_bb = self._make_node_bb()
    linker_bb = self._make_linker_bb()

    self._stk_cage = stk.ConstructedMolecule(
        topology_graph=self._build_stk_topology(
            node_bb=node_bb,
            linker_bb=linker_bb,
        )
    )

    if self.expansion_factor > 1.0:
        self._expand_cage_rigid()

    return self._stk_cage

count_unique_isomers_burnside ¶

count_unique_isomers_burnside() -> int

Count unique isomers via Burnside's lemma without full enumeration. Useful for large cages where iterating over 4^n_linkers is impractical.

Returns:

Name	Type	Description
`n_unique`	`int`

Source code in cage_isomer_builder/cage/base.py

def count_unique_isomers_burnside(self) -> int:
    """
    Count unique isomers via Burnside's lemma without full enumeration.
    Useful for large cages where iterating over 4^n_linkers is impractical.

    Returns
    -------
    n_unique : int
    """
    if self._fg_anchors is None:
        self.functionalise()
    n_fg_slots = len(self._fg_anchors)
    fg_per_linker = self._fg_per_linker if self._fg_per_linker is not None else 4
    n_linkers = n_fg_slots // fg_per_linker
    tl = self._get_transformation_library(self._fg_anchors.copy())
    raw_tl = list(zip(*tl))
    return count_unique_isomers(raw_tl, n_linkers, n_fg_slots, fg_per_linker)

enumerate_isomers ¶

enumerate_isomers(output_path: str = None, generate_files: bool = True, limit: int = None) -> list

Enumerate symmetry-unique isomers using the transformation library.

Walks candidate isomers with iter_unique_isomers, which checks each candidate directly against its own symmetry orbit instead of allocating a fg_per_linker**n_linkers flag array - the array is what used to raise MemoryError on a periodic structure with more than a couple dozen linkers (see docs/isomers.md). Memory use here is O(1) in n_linkers regardless of limit; runtime is unchanged (still one orbit check per raw combination visited) unless you pass limit, which stops early instead of walking the full space.

Parameters:

Name	Type	Description	Default
`output_path`	`str`	Directory to write XYZ structure files. No files are written if None or if generate_files is False.	`None`
`generate_files`	`bool`	Write one XYZ file per unique isomer (default True).	`True`
`limit`	`int`	Stop after this many unique isomers instead of walking the entire fg_per_linkern_linkers space. Use this for cages/MOFs where the full count (from count_unique_isomers_burnside) is too large to fully enumerate - canonical representatives occur roughly 1-in-\|symmetry group order\| along the walk, not 1-in-fg_per_linkern_linkers, so a modest limit is normally fast even when the full space isn't.	`None`

Returns:

Name	Type	Description
`unique_isomers`	`list of list of int`	Each element is a base-4 isomer descriptor.

Source code in cage_isomer_builder/cage/base.py

def enumerate_isomers(
    self,
    output_path: str = None,
    generate_files: bool = True,
    limit: int = None,
) -> list:
    """
    Enumerate symmetry-unique isomers using the transformation library.

    Walks candidate isomers with iter_unique_isomers, which checks each
    candidate directly against its own symmetry orbit instead of
    allocating a fg_per_linker**n_linkers flag array - the array is
    what used to raise MemoryError on a periodic structure with more
    than a couple dozen linkers (see docs/isomers.md). Memory use here
    is O(1) in n_linkers regardless of limit; runtime is unchanged
    (still one orbit check per raw combination visited) unless you
    pass limit, which stops early instead of walking the full space.

    Parameters
    ----------
    output_path : str, optional
        Directory to write XYZ structure files. No files are written
        if None or if generate_files is False.
    generate_files : bool
        Write one XYZ file per unique isomer (default True).
    limit : int, optional
        Stop after this many unique isomers instead of walking the
        entire fg_per_linker**n_linkers space. Use this for cages/MOFs
        where the full count (from count_unique_isomers_burnside) is
        too large to fully enumerate - canonical representatives occur
        roughly 1-in-|symmetry group order| along the walk, not
        1-in-fg_per_linker**n_linkers, so a modest limit is normally
        fast even when the full space isn't.

    Returns
    -------
    unique_isomers : list of list of int
        Each element is a base-4 isomer descriptor.
    """
    if self._fg_anchors is None:
        self.functionalise()

    n_fg_slots = len(self._fg_anchors)
    fg_per_linker = self._fg_per_linker if self._fg_per_linker is not None else 4
    n_linkers = n_fg_slots // fg_per_linker
    tl = self._get_transformation_library(self._fg_anchors.copy())

    if generate_files and output_path:
        os.makedirs(output_path, exist_ok=True)

    unique_isomers = []
    for _, descriptor in iter_unique_isomers(tl, n_linkers, fg_per_linker, limit=limit):
        unique_isomers.append(descriptor)

        if generate_files and output_path:
            generate_isomer_structure_file(
                self._ase_cage,
                self._fg_anchors,
                descriptor,
                self._fg_anchor_indices,
                output_path,
                fg_anchor_h_positions=self._fg_anchor_h_positions,
            )

    return unique_isomers

functionalise ¶

functionalise()

Functionalise the cage and identify functional-group anchor sites.

Returns:

Name	Type	Description
`fg_anchors`	`Atoms`	Functionalised ASE structure containing the generated anchor atoms.
`fg_anchor_indices`	`list[int]`	Indices of the anchor atoms used for isomer enumeration.

Source code in cage_isomer_builder/cage/base.py

def functionalise(self):
    """
    Functionalise the cage and identify functional-group anchor sites.

    Returns
    -------
    fg_anchors : ase.Atoms
        Functionalised ASE structure containing the generated anchor atoms.

    fg_anchor_indices : list[int]
        Indices of the anchor atoms used for isomer enumeration.
    """
    # Work on a copy so the raw (pre-At) cage is preserved if needed,
    # then store the At-decorated version as the canonical _ase_cage so
    # enumerate_isomers and generate_isomer_structure_file see the markers.
    functionalised = self.to_ase().copy()

    # When the STK cage is available, use its atom-info tracking to
    # identify subunits.  This is robust even when large building blocks
    # physically overlap (which breaks mofstructure's connectivity analysis
    # and the BFS node-center finder).
    if self._stk_cage is not None:
        linker_atom_indices, metal_atom_indices = self._get_sbus_from_stk(functionalised)
        self._node_centers = self._get_node_centers_from_stk(functionalised)
    else:
        linker_atom_indices, metal_atom_indices = get_sbus_indices(functionalised)
        self._node_centers = get_node_centers(functionalised, metal_atom_indices)

    self._linker_atom_indices = linker_atom_indices
    self._fg_anchors, self._fg_anchor_indices, self._fg_anchor_h_positions = full_functionalisation(
        functionalised, sum(linker_atom_indices, [])
    )
    self._fg_per_linker = len(self._fg_anchor_indices) // len(linker_atom_indices)
    self._ase_cage = functionalised

    return self._fg_anchors, self._fg_anchor_indices

get_bond_matrix ¶

get_bond_matrix() -> np.ndarray

This cage's bond-order matrix, taken directly from the STK construction rather than re-perceived from geometry - the same matrix write_run/write_gulp_gin use, for the same reason: a real bond order survives even where two bonded atoms sit unusually far apart, which a distance-based guess on the whole structure can miss (or, on a larger assembly like a host-guest complex, wrongly add between atoms that are merely close together but not bonded).

Needed by anything building on top of this cage's own bonding without re-perceiving it - e.g. utils.functionalise.place_guest_in_host's host_bond_matrix parameter, which keeps the host's real bonds and the union of each guest's own (separately, safely perceived) internal bonds, rather than running one neighbour-list pass over the merged host+guest structure.

Returns:

Type	Description
`(ndarray, shape(n_atoms, n_atoms))`

Raises:

Type	Description
`RuntimeError`	If build() hasn't been called - a structure loaded via load_from_file() has no STK-derived bond order information.

Source code in cage_isomer_builder/cage/base.py

def get_bond_matrix(self) -> np.ndarray:
    """
    This cage's bond-order matrix, taken directly from the STK
    construction rather than re-perceived from geometry - the same
    matrix write_run/write_gulp_gin use, for the same reason: a real
    bond order survives even where two bonded atoms sit unusually far
    apart, which a distance-based guess on the whole structure can
    miss (or, on a larger assembly like a host-guest complex, wrongly
    add between atoms that are merely close together but not bonded).

    Needed by anything building on top of this cage's own bonding
    without re-perceiving it - e.g.
    utils.functionalise.place_guest_in_host's host_bond_matrix
    parameter, which keeps the host's real bonds and the union of
    each guest's own (separately, safely perceived) internal bonds,
    rather than running one neighbour-list pass over the merged
    host+guest structure.

    Returns
    -------
    np.ndarray, shape (n_atoms, n_atoms)

    Raises
    ------
    RuntimeError
        If build() hasn't been called - a structure loaded via
        load_from_file() has no STK-derived bond order information.
    """
    if self._stk_cage is None:
        raise RuntimeError(
            "get_bond_matrix() requires the STK cage. Call build() "
            "first - a structure loaded via load_from_file() has no "
            "STK-derived bond order information."
        )
    _, bond_matrix = stk_2_ase_atoms(self._stk_cage)
    return bond_matrix

get_statistics ¶

get_statistics(n_kde_bins: int = 200) -> CageStatistics

Compute FG-to-FG distance statistics and a Gaussian KDE for this cage.

Parameters:

Name	Type	Description	Default
`n_kde_bins`	`int`	Number of evaluation points for the KDE curve (default 200).	`200`

Returns:

Name	Type	Description
`stats`	`CageStatistics`

Source code in cage_isomer_builder/cage/base.py

def get_statistics(self, n_kde_bins: int = 200) -> CageStatistics:
    """
    Compute FG-to-FG distance statistics and a Gaussian KDE for this cage.

    Parameters
    ----------
    n_kde_bins : int
        Number of evaluation points for the KDE curve (default 200).

    Returns
    -------
    stats : CageStatistics
    """
    if self._fg_anchors is None:
        self.functionalise()

    n_fg_slots = len(self._fg_anchors)
    fg_per_linker = self._fg_per_linker if self._fg_per_linker is not None else 4
    n_linkers = n_fg_slots // fg_per_linker
    n_unique = self.count_unique_isomers_burnside()

    distance_keys, total, endo_endo, endo_exo, exo_exo = fg2fg_distance_count(
        self._fg_anchors, self._fg_anchor_indices, fg_per_linker
    )
    kde_x, kde_y = kernel_density_estimation(distance_keys, total, n_kde_bins)

    return CageStatistics(
        n_linkers=n_linkers,
        n_fg_slots=n_fg_slots,
        n_unique_isomers=n_unique,
        distance_keys=distance_keys,
        total_fg_pairs=total,
        endo_endo_pairs=endo_endo,
        endo_exo_pairs=endo_exo,
        exo_exo_pairs=exo_exo,
        kde_x=kde_x,
        kde_y=kde_y,
    )

load_from_file `classmethod` ¶

load_from_file(path: str, **init_kwargs)

Create a builder from an existing molecular structure.

This method bypasses STK construction and instead loads an existing cage structure directly into the builder. The loaded structure may subsequently be functionalised, analysed, or enumerated.

Parameters:

Name	Type	Description	Default
`path`	`str or Path`	Path to a molecular structure readable by ASE (e.g. PDB, CIF, XYZ, MOL or SDF).	required
`**init_kwargs`		Additional keyword arguments forwarded to the constructor.	`{}`

Returns:

Type	Description
`CageBuilder`	Builder initialized from the supplied structure.

Source code in cage_isomer_builder/cage/base.py

@classmethod
def load_from_file(cls, path: str, **init_kwargs):
    """
    Create a builder from an existing molecular structure.

    This method bypasses STK construction and instead loads an existing
    cage structure directly into the builder. The loaded structure may
    subsequently be functionalised, analysed, or enumerated.

    Parameters
    ----------
    path : str or pathlib.Path
        Path to a molecular structure readable by ASE (e.g. PDB, CIF,
        XYZ, MOL or SDF).

    **init_kwargs
        Additional keyword arguments forwarded to the constructor.

    Returns
    -------
    CageBuilder
        Builder initialized from the supplied structure.
    """
    instance = cls.__new__(cls)
    CageBuilder.__init__(instance, **init_kwargs)
    instance._ase_cage = read(path)
    return instance

optimise ¶

optimise(fmax: float = 0.05, steps: int = 500, optimizer_cls=None, trajectory: Optional[str] = None, rigid_sbus: bool = True, include_charges: bool = True, charge: float = 0.0, method: str = 'UFF4MOF') -> Atoms

Run ASE geometry optimisation on the cage.

Two backends are available, selected via method:

"UFF4MOF" (default) - a classical force field evaluated from the STK bond matrix and explicit UFF atom types (see :mod:cage_isomer_builder.utils.optimise). Bond connectivity is taken from the STK-constructed cage and respected exactly, which matters when inter-atomic distances in the initial geometry are atypical (bulky multi-atom nodes, an expanded cage, ...) and would confuse a distance-only bond perception.
"GFN1-xTB" / "GFN2-xTB" - a tblite GFN-xTB electronic structure calculation. Unlike UFF4MOF, xTB needs no explicit bond-order matrix or atom typing - it derives its own electronic structure directly from atomic numbers and positions - but for the same reason it cannot be told to honour a specific bond matrix.

Call this after :meth:build and before :meth:functionalise. If called after :meth:functionalise, the At marker atoms are stripped before optimisation (they are not real chemistry) and :meth:functionalise should be re-run afterwards.

Parameters:

Name	Type	Description	Default
`fmax`	`float`	Force convergence threshold in eV/Å.	`0.05`
`steps`	`int`	Maximum number of optimisation steps.	`500`
`optimizer_cls`	`ASE optimizer class`	Any `ase.optimize` class (`LBFGS`, `BFGS`, `FIRE`, …). Defaults to `LBFGS`.	`None`
`trajectory`	`str`	Path for an ASE `.traj` trajectory file.	`None`
`rigid_sbus`	`bool`	Treat every node/linker SBU as a rigid body (see :class:`~cage_isomer_builder.utils.rigid_body.FixRigidBodies`), relaxing only inter-SBU placement. Metal nodes are kept fully rigid (their own cluster bonding is what full relaxation tends to get wrong), but each linker's donor atom(s) - whichever atom carries the bond crossing into a node - are excluded from the linker's rigid group. A node's binding vectors are fixed by its rigid shape, and a rigid linker's 6 whole-body DOF generally can't satisfy every M-L bond/angle around a node at once, which is what bends the M-L bond; freeing just the donor atom gives that bond the local freedom it needs without touching the rest of the linker (ring shape, aromaticity, ...). Set to False to relax every atom independently instead.	`False`
`include_charges`	`bool`	Store per-atom partial charges on the returned `Atoms` (via `set_initial_charges`). For `method="GFN1-xTB"`/`"GFN2-xTB"` these are a byproduct of the optimisation's own electronic structure. UFF4MOF itself has no notion of partial charges, so for `method="UFF4MOF"` a tblite GFN2-xTB single-point calculation is run on the pre-optimisation geometry to obtain them; those same fixed charges are then also fed into :class:`~cage_isomer_builder.utils.optimise.UFF4MOFCalculator` to add its Coulomb term (`keqiqj/r` between the same non-bonded pairs the van der Waals term already covers), rather than optimising with electrostatics entirely absent. Setting this to `False` skips both - no charges are computed or reported, and `method="UFF4MOF"` optimises without a Coulomb term at all.	`True`
`charge`	`float`	Total charge of the structure, passed to tblite - either as the xTB optimisation's own charge, or (for `method="UFF4MOF"`) as the charge of the pre-optimisation single-point calculation used to seed the Coulomb term.	`0.0`
`method`	`str`	`"UFF4MOF"`, or a tblite xTB Hamiltonian name such as `"GFN1-xTB"` / `"GFN2-xTB"`.	`"UFF4MOF"`

Returns:

Type	Description
`Atoms`	Optimised cage structure. `self._ase_cage` is updated in-place.

Source code in cage_isomer_builder/cage/base.py

def optimise(
    self,
    fmax: float = 0.05,
    steps: int = 500,
    optimizer_cls=None,
    trajectory: Optional[str] = None,
    rigid_sbus: bool = True,
    include_charges: bool = True,
    charge: float = 0.0,
    method: str = "UFF4MOF",
) -> Atoms:
    """
    Run ASE geometry optimisation on the cage.

    Two backends are available, selected via ``method``:

    * ``"UFF4MOF"`` (default) - a classical force field evaluated from
      the STK bond matrix and explicit UFF atom types (see
      :mod:`cage_isomer_builder.utils.optimise`). Bond connectivity is
      taken from the STK-constructed cage and respected exactly, which
      matters when inter-atomic distances in the initial geometry are
      atypical (bulky multi-atom nodes, an expanded cage, ...) and would
      confuse a distance-only bond perception.
    * ``"GFN1-xTB"`` / ``"GFN2-xTB"`` - a tblite GFN-xTB electronic
      structure calculation. Unlike UFF4MOF, xTB needs no explicit
      bond-order matrix or atom typing - it derives its own electronic
      structure directly from atomic numbers and positions - but for the
      same reason it cannot be told to honour a specific bond matrix.

    Call this after :meth:`build` and before :meth:`functionalise`.
    If called after :meth:`functionalise`, the At marker atoms are
    stripped before optimisation (they are not real chemistry) and
    :meth:`functionalise` should be re-run afterwards.

    Parameters
    ----------
    fmax : float, default 0.05
        Force convergence threshold in eV/Å.
    steps : int, default 500
        Maximum number of optimisation steps.
    optimizer_cls : ASE optimizer class, optional
        Any ``ase.optimize`` class (``LBFGS``, ``BFGS``, ``FIRE``, …).
        Defaults to ``LBFGS``.
    trajectory : str, optional
        Path for an ASE ``.traj`` trajectory file.
    rigid_sbus : bool, default False
        Treat every node/linker SBU as a rigid body (see
        :class:`~cage_isomer_builder.utils.rigid_body.FixRigidBodies`),
        relaxing only inter-SBU placement. Metal nodes are kept fully
        rigid (their own cluster bonding is what full relaxation tends
        to get wrong), but each linker's donor atom(s) - whichever atom
        carries the bond crossing into a node - are excluded from the
        linker's rigid group. A node's binding vectors are fixed by its
        rigid shape, and a rigid linker's 6 whole-body DOF generally
        can't satisfy every M-L bond/angle around a node at once, which
        is what bends the M-L bond; freeing just the donor atom gives
        that bond the local freedom it needs without touching the rest
        of the linker (ring shape, aromaticity, ...). Set to False to
        relax every atom independently instead.
    include_charges : bool, default True
        Store per-atom partial charges on the returned ``Atoms`` (via
        ``set_initial_charges``). For ``method="GFN1-xTB"``/``"GFN2-xTB"``
        these are a byproduct of the optimisation's own electronic
        structure. UFF4MOF itself has no notion of partial charges, so
        for ``method="UFF4MOF"`` a tblite GFN2-xTB single-point
        calculation is run on the *pre-optimisation* geometry to obtain
        them; those same fixed charges are then also fed into
        :class:`~cage_isomer_builder.utils.optimise.UFF4MOFCalculator`
        to add its Coulomb term (``ke*qi*qj/r`` between the same
        non-bonded pairs the van der Waals term already covers), rather
        than optimising with electrostatics entirely absent. Setting
        this to ``False`` skips both - no charges are computed or
        reported, and ``method="UFF4MOF"`` optimises without a Coulomb
        term at all.
    charge : float, default 0.0
        Total charge of the structure, passed to tblite - either as the
        xTB optimisation's own charge, or (for ``method="UFF4MOF"``) as
        the charge of the pre-optimisation single-point calculation used
        to seed the Coulomb term.
    method : str, default "UFF4MOF"
        ``"UFF4MOF"``, or a tblite xTB Hamiltonian name such as
        ``"GFN1-xTB"`` / ``"GFN2-xTB"``.

    Returns
    -------
    ase.Atoms
        Optimised cage structure. ``self._ase_cage`` is updated in-place.
    """
    from cage_isomer_builder.utils.rigid_body import FixRigidBodies

    if self._stk_cage is not None:
        # Always derive atoms from the original STK cage so that At
        # marker atoms from functionalise() are never included.
        ase_atoms, bond_matrix = stk_2_ase_atoms(self._stk_cage)
    else:
        ase_atoms = self.to_ase().copy()
        # Strip At atoms (FG markers) — not parameterised by any FF.
        non_at = [i for i, a in enumerate(ase_atoms) if a.symbol != 'At']
        ase_atoms = ase_atoms[non_at]
        bond_matrix = None

    if rigid_sbus:
        if self._stk_cage is None:
            raise RuntimeError(
                "rigid_sbus=True requires the STK cage to identify SBU "
                "groupings (unavailable after load_from_file()). Pass "
                "rigid_sbus=False to relax every atom independently instead."
            )
        bb_atoms, bb_is_node = self._classify_atoms_by_bb(ase_atoms)

        # Nodes stay fully rigid. Linkers drop their donor atom(s) - any
        # atom with a bond crossing into a node - out of the rigid group,
        # so each M-L bond keeps the local freedom a whole-rigid-body fit
        # can't give it (see the rigid_sbus docstring above).
        node_bb_ids = {bb_id for bb_id, is_node in bb_is_node.items() if is_node}
        node_atom_ids = {idx for bb_id in node_bb_ids for idx in bb_atoms[bb_id]}
        donor_atom_ids = {
            i
            for i in range(len(ase_atoms))
            if i not in node_atom_ids
            for j in np.nonzero(bond_matrix[i])[0]
            if j in node_atom_ids
        }

        groups = [bb_atoms[bb_id] for bb_id in node_bb_ids] + [
            [idx for idx in atoms if idx not in donor_atom_ids]
            for bb_id, atoms in bb_atoms.items()
            if bb_id not in node_bb_ids
        ]
        groups = [g for g in groups if len(g) > 1]
        ase_atoms.set_constraint(FixRigidBodies.from_atoms(ase_atoms, groups))

    if method.upper() == "UFF4MOF":
        if bond_matrix is None:
            raise RuntimeError(
                "method='UFF4MOF' requires the STK bond matrix (unavailable "
                "after load_from_file()). Pass method='GFN2-xTB' instead."
            )
        from cage_isomer_builder.utils.optimise import run_uff4mof_optimisation

        charges = None
        if include_charges:
            from cage_isomer_builder.utils.tblite_charges import predict_tblite_charges

            # Computed once, up front, from the pre-optimisation geometry,
            # and held fixed for the whole run (standard UFF/QEq
            # practice) so they can also drive UFF4MOF's Coulomb term -
            # which is otherwise entirely absent from the force field.
            charges = predict_tblite_charges(ase_atoms, charge=charge)

        optimised = run_uff4mof_optimisation(
            atoms=ase_atoms,
            bond_matrix=bond_matrix,
            fmax=fmax,
            steps=steps,
            optimizer_cls=optimizer_cls,
            trajectory=trajectory,
            charges=charges,
        )
    else:
        from cage_isomer_builder.utils.tblite_optimise import run_tblite_optimisation

        optimised, charges = run_tblite_optimisation(
            atoms=ase_atoms,
            fmax=fmax,
            steps=steps,
            optimizer_cls=optimizer_cls,
            trajectory=trajectory,
            charge=charge,
            method=method,
        )

    if include_charges and charges is not None:
        optimised.set_initial_charges(charges)

    self._ase_cage = optimised
    return optimised

save ¶

save(path: str) -> None

Save the cage structure to file.

If the cage was built with STK, STK's writer is used for PDB files. Otherwise, the ASE structure is written using ase.io.write.

Parameters:

Name	Type	Description	Default
`path`	`str`	Output file path, e.g. `"base_structure/Tet2Di4.pdb"`.	required

Source code in cage_isomer_builder/cage/base.py

def save(self, path: str) -> None:
    """
    Save the cage structure to file.

    If the cage was built with STK, STK's writer is used for PDB files.
    Otherwise, the ASE structure is written using ``ase.io.write``.

    Parameters
    ----------
    path : str
        Output file path, e.g. ``"base_structure/Tet2Di4.pdb"``.
    """
    path = str(path)

    if self._stk_cage is not None and path.lower().endswith(".pdb"):
        stk.PdbWriter().write(molecule=self._stk_cage, path=path)
        return

    ase_cage = self.to_ase()
    write(path, ase_cage)

to_ase ¶

to_ase() -> Atoms

Return the cage as an ASE Atoms object. Converts from the STK cage if not already done.

Source code in cage_isomer_builder/cage/base.py

def to_ase(self) -> Atoms:
    """
    Return the cage as an ASE Atoms object.
    Converts from the STK cage if not already done.
    """
    if self._ase_cage is not None:
        return self._ase_cage
    if self._stk_cage is None:
        raise RuntimeError(
            "No structure available. Call build() or load_from_file() first."
        )
    self._ase_cage, _ = stk_2_ase_atoms(self._stk_cage)
    return self._ase_cage

to_smiles ¶

to_smiles() -> str

Convert the cage to a SMILES string using mofstructure (OpenBabel).

Returns:

Name	Type	Description
`smiles`	`str`

Source code in cage_isomer_builder/cage/base.py

def to_smiles(self) -> str:
    """
    Convert the cage to a SMILES string using mofstructure (OpenBabel).

    Returns
    -------
    smiles : str
    """
    return ase_2_pybel(self.to_ase()).write('smi').strip()

write_ams_run ¶

write_ams_run(path: Union[str, Path]) -> None

Write an AMS geometry-optimisation run script for this cage.

The generated .run file contains:

Task GeometryOptimization
Cartesian coordinates under System / Atoms
Bond orders under System / BondOrders (from the STK bond matrix)
Lattice vectors if the structure is periodic
Engine DFTB block

Parameters:

Name	Type	Description	Default
`path`	`str or Path`	Output file path, e.g. `"Tet2Di4.run"`.	required

Raises:

Type	Description
`RuntimeError`	If :meth:`build` has not been called — the bond matrix is only available from the STK constructed molecule.

Source code in cage_isomer_builder/cage/base.py

def write_ams_run(self, path: Union[str, Path]) -> None:
    """
    Write an AMS geometry-optimisation run script for this cage.

    The generated ``.run`` file contains:

    * ``Task GeometryOptimization``
    * Cartesian coordinates under ``System / Atoms``
    * Bond orders under ``System / BondOrders`` (from the STK bond matrix)
    * Lattice vectors if the structure is periodic
    * ``Engine DFTB`` block

    Parameters
    ----------
    path : str or Path
        Output file path, e.g. ``"Tet2Di4.run"``.

    Raises
    ------
    RuntimeError
        If :meth:`build` has not been called — the bond matrix is
        only available from the STK constructed molecule.
    """
    if self._stk_cage is None:
        raise RuntimeError(
            "write_run() requires the STK cage. Call build() first."
        )

    atoms, bond_matrix = stk_2_ase_atoms(self._stk_cage)
    write_run(str(path), atoms, bond_matrix)

write_gulp_gin ¶

write_gulp_gin(path: Union[str, Path], lattice: str = 'conv', mechanical: bool = False) -> None

Write a GULP molecular-mechanics input file for this cage.

Uses gulp_setup - the same UFF4MOF atom-typing backend behind optimise(method="UFF4MOF") - to assign UFF4MOF atom types and write a ready-to-run .gin file: cartesian coordinates, explicit connect bond records taken from the STK bond matrix (so bond order is respected exactly rather than re-perceived from geometry), and a uff4mof species/library block.

Parameters:

Name	Type	Description	Default
`path`	`str or Path`	Output file path, e.g. `"Tet2Di4.gin"`.	required
`lattice`	`str`	GULP lattice-optimisation keyword: `"conv"` for constant volume (no cell optimisation) or `"conp"` for constant pressure (cell optimisation). Only relevant if the structure is periodic - a freshly built cage isn't, so this has no effect unless you've set a cell on it yourself.	`"conv"`
`mechanical`	`bool`	Use `gulp_setup`'s "mechanical properties" input variant (looser convergence, different optimiser settings) instead of the standard geometry-optimisation one.	`False`

Raises:

Type	Description
`RuntimeError`	If :meth:`build` has not been called — the bond matrix is only available from the STK constructed molecule.

Notes

gulp_setup.mmanalysis.write_gin maps each bond order to a GULP bond-type keyword via a fixed lookup ({4, 3, 2, 1.5, 1, 0.5, 0.25}); a bond order outside that set raises KeyError. The STK bond matrix normally only contains integer orders (plus 1.5 for perceived aromaticity), so this doesn't come up in practice.

Source code in cage_isomer_builder/cage/base.py

def write_gulp_gin(
    self, path: Union[str, Path], lattice: str = "conv", mechanical: bool = False
) -> None:
    """
    Write a GULP molecular-mechanics input file for this cage.

    Uses ``gulp_setup`` - the same UFF4MOF atom-typing backend behind
    ``optimise(method="UFF4MOF")`` - to assign UFF4MOF atom types and
    write a ready-to-run ``.gin`` file: cartesian coordinates, explicit
    ``connect`` bond records taken from the STK bond matrix (so bond
    order is respected exactly rather than re-perceived from geometry),
    and a ``uff4mof`` species/library block.

    Parameters
    ----------
    path : str or Path
        Output file path, e.g. ``"Tet2Di4.gin"``.
    lattice : str, default "conv"
        GULP lattice-optimisation keyword: ``"conv"`` for constant
        volume (no cell optimisation) or ``"conp"`` for constant
        pressure (cell optimisation). Only relevant if the structure
        is periodic - a freshly built cage isn't, so this has no
        effect unless you've set a cell on it yourself.
    mechanical : bool, default False
        Use ``gulp_setup``'s "mechanical properties" input variant
        (looser convergence, different optimiser settings) instead of
        the standard geometry-optimisation one.

    Raises
    ------
    RuntimeError
        If :meth:`build` has not been called — the bond matrix is
        only available from the STK constructed molecule.

    Notes
    -----
    ``gulp_setup.mmanalysis.write_gin`` maps each bond order to a GULP
    bond-type keyword via a fixed lookup ({4, 3, 2, 1.5, 1, 0.5, 0.25});
    a bond order outside that set raises ``KeyError``. The STK bond
    matrix normally only contains integer orders (plus 1.5 for
    perceived aromaticity), so this doesn't come up in practice.
    """
    if self._stk_cage is None:
        raise RuntimeError(
            "write_gulp_run() requires the STK cage. Call build() first."
        )

    from gulp_setup.mmanalysis import analyze_mm, write_gin

    atoms, bond_matrix = stk_2_ase_atoms(self._stk_cage)
    _, mmtypes = analyze_mm(atoms, bond_order=bond_matrix)
    write_gin(str(path), atoms, bond_matrix, mmtypes, lattice=lattice, mechanical=mechanical)

CageStatistics `dataclass` ¶

CageStatistics(n_linkers: int, n_fg_slots: int, n_unique_isomers: int, distance_keys: list, total_fg_pairs: list, endo_endo_pairs: list, endo_exo_pairs: list, exo_exo_pairs: list, kde_x: ndarray, kde_y: ndarray)

Statistics for a cage topology computed after functionalisation.

Attributes:

Name	Type	Description
`n_linkers`	`int`	Number of ditopic linkers in the cage.
`n_fg_slots`	`int`	Total number of functional group anchor slots (n_linkers × 4).
`n_unique_isomers`	`int`	Number of symmetry-unique isomers (Burnside's lemma result).
`distance_keys`	`list of float`	Sorted unique inter-linker FG-to-FG distances in Angstroms.
`total_fg_pairs`	`list of int`	Total active FG pair count at each distance.
`endo_endo_pairs`	`list of int`	Pair counts where both anchors point inward (endo-endo).
`endo_exo_pairs`	`list of int`	Pair counts with one inward and one outward anchor (endo-exo).
`exo_exo_pairs`	`list of int`	Pair counts where both anchors point outward (exo-exo).
`kde_x`	`ndarray`	Distance grid at which the KDE was evaluated.
`kde_y`	`ndarray`	KDE density values rescaled to total pair count.

Topology-specific builders¶

Tet2Di4CageBuilder ¶

Tet2Di4CageBuilder(linker: MoleculeInput = None, node: MoleculeInput = None, linker_fg=None, node_fg=None, node_fg_count: int = None, smarts_pattern: str = None, optimizer_step: float = 0.01, distance_threshold: float = 5.0, expansion_factor: float = 1.0, scale_multiplier: float = 1.0, linker_smiles: str = None, node_smiles: str = None)

Cage builders¶

CageBuilder ¶

build ¶

count_unique_isomers_burnside ¶

enumerate_isomers ¶

functionalise ¶

get_bond_matrix ¶

get_statistics ¶

load_from_file classmethod ¶

optimise ¶

save ¶

to_ase ¶

to_smiles ¶

write_ams_run ¶

write_gulp_gin ¶

CageStatistics dataclass ¶

Topology-specific builders¶

Tet2Di4CageBuilder ¶

Tet3_3Di3CageBuilder ¶

Tet4_4Di8CageBuilder ¶

Tet5Di10CageBuilder ¶

Tet6Di12CageBuilder ¶

Tet8Di16CageBuilder ¶

Tet16Di32CageBuilder ¶

Tet24Di48CageBuilder ¶

Tri2Di3CageBuilder ¶

Tri4Di6CageBuilder ¶

Tri4_2Di6CageBuilder ¶

Tri6Di9CageBuilder ¶

Tri8Di12CageBuilder ¶

Tri20Di30CageBuilder ¶

Tet6Tri8CageBuilder ¶

Tri4Tri4CageBuilder ¶

Tri1Tri1CageBuilder ¶

Tri2_2Tri2CageBuilder ¶

M2L4LanternCageBuilder ¶

M3L6CageBuilder ¶

M4L6TetrahedronSpacerCageBuilder ¶

M4L8CageBuilder ¶

M4L82CageBuilder ¶

M4L8xCageBuilder ¶

M6L12CubeCageBuilder ¶

M8L6CubeCageBuilder ¶

M9L18CageBuilder ¶

TenPlusTwentyCageBuilder ¶

TwelvePlusThirtyCageBuilder ¶

load_from_file `classmethod` ¶

CageStatistics `dataclass` ¶