Re: [gmx-users] Re: How to set the sigma and epsilon for Cu2+ in OPLS-AA/L

[fantasticqhl]

Dear Dr. Vitaly Chaban,

Thanks very much for your reply. My question is the relationship between 
van der Waals radius and sigma in the OPLS-AA/L force filed files of 
Gromacs.


Of course I did ab initio optimizations of my system, but I do not know 
there is some relation between the optimal bond length (copper--atom of 
the ligand) and sigma.
Could you please be more clear and give a little detailed explanation? 
Thanks very much!


All the best,
Qinghua

On 04/06/2013 06:07 PM, Dr. Vitaly Chaban wrote:

In systems of such kind, everything will depend on the atom of the ligand,
which coordinated by copper ion.

Perform ab initio geometry optimization and find the optimal distance. Then
adjust sigma(s).

Dr. Vitaly Chaban







There is a copper ion with four ligands in my system. I am going to

study this system using MD simulations.
For the vdW parameters, R*=1.74 angstrom and epsilon=1.14 kcal.mol from
one paper will be used in our
simulations. I already found the parameters of copper ion (Cu2+) in the
OPLS-AA/L force field files:
sigma= 2.08470e-01, epsilon=4.76976e+00, which are for Cu2+ without
ligands. The two epsilon are the same,
just with different units.

My question is that I do not know how to convert the vdW radius to
sigma. I found that the vdw radius of copper is
1.4 angstrom, and the sigma in the force field file is 2.08470e-01.
Could someone tell me how to do the converting?

Thanks very much!




--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] Re: How to set the sigma and epsilon for Cu2+ in OPLS-AA/L

[Dr. Vitaly Chaban]

Dear Qinghua -

The formal relation is diameter = pow (2, 1/6) * sigma, provided that you
have only LJ potential in your interacting subsystem.

If this is not the case, an optimal sigma can only be found iteratively.


Dr. Vitaly Chaban





On Sun, Apr 7, 2013 at 10:36 AM, fantasticqhl wrote:

> Dear Dr. Vitaly Chaban,
>
> Thanks very much for your reply. My question is the relationship between
> van der Waals radius and sigma in the OPLS-AA/L force filed files of
> Gromacs.
>
> Of course I did ab initio optimizations of my system, but I do not know
> there is some relation between the optimal bond length (copper--atom of the
> ligand) and sigma.
> Could you please be more clear and give a little detailed explanation?
> Thanks very much!
>
> All the best,
> Qinghua
>
> On 04/06/2013 06:07 PM, Dr. Vitaly Chaban wrote:
>
>> In systems of such kind, everything will depend on the atom of the ligand,
>> which coordinated by copper ion.
>>
>> Perform ab initio geometry optimization and find the optimal distance.
>> Then
>> adjust sigma(s).
>>
>> Dr. Vitaly Chaban
>>
>>
>>
>>
>>
>>
>>
>> There is a copper ion with four ligands in my system. I am going to
>>
>>> study this system using MD simulations.
>>> For the vdW parameters, R*=1.74 angstrom and epsilon=1.14 kcal.mol from
>>> one paper will be used in our
>>> simulations. I already found the parameters of copper ion (Cu2+) in the
>>> OPLS-AA/L force field files:
>>> sigma= 2.08470e-01, epsilon=4.76976e+00, which are for Cu2+ without
>>> ligands. The two epsilon are the same,
>>> just with different units.
>>>
>>> My question is that I do not know how to convert the vdW radius to
>>> sigma. I found that the vdw radius of copper is
>>> 1.4 angstrom, and the sigma in the force field file is 2.08470e-01.
>>> Could someone tell me how to do the converting?
>>>
>>> Thanks very much!
>>>
>>>
>>>
>
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

[gmx-users] Re: FEP + Ham. REMD

[Yuri Garmay]

Have anybody been constructing virtual sites for calculating PMF via FEP?


2013/4/6 Yuri Garmay 

