Trouble creating RNA system with espaloma

[kntkb]

Problem

I'm having trouble creating a RNA system with espaloma. I want to apply the espaloma force field to the RNA and treat all other components (water and ions) with the traditional force field (amber/tip3p_standard.xml). However, I'm getting the following error.

I would appreciate if someone can help me. All input files and scripts are attached for reproducibility.
Archive.zip

Error

ValueError: No template found for residue 1 (U). This might mean your input topology is missing some atoms or bonds, or possibly that you are using the wrong force field.

Script (brief summary)

#!/usr/bin/env python
# coding: utf-8

### create rna system using espaloma generator
import glob as glob
import numpy as np
import re
import warnings
from openmm import *
from openmm.app import *
from openmm.unit import *
from pdbfixer import PDBFixer
from openff.toolkit.topology import Molecule, Topology
from openmmforcefields.generators import EspalomaTemplateGenerator, GAFFTemplateGenerator


box_padding = 12.0 * angstrom
salt_conc = 0.15 * molar
nb_cutoff = 10 * angstrom
hmass = 3.5 * amu #Might need to be tuned to 3.5 amu 

water_model='tip3p'


#### espaloma for RNA system
molecule = Molecule.from_file("rna.pdb", file_format="pdb")
#molecule.visualize()


generator = EspalomaTemplateGenerator(molecules=molecule, forcefield='espaloma-0.2.2')
#EspalomaTemplateGenerator.INSTALLED_FORCEFIELDS

ff = ForceField('amber/protein.ff14SB.xml', 'amber/tip3p_standard.xml', 'amber/tip3p_HFE_multivalent.xml')
ff.registerTemplateGenerator(generator.generator)
pdb = PDBFile('amber_min.pdb')


system = ff.createSystem(pdb.topology, nonbondedMethod=PME, nonbondedCutoff=nb_cutoff, constraints=HBonds, rigidWater=True, hydrogenMass=hmass)

[jchodera]

The main issue is that OpenMM passes only one Residue at a time to the residue template generator, but your RNA molecule spans multiple residues. However, you don't want to combine the residue indices before loading with PDBFile or else OpenMM won't automatically infer the bond connectivity for your RNA atoms.

The easiest fix right now is to create a new Topology object following this example that simply adds all the RNA atoms to a single residue, so that the entire molecule will be passed to the espaloma residue template generator.

[jchodera]

Alternatively, your PDB file could contain the RNA as a small molecule (all atoms in a single residue, HETATM records, CONECT records for all bonds) and this should just work.

[kntkb]

Closing issue because the problem is solved.