FEP equilibration

[dlukauskis]

Hi, we are running some FEP validation tests using BSS from this paper. While some ligand pairs run fine, some can't be set up without allow_ring_breaking=True flag, and those fail to start as well, raising either parameter overflow-related errors or Traceback (most recent call last): File "/home/e620070/anaconda3/share/Sire/scripts/somd-freenrg.py", line 146, in <module> OpenMMMD.runFreeNrg(params) File "/home/e620070/anaconda3/lib/python3.7/site-packages/Sire/Tools/__init__.py", line 176, in inner retval = func() File "/home/e620070/anaconda3/lib/python3.7/site-packages/Sire/Tools/OpenMMMD.py", line 1637, in runFreeNrg system = integrator.minimiseEnergy(system, minimise_tol.val, minimise_max_iter.val) RuntimeError: Particle coordinate is nan

I ran the mapping and merge using the master and dev branches. I understand that this is a common issue, but if anyone could suggest possible solutions to the mapping/merging problem, that would be useful.

Attaching an input folder for replication.

Bace.zip

[lohedges]

Hi, can you please state the exact BioSimSpace version that you are using? You can get this as follows:

import BioSimSpace as BSS
print(BSS.__version__)

I can't even merge some of your problem ligand pairs with the latest devel commit of BioSimSpace, with it throwing the following error, e.g. for CAT-13a <-->CAT-17g:

IncompatibleError: Inconsistent number of bonds in merged molecule!

Others merge fine, but I've not run any simulations with them.

As for the parameter overflow error: For those sets, could you identify the number of different parameter records (bonds, angles, etc.) in the SOMD pertubation file and tell us which is triggering the error. (The error message might mention this.) As mentioned here, there is a limit to array sizes on the GPU which will be hardware dependent. I'm not sure how we would go about detecting this in order to reject any SOMD simulations which would exceed the limit, i.e. those with a large perturbation. This problem isn't so much a bug with BioSimSpace, rather a limitation of the way that SOMD / OpenMM manages the underlying simulation. (We could certainly fix SOMD if it's doing something in a silly way.)

The NaN errors from SOMD are likely due to steric clashes. This can be caused by poor equilibration or non-physical perturbations. Often, this problem is non-reproducible, so you might just have to run multiple simulations and take the ones that work. (The BioSimSpace free energy runner will try up to a maximum of 5 times per lambda leg before aborting the whole simulation.)

Do you have access to the exact mappings used in the paper? For the problem pairs, do the ones that are generated by BioSimSpace look sensible to you?

[lohedges]

I also note that the original paper uses a different force field, which will involve different terms in the perturbation, even if the mapping is the same. It could be interesting to find situations where a merge isn't possible, or gives bad results, with GAFF/GAFF2, but works with OPLS. (Or vice versa.)

[dlukauskis]

I'm running 2019.1.0+198.g54aa1b0 and I too can't merge the 13a_17g pair, attaching a list of problematic pairs for which merging failed/worked.

pairs_classified.txt

Indeed, there is a good correlation between parameter overflow and size of pert files. Most are between 700 and ~2000 lines, while the ones with parameter overflow are over 3000 lines long. I've looked at the error files and they all fail at 'torsionParams1', such as this one.

Is there an easy way to equilibrate with just somd-freenrg, or would you have to run 'normal' a simulation first for equilibration?

[lohedges]

Thanks for the additional information, that is helpful.

Yes, it's usual that dihedral/improper (torsion) terms are the problem with the overflow. This arises because of the modifications that the SOMD protocol makes to torsion terms involving dummy atoms for stability during the perturbation. If there are many dummy atoms, then the number of terms becomes large.

Is there an easy way to equilibrate with just somd-freenrg, or would you have to run 'normal' a simulation first for equilibration?

somd-freenrgy actually has a built in equilibration stage that can be activated. Because it's known to be unstable we decided to disable it in BioSimSpace. We normally just do a minimisation, then discard a chunk of the initial free energy data as an equilibration. However, you could simply use BioSimSpace to run minimisation and equilibration protocols before the free energy run. Both the BioSimSpace.Process.Gromacs and BioSimSpace.Process.Somd molecular dynamics drivers can handle systems containing a perturbable molecule. By default, they use the lambda = 0 state for their input and copy the updated coordinates back into your original system whenever you call getSystem on the object. At present its not possible to minimise or equilibrate an intermediate lambda state. (The somd-freenrgy equilibration equilibrates each lambda leg to the desired lambda value.) As an example, you could look at the script here. (The same concept should work with Gromacs or SOMD, but you might see stability issues.)

[jmichel80]

Could you share some details about the types of perturbations being attempted? I wouldn't expect issues if growing at most 1 ring per perturbation. Also we have processed the FEP+ dataset with a mix of BSS, FESetup and SOMD without running into those issues as far as I am aware so I'm wondering if you have run into mapping issues and are running calculations with sub-optimal mappings.

---

Dr. Julien Michel, Senior Lecturer Room 263, School of Chemistry University of Edinburgh David Brewster road Edinburgh, EH9 3FJ United Kingdom phone: +44 (0)131 650 4797 http://www.julienmichel.net/

On Mon, Jul 22, 2019 at 5:22 PM Lester Hedges notifications@github.com<mailto:notifications@github.com> wrote:

Thanks for the additional information, that is helpful.

Yes, it's usual that dihedral/improper (torsion) terms are the problem with the overflow. This arises because of the modifications that the SOMD protocol makes to torsion terms involving dummy atoms for stability during the perturbation. If there are many dummy atoms, then the number of terms becomes large.

Is there an easy way to equilibrate with just somd-freenrg, or would you have to run 'normal' a simulation first for equilibration?

somd-freenrgy actually has a built in equilibration stage that can be activated. Because it's known to be unstable we decided to disable it in BioSimSpace. We normally just do a minimisation, then discard a chunk of the initial free energy data as an equilibration. However, you could simply use BioSimSpace to run minimisation and equilibration protocols before the free energy run. Both the BioSimSpace.Process.Gromacs and BioSimSpace.Process.Somd molecular dynamics drivers can handle systems containing a perturbable molecule. By default, they use the lambda = 0 state for their input and copy the updated coordinates back into your original system whenever you call getSystem on the object. At present its not possible to minimise or equilibrate an intermediate lambda state. (The somd-freenrgy equilibration equilibrates each lambda leg to the desired lambda value.) As an example, you could look at the script herehttps://github.com/michellab/D3R2018/blob/master/CatS/BSS/binding_freenrg_gmx.py. (The same concept should work with Gromacs or SOMD, but you might see stability issues.)

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHubhttps://github.com/michellab/BioSimSpace/issues/97?email_source=notifications&email_token=ACZN3ZGYLYN36HJAXTW23OLQAXNCLA5CNFSM4IFY6BTKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD2QMQMI#issuecomment-513853489, or mute the threadhttps://github.com/notifications/unsubscribe-auth/ACZN3ZG55JK3PBPLUBO6KLDQAXNCLANCNFSM4IFY6BTA.

The University of Edinburgh is a charitable body, registered in Scotland, with registration number SC005336.

[dlukauskis]

I'm currently attempting to include equilibration into the protocol, analogous to the script you pointed to for Gromacs. However, I'm getting the following error:

Traceback (most recent call last): File "equilibration_binding_freenrg.py", line 66, in <module> minimised = BSS.Process.Somd(solvated, BSS.Protocol.Minimisation()) \ File "/home/e620070/anaconda3/lib/python3.7/site-packages/BioSimSpace/Process/_somd.py", line 152, in __init__ self._setup() File "/home/e620070/anaconda3/lib/python3.7/site-packages/BioSimSpace/Process/_somd.py", line 218, in _setup raise IOError("Failed to write system to 'PRM7' format.") from None OSError: Failed to write system to 'PRM7' format.

Here is the script that I used. I'm running on the dev branch, as described previously.

[lohedges]

Ah, this is because the GROMACS topology writer can handle perturbable molecules by default, i.e. the file contains both end states and we can write all the information to file. For SOMD, the topology file is simply AMBER format files for the lambda = 0 state. This is done for you when setting up a free energy simulation, along with some atom/residue renaming to match the required AMBER format (SOMD requires exact AMBER format waters, which means we have to rename those from the GROMACS solvated system.) This isn't done automatically when setting up minimisation or equilibration simulations.

As a workaround, you could try the following:

1. Simply run the minimisation and equilibration stages with GROMACS.
2. Manually extract and setup the lambda = 0 state for use with SOMD. This could be done as follows:

import BioSimSpace as BSS
import Sire.IO as IO
import Sire.Mol as Mol

# Here I assume we're at the point where you have a solvated system containing
# a merged molecule.

# Extract the waters from the existing system and convert them to AMBER format.
waters = IO.setAmberWater(solvated._sire_object.search("water"), "TIP3P")

# Remove all the original water molecules from the system, then add the renamed ones. 
solvated.removeWaterMolecules()
for wat in waters:
    solvated._sire_object.add(wat, Mol.MGName("all"))

# Now create a property map that maps properties in the system to those at lambda = 0,
# e.g. {"charge" : "charge0", "LJ" : "LJ0", ...}
property_map = {}
for prop in solvated._sire_object.propertyKeys():
    if prop[-1] == "0":
        property_map[prop[:-1]] = prop

# Now pass in the updated system and property map when initialising your process.
minimised = BSS.Process.Somd(solvated, BSS.Protocol.Minimisation(), property_map=property_map)

Going forward, I don't know whether we want to do all of this automatically with in BioSimSpace.Proces.Somd.

[dlukauskis]

Neither 1. or 2. seem to work for me.

For 2., I'm getting the same as before, unable to write 'PRM7'.

For 1., getting 'MDTraj failed to read ... ' trajectory reading error, which then fails to be passed further.

Including a folder for both.

[lohedges]

Are you positive that you are using the version of BioSImSpace that you stated? I've downgraded from the latest development version to 2019.1.0+198.g54aa1b0 and I am unable to merge the two ligands from your example folder:

Traceback (most recent call last):
  File "gmx_equil_binding_freenrg.py", line 58, in <module>
    merged = BSS.Align.merge(lig0, lig1, mapping, allow_ring_breaking=True)
  File "/home/lester/.conda/envs/biosimspace-dev/lib/python3.7/site-packages/BioSimSpace/Align/_align.py", line 576, in merge
    property_map0=property_map0, property_map1=property_map1)
  File "/home/lester/.conda/envs/biosimspace-dev/lib/python3.7/site-packages/BioSimSpace/_SireWrappers/_molecule.py", line 2443, in _merge
    raise _IncompatibleError("Inconsistent number of bonds in merged molecule!")