> Hi all!
>
> I have a problem. So as it can be seen in the theme I am trying to improve
> sampling of free energy calculating simulation by using replica exchange.
> In the gmx4.6 it is simple while using FEP technique, but does not
> implemented for umbrella sampling. But I want evaluate potential of mean
> force between two flexible molecules, so I have decided to do this with
> FEP+REMD.
>
> It is desirable to select distance between the mass centres as a reaction
> coordinate, but FEP has to do with bonds or constraints between single
> atoms. Luckily there virtual sites are possible to define in gromacs. So I
> have a thought of using it.
>
> But I have never been using virtual sites and so have some questions.
> I  introduced two virtual sites in topology with these lines:
>
> 
> *.top:
> [ atomtypes ]
> ;name mass charge ptype sigma epsilon
> DD 0.0 0.0 V 0 0 ; nb interaction off?
>
> [ implicit_genborn_params ]
> ; Atomtype sar st pi gbr hct
> DD 0.0 0 0.0 0.0 0.0 ; interaction with solvent off?
>
> [ atoms ]
> ...some lines .
>111 DD  9DDD DD111  0
>112 DD 10DDD DD112  0
>
> [ virtual_sitesn ]
> ; Site from funct a d
> 1112   12 ... ; virtual site at mass centre of some groupe
> 1122   82   83 ...
>
> [ bonds ]
> 
>   111   112 6   0.3   10001.3
> 1000 ; harmonic bond
> 
>
> *.gro:
> 9DDD DD  111   0.000   0.000   0.000 /// position have to be
> determined automatically?
>10DDD DD  112   0.000   0.000   0.000
> --
>
> The energy minimization showed correct constraint length at first sight.
> But does this method solve problem I have? I will be very much appreciated
> if someone experienced explains me whether this is correct or there are
> pitfalls at this way.
>
>
> --
> Best regards,
> Yuri
>
>


-- 
Best regards,
Yuri
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] Re: How to set the sigma and epsilon for Cu2+ in OPLS-AA/L

[fantasticqhl]

Dear Dr. Vitaly Chaban,

Thanks for the explanation. I know this equation. However, the van der 
Waals radius and its counterpart sigma in OPLS-AA/L force field files do 
not follow this equation.


For example, the vdw radius of copper ion is 1.4 angstrom, and its sigma 
is  2.08470e-01 (I guess the unit is nm). pow(2, 1/6) is more than 1, so 
obviously this equation
does not work with copper. So do other atoms. I guess that there might 
be an additional coefficient for this equation in gromacs. That's the 
purpose for asking. Thanks very much!



All the best,
Qinghua

On 04/07/2013 10:48 AM, Dr. Vitaly Chaban wrote:

Dear Qinghua -

The formal relation is diameter = pow (2, 1/6) * sigma, provided that 
you have only LJ potential in your interacting subsystem.


If this is not the case, an optimal sigma can only be found iteratively.


Dr. Vitaly Chaban





On Sun, Apr 7, 2013 at 10:36 AM, fantasticqhl > wrote:


Dear Dr. Vitaly Chaban,

Thanks very much for your reply. My question is the relationship
between van der Waals radius and sigma in the OPLS-AA/L force
filed files of Gromacs.

Of course I did ab initio optimizations of my system, but I do not
know there is some relation between the optimal bond length
(copper--atom of the ligand) and sigma.
Could you please be more clear and give a little detailed
explanation? Thanks very much!

All the best,
Qinghua

On 04/06/2013 06:07 PM, Dr. Vitaly Chaban wrote:

In systems of such kind, everything will depend on the atom of
the ligand,
which coordinated by copper ion.

Perform ab initio geometry optimization and find the optimal
distance. Then
adjust sigma(s).

Dr. Vitaly Chaban







There is a copper ion with four ligands in my system. I am
going to

study this system using MD simulations.
For the vdW parameters, R*=1.74 angstrom and epsilon=1.14
kcal.mol from
one paper will be used in our
simulations. I already found the parameters of copper ion
(Cu2+) in the
OPLS-AA/L force field files:
sigma= 2.08470e-01, epsilon=4.76976e+00, which are for
Cu2+ without
ligands. The two epsilon are the same,
just with different units.

My question is that I do not know how to convert the vdW
radius to
sigma. I found that the vdw radius of copper is
1.4 angstrom, and the sigma in the force field file is
2.08470e-01.
Could someone tell me how to do the converting?

Thanks very much!






--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] Re: How to set the sigma and epsilon for Cu2+ in OPLS-AA/L

[Dr. Vitaly Chaban]

The equation is a direct consequence of LJ-12-6 equation. This equation is
used in OPLS and most other force fields.

The difference you found originate from the fact that, besides LJ
potential, there is much stronger Coulomb potential in the copper-ion case.
If you run simulations, you will see that copper-ligand distance is smaller
than the sum of their sigmas multiplied by pow (2, 1/6).


Dr. Vitaly Chaban






On Sun, Apr 7, 2013 at 11:28 AM, fantasticqhl wrote:

