SEDenmarkLab / molli

Molecular Library Toolbox
https://molli.readthedocs.io
MIT License
21 stars 1 forks source link

Conformer Ensemble Initialization from Single Molecules #44

Closed blakeo2 closed 9 months ago

blakeo2 commented 1 year ago

Unable to correctly initialize conformer ensembles from single molecules.

Desired Code: conf_ens = ml.ConformerEnsemble(ml_mol)

Current Workaround: conf_ens = ml.ConformerEnsemble(ml_mol, n_conformers=1, coords=[ml_mol.coords])

austind4 commented 1 year ago

Hi I have investigated this issue and here are some things I found to help potentially resolve this issue.

The reason your work around works is because of how the code in ensemble.py is written.

The initialization of Conformer Ensemble has a default value of n_conformers = 0. Additionally, there is currently (from what I can tell) not a way to track the number of conformers from any class that ConformerEnsemble inherits from. Therefore, your work around works because you are manually assigning n_conformers. This allows line 54 of conformer ensemble to work.

Line 54:

        self._coords = np.full((n_conformers, self.n_atoms, 3), np.nan)

Now, if you add only n_conformers and not coordinates to your workaround solution, you get a matrix filled with np.nan. The reason you need to add 'coords = [ml_mol.coords]' is because of these lines 58-64:

if isinstance(other, ConformerEnsemble):
            self.atomic_charges = atomic_charges
            self.coords = other.coords
            self.weights = other.weights

The if statement is true in your example because this is an instance of Conformer Ensemble. Your workaround was setting other.coords to your manually designated coordinates. Then, the instance saved the coordinates as well.

        super().__init__(
            other,
            n_atoms=n_atoms,
            coords=coords, # Adding this line
            name=name,
            copy_atoms=copy_atoms,
            charge=charge,
            mult=mult,
            **kwds,

        coords = other.coords # Adding this line
        self._coords = np.full((n_conformers, self.n_atoms, 3), np.nan)
        self._atomic_charges = np.zeros((self.n_atoms,))
        self._weights = np.ones((n_conformers,))

This tested solution solves the issue, but I am not sure why.

blakeo2 commented 1 year ago

This helped me realize the reason we were running into the issues was because the self._coords is what is indicating both the number of conformers and construction of the coordinates, Since the default for coords is None and n_conformers is 0, it wasn't able to inherit it from a molecule object without it. Since the super().init is the same across most of the other classes in Molli, I figure we should probably avoid changing it. Here's the logic I used to fix it in mine (and I'll make a pull request for this shortly:

        super().__init__(
            other,
            n_atoms=n_atoms,
            name=name,
            copy_atoms=copy_atoms,
            charge=charge,
            mult=mult,
            **kwds,
        )

        if isinstance(other, Molecule):
            self._coords = np.full((1, self.n_atoms, 3), np.nan)
            self.coords = other.coords
        else:
            self._coords = np.full((n_conformers, self.n_atoms, 3), np.nan)