BioSimSpace._Exceptions._exceptions.IncompatibleError: Inconsistent number of bonds in merged molecule!

This is triggered because the number of bonds in the two lambda states differs and is not the result of a ring opening/closing (since allow_ring_breaking=True is set). This would suggest that the mapping is problematic and should be examined further.

However, for the errors that you see.

1. Are you running this as a batch script, or interactively like your SOMD example. If interactive, then you'll need to pass block=True to the getSystem call to make sure that you wait for the minimisation/equilibration to finish before attempting to grab the system. You might be failing to read the trajectory because it doesn't exist, or is partially written.

2. Apologies, my example above is wrong. You need a property map based on the properties of the merged molecule, not the solvated system. Try with the following replacement.

# Now create a property map that maps properties in the molecule to those at lambda = 0,
# e.g. {"charge" : "charge0", "LJ" : "LJ0", ...}
property_map = {}
for prop in merged._sire_object.propertyKeys():
    if prop[-1] == "0":
        property_map[prop[:-1]] = prop

[lohedges]

The merge error I see is for the GROMACS example script, where the batch command is:

/home/e620070/anaconda3/bin/sire_python --ppn=4 gmx_equil_binding_freenrg.py CAT-13a CAT-17g

I didn't notice that the SOMD notebook uses a different ligand pair:

# Extract the ligand numbers.
# CAT-17b_CAT-17e
num0 = 'CAT-17b'
num1 = 'CAT-17e'

[dlukauskis]

The merge error I see is for the GROMACS example script, where the batch command is: /home/e620070/anaconda3/bin/sire_python --ppn=4 gmx_equil_binding_freenrg.py CAT-13a CAT-17g I didn't notice that the SOMD notebook uses a different ligand pair: Extract the ligand numbers. CAT-17b_CAT-17e num0 = 'CAT-17b' num1 = 'CAT-17e'

Yeah, sorry, that's my mistake. Let's isolate this down to a pair that maps/merges okay, so 17b_17e.