>  Dear Dr. Vitaly Chaban,
>
> Thanks for the explanation. I know this equation. However, the van der
> Waals radius and its counterpart sigma in OPLS-AA/L force field files do
> not follow this equation.
>
> For example, the vdw radius of copper ion is 1.4 angstrom, and its sigma
> is  2.08470e-01 (I guess the unit is nm). pow(2, 1/6) is more than 1, so
> obviously this equation
> does not work with copper. So do other atoms. I guess that there might be
> an additional coefficient for this equation in gromacs. That's the purpose
> for asking. Thanks very much!
>
>
> All the best,
> Qinghua
>
> On 04/07/2013 10:48 AM, Dr. Vitaly Chaban wrote:
>
> Dear Qinghua -
>
>  The formal relation is diameter = pow (2, 1/6) * sigma, provided that
> you have only LJ potential in your interacting subsystem.
>
>  If this is not the case, an optimal sigma can only be found iteratively.
>
>
>  Dr. Vitaly Chaban
>
>
>
>
>
> On Sun, Apr 7, 2013 at 10:36 AM, fantasticqhl wrote:
>
>> Dear Dr. Vitaly Chaban,
>>
>> Thanks very much for your reply. My question is the relationship between
>> van der Waals radius and sigma in the OPLS-AA/L force filed files of
>> Gromacs.
>>
>> Of course I did ab initio optimizations of my system, but I do not know
>> there is some relation between the optimal bond length (copper--atom of the
>> ligand) and sigma.
>> Could you please be more clear and give a little detailed explanation?
>> Thanks very much!
>>
>> All the best,
>> Qinghua
>>
>> On 04/06/2013 06:07 PM, Dr. Vitaly Chaban wrote:
>>
>>> In systems of such kind, everything will depend on the atom of the
>>> ligand,
>>> which coordinated by copper ion.
>>>
>>> Perform ab initio geometry optimization and find the optimal distance.
>>> Then
>>> adjust sigma(s).
>>>
>>> Dr. Vitaly Chaban
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>> There is a copper ion with four ligands in my system. I am going to
>>>
 study this system using MD simulations.
 For the vdW parameters, R*=1.74 angstrom and epsilon=1.14 kcal.mol from
 one paper will be used in our
 simulations. I already found the parameters of copper ion (Cu2+) in the
 OPLS-AA/L force field files:
 sigma= 2.08470e-01, epsilon=4.76976e+00, which are for Cu2+ without
 ligands. The two epsilon are the same,
 just with different units.

 My question is that I do not know how to convert the vdW radius to
 sigma. I found that the vdw radius of copper is
 1.4 angstrom, and the sigma in the force field file is 2.08470e-01.
 Could someone tell me how to do the converting?

 Thanks very much!



>>
>
>
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

[gmx-users] Gromacs 4.6 GPU on Nvidia 250 GTS?

[Hovakim Grabski]

Dear Gromacs users,
Is it possible to use GROMACS 4.6 GPU Mode on Nvidia 250 GTS?
Thanks in advance,
Hovakim 
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] Gromacs 4.6 GPU on Nvidia 250 GTS?

[Justin Lemkul]

On Sun, Apr 7, 2013 at 6:39 AM, Hovakim Grabski
wrote:

> Dear Gromacs users,
> Is it possible to use GROMACS 4.6 GPU Mode on Nvidia 250 GTS?
>

No, the minimum compute capability is 2.0.  Your card has 1.1.

https://developer.nvidia.com/cuda-gpus

-Justin

-- 



Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540)
231-9080http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] unstable system

[Justin Lemkul]

On Sun, Apr 7, 2013 at 1:41 AM, Shima Arasteh
wrote:

> Hi all,
>
>  I have a system of peptide/POPC/water/ions. The energy minimization and
> NVT steps has passed successfully. I ran NPT step for around 10 ns with
> restraints of protein and P atoms at first nano seconds and then removing
> them gradually.
> I tried to go on MDRUN. I did not remove restraint of protein atoms
> completely and they are still restrained. When I run the mdrun command, I
> get error of "X particles communicated to PME node Y are more than a cell
> length out of the domain decomposition cell of their charge group" .
> I know this error means an unstable system. When I visualized the written
> pdb files, I see some popc hydrogen atoms are broken and located between
> two leaflets which are separated by a gap. The protein seems ok, however I
> don't get many pdb files to see.
>
>
> As what I see in Diagnosing unstable system web page,
> 1. it would be beneficial if one see what part of the system is unstable
> in first steps. As I saw, the unstable "POPC hydrogen atoms" are not fine.
> 2. The single molecules are supposed to examine in water or vacuum too. I
> have passed this step successfully.
> 3. I have not ignored any warning during the last steps.
> 4. And my mdp files to run md is as follow:
>
> integrator= md
> dt= 0.002
> nsteps= 500
>
>
> ns_type= grid
> nstlist= 5
> rlist= 1.2
> rlistlong   = 1.4
> rcoulomb= 1.2
> rvdw= 1.2
> pbc= xyz
> vdwtype = switch
> rvdw_switch = 0.1
>

What force field are you using? CHARMM? In any case, the value of
rvdw_switch does not make any sense. If you're using CHARMM, it should be
1.0.


