Dear scientific colleges,

I try to simulate a protein containing phosphorylation with OpenMM 8.1.1. To my understanding, this capability was added with #151
Let's have a look at this example:

from openmm.app import *
from openmm import *
from openmm.unit import *
from sys import stdout


forcefield = ForceField('./ff14SB.xml', './phosaa14SB.xml')
pdb = PDBFile('phospho14/5dfp_TPO.pdb')

# alternative force field:
#forcefield = ForceField('./ff10.xml')
#pdb = PDBFile('phospho10/5dfp_TPO.pdb')

system = forcefield.createSystem(pdb.topology, nonbondedMethod=PME, nonbondedCutoff=1*nanometer, constraints=HBonds)
integrator = LangevinIntegrator(300*kelvin, 1/picosecond, 0.002*picoseconds)
simulation = Simulation(pdb.topology, system, integrator)
simulation.context.setPositions(pdb.positions)
simulation.minimizeEnergy()

The pdbs are obtained from this repository (openmmforcefields-main/amber/files/phospho14 and phospho10).

Running this example yields for both force field choices the following error:
ValueError: No template found for residue 174 (SER). The set of atoms matches SER, but the bonds are different.

174 is the residue before the phosphorylation (residue 175). I'm pretty sure the error is triggered by residue 175 and not 174 as mentioned in the error message.

These might be related: #99 #203

Is there any way (at the moment) to simulate a protein containing phosphorylation? I feel like I'm missing out something obvious...

Thank you for providing OpenMM to the scientific community!

Best regards,
Johannes