Running the somd notebook still fails, even after correcting for property_map, BSS is unable to write PRM7.

[lohedges]

The AMBER PRM7 writing issue is because we only want to map the properties for the merged molecule, not for the protein and solvent. Sorry about that.

What you can do is as follows:

# Convert the merged molecule to a regular molecule at lambda = 0.
regular_molecule = merged._toRegularMolecule()

# Now create the composite system to solvate.
system = regular_molecule + protein

You should now be able to solvate and simulate with SOMD. When you get the output of the equilibration back you will need to reconstruct the solvated system with the merged molecule, then update with the coordinates from the equilibrated system.

# Make a new solvated system containing the merged molecule. We need to get the
# water molecules and ions from the original solvated system. Note that I have just manually
# worked out the indices of the ions, which are always at the end.
new_system = regular_molecule + protein + solvated.getWaterMolecules() + solvated.getMolecules()[7009:]

# Now copy the coordinates from the equilibrated system. You would then pass
# new_system to the free energy object.
new_system._updateCoordinates(equilibrated)

I've tested this up to the creation of new_system and I can successfully write it to PRM7 format.

Apologies for the clunkiness, the way that SOMD works makes this a pain. Since the GROMACS topology file contains all of the information for the perturbed molecule and is also able to handle any water naming convention, then these problems don't exist.

(Perhaps I need to write a getIons method! Could have also just done solvated.getMolecules()[2:] since we know that the protein and merged molecule are first.)

[lohedges]

Note that, by default, SOMD simulations run on the CPU. This is very slow, so you might want to pass platform="CUDA" to the BioSimSpace.Process.Somd constructor when setting up minimisation or equilibration simulations. Free energy simulations that are configured using the BioSimSpace.FreeEnergy package use platform="CUDA" by default when CUDA_VISIBLE_DEVICES is set. (SOMD free energy simulations are painfully slow on a CPU.)

[lohedges]

I've just realised that my example above won't work, since you also need to copy across the properties of the original system, such as the periodic box information. I'm currently incorporating all this into BioSimSpace.Process.Somd since it's such a pain to do things manually.

[dlukauskis]

Great, since using Gromacs for the equilibration and Somd for setting up FEP doesn't work consistently. Do keep me updated.

[lohedges]

This is now implemented in devel and an updated Conda package will be available soon. You should now be able to pass a system containing a perturbable molecule to BioSimSpace.Process.Somd and run any protocol. (It will raise a warning to say that it is using the lambda = 0 state.) When you call getSystem on the process it will copy the updated coordinates back into your original system meaning that all of the existing information is preserved. (You don't loose the merged molecule.)

[dlukauskis]

Does this work for you? It successfully minimises and equilibrates but fails to make the folder tree for FEP, running out of memory.

Specifically:

Traceback (most recent call last):
  File "test_somd.py", line 86, in <module>
    engine='SOMD', box=3*[35*BSS.Units.Length.angstrom])
  File "/home/e620070/anaconda3/envs/BSS_dev/lib/python3.7/site-packages/BioSimSpace/FreeEnergy/_binding.py", line 123, in __init__
    self._system1 = _Solvent.solvate(water_model, molecule=molecule, box=box)
  File "/home/e620070/anaconda3/envs/BSS_dev/lib/python3.7/site-packages/BioSimSpace/Solvent/_solvent.py", line 103, in solvate
    return _model_dict[model](molecule, box, shell, ion_conc, is_neutral, work_dir, property_map)
  File "/home/e620070/anaconda3/envs/BSS_dev/lib/python3.7/site-packages/BioSimSpace/Solvent/_solvent.py", line 261, in tip3p
    is_neutral, work_dir=work_dir, property_map=property_map)
  File "/home/e620070/anaconda3/envs/BSS_dev/lib/python3.7/site-packages/BioSimSpace/Solvent/_solvent.py", line 629, in _solvate
    proc = _subprocess.run(command, shell=True, stdout=stdout, stderr=stderr)
  File "/home/e620070/anaconda3/envs/BSS_dev/lib/python3.7/subprocess.py", line 472, in run
    with Popen(*popenargs, **kwargs) as process:
  File "/home/e620070/anaconda3/envs/BSS_dev/lib/python3.7/subprocess.py", line 775, in __init__
    restore_signals, start_new_session)
  File "/home/e620070/anaconda3/envs/BSS_dev/lib/python3.7/subprocess.py", line 1453, in _execute_child
    restore_signals, start_new_session, preexec_fn)
OSError: [Errno 12] Cannot allocate memory

[lohedges]

I think the problem here is that the free energy simulation is handled by a BioSimSpace.Process.ProcessRunner object, which is creating a BioSimSpace.Process.Somd object for each of the lambda windows for both the bound and free legs of the simulation. The ProcessRunner is designed to run arbitrary processes, either sequentially (as is done here) or distrubuted amongst the available hardware resources that are available (not yet implemented). You just give it a list of processes, which could be applying different protocols to the same system (as is done here), or using completely different systems.

When a Process object is instantiated, it creates a copy of the System internally. This is required since the job may be running interactively and the user could modify the system elsewhere. We need a record of the System so that getSystem can return information in a consistent format to what was originally passed in, i.e. just updating the coordinates. (Some processes modify naming, etc.) For ProcessRunner used here, you end up creating memory for the System object 24 times (plus the original object), which is causing the problem. (In reality, there should be one System for the bound leg.)

A free energy simulation is a special case, since we know that it will use the same System for each simulation, i.e. only the Protocol is changing. I think I'll need to implement a specific ProcessRunner for this type of simulation, where the memory for the System object is shared among all of the processes. (This would require some reworking of the BioSimSpace.Process package.) Alternatively, I could just create the individual processes one-by-one as the simulation progresses, although this doesn't provide the flexibility of distributing the processes in different ways. Perhaps this could be a different kind of ProcessRunner, e.g. SerialProcessRunner. I'll have a think about this to see what fits best with our plans.

[lohedges]