> ; Parameters for treating bonded interactions
> continuation= yes
> constraint_algorithm = LINCSNCS / SHAKE)
> constraints= all-bonds)
> lincs_iter= 1
> lincs_order= 4
>
> ; Parameters for treating electrostatic interactions
> coulombtype= PME; Long range electrostatic interactions
> treatment (cut-off, Ewald, PME)
> pme_order= 4; Interpolation order for PME (cubic interpolation
> is represented by 4)
> fourierspacing= 0.16; Maximum grid spacing for FFT grid using
> PME (nm)
>
> ; Temperature coupling parameters
> tcoupl= Nose-Hoover
> tc-grps= Protein_POPC Water_and_ions; Define groups to be
> coupled separately to temperature bath
> tau_t= 0.50.5 ; Group-wise coupling time constant
> (ps)
> ref_t= 310 310; Group-wise reference temperature (K)
>
> ; Pressure coupling parameters
> pcoupl= Parrinello-Rahman
> pcoupltype= semiisotropic
> tau_p= 2.0
> ref_p= 1.01325 1.01325
> compressibility = 4.5e-54.5e-5
>
> ; Miscellaneous control parameters
> ; Dispersion correction
> DispCorr= EnerPres
> ; Initial Velocity Generation
> gen_vel= no
> ; Centre of mass (COM) motion removal relative to the specified groups
> nstcomm= 1   y (steps)
> comm_mode= Linear
> comm_grps=Protein_POPC Water_and_ions; COM removal relative to the
> specified groups
>
>  Would you please let me know if these happen due to an improper
> equilibration? Do I need to extend the NPT step? Would that fix it?
>
>
Aside from the above comment, there is nothing particularly wrong about the
.mdp file aside from some odd characters here and there, which I will
assume are nothing more than quirks of transferring to an email, as they
otherwise would have triggered fatal errors in grompp.

-Justin

-- 



Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540)
231-9080http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

[gmx-users] (4/7/2013 3:33:29 PM)

[Acoot Brett]

  
http://www.rikgim.com/sq/ebetmmqleqljqo.ff  


 




















Acoot Brett
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] fail to pull

[Albert]

On 04/06/2013 08:52 PM, Justin Lemkul wrote:

Hard to tell. Does your ligand have a suitable exit pathway exactly aligned
along the x-axis? Have you tried increasing the pull rate? How long is the
simulation? I don't even see nsteps in the above .mdp file. How about
increasing the force constant? Is the vector connecting the COM of the
entire protein and the COM of the ligand suitable for describing the exit
pathway?

-Justin


Hello Justin:

 thanks a lot for kind rely.
 Yes, I adjust the conformation of whole protein/ligand so that it can 
exist from X-axies.  I only show part of the .mdp file so some of then 
are not shown.



; Run parameters
integrator  = md
dt  = 0.002
tinit   = 0
nsteps  = 50; 500 ps
nstcomm = 10
; Output parameters
nstxout = 5000  ; every 10 ps
nstvout = 5000
nstfout = 1000
nstxtcout   = 1000  ; every 1 ps
nstenergy   = 1000


probably I should consider use part of the protein such as residues 
around binding pocket as COM for protein instead of whole protein? I 
applied for 1ns with rate pull_rate1= 0.001, so at then end of pulling, 
the distance for COMprotein and COM ligand should be 10A. Probably this 
is too short for whole protein as COM?


thank you very much
best
Albert



--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] fail to pull

[Justin Lemkul]

On Sun, Apr 7, 2013 at 11:17 AM, Albert  wrote:

> On 04/06/2013 08:52 PM, Justin Lemkul wrote:
>
>> Hard to tell. Does your ligand have a suitable exit pathway exactly
>> aligned
>> along the x-axis? Have you tried increasing the pull rate? How long is the
>> simulation? I don't even see nsteps in the above .mdp file. How about
>> increasing the force constant? Is the vector connecting the COM of the
>> entire protein and the COM of the ligand suitable for describing the exit
>> pathway?
>>
>> -Justin
>>
>
> Hello Justin:
>
>  thanks a lot for kind rely.
>  Yes, I adjust the conformation of whole protein/ligand so that it can
> exist from X-axies.  I only show part of the .mdp file so some of then are
> not shown.
>
>
> ; Run parameters
> integrator  = md
> dt  = 0.002
> tinit   = 0
> nsteps  = 50; 500 ps
> nstcomm = 10
> ; Output parameters
> nstxout = 5000  ; every 10 ps
> nstvout = 5000
> nstfout = 1000
> nstxtcout   = 1000  ; every 1 ps
> nstenergy   = 1000
>
>
> probably I should consider use part of the protein such as residues around
> binding pocket as COM for protein instead of whole protein? I applied for
> 1ns with rate pull_rate1= 0.001, so at then end of pulling, the distance
> for COMprotein and COM ligand should be 10A. Probably this is too short for
> whole protein as COM?
>
>
Let me clear one thing up first. 1 ns of pulling with a 0.001 nm/ps pull
rate will not necessarily cause the ligand to be displaced by 1 nm. The
particle pulling the virtual spring will be displaced by 1 nm, but the
ligand will only move as a function of this applied force and the restoring
forces (i.e. interactions between the ligand and protein).

