MELD DNA system indexing error

[tkdoesntgetcomputers]

Hello,

I am working on developing code to study DNA-Protein interactions. I am currently trying to find a workaround for the forcefield restrictions in order to simulate with Parmbsc1 or for other more complicated systems that don't have commonly used residues. I am effectively trying to build an amber system with tleap and then convert that to an openmm format to build the MELD system. I am experiencing an indexing issue that I can't seem to find the source of. Below I laid out the steps I have taken with the snippet of relevant code giving me an issue. At the bottom is the error I am getting.

1. Use tleap on minimized system to create a prmtop, inpcrd and new pdb.

source leaprc.protein.ff14SB

loadoff parmBSC1.lib loadamberparams parmBSC1.frcmod

WTP = loadpdb minimized_no_hydrogen.pdb

saveamberparm WTP minimized.prmtop minimized.inpcrd

savepdb WTP min_amber_output.pdb

quit #remove all hydrogens

2. Adjust atom and residue indices to be 0 based for MELD
3. Write MELD system with class meld.system.meld_system.System(solvation, openmm_system, openmm_topology, integrator ...)

I have written functions to see if there were notation mismatches between the pdb and prmtop as well as to check that the indexing correctly lines up. Both functions came back without errors saying that the top and pdb matched up.

SETUP CODE****

import numpy as np from meld.remd import ladder, adaptor, leader from meld import comm, vault from meld import system from meld import parse import meld as meld from meld.system.scalers import LinearRamp from meld.system.amber import import glob from restraints import from openmm import * from openmm import unit as u import mdtraj

N_REPLICAS = 30 N_STEPS = 20000 BLOCK_SIZE = 100

minimized_prmtop='minimized.prmtop' minimized_inpcrd='minimized.inpcrd'

def amber_to_openmm(minimized_prmtop,minimized_inpcrd):

loading in amber parameters

#Note the pdb is a new output from the tleap script to ensure that indexing is kept consistent
prmtop=app.AmberPrmtopFile(minimized_prmtop)
inpcrd=app.AmberInpcrdFile(minimized_inpcrd)

system=prmtop.createSystem()
topology=prmtop.topology
positions=inpcrd.positions
velocities=inpcrd.velocities
box_vectors=inpcrd.boxVectors

return system, topology, positions, velocities, box_vectors

def setup_system(): templates = glob.glob('Templates/*.pdb')

system, topology, positions, velocities, box_vectors = amber_to_openmm(minimized_prmtop,minimized_inpcrd)

positions_array = np.array(positions*u.nanometers)
positions_array = positions_array[np.newaxis, :, :]

integrator = LangevinIntegrator(300*u.kelvin, 1/u.picoseconds, 0.004*u.picoseconds)
barostat = MonteCarloBarostat(1*u.bar, 300*u.kelvin)

s = meld.system.meld_system.System(
    solvation='gbneck2',
    openmm_system=system,
    openmm_topology=topology,
    integrator=integrator,
    barostat=barostat,
    template_coordinates=positions_array,
    template_velocities=velocities,
    template_box_vectors=box_vectors,
    builder_info={}
)

ERROR****

File "/home/tek5268/Brini_research/MELD-DNA/DNA-1a74/setupMELD_testing.py", line 213, in setup_system() File "/home/tek5268/Brini_research/MELD-DNA/DNA-1a74/setupMELD_testing.py", line 103, in setup_system s = meld.system.meld_system.System( File "/.autofs/tools/spack/var/spack/environments/brini-23102001/.spack-env/view/lib/python3.10/site-packages/meld/system/meld_system.py", line 112, in init self._setup_indexing() File "/.autofs/tools/spack/var/spack/environments/brini-23102001/.spack-env/view/lib/python3.10/site-packages/meld/system/meld_system.py", line 325, in _setup_indexing assert len(self._atom_names) == self._n_atoms AssertionError

Please let me know if you can see any issues with this or have any tips on how to accomplish what I want to accomplish. Thank you!