Actually, I don't think this is the problem. I've run your example (minus minimisation and equilibration since I'm on my laptop) and I can create the BindingFreeEnergy object passing in the solvated system. I checked the memory footprint of my iPython session both with and without copying the system inside each Process object and it barely exceeded 500 MiB in both cases, i.e. less resources than my browser is using while writing this reply.

I wonder if you're experiencing some left-over memory footprint from the equilibration. SOMD doesn't write restart files for individual molecular configurations so getSystem reads an entire trajectory (which might be large for a long simulation) and takes the final frame. However, this memory should be freed once the system is returned.

Could you let me know if you get the memory error immediately when instantiating the BindingFreeEnergy object, or whether you see several folders created before it occurs? Could you monitor the memory footprint while your script is running to see what happens when the error is raised? How much memory do you have on the machine where your script was running?

From your script I see that you would like a way to determine the box size that should be passed to the solvate function. You can do this with some of the hidden Sire functions within BioSimSpace, i.e. to get an AABox object that represents the axis-aligned bounding box of the system:

# Here the system is your protein plus the merged molecule.
aabox = system._getAABox()

# A buffer around the system.
buffer = BSS.Units.Length.angstrom

# Create a box size list from the minimum enclosing box plus a buffer.
box_size = [2*BSS.Units.Length.angstrom + buffer for x in aabox.halfExtents()]

[lohedges]

I've now exposed this directly to the BioSimSpace objects so you can do the following:

# Here the system is your protein plus the merged molecule.
box_min, box_max = system.getAxisAlignedBoundingBox()

# Work out the box size from the difference in the coordinates.
box_size = [y - x for x, y in zip(box_min, box_max)]

[dlukauskis]

Thanks a lot for the box tool, will check it out.

The script did seem to successfully get the equilibrated system, topping at 56% of MEM (on a machine with 64GB). It made the main work_dir for FEP, but no other directories inside it. Very strange.

[lohedges]

I've tested here with a short equilibration and I'm pretty confident that, as suggested above, the issue is the size of the trajectory that is loaded following the equilibration stage. I'll look into a way of working around this, since there should be no need to read the entire trajectory. The only issue is that it's sometimes hard to know the exact number of frames in the file. I'll let you know once I've come up with a solution.

[lohedges]

Basically, the problem is that SOMD only writes coordinates to a trajectory and since we need to use getSystem frequently we write to the file fairly often. I'll look at updating SOMD so that it can also write a single frame as a restart. That way we only ever need to load the trajectory if we actually want a trajectory, and can save frames at a less frequent interval.

[jmichel80]

Hi Lester,

Maybe I’m missing something but SOMD now updates sim_restart.s3 with coordinates every cycle so you don’t need to load a full trajectory just to get the last frame.

Best wishes

Julien

Sent from my iPhone

On 25 Jul 2019, at 14:47, Lester Hedges notifications@github.com wrote:

Basically, the problem is that SOMD only writes coordinates to a trajectory and since we need to use getSystem frequently we write to the file fairly often. I'll look at updating SOMD so that it can also write a single frame as a restart. That way we only ever need to load the trajectory if we actually want a trajectory, and can save frames at a less frequent interval.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub, or mute the thread.

[lohedges]

Ah great, I'll look into that.

[lohedges]

@jmichel80 How do you actually load and use the streamed Sire system? I've followed what's done in OpenMMMD.py, but get errors whenever I try to do anything to the system, e.g.:

from Sire.Stream import load

# Load a binary restart file.
system, _ = load("/tmp/tmpo5lhzy3c/sim_restart.s3")

# Try to get the molecule numbers.
system.molNums()

This results in the following error:

UserWarning: Exception 'SireError::program_bug' thrown by the thread 'master'.
There is a problem - the parent forcefield of molecules:molecules (7) is null!
Thrown from FILE: /home/lester/Code/Sire/corelib/src/libs/SireFF/ffmolgroup.cpp, LINE: 144, FUNCTION: void SireFF::detail::FFMolGroupPvt::assertNotNull() const
__Backtrace__
(  0) /home/lester/sire.app/lib/libSireError.so.2019 ([0x7f178da3aec2] ++0x42)
  -- SireError::getBackTrace()

(  1) /home/lester/sire.app/lib/libSireError.so.2019 ([0x7f178da379d6] ++0xa6)
  -- SireError::exception::exception(QString, QString)

/home/lester/sire.app/lib/libSireFF.so.2019(+0x758b0) [0x7f17357098b0]
(  3) /home/lester/sire.app/lib/libSireFF.so.2019 ([0x7f17357175a3] ++0x1e3)
  -- SireFF::detail::FFMolGroupPvt::assertNotNull() const

(  4) /home/lester/sire.app/lib/libSireFF.so.2019 ([0x7f1735717738] ++0x28)
  -- SireFF::FFMolGroup::FFMolGroup(SireFF::detail::FFMolGroupPvt const&)

(  5) /home/lester/sire.app/lib/libSireFF.so.2019 ([0x7f17357177c3] ++0x23)
  -- SireFF::detail::FFMolGroupPvt::clone() const

(  6) /home/lester/sire.app/lib/libSireBase.so.2019 ([0x7f1760b06612] ++0x32)
  -- SireBase::PropPtrBase::PropPtrBase(SireBase::Property const&)

(  7) /home/lester/sire.app/lib/libSireMol.so.2019 ([0x7f1735b05fa3] ++0xb3)
  -- SireMol::MolGroupsBase::selectAll() const

(  8) /home/lester/sire.app/lib/libSireMol.so.2019 ([0x7f1735b0611e] ++0x1e)
  -- SireMol::MolGroupsBase::groups() const

(  9) /home/lester/sire.app/lib/libSireMol.so.2019 ([0x7f1735b0ccdc] ++0x4c)
  -- SireMol::MolGroupsBase::getMoleculeNumbers() const

( 10) /home/lester/sire.app/lib/libSireMol.so.2019 ([0x7f1735b0d16e] ++0x1e)
  -- SireMol::MolGroupsBase::molNums() const

/home/lester/sire.app/lib/python3.7/site-packages/Sire/Mol/_Mol.so(+0x82de85) [0x7f1736b86e85]
( 12) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427442d] ++0x29d)
  -- boost::python::objects::function::call(_object*, _object*) const

