Error using harmonic restratins

[Michael Schauperl]

Hi, 
i have troubles using harmonic restraints with my calculation. 
If I set the flag 'use restraints = True' and try to run the calcuation as described in https://siremol.org/tutorials/somd/Binding_free_energy/FESetup.html

I am getting the following error message:

Traceback (most recent call last):

  File "/home/mschauperl/.conda/envs/sire/share/Sire/scripts/somd-freenrg.py", line 146, in <module>

    OpenMMMD.runFreeNrg(params)

  File "/home/mschauperl/.conda/envs/sire/lib/python3.7/site-packages/Sire/Tools/__init__.py", line 176, in inner

    retval = func()

  File "/home/mschauperl/.conda/envs/sire/lib/python3.7/site-packages/Sire/Tools/OpenMMMD.py", line 1576, in runFreeNrg

    moves = setupMovesFreeEnergy(system, ranseed, gpu.val, lambda_val.val)

  File "/home/mschauperl/.conda/envs/sire/lib/python3.7/site-packages/Sire/Tools/OpenMMMD.py", line 1185, in setupMovesFreeEnergy

    Integrator_OpenMM.initialise()

UserWarning: Exception 'SireBase::missing_property' thrown by the thread 'master:main'.

There is no property with name "d(0)". Available properties are [ AtomNum(7), k(3), x(7), nrestrainedatoms, z(0), k(4), z(3), x(0), x(6), x(3), k(1), AtomNum(6), AtomNum(2), z(1), z(5), x(2), z(2), AtomNum(0), AtomNum(4), k(6), z(6), y(7), y(4), y(1), z(7), AtomNum(1), y(6), k(7), y(3), AtomNum(3), y(0), x(1), y(2), y(5), k(0), x(4), AtomNum(5), z(4), k(2), k(5), x(5) ].

Thrown from FILE: /home/sireuser/Sire/corelib/src/libs/SireBase/properties.cpp, LINE: 469, FUNCTION: const SireBase::Property& SireBase::Properties::operator[](const SireBase::PropertyName&) const

__Backtrace__

.....

Not sure what I am doing wrong.

Best,

Michael

[Julien Michel]

Hi Michael,

Can you post your simulation config file ?

Is this run with Sire 2019 or newer?

Best wishes,

Julien

Sent from my iPhone

On 26 Aug 2020, at 09:37, Michael Schauperl <schauper...@gmail.com> wrote:

Hi, 

--
You received this message because you are subscribed to the Google Groups "Sire Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sire-users+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sire-users/ee34f13b-5a50-4bc5-8859-dcd2c4c219d3n%40googlegroups.com.

[Michael Schauperl]

Hi Julien,

Thanks for the fast response. 

I just used the standard input file from the tutorial and added the line. I tried already multiple different variations of it but always with the same outcome.

#Input file

morphfile = ../../input/MORPH.pert

topfile= ../../input/SYSTEM.top

crdfile= ../../input/SYSTEM.crd

nmoves = 20000

ncycles = 100

buffered coordinates frequency = 5000

save coordinates = True

timestep = 2 * femtosecond

constraint = hbonds-notperturbed

hydrogen mass repartitioning factor = 1.0

cutoff type = cutoffperiodic

cutoff distance = 10*angstrom

barostat = True

andersen = True

energy frequency = 250

precision = mixed

minimise = True

equilibrate = False

equilibration iterations = 5000

center solute = True

reaction field dielectric = 82.0

minimal coordinate saving = True

lambda array =  0.000, 0.125, 0.250, 0.375, 0.500, 0.625, 0.750, 0.875, 1.000

use restraints = True

Sire version: 2019.3.0 

Best,

Michael

[julie...@googlemail.com]

Hi Michael,

Sorry for the slow reply, I am away from work for one more week. It looks like there is something buggy indeed with the positional restraints code in this version of Sire.  Could you elaborate why you wish to restrain atoms during your FEP calculation ? You could alternatively try the keyword ''freeze residues''  if that fits your use case. 

Best wishes,

Julien

[Michael Schauperl]

Hi Julien,

Thanks for the answer. I thought it is a user error :) 

My ligand tends to move out of the pocket, that's why I wanted the restraints. Not sure what is the best solution for this case but I thought restraints might help.

Thanks for your help. 

Michael

[Julien Michel]

Hi Michael,

The distance restraints keyword together with a list of specific distances could help. See line 172-180 here https://github.com/michellab/Sire/blob/devel/wrapper/Tools/OpenMMMD.py for a description of the syntax.

I suggest selecting a protein atom that doesn’t move too much and a ligand atom close to the centre of mass of the ligand. The code works better if the distance is not too large (e.g. 6 Angstrom with a tolerance for deviation of +/-  1 or 2 Angstrom). A force constant of 10 kcal/mol/A**2 is usually sufficient.

BW,

Julien

Sent from my iPhone

You received this message because you are subscribed to a topic in the Google Groups "Sire Users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/sire-users/-MNiK9Tp6-M/unsubscribe.
To unsubscribe from this group and all its topics, send an email to sire-users+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sire-users/4cd154ae-617f-420f-814d-985742718c9an%40googlegroups.com.