Choosing a more suitable reference group and running the simulation for
longer will produce the desired result.

-Justin

-- 



Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540)
231-9080http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] fail to pull

[Albert]

hello Justine:

 thanks a lot for such kind comments.

 I may find the problem, the COM probably is not suitable to use whole 
protein since I found by g_dist that all the distance are between 
0.9-1.0 nm throughout the whole pulling process.


 BTW, I notice that we are using command:

grompp -f md_pull.mdp -c npt.gro -p topol.top -n index.ndx -t npt.cpt -o 
pull.tpr


for pulling .tpr generating. May I ask whey shall we use option

"-t npt.cpt"

in  this  step? Usually we only need to specify .mdp, .gro and .top file 
for .tpr file



thank you very much
best
Albert

On 04/07/2013 05:31 PM, Justin Lemkul wrote:

Let me clear one thing up first. 1 ns of pulling with a 0.001 nm/ps pull
rate will not necessarily cause the ligand to be displaced by 1 nm. The
particle pulling the virtual spring will be displaced by 1 nm, but the
ligand will only move as a function of this applied force and the restoring
forces (i.e. interactions between the ligand and protein).

Choosing a more suitable reference group and running the simulation for
longer will produce the desired result.

-Justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] fail to pull

[Justin Lemkul]

On Sun, Apr 7, 2013 at 12:05 PM, Albert  wrote:

> hello Justine:
>
>  thanks a lot for such kind comments.
>
>  I may find the problem, the COM probably is not suitable to use whole
> protein since I found by g_dist that all the distance are between 0.9-1.0
> nm throughout the whole pulling process.
>
>  BTW, I notice that we are using command:
>
> grompp -f md_pull.mdp -c npt.gro -p topol.top -n index.ndx -t npt.cpt -o
> pull.tpr
>
>
> for pulling .tpr generating. May I ask whey shall we use option
>
> "-t npt.cpt"
>
> in  this  step? Usually we only need to specify .mdp, .gro and .top file
> for .tpr file
>
>
I'm assuming you're getting that line from my tutorial. You pass the .cpt
file to grompp to preserve velocities from the previous equilibration
phase. If you don't, what was the point of equilibrating? Coordinates,
topology, and .mdp parameters are all that are strictly required to produce
a .tpr file, but other files may be necessary to produce a sensible .tpr
file based on previous steps.

-Justin

-- 



Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540)
231-9080http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] fail to pull

[Albert]

IC.

thanks a lot for explanations.

Albert


On 04/07/2013 06:08 PM, Justin Lemkul wrote:

I'm assuming you're getting that line from my tutorial. You pass the .cpt
file to grompp to preserve velocities from the previous equilibration
phase. If you don't, what was the point of equilibrating? Coordinates,
topology, and .mdp parameters are all that are strictly required to produce
a .tpr file, but other files may be necessary to produce a sensible .tpr
file based on previous steps.

-Justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] unstable system

[Shima Arasteh]

Yes, some letters were transferred by mistakes to email.
Yes, I am using charmm36. Also the rvdw_switch was correct the same as you 
wrote 1.0  in last email.

Now, I am wondering what the solution is in such cases of facing such an error? 
Would you please let me know your suggestions? 

Many thanks for your replies.
   

Sincerely,
Shima



 From: Justin Lemkul 
To: Shima Arasteh ; Discussion list for GROMACS 
users  
Sent: Sunday, April 7, 2013 4:46 PM
Subject: Re: [gmx-users] unstable system
 






On Sun, Apr 7, 2013 at 1:41 AM, Shima Arasteh  
wrote:

Hi all,
>
> I have a system of peptide/POPC/water/ions. The energy minimization and NVT 
>steps has passed successfully. I ran NPT step for around 10 ns with restraints 
>of protein and P atoms at first nano seconds and then removing them gradually.
>I tried to go on MDRUN. I did not remove restraint of protein atoms completely 
>and they are still restrained. When I run the mdrun command, I get error of "X 
>particles communicated to PME node Y are more than a cell length out of the 
>domain decomposition cell of their charge group" .
>I know this error means an unstable system. When I visualized the written pdb 
>files, I see some popc hydrogen atoms are broken and located between two 
>leaflets which are separated by a gap. The protein seems ok, however I  don't 
>get many pdb files to see.
>
>
>As what I see in Diagnosing unstable system web page,
>1. it would be beneficial if one see what part of the system is unstable in 
>first steps. As I saw, the unstable "POPC hydrogen atoms" are not fine.
>2. The single molecules are supposed to examine in water or vacuum too. I have 
>passed this step successfully.
>3. I have not ignored any warning during the last steps.
>4. And my mdp files to run md is as follow:
>
>integrator    = md       
>dt        = 0.002     
>nsteps        = 500   
>
>
>ns_type        = grid       
>nstlist        = 5       
>rlist        = 1.2       
>rlistlong   = 1.4
>rcoulomb    = 1.2       
>rvdw        = 1.2       
>pbc        = xyz       
>vdwtype = switch
>rvdw_switch = 0.1
>

What force field are you using? CHARMM? In any case, the value of rvdw_switch 
does not make any sense. If you're using CHARMM, it should be 1.0.
 
; Parameters for treating bonded interactions
>continuation    = yes      
>constraint_algorithm = LINCS    NCS / SHAKE)
>constraints    = all-bonds    )
>lincs_iter    = 1       
>lincs_order    = 4      
>
>; Parameters for treating electrostatic interactions
>coulombtype    = PME        ; Long range electrostatic interactions treatment 
>(cut-off, Ewald, PME)
>pme_order    = 4        ; Interpolation order for PME (cubic interpolation is 
>represented by 4)
>fourierspacing    = 0.16        ; Maximum grid spacing for FFT grid using PME 
>(nm)
>
>; Temperature coupling parameters
>tcoupl        = Nose-Hoover       
>tc-grps        = Protein_POPC Water_and_ions        ; Define groups to be 
>coupled separately to temperature bath
>tau_t        = 0.5    0.5         ; Group-wise coupling time constant (ps)
>ref_t        = 310     310        ; Group-wise reference temperature (K)
>
>; Pressure coupling parameters
>pcoupl        = Parrinello-Rahman     
>pcoupltype    = semiisotropic           
>tau_p        = 2.0               
>ref_p        = 1.01325 1.01325       
>compressibility = 4.5e-5    4.5e-5
>
>; Miscellaneous control parameters
>; Dispersion correction
>DispCorr    = EnerPres      
>; Initial Velocity Generation
>gen_vel        = no         
>; Centre of mass (COM) motion removal relative to the specified groups
>nstcomm        = 1           y (steps)
>comm_mode    = Linear      
>comm_grps    =Protein_POPC Water_and_ions    ; COM removal relative to the 
>specified groups
>
> Would you please let me know if these happen due to an improper 
>equilibration? Do I need to extend the NPT step? Would that fix it?
>
>

Aside from the above comment, there is nothing particularly wrong about the 
.mdp file aside from some odd characters here and there, which I will assume 
are nothing more than quirks of transferring to an email, as they otherwise 
would have triggered fatal errors in grompp.

-Justin
-- 

 Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080 
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin 

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] unstable system

[Justin Lemkul]

On Sun, Apr 7, 2013 at 1:58 PM, Shima Arasteh
wrote:

> Yes, some letters were transferred by mistakes to email.
> Yes, I am using charmm36. Also the rvdw_switch was correct the same as you
> wrote 1.0  in last email.
>
>
I would suggest you take great care to provide accurate information,
especially when it is requested of you. I thought I had identified a
mistake that would have lead to instability. It turns out that you posted
something fictitious and it feels a bit like I wasted my time. If you want
free help, make it easy for those who want to help you to do so.


> Now, I am wondering what the solution is in such cases of facing such an
> error? Would you please let me know your suggestions?
>
>
I see no reason why a system that was previously stable for 10 ns to
suddenly become unstable. Monitor energy terms very closely, and try
reducing the time step to 1 fs for a period of time.

-Justin