/home/lester/sire.app/lib/libboost_python37.so.1.67.0(+0x20599) [0x7f1784274599]
( 14) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427aceb] ++0x6b)
  -- boost::python::detail::exception_handler::operator()(boost::function0<void> const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/Chem/rdmolfiles.so(+0x53c94) [0x7f17462b9c94]
( 16) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a)
  -- boost::python::detail::exception_handler::operator()(boost::function0<void> const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/Chem/rdchem.so(+0x10fd54) [0x7f17488e6d54]
( 18) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a)
  -- boost::python::detail::exception_handler::operator()(boost::function0<void> const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/Chem/rdchem.so(+0xbf0c4) [0x7f17488960c4]
( 20) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a)
  -- boost::python::detail::exception_handler::operator()(boost::function0<void> const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/rdBase.so(+0x54d34) [0x7f1761f92d34]
( 22) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a)
  -- boost::python::detail::exception_handler::operator()(boost::function0<void> const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/rdBase.so(+0x54c54) [0x7f1761f92c54]
( 24) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a)
  -- boost::python::detail::exception_handler::operator()(boost::function0<void> const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/rdBase.so(+0x54c04) [0x7f1761f92c04]
( 26) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a)
  -- boost::python::detail::exception_handler::operator()(boost::function0<void> const&) const

/home/lester/sire.app/lib/python3.7/site-packages/Sire/Error/_Error.so(+0x2aa3) [0x7f178da24aa3]
( 28) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a)
  -- boost::python::detail::exception_handler::operator()(boost::function0<void> const&) const

/home/lester/sire.app/lib/python3.7/site-packages/Sire/Error/_Error.so(+0x2a61) [0x7f178da24a61]
( 30) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427aa6f] ++0x3f)
  -- boost::python::handle_exception_impl(boost::function0<void>)

/home/lester/sire.app/lib/libboost_python37.so.1.67.0(+0x1dd2a) [0x7f1784271d2a]
( 32) /home/lester/sire.app/bin/python ([0x55772b79f46b] ++0x49b)
  -- _PyObject_FastCallKeywords

( 33) /home/lester/sire.app/bin/python ([0x55772b7f9d26] ++0x5266)
  -- _PyEval_EvalFrameDefault

( 34) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9)
  -- _PyEval_EvalCodeWithName

( 35) /home/lester/sire.app/bin/python ([0x55772b73afa4] ++0x44)
  -- PyEval_EvalCodeEx

( 36) /home/lester/sire.app/bin/python ([0x55772b73afcc] ++0x1c)
  -- PyEval_EvalCode

/home/lester/sire.app/bin/python(+0x1df8ad) [0x55772b8048ad]
( 38) /home/lester/sire.app/bin/python ([0x55772b79eca9] ++0xe9)
  -- _PyMethodDef_RawFastCallKeywords

( 39) /home/lester/sire.app/bin/python ([0x55772b79ef41] ++0x21)
  -- _PyCFunction_FastCallKeywords

( 40) /home/lester/sire.app/bin/python ([0x55772b7f91f4] ++0x4734)
  -- _PyEval_EvalFrameDefault

( 41) /home/lester/sire.app/bin/python ([0x55772b797e52] ++0x2a2)
  -- _PyGen_Send

( 42) /home/lester/sire.app/bin/python ([0x55772b7f64b6] ++0x19f6)
  -- _PyEval_EvalFrameDefault

( 43) /home/lester/sire.app/bin/python ([0x55772b797e52] ++0x2a2)
  -- _PyGen_Send

( 44) /home/lester/sire.app/bin/python ([0x55772b7f64b6] ++0x19f6)
  -- _PyEval_EvalFrameDefault

( 45) /home/lester/sire.app/bin/python ([0x55772b797e52] ++0x2a2)
  -- _PyGen_Send

( 46) /home/lester/sire.app/bin/python ([0x55772b79ec4c] ++0x8c)
  -- _PyMethodDef_RawFastCallKeywords

( 47) /home/lester/sire.app/bin/python ([0x55772b79efaf] ++0x4f)
  -- _PyMethodDescr_FastCallKeywords

( 48) /home/lester/sire.app/bin/python ([0x55772b7f96a1] ++0x4be1)
  -- _PyEval_EvalFrameDefault

( 49) /home/lester/sire.app/bin/python ([0x55772b79decb] ++0xfb)
  -- _PyFunction_FastCallKeywords

( 50) /home/lester/sire.app/bin/python ([0x55772b7f51a6] ++0x6e6)
  -- _PyEval_EvalFrameDefault

( 51) /home/lester/sire.app/bin/python ([0x55772b79decb] ++0xfb)
  -- _PyFunction_FastCallKeywords

( 52) /home/lester/sire.app/bin/python ([0x55772b7f4ff0] ++0x530)
  -- _PyEval_EvalFrameDefault

( 53) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9)
  -- _PyEval_EvalCodeWithName

( 54) /home/lester/sire.app/bin/python ([0x55772b79e157] ++0x387)
  -- _PyFunction_FastCallKeywords

( 55) /home/lester/sire.app/bin/python ([0x55772b7f5fa1] ++0x14e1)
  -- _PyEval_EvalFrameDefault

 -- _PyEval_EvalFrameDefault

( 56) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9)
  -- _PyEval_EvalCodeWithName

( 57) /home/lester/sire.app/bin/python ([0x55772b79e0f5] ++0x325)
  -- _PyFunction_FastCallKeywords

( 58) /home/lester/sire.app/bin/python ([0x55772b7f4ff0] ++0x530)
  -- _PyEval_EvalFrameDefault

( 59) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9)
  -- _PyEval_EvalCodeWithName

( 60) /home/lester/sire.app/bin/python ([0x55772b79e0f5] ++0x325)
  -- _PyFunction_FastCallKeywords

( 61) /home/lester/sire.app/bin/python ([0x55772b7f4ff0] ++0x530)
  -- _PyEval_EvalFrameDefault

( 62) /home/lester/sire.app/bin/python ([0x55772b79decb] ++0xfb)
  -- _PyFunction_FastCallKeywords

( 63) /home/lester/sire.app/bin/python ([0x55772b7f4ff0] ++0x530)
  -- _PyEval_EvalFrameDefault

( 64) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9)
  -- _PyEval_EvalCodeWithName

( 65) /home/lester/sire.app/bin/python ([0x55772b73b3bc] ++0x3dc)
  -- _PyFunction_FastCallDict

( 66) /home/lester/sire.app/bin/python ([0x55772b752473] ++0x63)
  -- _PyObject_Call_Prepend

( 67) /home/lester/sire.app/bin/python ([0x55772b74719e] ++0x6e)
  -- PyObject_Call

( 68) /home/lester/sire.app/bin/python ([0x55772b7f69ff] ++0x1f3f)
  -- _PyEval_EvalFrameDefault

