Open Force Field Initiative
pdb loading #29
-
|
Hi All, I am new to learning about openFF, and I'm currently facing difficulties loading the pdb file (attached as pdb-file.txt) of the TBP molecule using the provided commands. Additionally, I am also attempting to load a data file from lammps instead of a pdb file. It would be greatly appreciated if someone could assist me in resolving this problem? " Error:UnassignedChemistryInPDBError Traceback (most recent call last) File ~/anaconda3/envs/openff/lib/python3.11/site-packages/openff/utilities/utilities.py:80, in requires_package..inner_decorator..wrapper(*args, **kwargs) File ~/anaconda3/envs/openff/lib/python3.11/site-packages/openff/toolkit/topology/topology.py:1731, in Topology.from_pdb(cls, file_path, unique_molecules, toolkit_registry, _additional_substructures) File ~/anaconda3/envs/openff/lib/python3.11/site-packages/openff/toolkit/utils/toolkit_registry.py:356, in ToolkitRegistry.call(self, method_name, raise_exception_types, *args, **kwargs) File ~/anaconda3/envs/openff/lib/python3.11/site-packages/openff/toolkit/utils/toolkit_registry.py:352, in ToolkitRegistry.call(self, method_name, raise_exception_types, *args, **kwargs) File ~/anaconda3/envs/openff/lib/python3.11/site-packages/openff/toolkit/utils/rdkit_wrapper.py:289, in RDKitToolkitWrapper._polymer_openmm_pdbfile_to_offtop(self, topology_class, pdbfile, substructure_dictionary, coords_angstrom) File ~/anaconda3/envs/openff/lib/python3.11/site-packages/openff/toolkit/utils/rdkit_wrapper.py:462, in RDKitToolkitWrapper._polymer_openmm_topology_to_rdmol(self, omm_top, substructure_library) UnassignedChemistryInPDBError: Some bonds or atoms in the input could not be identified. Hint: The following residue names with unassigned atoms were not found in the substructure library. While the OpenFF Toolkit identifies residues by matching chemical substructures rather than by residue name, it currently only supports the 20 'canonical' amino acids. Error: The following 44 atoms exist in the input but could not be assigned chemical information from the substructure library: |
Beta Was this translation helpful? Give feedback.
Replies: 1 comment 4 replies
-
|
Hi @faranak1991, The problem is that your PDB file is missing CONECT records indicating where the bonds should be. If possible, see if your original data source can produce a pdb file with CONECT information included. But if not, you can try letting RDKit guess where the bonds are: from rdkit import Chem
mol = Chem.MolFromPDBFile('file.pdb', removeHs=False)
Chem.MolToPDBFile(mol, 'file_CONECT.pdb')
from openff.toolkit import Molecule, Topology
TBP = Molecule.from_smiles("O=P(OCCCC)(OCCCC)OCCCC")
topology = Topology.from_pdb( "file_CONECT.pdb", unique_molecules=[TBP] )This leads to a correctly-loaded molecule for me - Could you let me know if it works for you as well? |
Beta Was this translation helpful? Give feedback.
-
|
Dear Jeff Appreciate your quick reply. It works for me as well :)) Following after conversion to LAMMPS file, shows " UnsupportedExportError: Unsupported electrostatics method Coulomb " . Is this error related to the topology? Thanks, |
Beta Was this translation helpful? Give feedback.
-
|
Is your topology periodic? (What is |
Beta Was this translation helpful? Give feedback.
-
|
The topology is periodicity. I appreciate you bringing up the command "topology.box_vectors," which I overlooked. By applying the following command, I found that it effectively resolved that issue. topology.box_vectors = [30, 30, 30] * unit.angstrom If I modify the line "pair_style lj/cut/coul/long 9.0 9.0" in the tmp.in file to " 15 A," would that be correct? Thanks, |
Beta Was this translation helpful? Give feedback.
-
|
Not necessarily - if you're using an OpenFF force field, it's been trained on a 9 A cutoff so that's the recommended use. Using a longer cutoff is sometimes desirable, but my understanding is that it's better to stick with what the force field prescribes because that probably produces more accurate physical properties, etc. Also, in case this is what you happened to be doing, you probably don't want to keep the cutoff at half the box length, especially with larger boxes, because the cutoff error becomes vanishingly small but the extra computational cost increases significantly. If you wanted to modify the force field and use a different cutoff, and I'm not recommending you do, you can do that by setting it on the >>> from openff.toolkit import ForceField
>>> sage = ForceField("openff-2.1.0.offxml")
>>> sage['vdW'].cutoff
<Quantity(9.0, 'angstrom')>
>>> from openff.units import unit
>>> sage['vdW'].cutoff = 15 * unit.angstrom
>>> sage['vdW'].cutoff
<Quantity(15, 'angstrom')>before you do anything like |
Beta Was this translation helpful? Give feedback.
Not necessarily - if you're using an OpenFF force field, it's been trained on a 9 A cutoff so that's the recommended use. Using a longer cutoff is sometimes desirable, but my understanding is that it's better to stick with what the force field prescribes because that probably produces more accurate physical properties, etc. Also, in case this is what you happened to be doing, you probably don't want to keep the cutoff at half the box length, especially with larger boxes, because the cutoff error becomes vanishingly small but the extra computational cost increases significantly.
If you wanted to modify the force field and use a different cutoff, and I'm not recommending you do, you can do …