> Many thanks for your replies.
>
>
> Sincerely,
> Shima
>
>
> 
>  From: Justin Lemkul 
> To: Shima Arasteh ; Discussion list for
> GROMACS users 
> Sent: Sunday, April 7, 2013 4:46 PM
> Subject: Re: [gmx-users] unstable system
>
>
>
>
>
>
>
> On Sun, Apr 7, 2013 at 1:41 AM, Shima Arasteh 
> wrote:
>
> Hi all,
> >
> > I have a system of peptide/POPC/water/ions. The energy minimization and
> NVT steps has passed successfully. I ran NPT step for around 10 ns with
> restraints of protein and P atoms at first nano seconds and then removing
> them gradually.
> >I tried to go on MDRUN. I did not remove restraint of protein atoms
> completely and they are still restrained. When I run the mdrun command, I
> get error of "X particles communicated to PME node Y are more than a cell
> length out of the domain decomposition cell of their charge group" .
> >I know this error means an unstable system. When I visualized the written
> pdb files, I see some popc hydrogen atoms are broken and located between
> two leaflets which are separated by a gap. The protein seems ok, however I
> don't get many pdb files to see.
> >
> >
> >As what I see in Diagnosing unstable system web page,
> >1. it would be beneficial if one see what part of the system is unstable
> in first steps. As I saw, the unstable "POPC hydrogen atoms" are not fine.
> >2. The single molecules are supposed to examine in water or vacuum too. I
> have passed this step successfully.
> >3. I have not ignored any warning during the last steps.
> >4. And my mdp files to run md is as follow:
> >
> >integrator= md
> >dt= 0.002
> >nsteps= 500
> >
> >
> >ns_type= grid
> >nstlist= 5
> >rlist= 1.2
> >rlistlong   = 1.4
> >rcoulomb= 1.2
> >rvdw= 1.2
> >pbc= xyz
> >vdwtype = switch
> >rvdw_switch = 0.1
> >
>
> What force field are you using? CHARMM? In any case, the value of
> rvdw_switch does not make any sense. If you're using CHARMM, it should be
> 1.0.
>
> ; Parameters for treating bonded interactions
> >continuation= yes
> >constraint_algorithm = LINCSNCS / SHAKE)
> >constraints= all-bonds)
> >lincs_iter= 1
> >lincs_order= 4
> >
> >; Parameters for treating electrostatic interactions
> >coulombtype= PME; Long range electrostatic interactions
> treatment (cut-off, Ewald, PME)
> >pme_order= 4; Interpolation order for PME (cubic
> interpolation is represented by 4)
> >fourierspacing= 0.16; Maximum grid spacing for FFT grid using
> PME (nm)
> >
> >; Temperature coupling parameters
> >tcoupl= Nose-Hoover
> >tc-grps= Protein_POPC Water_and_ions; Define groups to be
> coupled separately to temperature bath
> >tau_t= 0.50.5 ; Group-wise coupling time constant
> (ps)
> >ref_t= 310 310; Group-wise reference temperature (K)
> >
> >; Pressure coupling parameters
> >pcoupl= Parrinello-Rahman
> >pcoupltype= semiisotropic
> >tau_p= 2.0
> >ref_p= 1.01325 1.01325
> >compressibility = 4.5e-54.5e-5
> >
> >; Miscellaneous control parameters
> >; Dispersion correction
> >DispCorr= EnerPres
> >; Initial Velocity Generation
> >gen_vel= no
> >; Centre of mass (COM) motion removal relative to the specified groups
> >nstcomm= 1   y (steps)
> >comm_mode= Linear
> >comm_grps=Protein_POPC Water_and_ions; COM removal relative to
> the specified groups
> >
> > Would you please let me know if these happen due to an improper
> equilibration? Do I need to extend the NPT step? Would that fix it?
> >
> >
>
> Aside from the above comment, there is nothing particularly wrong about
> the .mdp file aside from some odd characters here and there, which I will
> assume are nothing more than quirks of transferring to an email, as they
> otherwise would have triggered fatal errors in grompp.
>
> -Justin
> --
>
>  Justin A. Lemkul, Ph.D.
> Research Scientist
> Departmen

[gmx-users] Re: fail to pull

[Thomas Schlesier]

Must get the bus, so only short answer.
You could try to use constraint pulling instead of an umbrella 
potential. Then the ligand should move 1nm in 1ns. And you could see is 
the setup is ok, or if it would be better to pull into another direction...

greetings
thomas


Am 07.04.2013 18:05, schrieb gmx-users-requ...@gromacs.org:

Hard to tell. Does your ligand have a suitable exit pathway exactly
>>aligned
>>along the x-axis? Have you tried increasing the pull rate? How long is the
>>simulation? I don't even see nsteps in the above .mdp file. How about
>>increasing the force constant? Is the vector connecting the COM of the
>>entire protein and the COM of the ligand suitable for describing the exit
>>pathway?
>>
>>-Justin
>>

>
>Hello Justin:
>
>  thanks a lot for kind rely.
>  Yes, I adjust the conformation of whole protein/ligand so that it can
>exist from X-axies.  I only show part of the .mdp file so some of then are
>not shown.
>
>
>; Run parameters
>integrator  = md
>dt  = 0.002
>tinit   = 0
>nsteps  = 50; 500 ps
>nstcomm = 10
>; Output parameters
>nstxout = 5000  ; every 10 ps
>nstvout = 5000
>nstfout = 1000
>nstxtcout   = 1000  ; every 1 ps
>nstenergy   = 1000
>
>
>probably I should consider use part of the protein such as residues around
>binding pocket as COM for protein instead of whole protein? I applied for
>1ns with rate pull_rate1= 0.001, so at then end of pulling, the distance
>for COMprotein and COM ligand should be 10A. Probably this is too short for
>whole protein as COM?
>
>