( 69) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9)
  -- _PyEval_EvalCodeWithName

( 70) /home/lester/sire.app/bin/python ([0x55772b79e157] ++0x387)
  -- _PyFunction_FastCallKeywords

( 71) /home/lester/sire.app/bin/python ([0x55772b7f51a6] ++0x6e6)
  -- _PyEval_EvalFrameDefault

( 72) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9)
  -- _PyEval_EvalCodeWithName

( 73) /home/lester/sire.app/bin/python ([0x55772b73afa4] ++0x44)
  -- PyEval_EvalCodeEx

( 74) /home/lester/sire.app/bin/python ([0x55772b73afcc] ++0x1c)
  -- PyEval_EvalCode

/home/lester/sire.app/bin/python(+0x22f664) [0x55772b854664]
( 76) /home/lester/sire.app/bin/python ([0x55772b85e881] ++0xa1)
  -- PyRun_FileExFlags

( 77) /home/lester/sire.app/bin/python ([0x55772b85ea73] ++0x1c3)
  -- PyRun_SimpleFileExFlags

/home/lester/sire.app/bin/python(+0x23ab2f) [0x55772b85fb2f]
( 79) /home/lester/sire.app/bin/python ([0x55772b85fc4c] ++0x3c)
  -- _Py_UnixMain

( 80) /usr/lib/libc.so.6 ([0x7f179026bee3] ++0xf3)
  -- __libc_start_main

/home/lester/sire.app/bin/python(+0x1df982) [0x55772b804982]
__EndTrace__
Exception 'SireError::program_bug' thrown by the thread 'master'.
There is a problem - the parent forcefield of molecules:molecules (7) is null!
Thrown from FILE: /home/lester/Code/Sire/corelib/src/libs/SireFF/ffmolgroup.cpp, LINE: 144, FUNCTION: void SireFF::detail::FFMolGroupPvt::assertNotNull() const

I think it is easier for me to just modify OpenMMMD.py so that it writes a separate DCD file for coordinates only. This is trivial to do and I've tested that it works properly, unlike the restart file above.

[jmichel80]

Odd...check that somd actually restarts correctly when run from the command line in a folder that already contains a restart file (current behavior). Could be an import problem, i.e. you need to load more modules.

Sent from my iPhone

On 25 Jul 2019, at 16:03, Lester Hedges notifications@github.com wrote:

@jmichel80 How do you actually load and use the streamed Sire system? I've followed what's done in OpenMMMD.py, but get errors whenever I try to do anything to the system, e.g.:

from Sire.Stream import load

Load a binary restart file.

system, _ = load("/tmp/tmpo5lhzy3c/sim_restart.s3")

Try to get the molecule numbers.

system.molNums() This results in the following error:

UserWarning: Exception 'SireError::program_bug' thrown by the thread 'master'. There is a problem - the parent forcefield of molecules:molecules (7) is null! Thrown from FILE: /home/lester/Code/Sire/corelib/src/libs/SireFF/ffmolgroup.cpp, LINE: 144, FUNCTION: void SireFF::detail::FFMolGroupPvt::assertNotNull() const Backtrace ( 0) /home/lester/sire.app/lib/libSireError.so.2019 ([0x7f178da3aec2] ++0x42) -- SireError::getBackTrace()

( 1) /home/lester/sire.app/lib/libSireError.so.2019 ([0x7f178da379d6] ++0xa6) -- SireError::exception::exception(QString, QString)

/home/lester/sire.app/lib/libSireFF.so.2019(+0x758b0) [0x7f17357098b0] ( 3) /home/lester/sire.app/lib/libSireFF.so.2019 ([0x7f17357175a3] ++0x1e3) -- SireFF::detail::FFMolGroupPvt::assertNotNull() const

( 4) /home/lester/sire.app/lib/libSireFF.so.2019 ([0x7f1735717738] ++0x28) -- SireFF::FFMolGroup::FFMolGroup(SireFF::detail::FFMolGroupPvt const&)

( 5) /home/lester/sire.app/lib/libSireFF.so.2019 ([0x7f17357177c3] ++0x23) -- SireFF::detail::FFMolGroupPvt::clone() const

( 6) /home/lester/sire.app/lib/libSireBase.so.2019 ([0x7f1760b06612] ++0x32) -- SireBase::PropPtrBase::PropPtrBase(SireBase::Property const&)

( 7) /home/lester/sire.app/lib/libSireMol.so.2019 ([0x7f1735b05fa3] ++0xb3) -- SireMol::MolGroupsBase::selectAll() const

( 8) /home/lester/sire.app/lib/libSireMol.so.2019 ([0x7f1735b0611e] ++0x1e) -- SireMol::MolGroupsBase::groups() const

( 9) /home/lester/sire.app/lib/libSireMol.so.2019 ([0x7f1735b0ccdc] ++0x4c) -- SireMol::MolGroupsBase::getMoleculeNumbers() const

( 10) /home/lester/sire.app/lib/libSireMol.so.2019 ([0x7f1735b0d16e] ++0x1e) -- SireMol::MolGroupsBase::molNums() const

/home/lester/sire.app/lib/python3.7/site-packages/Sire/Mol/_Mol.so(+0x82de85) [0x7f1736b86e85] ( 12) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427442d] ++0x29d) -- boost::python::objects::function::call(_object, _object) const

/home/lester/sire.app/lib/libboost_python37.so.1.67.0(+0x20599) [0x7f1784274599] ( 14) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427aceb] ++0x6b) -- boost::python::detail::exception_handler::operator()(boost::function0 const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/Chem/rdmolfiles.so(+0x53c94) [0x7f17462b9c94] ( 16) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a) -- boost::python::detail::exception_handler::operator()(boost::function0 const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/Chem/rdchem.so(+0x10fd54) [0x7f17488e6d54] ( 18) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a) -- boost::python::detail::exception_handler::operator()(boost::function0 const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/Chem/rdchem.so(+0xbf0c4) [0x7f17488960c4] ( 20) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a) -- boost::python::detail::exception_handler::operator()(boost::function0 const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/rdBase.so(+0x54d34) [0x7f1761f92d34] ( 22) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a) -- boost::python::detail::exception_handler::operator()(boost::function0 const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/rdBase.so(+0x54c54) [0x7f1761f92c54] ( 24) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a) -- boost::python::detail::exception_handler::operator()(boost::function0 const&) const

/home/lester/sire.app/lib/python3.7/site-packages/rdkit/rdBase.so(+0x54c04) [0x7f1761f92c04] ( 26) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a) -- boost::python::detail::exception_handler::operator()(boost::function0 const&) const

/home/lester/sire.app/lib/python3.7/site-packages/Sire/Error/_Error.so(+0x2aa3) [0x7f178da24aa3] ( 28) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427acba] ++0x3a) -- boost::python::detail::exception_handler::operator()(boost::function0 const&) const

/home/lester/sire.app/lib/python3.7/site-packages/Sire/Error/_Error.so(+0x2a61) [0x7f178da24a61] ( 30) /home/lester/sire.app/lib/libboost_python37.so.1.67.0 ([0x7f178427aa6f] ++0x3f) -- boost::python::handle_exception_impl(boost::function0)

/home/lester/sire.app/lib/libboost_python37.so.1.67.0(+0x1dd2a) [0x7f1784271d2a] ( 32) /home/lester/sire.app/bin/python ([0x55772b79f46b] ++0x49b) -- _PyObject_FastCallKeywords

( 33) /home/lester/sire.app/bin/python ([0x55772b7f9d26] ++0x5266) -- _PyEval_EvalFrameDefault

( 34) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9) -- _PyEval_EvalCodeWithName

( 35) /home/lester/sire.app/bin/python ([0x55772b73afa4] ++0x44) -- PyEval_EvalCodeEx

( 36) /home/lester/sire.app/bin/python ([0x55772b73afcc] ++0x1c) -- PyEval_EvalCode

/home/lester/sire.app/bin/python(+0x1df8ad) [0x55772b8048ad] ( 38) /home/lester/sire.app/bin/python ([0x55772b79eca9] ++0xe9) -- _PyMethodDef_RawFastCallKeywords

( 39) /home/lester/sire.app/bin/python ([0x55772b79ef41] ++0x21) -- _PyCFunction_FastCallKeywords

( 40) /home/lester/sire.app/bin/python ([0x55772b7f91f4] ++0x4734) -- _PyEval_EvalFrameDefault

( 41) /home/lester/sire.app/bin/python ([0x55772b797e52] ++0x2a2) -- _PyGen_Send

( 42) /home/lester/sire.app/bin/python ([0x55772b7f64b6] ++0x19f6) -- _PyEval_EvalFrameDefault

( 43) /home/lester/sire.app/bin/python ([0x55772b797e52] ++0x2a2) -- _PyGen_Send

( 44) /home/lester/sire.app/bin/python ([0x55772b7f64b6] ++0x19f6) -- _PyEval_EvalFrameDefault

( 45) /home/lester/sire.app/bin/python ([0x55772b797e52] ++0x2a2) -- _PyGen_Send

( 46) /home/lester/sire.app/bin/python ([0x55772b79ec4c] ++0x8c) -- _PyMethodDef_RawFastCallKeywords

( 47) /home/lester/sire.app/bin/python ([0x55772b79efaf] ++0x4f) -- _PyMethodDescr_FastCallKeywords

( 48) /home/lester/sire.app/bin/python ([0x55772b7f96a1] ++0x4be1) -- _PyEval_EvalFrameDefault

( 49) /home/lester/sire.app/bin/python ([0x55772b79decb] ++0xfb) -- _PyFunction_FastCallKeywords

( 50) /home/lester/sire.app/bin/python ([0x55772b7f51a6] ++0x6e6) -- _PyEval_EvalFrameDefault

( 51) /home/lester/sire.app/bin/python ([0x55772b79decb] ++0xfb) -- _PyFunction_FastCallKeywords

( 52) /home/lester/sire.app/bin/python ([0x55772b7f4ff0] ++0x530) -- _PyEval_EvalFrameDefault

( 53) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9) -- _PyEval_EvalCodeWithName

( 54) /home/lester/sire.app/bin/python ([0x55772b79e157] ++0x387) -- _PyFunction_FastCallKeywords

( 55) /home/lester/sire.app/bin/python ([0x55772b7f5fa1] ++0x14e1) -- _PyEval_EvalFrameDefault

-- _PyEval_EvalFrameDefault

( 56) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9) -- _PyEval_EvalCodeWithName

( 57) /home/lester/sire.app/bin/python ([0x55772b79e0f5] ++0x325) -- _PyFunction_FastCallKeywords

( 58) /home/lester/sire.app/bin/python ([0x55772b7f4ff0] ++0x530) -- _PyEval_EvalFrameDefault

( 59) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9) -- _PyEval_EvalCodeWithName

( 60) /home/lester/sire.app/bin/python ([0x55772b79e0f5] ++0x325) -- _PyFunction_FastCallKeywords

( 61) /home/lester/sire.app/bin/python ([0x55772b7f4ff0] ++0x530) -- _PyEval_EvalFrameDefault

( 62) /home/lester/sire.app/bin/python ([0x55772b79decb] ++0xfb) -- _PyFunction_FastCallKeywords

( 63) /home/lester/sire.app/bin/python ([0x55772b7f4ff0] ++0x530) -- _PyEval_EvalFrameDefault

( 64) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9) -- _PyEval_EvalCodeWithName

( 65) /home/lester/sire.app/bin/python ([0x55772b73b3bc] ++0x3dc) -- _PyFunction_FastCallDict

( 66) /home/lester/sire.app/bin/python ([0x55772b752473] ++0x63) -- _PyObject_Call_Prepend

( 67) /home/lester/sire.app/bin/python ([0x55772b74719e] ++0x6e) -- PyObject_Call

( 68) /home/lester/sire.app/bin/python ([0x55772b7f69ff] ++0x1f3f) -- _PyEval_EvalFrameDefault

( 69) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9) -- _PyEval_EvalCodeWithName

( 70) /home/lester/sire.app/bin/python ([0x55772b79e157] ++0x387) -- _PyFunction_FastCallKeywords

( 71) /home/lester/sire.app/bin/python ([0x55772b7f51a6] ++0x6e6) -- _PyEval_EvalFrameDefault

( 72) /home/lester/sire.app/bin/python ([0x55772b73a0d9] ++0x2f9) -- _PyEval_EvalCodeWithName

( 73) /home/lester/sire.app/bin/python ([0x55772b73afa4] ++0x44) -- PyEval_EvalCodeEx

( 74) /home/lester/sire.app/bin/python ([0x55772b73afcc] ++0x1c) -- PyEval_EvalCode

/home/lester/sire.app/bin/python(+0x22f664) [0x55772b854664] ( 76) /home/lester/sire.app/bin/python ([0x55772b85e881] ++0xa1) -- PyRun_FileExFlags

( 77) /home/lester/sire.app/bin/python ([0x55772b85ea73] ++0x1c3) -- PyRun_SimpleFileExFlags

/home/lester/sire.app/bin/python(+0x23ab2f) [0x55772b85fb2f] ( 79) /home/lester/sire.app/bin/python ([0x55772b85fc4c] ++0x3c) -- _Py_UnixMain

( 80) /usr/lib/libc.so.6 ([0x7f179026bee3] ++0xf3) -- __libc_start_main

/home/lester/sire.app/bin/python(+0x1df982) [0x55772b804982] EndTrace Exception 'SireError::program_bug' thrown by the thread 'master'. There is a problem - the parent forcefield of molecules:molecules (7) is null! Thrown from FILE: /home/lester/Code/Sire/corelib/src/libs/SireFF/ffmolgroup.cpp, LINE: 144, FUNCTION: void SireFF::detail::FFMolGroupPvt::assertNotNull() const I think it is easier for me to just modify OpenMMMD.py so that it writes a separate DCD file for coordinates only. This is trivial to do and I've tested that it works properly, unlike the restart file above.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

[lohedges]

It seems to restart fine, although I get the same error even after reproducing all of the imports from OpenMMMD.py. I think I'll stick with my single frame DCD for now.

[jmichel80]

Ok. Could you post an issue on the Sire git page? This sounds like a bug in SireFF, maybe @chryswoods can advise.

Julien

Sent from my iPhone

On 25 Jul 2019, at 16:38, Lester Hedges notifications@github.com wrote:

It seems to restart fine, although I get the same error even after reproducing all of the imports from OpenMMMD.py. I think I'll stick with my single frame DCD for now.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

[lohedges]

@dlukauskis I've updated OpenMMMD.py in Sire so that it writes an AMBER RST coordinate file each SOMD cycle. This can be used to get the latest coordinates and generate the updated system, rather than having to load the whole trajectory.

Sire is currently rebuilding, which will take an hour or so. After that I'll push another BioSimSpace commit to generate an updated Conda package. This should be available around lunch time.

[lohedges]

A new BioSimSpace package won't be ready for a while. There's a corrupt package in the Conda archives that is causing the Linux Sire build to fail, which means I can't rebuild BioSimSpace on top of the Sire updates.

(Does a week go by without yet another ridiculous Conda issue!)

[dlukauskis]

No worries, keep me updated.

[lohedges]

This should now be fixed in the latest Conda package. When you call getSystem on a BioSimSpace.Process.Somd object it reads the configuration from a single restart file (not a trajectory) and updates coordinates in the original system. Hopefully you'll no longer have issues with memory. Because of this change, trajectory frames aren't saved so frequently so we also end up with much smaller file sizes.

I'll work on trying to implement something similar for GROMACS too, since that also only writes trajectory data.

[lohedges]

Actually, hold fire. It looks like the latest Linux build failed, because of (you guessed it!) yet another Conda issue.

[lohedges]

Okay, we now have a working Conda package. (Conda 4.7 breaks the RDKit install within our Conda environment. Rolling back to 4.6.14 fixes things.)

I've also updated BioSimSpace.Process.Gromacs to make use of gmx trjconv to extract single trajectory frames when calling getSystem.

[dlukauskis]

I've tried testing the equilibration feature, and I'm getting:

TypeError: 'system' must be of type 'BioSimSpace._SireWrappers.System'

Tested on 2019.1.0+220.gfcc2c13.

[lohedges]

In future, could you please post the full error trace. This should include additional context about the error, such as the file and line number, which are super helpful for debugging.

What I think has happened is that your BioSimSpace Conda package has updated, but Sire hasn't updated too. (This contains the OpenMMMD.py fix to save restarts every cycle.) Could you try installing into a clean Conda environment, rather than doing an upgrade. (You could also try upgrading Sire separately.) This shouldn't have happened, but it can be difficult to get Conda's dependency pinning rules to work correctly, particularly when packages, e.g. Sire, don't use semantic versioning.

[dlukauskis]

I've done a clean conda env reinstall of BSS (version is same as mentioned before) and Sire (2019.1.0=py3.7hf484d3e_23). Same issue persists. Sorry I didn't attach the full trace, here it is. It looks like BSS doesn't let the process run and finish before getting the system.

[lohedges]

I've figured it out. The getSystem method in BioSimSpace.Process.Somd was missing the code needed to block until the process finishes before trying to get the system. (All other methods were fine, so I'm not sure how I missed it there.) I've updated the code and a new package should be available in an hour or so.

[lohedges]

I've just run your example using the updated code and the equilibration fails with the following error:

Traceback (most recent call last):
  File "/home/lester/sire.app/share/Sire/scripts/somd.py", line 125, in <module>
    OpenMMMD.run(params)
  File "/home/lester/sire.app/lib/python3.7/site-packages/Sire/Tools/__init__.py", line 176, in inner
    retval = func()
  File "/home/lester/sire.app/lib/python3.7/site-packages/Sire/Tools/OpenMMMD.py", line 1487, in run
    system = moves.move(system, nmoves.val, True)
RuntimeError: Particle coordinate is nan

At this stage, I think it's probably best if you take a closer look at the merged system for this ligand pair to see if there are any issues with the mapping.

[dlukauskis]

Yeah I'm getting the same error. Thanks.

[lohedges]

Same issue if minimising/equilibrating with GROMACS:

fatal error:
Step 0: The total potential energy is 4.8604e+14, which is extremely high. The
LJ and electrostatic contributions to the energy are 4.8385e+14 and -474590,
respectively. A very high potential energy can be caused by overlapping
interactions in bonded interactions or very large coordinate values. Usually
this is caused by a badly- or non-equilibrated initial configuration,
incorrect interactions or parameters in the topology.

I'd look at your input structure to see if there are any steric clashes.

[dlukauskis]

The issue was poorly parameterised ligand, not assigning a planar dihedral due to a failed conversion from an sdf to a pdb format. That issue is now fixed.