Let me clear one thing up first. 1 ns of pulling with a 0.001 nm/ps pull
rate will not necessarily cause the ligand to be displaced by 1 nm. The
particle pulling the virtual spring will be displaced by 1 nm, but the
ligand will only move as a function of this applied force and the restoring
forces (i.e. interactions between the ligand and protein).

Choosing a more suitable reference group and running the simulation for
longer will produce the desired result.

-Justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

[gmx-users] Simulating a large system

[Juan Antonio Raygoza Garay]

Hi, is it possible to instruct gromacs to only perform the dynamics on half of 
the system or protein while ignoring the rest?

thanks--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] Simulating a large system

[Justin Lemkul]

On Sun, Apr 7, 2013 at 9:55 PM, Juan Antonio Raygoza Garay  wrote:

> Hi, is it possible to instruct gromacs to only perform the dynamics on
> half of the system or protein while ignoring the rest?
>

You can use freezegrps to immobilize whatever you like, but you won't gain
any performance by doing so.

-Justin

-- 



Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540)
231-9080http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] fail to pull

[Jianguo Li]

You switched on the position restraint in your mdp file, is that the reason?



 From: Albert 
To: gromacs maillist  
Sent: Sunday, 7 April 2013, 2:47
Subject: [gmx-users] fail to pull
 
Dear:

I am trying to pull my ligand outside of the binding pocket with following 
configurations:

title       = Umbrella pulling simulation
define      = -DPOSRES
; Pull code
pull            = umbrella
pull_geometry   = distance  ; simple distance increase
pull_dim        = Y N N
pull_start      = yes       ; define initial COM distance > 0
pull_ngroups    = 1
pull_group0     = Protein
pull_group1     = LIG
pull_rate1      = 0.001      ; 0.01 nm per ps = 10 nm per ns
pull_k1         = 1000      ; kJ mol^-1 nm^-2

Tcoupl      = v-rescale
tc_grps     = Protein_LIG   Water_and_ions
tau_t       = 0.5       0.5
ref_t       = 310       310
; Pressure coupling is on
Pcoupl          = Parrinello-Rahman
pcoupltype      = isotropic
tau_p           = 1.0   1.0
compressibility = 4.5e-5
ref_p           = 1.0     1.0
refcoord_scaling = com
; Generate velocities is off
gen_vel     = no
; Periodic boundary conditions are on in all directions
pbc     = xyz
; Long-range dispersion correction
DispCorr    = EnerPres



It is quite strange, the ligand is still in place and not outside the pocket at 
the end of simulations. I am just wondering where is the problem?

thank you very much
best
Albert

-- gmx-users mailing list    gmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the www interface 
or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

[gmx-users] GROMACS 4.6v - Myrinet2000

[Hrachya Astsatryan]

Dear all,

We have installed the latest version of Gromacs (version 4.6) on our 
cluster by the following step:


 * cmake .. -DGMX_MPI=ON -DCMAKE_INSTALL_PREFIX=/backup/sicnas/gromacs
   -DGMX_BUILD_OWN_FFTW=ON -DGMX_PREFER_STATIC_LIBS=ON

The interconnection of the cluster is Myrinet2000 and the driver is MX 
(and mpich 1.3.1).


The job scripts (see below) and different simulations show that there is 
no adequate performance (usually it should be 2-3 time more).


Could you, pls. help us to overcome the problem in order to get better 
performance?


With regards,
Hrach

#PBS -l nodes=8:ppn=2
##PBS -l walltime=360:00
#PBS -q armcluster
#PBS -e 33mM_16_new.err
#PBS -o 33mM_16_new.log
## Specify the shell to be bash
#PBS -S /bin/bash

export 
LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/backup/sicnas/gromacs/lib:/opt/mpi/lib


/opt/mpi/bin/mpirun -machinefile /backup/sicnas/npt101/machine -np 16 
/backup/sicnas/gromacs/bin/mdrun_mpi -s /backup/sicnas/npt101/topol.tpr -v






--
Hrachya Astsatryan
Head of HPC Laboratory,
Institute for Informatics and Automation Problems,
National Academy of Sciences of the Republic of Armenia
1, P. Sevak str., Yerevan 0014, Armenia
t: 374 10 284780
f: 374 10 285812
e: hr...@sci.am
skype: tighra
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists