Re: [gmx-users] Conceptual question: Why are MD production simulations usually run in the NVT ensemble?

2012-07-14 Thread Justin A. Lemkul



On 7/14/12 11:06 AM, Andrew DeYoung wrote:

Hi,

I'm not sure if this is an appropriate question for the Gromacs users'
mailing list.  If it isn't, please forgive me and disregard this message.

I asked this question on Physics Stack Exchange: Why is the canonical (NVT)
ensemble often used for molecular dynamics (MD) simulations?  My question
and an answer are posted here:

http://physics.stackexchange.com/questions/31997/why-is-the-canonical-nvt-en
semble-often-used-for-classical-molecular-dynam

My question is, do you agree with the answer?  I'm not sure that I do.  The
responder seems to imply that it is practically impossible to use periodic
boundary conditions in an NPT simulation.  But, I think that the algorithms
in Gromacs do just this; am I correct in this belief?



I also disagree with the answer you got, on several grounds.  I will admit that 
I have never gone into the code for how boxes are scaled and how this is all 
done in terms of bookkeeping, but I think there is plenty of evidence that we 
have decent barostat algorithms.  The response you got seemed to just say "it's 
all too inconvenient to try to think about, so we shouldn't do it."  Part of the 
response was simply that if you increase the system size, the ensembles become 
the same, which may be true in a theoretical sense (since P fluctuations 
decrease proportionally with system size) but it becomes totally impractical for 
doing simulations.  How then, does one account for heterogeneous systems like 
I've got with a membrane, protein, water, and maybe small molecules?  I can't 
just keep scaling that up and assume that the relative majority of water 
molecules is going to make everything happy.  That, and the lipid model was 
parameterized for NPT...



People do often run _equilibration_ simulations in the NPT ensemble, I
think.  They do this usually, I think, to obtain the proper density -- the
code changes the box dimensions until the proper average pressure is
reached.  Then, for _production_ simulations, people usually use the NVT
ensemble, where the dimensions of the simulation box are held fixed.  My
question here is, why is the NVT ensemble, rather than the NPT ensemble,
typically used for production MD simulations?



I disagree with this premise as well.  The vast majority of all biological 
simulations I see are done under NPT conditions.  I think the heart of the issue 
is the discipline you're studying and what is commonly done.  Maybe in pure 
physics or chemistry simulations, NVT is more common because the experiments are 
done that way.  For us in biological science, the NVT ensemble is less relevant 
than NPT because we are studying processes that go on in cells, within 
membranes, or even just in solution an open-topped test tube or flask that is 
under atmospheric pressure.  The general advice from many on this list has 
always been, "pick an ensemble that models the reality of your system best."  If 
it makes sense to switch back to NVT from NPT for whatever you're doing, by all 
means do it.



Do you have any suggestions for helpful reading that I can do on this topic?
Thank you for your time!



I think it's a philosophical discussion, but comes down to what it is you're 
modeling and what the most appropriate conditions are.  Knowing the ins and outs 
of thermostat and barostat algorithms would be very important in debating these 
issues.


-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
* Only plain text messages are allowed!
* 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] (no subject)

2012-07-14 Thread Justin A. Lemkul



On 7/14/12 3:34 AM, sara elham wrote:

Dear Gromacs users

I want to find the solution for my problem from mailing list, but it
give me this massage:
"Timeout expired. The timeout period elapsed prior to completion of
the operation or the server is not responding."


The website experiences intermittent issues.  Hopefully they will be sorted out 
soon.



I tried to registration again, but in the section " Type the
characters you see in the image below " , it isn't shown anything!!!



You don't need to register to search the mailing list archive.

-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
* Only plain text messages are allowed!
* 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: Melting simulations - regd

2012-07-14 Thread Justin A. Lemkul



On 7/14/12 1:02 AM, rapolu sukumar wrote:

Dear Justin,

Thank you for your reply,  here I have one more doubt, I am
using 'pcoupltype= anisotropic' in my simulationa as the polymers are
aniostropic, I have searched for compressibility values of system but i
couldn't find is there any generalized way to set up compressibility values.




Most simulations probably just use the compressibility of water.  As for the 
accuracy of your simulation if you do this, I cannot say.


-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
* Only plain text messages are allowed!
* 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] g_saltbr speed

2012-07-13 Thread Justin A. Lemkul



On 7/13/12 12:05 PM, Kavyashree M wrote:

Its 50ns 25000 frames. the xtc file is 695MB. it has 16GB RAM. So will
that be insufficient? I have previously run other analysis which used to take
huge memory, for eg. covariance analysis, in a system with much lesser memory
even though CPU usage was low the job used to finish. But in this case its not
so.



I can't see a reason why the file itself would present a problem.  I have run 
g_saltbr on similarly sized systems.  Sorry, I can't offer an explanation as to 
why it's stopping prematurely.


-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
* Only plain text messages are allowed!
* 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] g_saltbr speed

2012-07-13 Thread Justin A. Lemkul



On 7/13/12 11:50 AM, Kavyashree M wrote:

Dear users,

Its the continuation of the question I asked yesterday, Inorder to reduce
the memory usage during g_saltbr calculations i got the trajectory of only
protein, and tpr file without water and was able to successfully run it. But
unfortunately this again got stopped at 36ns as it had stopped when i was
using the whole trajectory. I tried with -dt 2, still the same problem exists.
Kindly suggest a way out of this situation.



How long is the trajectory?  How many frames?  What is the size of the file on 
disk?  It sounds to me like you're simply exhausting available memory, so the 
only advice is in the link I posted before - use fewer frames or use a machine 
that has more memory to do the analysis.


-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
* Only plain text messages are allowed!
* 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: Error in Membrane simulations with POPC bilayer

2012-07-13 Thread Justin A. Lemkul



On 7/13/12 6:44 AM, J Peterson wrote:

Thanks Justin,

The explanations are very very useful during my course of simulating a
protein with POPC.

I also would like to get explanation on how to simulate a protein which has
only its N-terminal region embedded in the membrane but the rest in solvent.

What is the easy and accurate way to do it?



Position the protein with editconf -center and follow all the other steps in the 
same way.  The dimensions of the box become quite important in such a case. 
There are some general tips on such placement in another of my tutorials:


http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/biphasic/03_tricks.html

-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
* Only plain text messages are allowed!
* 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] Melting simulations - regd

2012-07-13 Thread Justin A. Lemkul



On 7/13/12 2:27 AM, sukumar rapolu wrote:

Hello gmx users,

I am doing melting simulations of a polymer in gromacs, for this am
giving output from a  production run of  particular temperature as the
input for equilibration at  next temperature and generating velocities
in equilibration at each temperature by using gen_vel = yes.
  Here my doubt is, as I am using output from production  as input for
equilibration so do I have to use   continuation= yes in my
equilibration mdp file or as I am generating velocities at each
temperature during equilibration do I have to use   continuation=
no ?



If you're using the structure from a previous run that used constraints, the 
difference between these two settings should be minimal.  But since you're 
basically starting totally new simulations, you can't go wrong with 
"continuation = no."


-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
* Only plain text messages are allowed!
* 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: Error in Membrane simulations with POPC bilayer

2012-07-13 Thread Justin A. Lemkul



On 7/13/12 4:14 AM, J Peterson wrote:

Hi Justin,

I have another doubt on the strong posres that was included in the topology
file. When do we need to remove that position restraint? Does it really
affect at point of time the system?



The strong restraints are only needed for InflateGRO.  During equilibration, 
position restraints with a normal magnitude are appropriate.


-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
* Only plain text messages are allowed!
* 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: Error in Membrane simulations with POPC bilayer

2012-07-13 Thread Justin A. Lemkul



On 7/13/12 1:44 AM, J Peterson wrote:

Thanks for the comment. By the way how to make a bigger box at this time of
the tutorial without affecting any part of the system. Can I use editconf
with slightly bigger number for z-axis (something like 6.7 which was 5.7
before)?



The box should be sufficiently large so as to avoid spurious interactions across 
periodic boundaries.  The size is thus motivated by the length of the cutoffs 
used.  Increasing the box by 1.0 nm might be just enough, but you should 
calculate this very carefully.


-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
* Only plain text messages are allowed!
* 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: Final state not reached in pulling simulation

2012-07-13 Thread Justin A. Lemkul



On 7/13/12 2:32 AM, Neeru Sharma wrote:

Dear Justin,

Thanks for the suggestion regarding the pull force.
1)I will try by reducing the pull force to 1000 and 100 KJmol/nm2 now.
Regarding the force being applied in z-direction, after visualizing my
trajectory i thought it would work by providing force in z-direction.
2)Regarding the protein rotation, it does not rotate over the time. There
are just some localized changes in it.
3)Regarding the distance measurement, I am measuring the distance between
the specific atoms of the residue and the atoms of the GTP. Till now, these
distances as well as the overall protein geometry is maintained well in the
range too.
4)If this kind of pulling does not work out in my case, I will again try it
with using "position" geometry too with pull_vec1.

Can you suggest why have I not reached the final distance at the end of the
simulation? Is it because of the geometry that I have used, because the
force constant is already too high in this case?



I don't know why this is happening.  That's why you need to try more things to 
see if you can root out the issue.  There are a whole host of factors that can 
act against the pulling force.  You're not using any sort of position 
restraints, are you?


The other thing to keep in mind is that if your initial distance is 0.7 nm, and 
you pull for 10 ns at 0.05 nm/ns, you should in theory end up with a distance of 
0.2 nm, not 0.3 nm.  That is, assuming that the pulled group is free to move and 
thus smoothly goes along with the spring.  Other forces within the structure may 
act in opposition.


-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
* Only plain text messages are allowed!
* 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: Error in Membrane simulations with POPC bilayer

2012-07-12 Thread Justin A. Lemkul



On 7/13/12 12:24 AM, J Peterson wrote:

Hi Justin,

I followed your comments and now at the stage of adding solvents.
I wonder to see the protein after shrinking step to have no SOL molecules as
there were SOL molecules in the source popc128b.pdb. Had we removed all the
original SOL molecules anywhere during the course of tutorial?



InflateGRO removes all solvent molecules.


I also see one of the POPC molecules standing out (upside down) of the
bilayer area. How would that be adjusted?

http://gromacs.5086.n6.nabble.com/file/n4999393/POPC..jpg



You need a bigger box.  The lipids are moving into a neighboring periodic cell, 
which in itself is not a problem, but I suspect you'll have inadequate hydration 
and thus you'll be calculating headgroup-headgroup interactions across periodic 
boundaries, which you probably want to avoid.


-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
* Only plain text messages are allowed!
* 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: Error in Membrane simulations with POPC bilayer

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 10:23 PM, J Peterson wrote:

Hi Justin,

Thanks for the suggestion I got it solved somehow. The main problem was in
the popc128b.pdb itself, it has first 64 lipids and half of the SOl
molecules followed by rest of the POPC and SOL molecules. When I rearranged
them the error was solved.

But now another thing I would like to confirm with you. During shrinking
stage, I get the following notes, please tell me are they OK at this stage,

NOTE 1 [file topol.top, line 3070]:
   System has non-zero total charge: 5.00e+00



Assuming this is the net charge on your protein, there's nothing wrong here.


Analysing residue names:
There are:46Protein residues
There are:   126  Other residues
Analysing Protein...
Analysing residues not classified as Protein/DNA/RNA/Water and splitting
into groups...
Number of degrees of freedom in T-Coupling group rest is 21063.00
Largest charge group radii for Van der Waals: 0.249, 0.247 nm
Largest charge group radii for Coulomb:   0.249, 0.247 nm
Calculating fourier grid dimensions for X Y Z
Using a fourier grid of 160x160x48, spacing 0.110 0.114 0.119
Estimate for the relative computational load of the PME mesh part: 0.92

NOTE 2 [file em_st.mdp]:
   The optimal PME mesh load for parallel simulations is below 0.5
   and for highly parallel simulations between 0.25 and 0.33,
   for higher performance, increase the cut-off and the PME grid spacing



PME load for EM runs is usually irrelevant.  During actual MD, you'd strive for 
better balance, but since an inflated lipid system in vacuo is rather weird, you 
can ignore this as well.


-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
* Only plain text messages are allowed!
* 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] Final state not reached in pulling simulation

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 3:43 PM, lloyd riggs wrote:

your pull force looks insanly high especially if your pulling a small piece
of residue?  But for a whole protein of averidge 40 KDa, or 350 amino acids
its around 2000 to 3000 from liturature (only about 6 that I could find
anyways).  I might thus be wrong, but wounder if you have a pull rate and
force what one wins?



The force constant and pull rate do not compete.  The pull rate governs how fast 
the virtual particle is pulled (at constant velocity) and the stiffness of the 
spring dictates the response of the pulled group to which it is attached.


I agree that the force constant seems excessive.  I also wonder why the pulling 
is only being applied in the z-direction.  To the OP - are the GTP and chosen 
residue completely linear in the z-direction?  Does the protein structure rotate 
over time?  If it does, then you're doing nonproductive pulling.  How are you 
measuring the distance (achieved and desired)?  Is it the z-component of the 
distance only, or the complete distance?


I would also note that one-dimensional pulling using the "distance" geometry is 
relatively inflexible.  Using "direction" or "position" geometries in concert 
with pull_vec1 and/or more sophisticated index groups may be more appropriate.


-Justin



 Original-Nachricht 

Datum: Thu, 12 Jul 2012 18:15:53 +0200 Von: Thomas Schlesier
 An: gmx-users@gromacs.org Betreff: [gmx-users] Final
state not reached in pulling simulation



It could be possible tht you do not pull into the 'right' direction. if
there is another group between 'GTP' and 'Residue' you will get clashes and
'Residue' won't move further (could be a water molecule, or some other part
of 'GTP'). If this happens you should observe an increase in the force due
to the umbrella potential. If the problems are due to waters molecules
which block the 'pathway' you could just delete them. If another group is
in the way, you might want to change the pull-vector (and if lucky find the
right one). But don't know what would be the best strategy in this case.
Maybe you can look into what docking-people do, seems to me that your
simulation is related to what they do (but myself has absolute no knowledge
about docking simulations).

greetings thomas



Am 12.07.2012 17:26, schrieb gmx-users-requ...@gromacs.org:

Hello all,

I m performing a pulling simulation on my Protein-Mg-GTP complex. I have
considered pulling between the GTP and a residue of protein. The pull
code in the .mdp file im using is as follows:


; Pull code pull= umbrella pull_geometry   = distance  ;
simple distance increase pull_dim= N N Y pull_start  = yes
; define initial COM distance>  0 pull_ngroups= 1 pull_group0 =
GTP pull_group1 = Residue pull_rate1  = -0.5  ; 0.5
nm per ps = .05 nm per ns pull_k1 = 1  ; kJ mol^-1
nm^-2



The initial distance between GTP and the residue was 7 A and the desired
one was 3A. After the completion of run (10ns), I could get a trajectory
where the final distance was still 4.25 A.

I tried to continue the simulation for another 10ns with the same value
for pull_k1 parameter and one by increasing the value to 100,000 also. In
both of the case, the  trajectories showed the distance stabilized near
_4.25 A only. Can anyone please tell me the reason behind it? What should
I do, so that I could get the desired distance ?

Any suggestion and help is welcome !!!


Thanks,

Neeru Sharma


-- gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users * Only plain text
messages are allowed! * 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


--


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
* Only plain text messages are allowed!
* 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] more than one protonation per residue with pdb2gmx

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 9:00 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Can you please give me an example of such a modeling software?
I tried it with PYMOL but the results are not very good.
And I also found several programs but the are all not free.



http://www.gromacs.org/Documentation/File_Formats/Coordinate_File#Sources

-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
* Only plain text messages are allowed!
* 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] g_saltbr speed

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 8:31 AM, Kavyashree M wrote:

Ok may be i need to specify an index file. I will try that.
And regarding the WARNING: this .tpx file is not fully functional.
I hope it will work fine enough to finish g_saltbr calculation?



In principle, you should be prompted to choose a default group, but you can also 
use a custom index group as needed.


The warning message is intended to note that the .tpr file you produce will not 
likely work for running an actual simulation.  It should be fine for analysis.


-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
* Only plain text messages are allowed!
* 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] g_saltbr speed

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 8:15 AM, Kavyashree M wrote:

Dear Sir,

I had a problem again during g_saltbr calculation it needs
.xtc and .tpr file, I can reduce the .xtc file to have only
protein but .tpr file will have water also in it. inorder to generate
new tpr file without water using grompp, i need topology file
without water .. so do you suggest me to go this way.. it
appears quite complicated!



In fact, it's quite easy with tpbconv.  From tpbconv -h:

"3. by creating a .tpx file for a subset of your original tpx file, which is
useful when you want to remove the solvent from your .tpx file, or when you
want to make e.g. a pure Calpha .tpx file. Note that you may need to use
-nsteps -1 (or similar) to get this to work. WARNING: this .tpx file is not
fully functional."

-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
* Only plain text messages are allowed!
* 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] water bridge

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 8:02 AM, Raj wrote:

Hi all,

I have an protein and drg complex. I've done the MD in SPC water
environment. Now I would like to check number of hydrogen bond formed
between the drug and protein, water as well. What will the right command for
me to analyse. Thanks in advance



g_hbond with suitable index groups.  If you're looking for water-mediated 
hydrogen bonds, there are numerous discussions about ways to accomplish that if 
you search the mailing list archive.


-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
* Only plain text messages are allowed!
* 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] g_saltbr speed

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 6:38 AM, Kavyashree M wrote:

Dear Sir,

That is true as the number of the frames increased the
memory had almost reached 95% but still it has been in
95% since long and CPU usage drops down to 1.5 -2 %
but in many cases i have seen that still it will run (off course
slowly) and finish. But this was too long. So any suggestions?



http://www.gromacs.org/Documentation/How-tos/Reducing_Trajectory_Storage_Volume

-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
* Only plain text messages are allowed!
* 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] more than one protonation per residue with pdb2gmx

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 6:12 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi Justin,
so you mean that I first have to add the capping groups to my structure
and then run pdb2gmx with the -ter function?



Yes.  The capping groups (see your force field's .rtp file for available 
entries) are treated as any other residue.  Gromacs will not magically build 
them; they need to be present.



But how can I add the capping groups to the structure?




You'll have to use some external modeling software for that.

-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
* Only plain text messages are allowed!
* 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: Error in Membrane simulations with POPC bilayer

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 4:59 AM, J Peterson wrote:

Hi Justin,

Thanks for the effort to help me.

I still no out of the error. The following is the content of my
topol_popc.top

; Include chain topologies
#include "gromos53a6_lipid.ff/forcefield.itp"
#include "popc.itp"

; Include water topology
#include "gromos53a6_lipid.ff/spc.itp"

; Include ion topologies
#include "gromos53a6_lipid.ff/ions.itp"

; System specifications
[ system ]
128-Lipid POPC Bilayer

[ molecules ]
; molecule name nr.
POPC 128
SOL 2460

And the following is the first part of my popc128b.gro file

Alm on surf + relaxed popc
14036
 1POP C11   0.253   5.425   1.688
 1POP C22   0.428   5.314   1.792
 1POP C33   0.334   5.243   1.571
 1POP N44   0.378   5.352   1.660
 1POP C55   0.474   5.439   1.590
 1POP C66   0.606   5.390   1.531
 1POP O77   0.692   5.366   1.643
 1POP P88   0.834   5.316   1.587
 1POP O99   0.800   5.197   1.505
 1POPO10   10   0.903   5.435   1.533
 1POPO11   11   0.890   5.266   1.729
 1POPC12   12   0.823   5.142   1.753
 1POPC13   13   0.836   5.086   1.895
 1POPO14   14   0.766   4.964   1.924
 1POPC15   15   0.632   4.959   1.944
 1POPO16   16   0.555   5.019   1.869
 1POPC17   17   0.589   4.833   2.020
 1POPC18   18   0.615   4.703   1.944
 1POPC19   19   0.595   4.575   2.025
 1POPC20   20   0.625   4.467   1.920
 1POPC21   21   0.677   4.340   1.987
 1POPC22   22   0.813   4.357   2.055
 1POPC23   23   0.848   4.223   2.120
 1POPC24   24   0.884   4.121   2.012
 1POPC25   25   0.957   4.013   2.043
 1POPC26   26   0.994   3.986   2.189
 1POPC27   27   1.069   3.853   2.202
 1POPC28   28   0.967   3.739   2.204
 1POPC29   29   1.012   3.593   2.202
 1POPC30   30   0.916   3.478   2.169
 1POPC31   31   0.869   3.518   2.029
 1POPC32   32   0.801   5.205   1.983
.
.
.
.
.

What are the parts mismatching in these files?


From what's shown, it's impossible to tell.  Based on this information, 
everything lines up.  I'm assuming the original coordinate file and the popc.itp 
topology came from Tieleman's site?



Since it says the atom names from the top file will be used, is it safe to
ignore this warning too?



No.  The naming mismatches imply that grompp is trying to map parameters for 
lipids onto water molecules.  The resulting simulation will be nonsense and will 
probably crash immediately.


-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
* Only plain text messages are allowed!
* 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] g_saltbr speed

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 4:51 AM, Kavyashree M wrote:

Dear Gromacs users,

I was running the saltbridge calculations for a dimeric
protein simulation using g_saltbr, But its taking very
long time, almost four days still its not completed.
Could anyone has suggestion regarding this issue? I am
using the same system -
Intel(R) Core(TM) i7-2600 CPU @ 3.40GHz
where i ran MD. Please give some suggestion as to how to
increase the speed of calculation.
Command i issued was:
g_saltbr -f ../traj.xtc -s md.tpr -t 0.4 -sep



g_saltbr calculates properties of all possible ionic pairs in the system, so if 
there are many, the calculation might take a long time.  Four days sounds 
ridiculous, and perhaps the program has frozen by exhausting the available memory.


-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
* Only plain text messages are allowed!
* 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] more than one protonation per residue with pdb2gmx

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 4:37 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

I have another question about the option -ter of the pdb2gmx command.

I choose it because I thought that this is a way that I can determine what
shell happen with the termini but I was not ask anything by the program.

My aim is it to block the termini with a neutral group. Is there a way to
do this with gromacs?



Yes.  You need to build the appropriate groups onto your protein's structure 
using normal capping groups present within the force field (ACE, NME, NH2, etc) 
and choose "None" for the termini when running pdb2gmx so that no additional 
protons are added or removed; the first and last amino acids are then treated as 
internal residues with normal peptide bonds.


-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
* Only plain text messages are allowed!
* 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] pH and protein

2012-07-12 Thread Justin A. Lemkul



On 7/12/12 1:42 AM, tarak karmakar wrote:

Dear All,


I am simulating a protein in gromacs with amber force field. The
protein shows maximum biological activity at pH 5.0 and at pH 7.4 it
shows no activity. So whichever biological process I am going to model
should be at the biologically active pH . So can anyone suggest me

1) how to get the protonation state of each and every residues present
in the protein at that particular pH [i.e. at pH 5.00] ? Should I look
into the pKa values of each and every residues from standard table to
or is there any software to find the same ?


pKa values can vary greatly depending on the local environment of the residue. 
You should find a method to calculate the pKa values.  One possibility:


http://biophysics.cs.vt.edu/


2) how to get the net charge of the protein at that particular pH ? [
Protein Calculator !! ]



The net charge is printed in the topology produced by pdb2gmx.

-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
* Only plain text messages are allowed!
* 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] More accurate potential energy in output?

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 3:14 PM, Markus Kaukonen wrote:

Dear Gromacs,

I'm trying to build a QM/MM with gromacs as MM.
The thing would be run by ASE-simulation environment
(https://wiki.fysik.dtu.dk/ase/).

Question: I'm trying to get the potential energy of a single atomic
configutation.
Is it possible to get the potential energy of a given configuration
with more than 5 digits?



Compile and run in double precision.

-Justin


Both the log file and gmxdump -e give me only 5 digits. The .mdp file
I used is below.

terveisin, Markus

;===
;Gromacs input file
;Created using the Atomic Simulation Environment (ASE)
;===
rlist   = 0.0 ;
nstlist = 1 ;
pbc = no ;
integrator  = md ; md: molecular dynamics(Leapfrog),
rvdw= 0.0 ;
nstlog  = 1 ;
nstenergy   = 1 ;
rcoulomb= 0.0 ;
energygrps  = System ;
ns_type = grid ;
nsteps  = 0 ;
nstfout = 1 ;
define  = -DFLEXIBLE ; flexible/ rigid water






--
====

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
* Only plain text messages are allowed!
* 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] Time for NPT or NVT equilibration

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 2:08 PM, Shima Arasteh wrote:

Dear gmx friends,How much time should have been spent in NPT equilibrium for
a system composed of  protein and water? In Justin's tutorial, I saw it is
100 ps for a system composed of Lysozym and water.



Long enough for the observables of interest to converge.

-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
* Only plain text messages are allowed!
* 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] Forcefield parameters for AcetylCoA

2012-07-11 Thread Justin A. Lemkul


Please keep the discussion on the gmx-users mailing list.

On 7/11/12 10:43 AM, panbazha wrote:

Dear Justin,
That file was generated by me in ATB, The issue is when i looked into
the united atom itp file, I found many bonds and angles has not determined. So,
just wondering if there are any published reports to edit it manually.



A quick literature search should turn something up.  Google returns lots of 
results if you search "Gromos96 acetyl CoA" (without the quotes).


I looked at the missing bonded parameters and they should all be quite easy to 
define based on known values in ffbonded.itp - it is not clear to me why they 
are missing.  Most, if not all, of the missing bonds are between OA and H atoms, 
which are simple alcohol groups (bond type gb_1).


-Justin


Best regards





Quoting "Justin A. Lemkul" :




On 7/11/12 8:18 AM, Padmanabhan Anbazhagan wrote:

Dear All,
Could anyone please suggest me some article or information that contains
force field parameters (gromacs53A6) for acetyl Coenzyme A



Google turned up:

http://compbio.biosci.uq.edu.au/atb/download.py?molid=5470

-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
* Only plain text messages are allowed!
* 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








--
========

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
* Only plain text messages are allowed!
* 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] more than one protonation per residue with pdb2gmx

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 10:24 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,

I wanted to ask if there is a possibility to tell pdb2gmx which residues I
want to be protonated or not.
So that I can for example say HIS at position 29 shell be protonated twice.

I need this for further electrostatic analysis.



Please read pdb2gmx -h.  There are tons of available options to control 
protonation of any and all titratable residues.


-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
* Only plain text messages are allowed!
* 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] Forcefield parameters for AcetylCoA

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 8:18 AM, Padmanabhan Anbazhagan wrote:

Dear All,
Could anyone please suggest me some article or information that contains
force field parameters (gromacs53A6) for acetyl Coenzyme A



Google turned up:

http://compbio.biosci.uq.edu.au/atb/download.py?molid=5470

-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
* Only plain text messages are allowed!
* 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] question about the output of minimization

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 7:34 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Yes I tried the command you wrote but it is still in the corner of the box.
Is there anything else I can do?
The coordinates I used to center the protein are the ones which are in the
last line of the .gro file.



The last line of the .gro file are the box vectors, not the center of the 
system, so you are not, in fact, centering the protein.  You are placing it in 
the corner of the box yourself.  The -center option overrides the -c option.  If 
you want the protein centered, omit -center in your editconf command.


-Justin


Thank you




On 7/11/12 7:12 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi Justin,
ah okey. Thank you.
And I have another question. Why is the protein in the corner of the box
and not in the middle?
I thought I centered it with the command:

editconf -f wholeProtein.gro -o 3m71_box.gro -center 4.59340   4.59470
   5.17330 -c -bt dodecahedron -d 1.0 2>>logErr 1>>logOut

Or not?



Have you tried trjconv -pbc mol -ur compact?  The "center" of an infinite
system
is an arbitrary location; you have to re-wrap the unit cell to achieve the
expected outcome.  It's centered, it just may not look like it in the
current
representation.

-Justin




On 7/11/12 6:13 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
I did a minimization of my structure. But the output seems a bit
strange
for me, since my input was the protein with its membrane in a box like
this:
# = box  P = Protein  M=Membrane


#
#   #
# MM#
#   #
#   #
#   #
#

But the output of the minimization looks like this

#
##
#M M #
 #    #
  ##
   ##
##

How can that be?


It's a triclinic representation of the unit cell.  Nothing is wrong.


My box was created with:

editconf -f wholeProtein.gro -o 3m71_box.gro -center 4.59340   4.59470
5.17330 -c -bt dodecahedron -d 1.0 2>>logErr 1>>logOut



Is a dodecahedral box appropriate for a protein in a membrane?  The
inherent
symmetry of a slab of any sort dictates that you should be using a
cubic
or
rectangular box.  Dodecahedral and octahedral boxes are better for
systems
with
spherical symmetry, like globular proteins in water.



I put solvent in it with:

genbox -cp 3m71_box.gro -cs spc216.gro -p 3m71.top -o 3m71_water.gro
2>>logErr 1>>logOut


The mdp file for the minimization looks like this:

define  = -DPOSRES
integrator  = steep
emtol   = 10
nsteps  = 1500
nstenergy   = 1
energygrps  = System
coulombtype = PME
rcoulomb= 0.9
rvdw= 0.9
rlist   = 0.9
fourierspacing  = 0.12
pme_order   = 4
ewald_rtol  = 1e-5
pbc = xyz


Can you please tell me how I can prevent my box to shift?



There is no box shift; it's just a visual representation of the
dodecahedral
box.  You can "correct" it (for visualization purposes) by running
trjconv
-pbc
mol -ur compact.

-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
* Only plain text messages are allowed!
* 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







--


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
* Only plain text messages are allowed!
* 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







--


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-user

Re: [gmx-users] question about the output of minimization

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 7:12 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi Justin,
ah okey. Thank you.
And I have another question. Why is the protein in the corner of the box
and not in the middle?
I thought I centered it with the command:

editconf -f wholeProtein.gro -o 3m71_box.gro -center 4.59340   4.59470
  5.17330 -c -bt dodecahedron -d 1.0 2>>logErr 1>>logOut

Or not?



Have you tried trjconv -pbc mol -ur compact?  The "center" of an infinite system 
is an arbitrary location; you have to re-wrap the unit cell to achieve the 
expected outcome.  It's centered, it just may not look like it in the current 
representation.


-Justin




On 7/11/12 6:13 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
I did a minimization of my structure. But the output seems a bit strange
for me, since my input was the protein with its membrane in a box like
this:
# = box  P = Protein  M=Membrane


#
#   #
# MM#
#   #
#   #
#   #
#

But the output of the minimization looks like this

#
##
   #M M #
#    #
 ##
  ##
   ##

How can that be?


It's a triclinic representation of the unit cell.  Nothing is wrong.


My box was created with:

editconf -f wholeProtein.gro -o 3m71_box.gro -center 4.59340   4.59470
5.17330 -c -bt dodecahedron -d 1.0 2>>logErr 1>>logOut



Is a dodecahedral box appropriate for a protein in a membrane?  The
inherent
symmetry of a slab of any sort dictates that you should be using a cubic
or
rectangular box.  Dodecahedral and octahedral boxes are better for systems
with
spherical symmetry, like globular proteins in water.



I put solvent in it with:

genbox -cp 3m71_box.gro -cs spc216.gro -p 3m71.top -o 3m71_water.gro
2>>logErr 1>>logOut


The mdp file for the minimization looks like this:

define  = -DPOSRES
integrator  = steep
emtol   = 10
nsteps  = 1500
nstenergy   = 1
energygrps  = System
coulombtype = PME
rcoulomb= 0.9
rvdw= 0.9
rlist   = 0.9
fourierspacing  = 0.12
pme_order   = 4
ewald_rtol  = 1e-5
pbc = xyz


Can you please tell me how I can prevent my box to shift?



There is no box shift; it's just a visual representation of the
dodecahedral
box.  You can "correct" it (for visualization purposes) by running trjconv
-pbc
mol -ur compact.

-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
* Only plain text messages are allowed!
* 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







--


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
* Only plain text messages are allowed!
* 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] (no subject)

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 6:19 AM, amir abbasi wrote:

Thanks Justin,
But I want to neutralize my system in implicit solvent.
In Amber I had use Debye screening but in gromacs I don't know what should I do.



From my understanding, this remains an unresolved issue.

-Justin



On Wed, Jul 11, 2012 at 2:45 PM, Justin A. Lemkul  wrote:



On 7/11/12 6:00 AM, amir abbasi wrote:


Hi All!
I want to use Implicit solvent to simulate a nucleic acid sequence.
How can I do it?
I use this command:
genion -s ions.tpr -o nucleic_ions.gro -p nucleic.top -pname K+ -nname
CL -neutral -conc 0.1

ions.tpr file is same as umbrella sampling tutorial.

I got this error message:
Fatal error:
Your solvent group size (2898) is not a multiple of 31
what should I do?



One does not typically add explicit ions in an implicit solvent system.
genion fails because it appears you are trying to replace parts of your
nucleic acid with ions.

-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
* Only plain text messages are allowed!
* 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




--
============

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
* Only plain text messages are allowed!
* 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] question about the output of minimization

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 6:13 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
I did a minimization of my structure. But the output seems a bit strange
for me, since my input was the protein with its membrane in a box like
this:
# = box  P = Protein  M=Membrane


#
#   #
# MM#
#   #
#   #
#   #
#

But the output of the minimization looks like this

#
##
  #M M #
   #    #
##
 ##
  ##

How can that be?


It's a triclinic representation of the unit cell.  Nothing is wrong.


My box was created with:

editconf -f wholeProtein.gro -o 3m71_box.gro -center 4.59340   4.59470
5.17330 -c -bt dodecahedron -d 1.0 2>>logErr 1>>logOut



Is a dodecahedral box appropriate for a protein in a membrane?  The inherent 
symmetry of a slab of any sort dictates that you should be using a cubic or 
rectangular box.  Dodecahedral and octahedral boxes are better for systems with 
spherical symmetry, like globular proteins in water.




I put solvent in it with:

genbox -cp 3m71_box.gro -cs spc216.gro -p 3m71.top -o 3m71_water.gro
2>>logErr 1>>logOut


The mdp file for the minimization looks like this:

define  = -DPOSRES
integrator  = steep
emtol   = 10
nsteps  = 1500
nstenergy   = 1
energygrps  = System
coulombtype = PME
rcoulomb= 0.9
rvdw= 0.9
rlist   = 0.9
fourierspacing  = 0.12
pme_order   = 4
ewald_rtol  = 1e-5
pbc = xyz


Can you please tell me how I can prevent my box to shift?



There is no box shift; it's just a visual representation of the dodecahedral 
box.  You can "correct" it (for visualization purposes) by running trjconv -pbc 
mol -ur compact.


-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
* Only plain text messages are allowed!
* 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] (no subject)

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 6:00 AM, amir abbasi wrote:

Hi All!
I want to use Implicit solvent to simulate a nucleic acid sequence.
How can I do it?
I use this command:
genion -s ions.tpr -o nucleic_ions.gro -p nucleic.top -pname K+ -nname
CL -neutral -conc 0.1

ions.tpr file is same as umbrella sampling tutorial.

I got this error message:
Fatal error:
Your solvent group size (2898) is not a multiple of 31
what should I do?



One does not typically add explicit ions in an implicit solvent system.  genion 
fails because it appears you are trying to replace parts of your nucleic acid 
with ions.


-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
* Only plain text messages are allowed!
* 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] sovation with tip3p

2012-07-11 Thread Justin A. Lemkul



On 7/11/12 5:36 AM, Shima Arasteh wrote:

Hi all,

I'm going to do the step 3 in Justin's tutorial (Lysozyme in water) . In this 
step, the solvation is accomplished through this command:
genbox -cp 1AKI_newbox.gro -cs spc216.gro -o 1AKI_solv.gro -p topol.top

in which that spc216.gro is used. In first step I had used the TIP3P water 
model, here I'd rather to call TIP3P too. But there is not tip3p.gro in 
share/top .

I would appreciate you if you give me any suggestion.


http://www.gromacs.org/Documentation/How-tos/TIP3P_coordinate_file

-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
* Only plain text messages are allowed!
* 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] refcoord_scaling

2012-07-10 Thread Justin A. Lemkul



On 7/10/12 5:48 PM, Katie Maerzke wrote:

Hi all -

I'm new to Gromacs (and MD).  I can't figure out whether or not
refcoord_scaling is important for an NpT simulation.  In an MC simulation,
when the volume changes, the coordinates of the molecules are also scaled.
Does refcoord_scaling control how to handle scaling the molecular coordinates
(e.g., scale all coordinates or scale based on COM), or does it do something
else?  Is it even important?  I've searched the manual and the user forum,
but I haven't found an answer to my question.



The important distinction is between reference coordinates and the actual 
coordinates of the atoms.  When employing position restraints, the restraint 
potential is based on the reference positions of the atoms.  It is these that 
are updated based on the refcoord_scaling option.  In most cases, the box 
shouldn't change drastically, and if it does, you have bigger problems.  Thus, I 
imagine the contribution to energy is likely small, but I have never tested this 
myself.  In theory, if the reference coordinates are not updated correctly, you 
will get incorrect contributions to the virial and thus some errors in pressure 
calculation.


This is my understanding of the implementation; if it is incorrect or 
incomplete, hopefully someone will chime in.


-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
* Only plain text messages are allowed!
* 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] POPC in water

2012-07-10 Thread Justin A. Lemkul



On 7/10/12 6:09 AM, Shima Arasteh wrote:

Dear gmx users,


I want to simulate a protein in bilayer. The chosen bilayer is POPC.
According to Justin's tutorial about KALP15 in DPPC, I would simulate the
protein in lipid bilayer and water. In this tutorial I didn't find the
simulation of bilayer in water seperately and it is just going through the
simulation with protein.



See part 1 of step 3:

http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/membrane_protein/03_solvate.html

It's not a full simulation of a membrane, but it presents all the essential 
elements - coordinate file, topology (very simple), and .mdp (though it is for 
EM, not an MD simulation).



I need to say that I got the POPC.itp and .top through the link sent me by
dear Peter.


Now I'm wondering if POPC is needed to simulate in water before applying in
simulation of protein in POPC?



The answer depends on how well equilibrated the starting structure is.  If 
you're going to be putting a protein in it, you'll have to re-equilibrate that 
new system for a considerable amount of time, so simulating the membrane 
beforehand may or may not be relevant, but a thoroughly equilibrated membrane 
will reduce the time needed in equilibrating the membrane protein system.


-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
* Only plain text messages are allowed!
* 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] node decomposition' problem

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 4:25 PM, Justin A. Lemkul wrote:



On 7/9/12 4:23 PM, Thales Kronenberger wrote:

I'm trying to run a kinase (what means that I had ATP - large charged
group) energy minimization and then MD.

But when I put my protein together with its ligands I gotcha the
follow error message:

"There is no domain decomposition for 6 nodes that is compatible with
the given box and a minimum cell size of 6.47943 nm
Change the number of nodes or mdrun option -rdd
Look in the log file for details on the domain decomposition"

I read about in gromacs forums and I can force the thing's running
with the option nt 1 (one node...).

My problem is that I still want to run in parallel. Is it still that
possible for my system or I'm doomed to the 1 core simulation



The minimum charge group size depends on a whole host of factors:



*Edit* "minimum domain decomposition cell" - I don't know where "charge group" 
came from...



http://www.gromacs.org/Documentation/Errors#There_is_no_domain_decomposition_for_n_nodes_that_is_compatible_with_the_given_box_and_a_minimum_cell_size_of_x_nm


The large size you have obtained indicates there are likely problems with the
.mdp file, topology, or both.

-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
* Only plain text messages are allowed!
* 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] node decomposition' problem

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 4:23 PM, Thales Kronenberger wrote:

I'm trying to run a kinase (what means that I had ATP - large charged
group) energy minimization and then MD.

But when I put my protein together with its ligands I gotcha the
follow error message:

"There is no domain decomposition for 6 nodes that is compatible with
the given box and a minimum cell size of 6.47943 nm
Change the number of nodes or mdrun option -rdd
Look in the log file for details on the domain decomposition"

I read about in gromacs forums and I can force the thing's running
with the option nt 1 (one node...).

My problem is that I still want to run in parallel. Is it still that
possible for my system or I'm doomed to the 1 core simulation



The minimum charge group size depends on a whole host of factors:

http://www.gromacs.org/Documentation/Errors#There_is_no_domain_decomposition_for_n_nodes_that_is_compatible_with_the_given_box_and_a_minimum_cell_size_of_x_nm

The large size you have obtained indicates there are likely problems with the 
.mdp file, topology, or both.


-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
* Only plain text messages are allowed!
* 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: DNA simulations

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 4:07 PM, SatyaK wrote:

  Hello,

  I followed below steps using VMD and GROMACS but something went wrong in
using GROMACS which I am not able to figure out. Appreciate your help.

  1. editconf -f initialfile.pdb -o initialfile.gro -d 0.2
  2. VMD: within 5 of nucleic
  $sel writepdb initialfile_updated.pdb
  $sel delete

   (initialfile_updated.pdb has only those molecules that are at a
distance of 5A)

  3. I used GROMACS to convert back to .gro:editconf -f
initialfile_updated.pdb -o initialfile_updated.gro -d 0.2

  The coordinates in ininitialfile.gro and initialfile_updared.gro are
different. I quite don't understand the reason for the same.


You are resetting the box by invoking the -d option with editconf.  If 
initialfile.gro had all the original atoms and initialfile_updated.gro has 
significantly fewer (due to the selection made in VMD), then the size associated 
with the system is much smaller.  Using -d re-centers the system within the 
defined box, thus shifting the coordinates.


-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
* Only plain text messages are allowed!
* 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] force field parametrs for Mn2+

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 1:21 PM, tarak karmakar wrote:

Oh !! nice work
Thanks a lot for the quick reply. But I'm very sorry to inform you
that whichever table [supplementary table S4] you are specifying in
the supporting info, I couldn't find anywhere. So, may be I'm finding
my way in wrong track. Can you please provide me the link and / the
table containing the parameters for the Manganese ?


It's on p. 11 of the supplement.

http://onlinelibrary.wiley.com/store/10.1002/anie.201202032/asset/supinfo/anie_201202032_sm_miscellaneous_information.pdf?v=1&s=cff1986017d1843da85eb75dbd174c8e11727dc8

-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
* Only plain text messages are allowed!
* 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] Reg dimers

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 9:48 AM, Ramya LN wrote:

Dear all,
I have done protein-ligand dynamics.I got the final gro file.When
converted to PDB, i observed that my active site has ligand but two
chains of teh protein got separated. What might be the reason for
this?should i consider this as an error in my simulation???kindly help
me in this regard.



http://www.gromacs.org/Documentation/Terminology/Periodic_Boundary_Conditions

-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
* Only plain text messages are allowed!
* 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] error in mdrun

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 9:46 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

With the mentioned below options I get the following error:

Fatal error:
1 particles communicated to PME node 1 are more than 2/3 times the cut-off
out of the domain decomposition cell of their charge group in dimension x.
This usually means that your system is not well equilibrated.


But it does not occur immediatly but only at step

Step 695, time 1.39 (ps)  LINCS WARNING
relative constraint deviation after LINCS:
rms 2.011725, max 78.611298 (between atoms 3778 and 3781)
bonds that rotated more than 30 degrees:
  atom 1 atom 2  angle  previous, current, constraint length
3782   3785   92.21.3061   0.2785  0.1090
3782   3784   93.00.1055   0.2552  0.1090
3782   3783   91.00.1062   0.3086  0.1090
3204   3206   89.60.1449   1.0256  0.1449
3204   3205   89.70.1010   0.9949  0.1010


Is there something wrong with my temperature coupling?



I doubt it.  Decrease your timestep (as I've said twice now) and try again with 
something sensible for dt.


-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
* Only plain text messages are allowed!
* 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] error in mdrun

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 9:40 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi Justin,
  okey then I will try it with this timestep.
No it is not my goal to do a NVE. I already had temperature coupling
options in my .mdp file but on the blowing up side was written "you are
using inappropriate temperature coupling" so I thought that that might be
the reason and deleted it from my .mdp file.



Temperature coupling in itself is not the source of the problem, but what is 
stated on the wiki is that inappropriate use of it can lead to instability.  The 
settings below look fine.  I think your biggest problem is the large timestep.


-Justin


I had the following temperature coupling options:

tcoupl  = V-rescale
tc-grps = Protein  Non-Protein
tau_t   = 0.1  0.1
ref_t   = 298  298
pcoupl  = no



Thank you for your answer.
Eva




On 7/9/12 9:25 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
I want to do a md for a protein with a membrane around it.
I already minimised the energy of the protein.



Output of the minimization:

^MStep=  812, Dmax= 1.4e-06 nm, Epot= -7.59283e+05 Fmax= 3.59781e+02,
atom= 1653
Stepsize too small, or no change in energy.
Converged to machine precision,
but not to the requested precision Fmax < 10

Double precision normally gives you higher accuracy.
You might need to increase your constraint accuracy, or turn
off constraints alltogether (set constraints = none in mdp file)

writing lowest energy coordinates.


Steepest Descents converged to machine precision in 813 steps,
but did not reach the requested Fmax < 10.
Potential Energy  = -7.5928356e+05
Maximum force =  3.6197971e+02 on atom 1653
Norm of force =  2.6517429e+00





In my eyes this output was okey so I went on with the md.
And here I get the error:

Fatal error:
A charge group moved too far between two domain decomposition steps
This usually means that your system is not well equilibrated




I already read the error page for this error and also the blowing up
page
but I still do not know what to do now.

My .mdp file for the md runs looks like this:

define  = -DPOSRES
integrator  = md
dt  = 0.005


This timestep is huge.  Even with constraints, you probably can't exceed 2
fs
stably (0.002 ps).


nsteps  = 2000
nstxout = 0
nstvout = 0
nstfout = 0
nstlog  = 1000
nstxtcout   = 0
nstenergy   = 5
energygrps  = Protein Non-Protein
nstcalcenergy   = 5
nstlist = 10
ns-type = Grid
pbc = xyz
rlist   = 0.9
coulombtype = PME
rcoulomb= 0.9
rvdw= 0.9
fourierspacing  = 0.12
pme_order   = 4
ewald_rtol  = 1e-5
gen_vel = yes
gen_temp= 200.0
gen_seed= 
constraints = all-bonds


Are here any optinos which can cause the error?



In the absence of temperature and/or pressure coupling, the ensemble
you're
trying to simulate is NVE, which is very tricky to get stabilized.

http://www.gromacs.org/Documentation/Terminology/NVE

If you're not going for an NVE ensemble, you need several adjustments in
the
.mdp file.  See any basic tutorial for examples of how to simulate other
ensembles, if this is your goal.

-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
* Only plain text messages are allowed!
* 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







--
============

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
* Only plain text messages are allowed!
* 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] error in mdrun

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 9:25 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
I want to do a md for a protein with a membrane around it.
I already minimised the energy of the protein.



Output of the minimization:

^MStep=  812, Dmax= 1.4e-06 nm, Epot= -7.59283e+05 Fmax= 3.59781e+02,
atom= 1653
Stepsize too small, or no change in energy.
Converged to machine precision,
but not to the requested precision Fmax < 10

Double precision normally gives you higher accuracy.
You might need to increase your constraint accuracy, or turn
off constraints alltogether (set constraints = none in mdp file)

writing lowest energy coordinates.


Steepest Descents converged to machine precision in 813 steps,
but did not reach the requested Fmax < 10.
Potential Energy  = -7.5928356e+05
Maximum force =  3.6197971e+02 on atom 1653
Norm of force =  2.6517429e+00





In my eyes this output was okey so I went on with the md.
And here I get the error:

Fatal error:
A charge group moved too far between two domain decomposition steps
This usually means that your system is not well equilibrated




I already read the error page for this error and also the blowing up page
but I still do not know what to do now.

My .mdp file for the md runs looks like this:

define  = -DPOSRES
integrator  = md
dt  = 0.005


This timestep is huge.  Even with constraints, you probably can't exceed 2 fs 
stably (0.002 ps).



nsteps  = 2000
nstxout = 0
nstvout = 0
nstfout = 0
nstlog  = 1000
nstxtcout   = 0
nstenergy   = 5
energygrps  = Protein Non-Protein
nstcalcenergy   = 5
nstlist = 10
ns-type = Grid
pbc = xyz
rlist   = 0.9
coulombtype = PME
rcoulomb= 0.9
rvdw= 0.9
fourierspacing  = 0.12
pme_order   = 4
ewald_rtol  = 1e-5
gen_vel = yes
gen_temp= 200.0
gen_seed= 
constraints = all-bonds


Are here any optinos which can cause the error?



In the absence of temperature and/or pressure coupling, the ensemble you're 
trying to simulate is NVE, which is very tricky to get stabilized.


http://www.gromacs.org/Documentation/Terminology/NVE

If you're not going for an NVE ensemble, you need several adjustments in the 
.mdp file.  See any basic tutorial for examples of how to simulate other 
ensembles, if this is your goal.


-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
* Only plain text messages are allowed!
* 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] large radius problem

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 7:09 AM, siddhant jain wrote:

I was performing simulations on urea unfolding of protein. I performed
one set of simulation for 50 ns with a velocity (set by gen_seed) and
it went fine.
Now when I am doing simulation using a different gen_seed from 10 to
20 ns many long bonds are formed. I could visualize them in vmd.
Subsequently rmsd and radius of gyration plots also show the large
change. But,  energy is almost constant during the whole time period.
Also after 20 ns this problem goes by itself and all bonds become
normal. Why it could be so and how could I correct it. These bonds are
so long that radius of gyration changes from 1.2 nm to 3.5 nm.



That sounds like a PBC issue to me.

http://www.gromacs.org/Documentation/Terminology/Periodic_Boundary_Conditions

-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
* Only plain text messages are allowed!
* 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] Peptide folding simulation

2012-07-09 Thread Justin A. Lemkul



On 7/9/12 12:02 AM, bharat gupta wrote:

Hi,


I have been trying to study folding of a peptide 24 residues long. I
did  a  simulation of 50 ns with explicit solvent, CHARMM FF, but I
was not able to find even a single folding event. Then I decided use
explicit solvent for simulation and I again simulated the peptide for
100 ns . This time again I ended with no folding events.


  I know that in case of explicit solvent , a 50ns simulation time is
not enough to observe anything. But I did it to see the initial
behavior of the peptide in water. In take many random like
conformation but doesnot fold into a desired beta-hairpin. For the
explicit solvent simulation, I followed the lysozyme tutorial
parameters.



You shouldn't.  The .mdp settings are appropriate for OPLS-AA, not CHARMM.


For implicit solvent simulation, I used the following parameters for
Energy minimization  :

  define  =  -DFLEXIBLE
  constraints =  none
  integrator  =  steep
  dt  =  0.001; ps
  nsteps  =  3
  vdwtype =  cut-off
  coulombtype =  cut-off
  pbc =  no
  nstlist =  0
  ns_type =  simple
  rlist   =  0   ; this means all-vs-all (no cut-off),
which   gets expensive for bigger systems
  rcoulomb=  0
  rvdw=  0
  comm-mode   =  angular
  comm-grps   =  Protein
  optimize_fft=  yes
  ;
  ;   Energy minimizing stuff
  ;
  emtol   =  5.0
  emstep  =  0.01
  ;
  ; Implicit solvent
  ;
  implicit_solvent=  GBSA
  gb_algorithm=  Still ; HCT ; OBC
  nstgbradii  =  1
  rgbradii=  0   ; [nm] Cut-off for the calculation of the
Born  radii. Currently must be equal to rlist
  gb_epsilon_solvent  =  80; Dielectric constant for the implicit  solvent
   ; gb_saltconc   =  0 ; Salt concentration for implicit
solvent  models, currently not used
  sa_algorithm=  Ace-approximation
  sa_surface_tension  = -1



For MD I used the following : -


define  =  -DPOSRESHELIX ; -DFLEXIBLE -DPOSRES
  constraints =  none
  integrator  =  md
  dt  =  0.001   ; ps
  nsteps  =  10 ; 10 ps = 100 ns
  nstcomm =  10
  nstcalcenergy   =  10
  nstxout =  1000 ; frequency to write coordinates to output
  trajectory
  nstvout =  0   ; frequency to write velocities to output
  trajectory; the last velocities are always written
  nstfout =  0   ; frequency to write forces to output
  trajectory
  nstlog  =  1000 ; frequency to write energies to log
  file
  nstenergy   =  1000 ; frequency to write energies to edr file

  vdwtype =  cut-off
  coulombtype =  cut-off

  pbc =  no

  nstlist =  0
  ns_type =  simple
  rlist   =  0   ; this means all-vs-all (no cut-off), which
  gets expensive for bigger systems
  rcoulomb=  0
  rvdw=  0

  comm-mode   =  angular
  comm-grps   =  system

  optimize_fft=  yes

  ; V-rescale temperature coupling is on
  Tcoupl  =  v-rescale
  tau_t   =  0.1
  tc_grps =  system
  ref_t   =  300
  ; Pressure coupling is off
  Pcoupl  =  no
  ; Generate velocites is on
  gen_vel =  yes
  gen_temp=  300
  gen_seed=  -1

  ;
  ; Implicit solvent
  ;
  implicit_solvent=  GBSA
  gb_algorithm=  Still ; HCT ; OBC
  nstgbradii  =  1
  rgbradii=  0   ; [nm] Cut-off for the calculation of the
Born radii. Currently must be equal to rlist
  gb_epsilon_solvent  =  80; Dielectric constant for the implicit   solvent
  ; gb_saltconc   =  0 ; Salt concentration for implicit
solvent   models, currently not used
  sa_algorithm=  Ace-approximation
  sa_surface_tension  = -1

So, finally I want to know where have I gone in my simulation
experiments, both implicit and explicit ?? ... Please reply .



What evidence do you have that you should expect to see a folding event in such 
a short time?  Most people will use more extensive sampling methods like REMD to 
observe peptide folding.


-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
* Only plain text messages are allowed!
* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests t

Re: [gmx-users] Re: Two questions about index files

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 4:14 PM, Andrew DeYoung wrote:

Hi,

Thanks, Justin! Do you by any chance have any experience with the "expect"
tool (http://www.nist.gov/el/msid/expect.cfm)? I have never used it. I guess
that this may not help in this case because, of course, the printing of
selections is due to Gromacs. But I wonder if there is any way that "expect"
can be used in tandem with modification of the source code to run the
Gromacs utility in the background.



No, I've never used it.  I think there has been some mention across the list, 
but I don't know how it would be necessary (or exceptionally useful) here. 
Maybe someone with more experience with it can comment.  I certainly doubt it 
will suppress Gromacs screen output.


-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
* Only plain text messages are allowed!
* 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] define a new residue

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 12:43 PM, Shima Arasteh wrote:

OK.
What about generating an output file through CGenFF by the first 3 residues of 
the protein, rather thn the first 2 (formyl+valine)?




Maybe.  Try it and see, rather than waiting a few hours for someone to get 
around to replying :)


-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
* Only plain text messages are allowed!
* 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] Two questions about index files

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 3:59 PM, Andrew DeYoung wrote:

Hi,

I hope that all is well. If you have time, I have two questions about index
files:

(1) Do you know if there is a limit on the number of entries an index file
can have? I am trying to write a shell script which would allow me to run
g_traj 20 times, feeding a different index entry to it on each
iteration. This is because I have used g_select to (dynamically) pick out
atoms that meet certain criteria. I have 20 time steps in my trajectory,
so I have 20 entries in my index file. I will try this, but do you have
any experience with this about whether such a large index file will work?



The only limitation here is disk space and available memory, both of which are 
external to Gromacs.



(2) Normally, when I call any Gromacs program with the -n switch, the
program prints to the screen all of the choices available in the index file.
I am planning to use the "<http://www.gromacs.org/Documentation/How-tos/Using_Commands_in_Scripts#Withi
n_Script

This works, but the problem is that the program still prints all of the
index file selections to the screen, even though the input is now automated
from the shell script. This means that the program insists upon printing out
all 20 index selections on each of my 20 calls to the program, and
this can take a lot of time.  Do you know if there is a way to suppress the
printing of the index files selections and somehow feed input to the program
in the background/invisibly? Simply putting & at the end of the command does
not seem to work.

In other words, is there a way to make commands non-interactive _in the
background_?  Thank you!



I don't think so.  There is a hidden option called -quiet that reduces some 
screen output, but since selections are dependent upon user choice, they can't 
be suppressed, at least in the quick trials I just did.  I suppose you could 
modify the source to prevent this printing, but whether or not that's worth the 
effort is up to you.


-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
* Only plain text messages are allowed!
* 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] define a new residue

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 11:29 AM, Shima Arasteh wrote:

Thanks all.
So if I find a protein which is parametrized by CHARMM and then find the valine 
residue there, I might use the parameters of side chains of it in my own rtp 
file. Right?




You can look this up in the force field's .rtp file.  For full parameterization 
procedures, refer to the primary literature for CHARMM.  It should detail how to 
derive parameters for new species.


-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
* Only plain text messages are allowed!
* 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: Re: g_sham problem

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 9:39 AM, siddhant jain wrote:

Thanks Justin
Also I wanted to ask that to know my structure which has minimum energy, can
I take the corresponding values of principle components, look at the time at
which they occur simultaneously, and hence note the structure corresponding
to that time. Will it be reasonable enough or its vague.



You can certainly identify structures within different regions of the free 
energy surface in this way.  Just be careful about saying any structure has a 
particular "energy" value, since this can mean a number of different things.


-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
* Only plain text messages are allowed!
* 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] g_sham problem

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 9:16 AM, siddhant jain wrote:

I used the g_sham command to get 2-d energy field where the two
components are eigenvector 1 and 2. The graph is fine but it has no
labelling along the axis.
I used
g_sham -f proj124.xvg -ls 124.xpm -notime
I want the labels along the axes so that i can comprehend it better.



Use the -xmin and -xmax options to set suitable values.  g_sham tries to guess, 
but in my experience, does a poor job of it and I usually have to set everything 
manually.


-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
* Only plain text messages are allowed!
* 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] Crashes during protein-ligand simulation

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 6:57 AM, James Starlight wrote:

Justin,

unfortunately my last system have also been crashed after 35ns of
simulation with the links warnings accompanied by the error

Fatal error:
A charge group moved too far between two domain decomposition steps
This usually means that your system is not well equilibrated
For more information and tips for troubleshooting, please check the GROMACS
website at http://www.gromacs.org/Documentation/Errors

This time I've devided all functional groups of my ADENOSINE ligand
into separate charge groups in topol.itp.



Note that these charges are not the same as the adenosine moiety of ATP, which 
is standard in the Gromos96 force fields.  Perhaps you should re-evaluate some 
charges (see aminoacids.rtp for the existing ATP implementation).



;  nr  type  resnr  resid  atom  cgnr  chargemasstotal_charge
 1NT1ADN N61   -0.844  14.0067
 2 H1ADNH1110.422   1.0080
 3 H1ADNH1210.422   1.0080  ;  0.000
 4 C1ADN C820.097  12.0110
 5HC1ADNH0120.177   1.0080
 6NR1ADN N32   -0.642  14.0067
 7 C1ADN C420.175  12.0110
 8 C1ADN C520.092  12.0110
 9NR1ADN N72   -0.556  14.0067
10 C1ADN C620.657  12.0110  ;  0.000
11 C1ADNC5'3   -0.677  12.0110
12 C1ADNC4'30.834  12.0110
13OE1ADNO4'3   -0.248  15.9994
14 C1ADNC1'3   -0.558  12.0110
15 C1ADNC2'40.603  12.0110
16 C1ADNC3'5   -0.212  12.0110
17NR1ADN N930.415  14.0067
18OA1ADNO2'4   -0.606  15.9994
19 H1ADNH0840.482   1.0080
20OA1ADNO3'5   -0.606  15.9994
21 H1ADNH0650.482   1.0080
22OA1ADNO5'3   -0.246  15.9994
23 H1ADNH0330.337   1.0080  ; -0.000
24 C1ADN C260.502  12.0110
25HC1ADNH1060.106   1.0080
26NR1ADN N16   -0.608  14.0067  ;  0.000


Also I've done proper equilibration in tho steps (I've used 1fs
integrator steps on both stages of equilibration)
1- I've made equilibration with posres ( 500ps) on protein backbone
atoms as well as ligand with x-ray water

2- The next step (5ns) was done without posres only with smaller
integrator steps.

What another possible sollutions could be ?



Try some diagnostics:

http://www.gromacs.org/Documentation/Terminology/Blowing_Up#Diagnosing_an_Unstable_System

If simulations of the same system without a ligand run fine, then the error is 
still related to the ligand, either its topology or placement.


-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
* Only plain text messages are allowed!
* 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] error in grompp

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 6:02 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi Justin,
thank you for your answer.

Now I tried it with two different restraint .itp files. One for the
protein and one for the dummy atoms.
But still it doesn't work. Now the error is:

[ file posre_memb.itp, line 5 ]:
Atom index (4942) in position_restraints out of bounds (1-1).
This probably means that you have inserted topology section
"position_restraints"
in a part belonging to a different molecule than you intended to.
In that case move the "position_restraints" section to the right molecule.


But I think I included it the right way:


; Include Position restraint file
#ifdef POSRES
#include "posre.itp"
#endif

; Include water topology
#include "amber03.ff/tip3p.itp"

#ifdef POSRES_WATER
; Position restraint for each water oxygen
[ position_restraints ]
;  i funct   fcxfcyfcz
11   1000   1000   1000
#endif

; Include topology for ions
#include "amber03.ff/ions.itp"

#include "amber03.ff/dum.itp"
#ifdef POSRES
#include "posre_memb.itp"
#endif



In my coordiate file the difference between them look like this:

313LEU   HD23 4938   3.813   4.505   3.308
   313LEU  C 4939   3.435   4.335   3.190
   313LEUOC1 4940   3.429   4.330   3.090
   313LEUOC2 4941   3.337   4.305   3.259
   314DUMDUM 4942   1.996   2.371   6.171
   314DUMDUM 4943   1.996   2.371   6.271
   314DUMDUM 4944   1.996   2.471   6.171
   314DUMDUM 4945   1.996   2.471   6.271


my restraint file for the protein looks like this:


; position restraints for Protein-H of GROup of MAchos and Cynical Suckers

[ position_restraints ]
;  i funct   fcxfcyfcz
11   1000   1000   1000
41   1000   1000   1000
71   1000   1000   1000
   101   1000   1000   1000
   131   1000   1000   1000


and my restraint file for the dummy atoms look like this:

; position restraints for Protein-H of GROup of MAchos and Cynical Suckers

[ position_restraints ]
;  i funct   fcxfcyfcz
49421   1000   1000   1000
49431   1000   1000   1000
49441   1000   1000   1000
49451   1000   1000   1000


What is wrong?



Atom numbering is done per [moleculetype] and has nothing to do with the atom 
numbers in the coordinate file.  If you have a one-atom dummy [moleculetype], 
then the only valid content of posre_memb.itp is:


[ position_restraints ]
;  i funct   fcxfcyfcz
   11   1000   1000   1000

-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
* Only plain text messages are allowed!
* 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] top file

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 3:12 AM, Shima Arasteh wrote:





  Dear gmx users,

I have gotten 2 topology files from a unique simulation ( a protein in water) 
with gmx and CHARMM36. I  know the charges of atoms derived by PRODRG and 
CGENFF are different, but I'm wondering if the total charge which is visible in 
top file, is supposed to be the same? What items are expected to be the same? 
or be in agreement?



The total charge is given in the qtot column of the .top file.  In theory, at a 
given protonation state, any protein should have the same net charge under any 
force field.


-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
* Only plain text messages are allowed!
* 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] define a new residue

2012-07-08 Thread Justin A. Lemkul



On 7/8/12 7:10 AM, Shima Arasteh wrote:

Because I saw that residues defined in aminoacids.rtp file don't have H in 
their carboxyl and atom N. So I decided to remove H !



Chemically, there is never an H atom on the C-terminal carbonyl of a group 
involved in a peptide bond.  If you parameterized a molecule that was protonated 
in such a way, it will give an erroneous output.  You can't just delete atoms, 
and shifting the charge by adding it back somewhere else is also very suspect, 
because you're changing the electronic nature of the compound.


I would say you need a better parameterization protocol using a more suitable 
model compound.


-Justin





Sincerely,
Shima



From: francesco oteri 
To: Shima Arasteh 
Cc: Discussion list for GROMACS users 
Sent: Sunday, July 8, 2012 3:28 PM
Subject: Re: [gmx-users] define a new residue


Why have you removed the hydrogen?


2012/7/8 Shima Arasteh 

Dear Francesco,

Thanks. Honestly I thought about this, but I don't know how much charges I need 
to increase or decrease of other atoms?

Is it possible to add the FVAL with the H atom( which I removed before)? I mean 
that I apply the complete formyl-valine and don't remove the H atom.

Thanks for your suggestions.


Cheers,
Shima




From: francesco oteri 
To: Shima Arasteh ; Discussion list for GROMACS users 

Sent: Sunday, July 8, 2012 3:04 PM
Subject: Re: [gmx-users] define a new residue



Hi Shima,
usually charge calculation are carried out on the system that is supposed to be 
used in the MD.
So in my opinion you shouldn't remove the hydrogen.
Anyway, if you wanna remove the hydrogen, either you increase the charge of 
some other atom or
calculate the charges on the new residue without hydrogen.

Francesco


2012/7/8 Shima Arasteh 

Hi dear gmx friends,


I got the parameters of formyl-valine through the CHARMM website. Now I need to 
define it as a new residue FVAL( as Justin suggested me earlier) in .rtp file. 
To set the correct charges for atoms, I used the CHARMM output. In order to 
define FVAL to rtp file, I added these lines as below:

CNC0.3490
 ONO-0.4941
 H1HC0.1002
 NNH1-0.4233
 HNH0.094
 CACT10.1445
 HAHB0.096
 CBCT1-0.0977
 HBHA0.098
 CG1CT3-0.2689
 HG11HA0.0910
 HG12HA0.0911
 HG13HA0.0912
 CG2CT3-0.26813
 HG21HA0.0914
 HG22HA0.0915
 HG23HA0.0916
 CC0.20917
 OO-0.39518

As you see I omitted the H connected to Carboxyl. Now the total charge of the 
new-defined residue is not zero (-0.333) . How can I correct it?


I would appreciate you for your suggestions.
Thanks in advance.

Sincerely,
Shima
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Only plain text messages are allowed!
* 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




--
Cordiali saluti, Dr.Oteri Francesco






--
============

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
* Only plain text messages are allowed!
* 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] Crashes during protein-ligand simulation

2012-07-07 Thread Justin A. Lemkul



On 7/7/12 11:08 AM, James Starlight wrote:

justin,

It seems that problem was in the big charge groups in the ligand.itp
file. In particularly I've devided largest group into several smaller
and there haven't any crashes been occured yet.

With my last system with the the default COM group I've obtained crash
on the 25ns with the error about ligand's big cngr exactly.


By the way could you tell me how I could analyse stability of the
protein-ligand system as well as contributions of the individual
non-covalent contacts with the .EDR file?

In my mdp file I've defined energygrps= Protein Ligand   as the
separate energy terms.

During analysis of the edr file I've observed some options like
LJ-SR:Protein-ADN ( this value was slightly negative (-100) during my trr )
Coul-14:Protein-ADN ( this value was constantly zero )
Coul-SR:ADN-rest ( this value was equal to the LJ-SR:Protein-ADN  -100 )

What conclusions in general could I do based on that values ? How I
could measure stability of the ligand ( beside dirrect RMSD
measurement) from such energy terms?



In isolation, these values are not indicative of anything.  Their absolute 
values are completely dependent upon the parameters set in the topology.  Force 
field parameterization methodology does not require that any topology reproduce 
protein-ligand energetics; only true free energy calculations might tell you this.


In conjunction with other metrics (hydrogen bonding, contacts, RMSD, SASA, etc) 
you may be able to make some meaningful observations about ligand stability in 
the protein binding site.


-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
* Only plain text messages are allowed!
* 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] error in grompp

2012-07-07 Thread Justin A. Lemkul



On 7/7/12 8:20 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
I want do to an energy minimization and simulation with position restraints.
Additionally to the protein I have a membrane around my protein which I
also ant to fix. If I won't fix it, it is in the whole box after the
minimization and simulation but not around my protein.
Since I only want the hydrogen atoms to be flexible I use the second
option "protein-h" when I was asked by genrestr. This makes that the whole
protein including the membrane of dummy atoms is fix.
But when I want to use grompp I get the error:

Atom index (4942) in position_restraints out of bounds (1-4941).
This probably means that you have inserted topology section
"position_restraints"
in a part belonging to a different molecule than you intended to.
In that case move the "position_restraints" section to the right molecule.


The first 4941 atoms are of the "real" protein without the membrane and
after this the membrane starts:

  313LEUOC1 4940   3.429   4.330   3.090
   313LEUOC2 4941   3.337   4.305   3.259
   314DUMDUM 4942   1.996   2.371   6.171
   315DUMDUM 4943   1.996   2.371   6.271
   316DUMDUM 4944   1.996   2.471   6.171


When I remove the restriction of the DUM atoms I don't get the error but
that is not what I want.

Can you please tell me how I can fix the membrane?



You need a separate position restraint .itp file.  As the error message states, 
position restraints are only applicable per [moleculetype].


http://www.gromacs.org/Documentation/Errors#Atom_index_n_in_position_restraints_out_of_bounds

-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
* Only plain text messages are allowed!
* 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: specifying the direction of Pull in US

2012-07-07 Thread Justin A. Lemkul



On 7/7/12 12:15 AM, Raj wrote:

Thanks for ur suggestion Justin,

I'm facing trouble in setting that vector, actually I cant figure out how
can i set up a vector. Is there any easier way with which i can set up a
vector. Thanks



It can be set simply to the difference in the (x,y,z) coordinates of the COMs of 
the two groups.  You can calculate the COM position of each group with g_traj 
-ox -com -x -y -z.


Alternatively, if you use "pull_geometry = distance" and set "pull_dim = Y Y Y" 
then the pulling direction is the vector between the two groups in all dimensions.


-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
* Only plain text messages are allowed!
* 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] Crashes during protein-ligand simulation

2012-07-06 Thread Justin A. Lemkul



On 7/6/12 3:56 PM, James Starlight wrote:

Justin,

I've experimented with 2 dirrerent COM groups

comm-grps   = SOL_NA_CL XW Protein_CCl4_ADN ; 3 groups

comm-grps   = SOL_NA_CL_XW Protein_CCl4_ADN; 2 groups

but the crashes were in both cases after 12- 15ns of simulation

this time I've changed to the

comm-grps   = System

and there have not been any crashes yet (to this time I've already
calculated addition 20ns after privious crhased simulation using
checkpoint file for the last simulation ). But it's posible that it
was lucky coincidence :)

Could you tell me how I could devide largest group in the above
axample into several smaller sub-groups ? Should I do that separation
randomly or there are most correct way for that ? (e.g within third
cgnr separate all nitrogens and oxygens with corresponded hydrogens in
the separate cgrp's from carbons)



Look at existing examples in the force field .rtp file.  In general, a charge 
group consists of a functional group (amine, carboxylate, etc) or small CHn 
units.  Generally there are 2-4 atoms per charge group.  There is some 
discussion on these topics in the manual.


-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
* Only plain text messages are allowed!
* 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] Crashes during protein-ligand simulation

2012-07-06 Thread Justin A. Lemkul



On 7/6/12 2:05 PM, James Starlight wrote:

Justin,


I've done all steps in accordance to your tutorial. I've already done
the same systems with another ligands but had no problem.

This time I've made topology of the ligand via ATB server. I've only
noticed that some cgnr are too big in that topology . This is the
example

ADN 3
[ atoms ]
;  nr  type  resnr  resid  atom  cgnr  chargemasstotal_charge
 1NT1ADN N61   -0.844  14.0067
 2 H1ADNH1110.422   1.0080
 3 H1ADNH1210.422   1.0080  ;  0.000
 4 C1ADN C820.097  12.0110
 5HC1ADNH0120.177   1.0080
 6NR1ADN N32   -0.642  14.0067
 7 C1ADN C420.175  12.0110
 8 C1ADN C520.092  12.0110
 9NR1ADN N72   -0.556  14.0067
10 C1ADN C620.657  12.0110  ;  0.000
11 C1ADNC5'3   -0.677  12.0110
12 C1ADNC4'30.834  12.0110
13OE1ADNO4'3   -0.248  15.9994
14 C1ADNC1'3   -0.558  12.0110
15 C1ADNC2'30.603  12.0110
16 C1ADNC3'3   -0.212  12.0110
17NR1ADN N930.415  14.0067
18OA1ADNO2'3   -0.606  15.9994
19 H1ADNH0830.482   1.0080
20OA1ADNO3'3   -0.606  15.9994
21 H1ADNH0630.482   1.0080
22OA1ADNO5'3   -0.246  15.9994
23 H1ADNH0330.337   1.0080  ; -0.000
24 C1ADN C240.502  12.0110
25HC1ADNH1040.106   1.0080
26NR1ADN N14   -0.608  14.0067  ;  0.000
; total charge of the molecule:  -0.000



Large charge groups could account for errors in neighbor searching, leading to 
clashes that cause the simulation to collapse.




2) To the binding pocket I've inserted this ligand manually by means
of superimposition with the reference x-ray structure wich include the
same protein in the same conformation with the same ligand. I've done
some systems already and that aproach was good :)



OK, just be ready for reviewers to ask why you didn't do docking ;)


3) It's strange that the simulation crashes without any reasons ( the
system is very stable during calculated 10-15ns trajectory)



There's always a reason, you just haven't found it yet.  The charge group size 
could indeed be the problem; neighbor searching can fail at any time when some 
atoms run into one another.



Also I suppose that such problems could be with the COM groups

this is the example from my mdp

comm-grps   = SOL_NA_CL XW Protein_CCl4_ADN

here XW is the water wich were coppied from X-ray structure .
Also in that system Ccl4 is the membrane mimicking layer so I've
merged it with protein and ligand in the same group.



I see no reason to add such complexity to the system.  Breaking the crystal 
waters into their own COM removal group does not make sense to me.  Physically, 
they are basically part of the protein.



On the current stage I've tried to make changes in the mdp on

comm-grps   = System

to check if the problem was with that COM motion



And what was the outcome?  I see no reason that two COM motion removal groups 
wouldn't be appropriate (as layers can slide with respect to one another, like a 
membrane) but three groups does not sound appropriate.


-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
* Only plain text messages are allowed!
* 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] Boundary

2012-07-06 Thread Justin A. Lemkul



On 7/6/12 3:06 PM, dariush wrote:

Hi all,

Why does GROMACS just provide PBC for boundary condition?
However, LAMPPS as an example provides four kind of boundary:
periodic, non-periodic and fixed, non-periodic and shrink-wrapped and
non-periodic and shrink-wrapped with a minimum value.



Gromacs can also do non-periodic systems by setting "pbc = no" in the .mdp file. 
 A variety of options exist for walls, as well.


If there are particular features that users want, they are welcome to implement 
them and submit them for review in the development code.  Otherwise, if no one 
files a feature request on redmine, the developers aren't going to invest time 
in the feature unless they need it themselves.  I would hazard a guess that 3-D 
periodic boundary conditions are the most commonly used in simulations of 
biomolecules.



PBC makes problem when you want to make movie in VMD, even you add some
mirror-wise molecule in each direction.

Is there anyway to figure out this problem?



Any trajectory can be "fixed" for visualization purposes with trjconv.  You may 
need several iterations to achieve the desired effect.


http://www.gromacs.org/Documentation/Terminology/Periodic_Boundary_Conditions#Suggested_trjconv_workflow

-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
* Only plain text messages are allowed!
* 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: bilayers move apart by nanometers upon implementing dihedral restraints on lipid tails

2012-07-06 Thread Justin A. Lemkul



On 7/6/12 10:40 AM, khandelia wrote:

Justin, Vitaly

The topology is fine, I double-checked

The simulation runs perfectly fine without the restraints.

It is not a PBC effect, since the box size along z is > 50 nm after a ns or
so.

Does one need yet another restraint to hold the bilayer together?



I have never had a need for any restraints to keep a bilayer intact.


There has been some discussion about problems with dihedral restrains in the
list earlier, but nothing like this.



The problem you're observing seems to indicate that your manipulation of the 
lipid chain causes physical instability.  How extensive are the restraints?  How 
many atoms do they involve?  You provided an "etc" in your previous message, so 
I'm trying to clarify what's going on.  Is it even physically possible to orient 
the lipid chain in such a way?  You've got basically all the consecutive 
dihedrals in a very specific orientation - is that compatible with your system? 
 Can you run a simulation of a single lipid in vacuo using these restraints?


-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
* Only plain text messages are allowed!
* 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] Error in Membrane simulations with POPC bilayer

2012-07-06 Thread Justin A. Lemkul



On 7/6/12 10:48 AM, J Peterson wrote:

Hi Justin and Anirban,

I started a membrane simulation with POPC bilayer after a training with the
given KALP peptide and DPPC bilayer. I am following both of your tutorials
(mainly the Justin's). I have problem at where I generate a .tpr file for a
DPPC (POPC here)-only system using grompp.

I see another warning on non-matching number of atoms along with the error
that you recommended a safe one.

My error is

Warning: atom name 3340 in topol_popc.top and popc_128b_H.pdb does not match
(C12 - H2)
Warning: atom name 3341 in topol_popc.top and popc_128b_H.pdb does not match
(C13 - O)
Warning: atom name 3342 in topol_popc.top and popc_128b_H.pdb does not match
(O14 - H1)
Warning: atom name 3343 in topol_popc.top and popc_128b_H.pdb does not match
(C15 - H2)
Warning: atom name 3344 in topol_popc.top and popc_128b_H.pdb does not match
(O16 - O)
Warning: atom name 3345 in topol_popc.top and popc_128b_H.pdb does not match
(C17 - H1)
Warning: atom name 3346 in topol_popc.top and popc_128b_H.pdb does not match
(C18 - H2)
Warning: atom name 3347 in topol_popc.top and popc_128b_H.pdb does not match
(C19 - O)
Warning: atom name 3348 in topol_popc.top and popc_128b_H.pdb does not match
(C20 - H1)
(more than 20 non-matching atom names)

WARNING 1 [file topol_popc.top, line 26]:
   10708 non-matching atom names
   atom names from topol_popc.top will be used
   atom names from popc_128b_H.pdb will be ignored

Analysing residue names:
There are:   128  Other residues
There are:  2460  Water residues
Analysing residues not classified as Protein/DNA/RNA/Water and splitting
into groups...
Number of degrees of freedom in T-Coupling group rest is 42105.00
Largest charge group radii for Van der Waals: 6.115, 5.932 nm
Largest charge group radii for Coulomb:   6.546, 6.115 nm

WARNING 2 [file em_st.mdp]:
   The sum of the two largest charge group radii (12.661407) is larger than
   rlist (0.90)

Calculating fourier grid dimensions for X Y Z
Using a fourier grid of 54x56x48, spacing 0.117 0.117 0.119
Estimate for the relative computational load of the PME mesh part: 0.45
This run will generate roughly 34 Mb of data

There were 2 warnings

---
Program grompp, VERSION 4.5.3
Source code file: grompp.c, line: 1563

Fatal error:
Too many warnings (2), grompp terminated.
If you are sure all warnings are harmless, use the -maxwarn option.

Can any one of you help me move from here?



The order of your [molecules] directive does not agree with the contents of the 
coordinate file.  What you're seeing is a mismatch between lipid atom names and 
water (O, H1, and H2).  The second warning may or may not be problematic.  If 
your charge groups are split across periodic boundaries, they will be 
reconstructed properly.  If your molecules are already whole, then you have a 
separate issue.


-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
* Only plain text messages are allowed!
* 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: specifying the direction of Pull in US

2012-07-06 Thread Justin A. Lemkul



On 7/6/12 9:52 AM, Raj wrote:

hi ,

When i applied distance I cant have control over the direction of the pull.
The ligand is not exactly pulled along the direction i meant to pull



You can specify the pull vector using the difference between the COM positions 
of the pull group and its reference.  Set those (x,y,z) values to pull_vec1.


-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
* Only plain text messages are allowed!
* 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] Crashes during protein-ligand simulation

2012-07-06 Thread Justin A. Lemkul



On 7/6/12 1:12 AM, James Starlight wrote:

Dear Gromacs users!

I have some problems with the simulation of protein-ligand complex
embedded in the ccl4-water environment. In addition there are some
crystallography waters (xw) embedded in the protein interiour of the
protein. I've done equilibration and minimisation of my system and run
it in NVT ensemble.

Finally I've already simulated this system in the apo form as well as
without XW and there were no any problems.

In the current case my system always crashed after 10-15 ns of
simulation with the errors like

Step 6651310, time 13302.6 (ps)  LINCS WARNING
relative constraint deviation after LINCS:
rms 0.060675, max 1.520945 (between atoms 3132 and 3130)
bonds that rotated more than 30 degrees:
  atom 1 atom 2  angle  previous, current, constraint length
3148   3147   90.00.1281   0.1483  0.1000
3150   3149   90.00.1084   0.1444  0.1000
3131   3130   90.00.1321   0.1325  0.1000
3132   3130   90.00.1067   0.2521  0.1000

---
Program mdrun_mpi.openmpi, VERSION 4.5.5
Source code file: /tmp/build/gromacs-4.5.5/src/mdlib/constr.c, line: 189


here both atoms 3132 and 3130 are from LIGAND.

During data analysing I didnt observed any serious artifacts in that
system. In addition RMSD both of protein and ligand were very stable.
Finally there are no fluctuations in energy or temperature. So I could
understand why this crasshes could occur. If I try to continue this
simulation from the crasshed checkpoint my simulation always goon but
within next 5-10ns I've always obtain second crash etc.

This is the last step from log file

DD  step 664  vol min/aver 0.758  load imb.: force  1.0%  pme
mesh/force 0.708

Step   Time Lambda
 66513300.00.0

Energies (kJ/mol)
   Angle   G96AngleProper Dih.  Improper Dih.  LJ-14
 5.78299e+011.24988e+042.10966e+031.81619e+038.94229e+01
  Coulomb-14LJ (SR)LJ (LR)  Disper. corr.   Coulomb (SR)
 4.65481e+048.27301e+04   -6.89535e+03   -2.15286e+03   -7.16722e+05
Coul. recip.  PotentialKinetic En.   Total Energy  Conserved En.
-1.72813e+05   -7.52733e+051.46708e+05   -6.06025e+05   -1.37045e+06
 Temperature Pres. DC (bar) Pressure (bar)   Constr. rmsd
 3.10968e+02   -9.74796e+012.52993e+021.55693e-05


Could you explain me what could be wrong with that system and what
addition data should I provide to help sheld light on that problem ?



If the addition of a ligand causes the simulation to crash (and the simulation 
runs normally in the apo form with and without crystal waters), then that sounds 
like a problem with the ligand topology or its initial placement.


What is the ligand?  How did you generate and validate its topology?  How did 
you place it?


-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
* Only plain text messages are allowed!
* 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: Problem with minimizing the energy of the solvated system

2012-07-05 Thread Justin A. Lemkul



On 7/5/12 5:51 PM, jonas87 wrote:

Protein-water.pdb has the following line:
CRYST1   75.324   75.324   75.324  60.00  60.00  90.00 P 1   1

But after adding the NA/CL ions with genion the CRYST1 line is missing
entirely from the output file protein-solvated.pdb
Could this be because of the bug u mentioned?


Sounds about right.


I'm using the latest version of gromacs, that means this bug was
reintroduced?


Likely it was not fully fixed; I recall a change to editconf, but maybe genion 
needs a look as well.



So it is ok to always output to .gro instead of .pdb?



Always.

-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
* Only plain text messages are allowed!
* 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: Problem with minimizing the energy of the solvated system

2012-07-05 Thread Justin A. Lemkul



On 7/5/12 5:05 PM, jonas87 wrote:

I'm following the tutorial exactly. Even have my files named the same way.
The contents of my minil.mdp (pcb is changed from no to xyz right before
running the energy minimazation of the solvated system):

title   = Energy Minimization   ; Title of run
cpp = /lib/cpp  ; Preprocessor
define  = -DFLEXIBLE
integrator  = steep ; Algorithm (steep = steepest descent 
minimization)
emtol   = 1.0   ; Stop minimization when the maximum force < 
1.0 kJ/mol
nsteps  = 500   ; Maximum number of (minimization) steps to 
perform
nstenergy   = 1 ; Write energies to disk every nstenergy steps
energygrps  = System; Which energy group(s) to write to disk
ns_type = simple; Method to determine neighbor list (simple, 
grid)
coulombtype = cut-off   ; Treatment of long range electrostatic 
interactions
rcoulomb= 1.0   ; long range electrostatic cut-off
rvdw= 1.0   ; long range Van der Waals cut-off
constraints = none  ; Bond types to replace by constraints
pbc = no; Periodic Boundary Conditions (yes/no)




You're probably losing the box information somewhere along the line while 
switching back and forth between .gro and .pdb.  The CRYST1 line in 
protein-solvated.pdb should agree with the previous settings.  If it doesn't, 
something has gone wrong.  I recall some previous version of Gromacs had issues 
reading and writing correct box vectors in .pdb files; I don't know which one. 
It's always safe to use .gro files for everything, though in principle, .pdb 
files should work as well.


-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
* Only plain text messages are allowed!
* 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] Problem with minimizing the energy of the solvated system

2012-07-05 Thread Justin A. Lemkul



On 7/5/12 3:15 PM, jonas87 wrote:

Chopping out a bunch of unnecessary stuff...


*editconf -f protein-EM-vacuum.pdb -o protein-PBC.gro -bt dodecahedron -d
1.0 *

Read 1202 atoms
No velocities found
 system size :  5.245  4.320  2.507 (nm)
 diameter:  5.532   (nm)
 center  :  0.010  0.047 -0.011 (nm)
 box vectors :  0.000  0.000  0.000 (nm)
 box angles  :   0.00   0.00   0.00 (degrees)
 box volume  :   0.00   (nm^3)
 shift   :  5.640  5.602  2.674 (nm)
new center  :  5.649  5.649  2.663 (nm)
new box vectors :  7.532  7.532  7.532 (nm)


Here, you're setting the box.  The dimensions look sensible.

(Chopping out more extraneous stuff)

Here's where something is going wrong.  Without seeing the contents of 
minim.mdp, it's hard to say exactly what, aside from a generic comment that the 
settings in the .mdp file are incompatible with the size of the box.  Typically 
tutorials work quite well, so you should ensure that you're following everything 
exactly and keeping track of all your files so nothing has gotten switched up.



grompp -v -f minim.mdp -c protein-solvated.pdb -p protein.top -o
protein-EM-solvated.tpr *

Ignoring obsolete mdp entry 'title'
Ignoring obsolete mdp entry 'cpp'

Back Off! I just backed up mdout.mdp to ./#mdout.mdp.3#
checking input for internal consistency...
processing topology...
Generated 141 of the 1176 non-bonded parameter combinations
Excluding 3 bonded neighbours molecule type 'Protein_chain_A'
Excluding 2 bonded neighbours molecule type 'SOL'
Excluding 1 bonded neighbours molecule type 'NA'
Excluding 1 bonded neighbours molecule type 'CL'
processing coordinates...
double-checking input for internal consistency...

ERROR 1 [file protein.top, line 7583]:
   ERROR: One of the box lengths is smaller than twice the cut-off length.
   Increase the box size or decrease rlist.


---
Program grompp, VERSION 4.5.5
Source code file: /build/buildd/gromacs-4.5.5/src/kernel/grompp.c, line:
1372

Fatal error:
There was 1 error in input file(s)
For more information and tips for troubleshooting, please check the GROMACS
website at http://www.gromacs.org/Documentation/Errors
---



Did you follow the link?  It's where you can find an answer to almost every 
error you'll encounter.  For instance:


http://www.gromacs.org/Documentation/Errors#The_cut-off_length_is_longer_than_half_the_shortest_box_vector_or_longer_than_the_smallest_box_diagonal_element._Increase_the_box_size_or_decrease_rlist

-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
* Only plain text messages are allowed!
* 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] question about energy minimization

2012-07-05 Thread Justin A. Lemkul



On 7/5/12 10:52 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
I want to minimize the energy of my structure.
I use the following mdp file:

define  = -DPOSRES
integrator  = steep
emtol   = 10
nsteps  = 5000
nstenergy   = 1
energygrps  = System
coulombtype = PME
rcoulomb= 0.9
rvdw= 0.9
rlist   = 0.9
fourierspacing  = 0.12
pme_order   = 4
ewald_rtol  = 1e-5
pbc = xyz

But it nevers runs till the end. It allways stops before the end and I get
no pdb file. When I use nsteps = 2000 there is no problem and I receive a
pdb file but it says that there were not enough steps to get under a
specific force value.

My question is: is there a possibility to make the minimization run
longer? Can I set a parameter that makes it run longer?



Set nsteps to either a very high value or to -1 (which says run infinitely 
long).  Note that several factors are working against you here, including the 
use of position restraints.  These disfavor the motion of whatever the 
restrained atoms are, which can impede energy minimization.  Also, assuming 
you're using single-precision Gromacs, an emtol of 10 is very low and unlikely 
to be achievable, especially using a single run of steepest descent EM.


-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
* Only plain text messages are allowed!
* 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] bilayers move apart by nanometers upon implementing dihedral restraints on lipid tails

2012-07-05 Thread Justin A. Lemkul



On 7/5/12 10:25 AM, himanshu khandelia wrote:

I am trying to implement dihedral restraints for lipids in a bilayer
using what is suggested here:

http://www.gromacs.org/Documentation/How-tos/Dihedral_Restraints

However, although the dihedral angles seem to be restrained fine, the
leaflets move apart by 10s of nanometers along +z over a nanosecond or
so, after which of course, the simulation crashes.

Can anyone suggest what I might be doing wrong?



Are the simulations stable in the absence of these restraints?  Are you sure 
this isn't just a matter of periodicity creating the illusion of large-scale 
diffusion?


-Justin


  version 4.5.4

In the mdp file:

;dihedral restraints
dihre   =  yes
dihre_fc   =  100

In the topology:

[ dihedral_restraints ]

17 18 19 20 11100012
18 19 20 21 11100012
19 20 21 22 11100012
20 21 22 23 11100012
21 22 23 24 11100012
...
etc.

Himanshu



--
====

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
* Only plain text messages are allowed!
* 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] COM RDF

2012-07-04 Thread Justin A. Lemkul



On 7/4/12 1:31 PM, Dr. Vitaly V. G. Chaban wrote:

Dear GROMACS people -

I am calculating radial distribution function between the centers of
mass of two large particles in a periodic box.

My command is --  g_rdf_455 -n 2drops.ndx -rdf mol_com -s topol.tpr -f
conf.gro

My index file contains two groups of atoms standing for the first and
the second supraparticles whose centers-of-mass I am interested in.
Since here I provide the "conf.gro" file rather than a trajectory, I
expect to get just one peak in my RDF, but instead I get multiple
peaks at different separation distances (please, see
http://i47.tinypic.com/24m826v.jpg). The groups of atoms in the index
file are separated spatially, if this may matter.

Would anyone kindly explain how the COM RDF function works? My
ultimate purpose is to depict how the distance between those two
centers-of-mass evolves in time.



You also need to use the -com flag to use the COM of the reference group.

-Justin




The -com flag changed the output (please, see
http://i48.tinypic.com/29uuo7b.jpg) but there is still a number of
embarassing peaks after 7nm. Given the positions of these
supraparticles, I would expect maximum about 10nm whereas the box side
is 26nm.


From the g_rdf help message --


The option -rdf sets the type of RDF to be computed. Default is for atoms or
particles, but one can also select center of mass or geometry of molecules or
residues. In all cases, only the atoms in the index groups are taken into
account.


I understand this in the following way. The centers-of-mass of the
groups in the index file are computed and further the calculation
proceeds as if these COMs were regular atoms. If I have only two
groups, I should get one maximum.

Am I not correct?



In theory, yes.  I have no idea how effective g_rdf is for this purpose.  From 
your initial description, if you're just trying to track COM distance over time, 
why not use g_dist?  That's precisely what it does.


-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
* Only plain text messages are allowed!
* 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] COM RDF

2012-07-04 Thread Justin A. Lemkul



On 7/4/12 1:07 PM, Dr. Vitaly V. G. Chaban wrote:

Dear GROMACS people -

I am calculating radial distribution function between the centers of
mass of two large particles in a periodic box.

My command is --  g_rdf_455 -n 2drops.ndx -rdf mol_com -s topol.tpr -f conf.gro

My index file contains two groups of atoms standing for the first and
the second supraparticles whose centers-of-mass I am interested in.
Since here I provide the "conf.gro" file rather than a trajectory, I
expect to get just one peak in my RDF, but instead I get multiple
peaks at different separation distances (please, see
http://i47.tinypic.com/24m826v.jpg). The groups of atoms in the index
file are separated spatially, if this may matter.

Would anyone kindly explain how the COM RDF function works? My
ultimate purpose is to depict how the distance between those two
centers-of-mass evolves in time.



You also need to use the -com flag to use the COM of the reference group.

-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
* Only plain text messages are allowed!
* 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] Nucleic acid simulation

2012-07-04 Thread Justin A. Lemkul



On 7/4/12 6:45 AM, Ravi Raja Merugu wrote:

Hello every one,

Im interesting in performing a MD for Protein - RNA complex , Can any
one suggest a good  tutorial.



Such systems do not differ significantly from simulations of simple proteins in 
water, since pdb2gmx can produce topologies for both proteins and nucleic acids. 
 The remaining workflow is basically the same.


-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
* Only plain text messages are allowed!
* 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] Regarding parameters for pull_groups in mdp file for the pull code

2012-07-04 Thread Justin A. Lemkul



On 7/4/12 2:31 AM, neeru sharma wrote:

Dear Gromacs Users,


I have some queries about the parameters in the .mdp file for the pull code.

If I want to pull my ligand, towards specific atom/group of atoms from
the protein, how am I supposed to mentioned these in the mdp file?

*
  pull = umbrella
  pull_geometry  = distance
  pull_start = yes
  pull_ngroups= 1
  pull_group0 = Ligand
  pull_group1 = Atom/group of atoms from the protein
*

Here, I want to specify the pull_group0 as "Ligand" and pull_group1 as
"Atoms of the protein".  My query is regarding these groups. Shall I
just write the name of the Ligand and Atoms (specifying the atom no)
or am I supposed to create a separate index file for each of them (one
for ligand and other for group of atoms) ?



All groups specified in the .mdp file must be either valid default groups or 
custom groups provided in an index file.


-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
* Only plain text messages are allowed!
* 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] PMF trails off to infinity.

2012-07-03 Thread Justin A. Lemkul



On 7/3/12 2:50 PM, Thomas Schlesier wrote:

think you encounter the problem, that you construct your pmf from a 3d
simulation and project it onto 1d, but do no correction.

For TI (if you constrain the distance in all three directions) the pmf is given 
by
V_pmf(r) = - \int [ F_c + 2/(beta*r) ] dr
with F_c the constraint force and \beta = 1/(k_B*T).

for umbrella sampling one should need the same factor
-2/beta \int 1/r dr
if one restraints the system in 3 directions.

Since 'g_wham' uses no *.tpr it can not know the value of 'pull_dim' one should
introduce the factor afterwards.



The input to g_wham (at least in modern Gromacs versions) is a list of .tpr 
files (passed to the -it flag) and a list of pullx or pullf files (with -ix or 
-if).  How is g_wham ignorant of these pull settings?


-Justin


There could also the problem that the dummy-atom interacts with the ligand, or
other clashes, but i don't think so, because the force looks too fine. If the
ligand would crash into something one should see a greatly increasing force.
like around 1000-2500 in your force-profile, when the protein still helds the
ligand in the binding pocket.



Additional comment:
'pull_geometry = distance'
is a really bad idea for pmf generation if one has an actual distance of 0 and
the system is not isotropic, since distance only knows distances and has no
ideas about directrions:

place a reference particle at the origin of a box
put a second particle on top of it
pull along x
pull along -x
in both cases you get the distance r=|x|=|-x|
if the system is isotropic, everything is fine
if system is anisotropic you get problems

on the nice side, this problem should affect simulations where one pulls a
ligand away from a binding pocket, only for the windows which is centered the
nearest to 0.
for mapping an ion-channel with a reference group in the middle of the channel
it's (/ it should be) far more worse, (if the system is not isotropic).
Somewhere in the archive is a longer explainition, i discussed the topic with
2-3 different people...






Am 03.07.2012 11:34, schrieb gmx-users-requ...@gromacs.org:

Basically what I'm doing is pulling a ligand out of a protein towards a
dummy atom, which has a mass of 1 and no charge. I've attached the a
portion .mdp files for both the smd portion and the umbrella sampling. I
know that the ligand gets very close possibly crashing into the dummy
atom. So from what you're saying, I'm thinking this might be the source
of the problem. - Laura ## smd.mdp## ; Generate velocites is
on at 300 K. gen_vel = yes gen_temp = 300.0 gen_seed = 123456 ;COM
PULLING ; Pull type: no, umbrella, constraint or constant_force pull =
umbrella ; Pull geometry: distance, direction, cylinder or position
pull_geometry = distance ; Select components for the pull vector.
default: Y Y Y pull_dim = Y Y Y ; Cylinder radius for dynamic reaction
force groups (nm) pull_r1 = 1 ; Switch from r1 to r0 in case of dynamic
reaction force pull_r0 = 1.5 pull_constr_tol = 1e-06 pull_start = yes
pull_nstxout = 10 pull_nstfout = 1 ; Number of pull groups pull_ngroups
= 1 ; Group name, weight (default all 1), vector, init, rate (nm/ps),
kJ/(mol*nm^2) pull_group0 =JJJ pull_weights0 = pull_pbcatom0 = 0
pull_group1 = CPZ pull_weights1 = pull_pbcatom1 = 0 pull_vec1 =
pull_init1 = pull_rate1 = -0.0006 pull_k1 = 800 pull_kB1 = ###
um.mdp ## ;COM PULLING ; Pull type: no, umbrella, constraint
or constant_force pull = umbrella ; Pull geometry: distance, direction,
cylinder or position pull_geometry = distance ; Select components for
the pull vector. default: Y Y Y pull_dim = Y Y Y ; Cylinder radius for
dynamic reaction force groups (nm) pull_r1 = 1 ; Switch from r1 to r0 in
case of dynamic reaction force pull_r0 = 1.5 pull_constr_tol = 1e-06
pull_start = yes pull_nstxout = 100 pull_nstfout = 100 ; Number of pull
groups pull_ngroups = 1 ; Group name, weight (default all 1), vector,
init, rate (nm/ps), kJ/(mol*nm^2) pull_group0 =JJJ pull_weights0 =
pull_pbcatom0 = 0 pull_group1 = CPZ pull_weights1 = pull_pbcatom1 = 0
pull_vec1 = pull_init1 = 0 pull_rate1 = 0 pull_k1 = 1000 pull_kB1 = On
07/02/2012 04:38 PM, Justin A. Lemkul wrote:

>
>  On 7/2/12 4:30 PM, Laura Kingsley wrote:

>>  Just for clarification, the PMF is read from right to left and the force
profile
>>  is read from left to right.
>>

>  The dramatic change in the magnitude and sign of the force, coupled with the
>  steady increase in PMF, indicates to me that some elements of your system are
>  crashing into one another.  In the absence of an accompanying explanation of
>  what you're doing (description of system, procedure with .mdp parameters,
etc)
>  that's the best I can offer.
>
>  -Justin
>

>>  On 07/02/2012 04:27 PM, Kingsley, Laura J wrote:

>>>  Here is a link to both the PMF profile and the force profil

Re: [gmx-users] inquiry about SAS

2012-07-03 Thread Justin A. Lemkul



On 7/3/12 10:49 AM, Turgay Cakmak wrote:

Hi all,

I am calculating SAS using g_sas of my system (several peptides in
water and ions, Na and Cl). I choose:
  for calculation group: non-water
  for output group: protein
  (400 out of 750 atoms were classified as hydrophobic)

When I plot the Area vs time graphs, both the hydrophobic-SAS and
hydrophilic-SAS decrease over the time. I think, the reason of the
hydrophobic-SAS decrease from 65 nm2 to 30nm2 is due to aggregation of
the peptides. But, I couldn't understant why hydrophilic-SAS is
decreasing from 40 nm2 to 20nm2 over time. Do the results make sense?



Sure, you've got both hydrophobic and hydrophilic contacts happening.


May be I mis-understand the meanings of the hydrophobic and
hydrophilic SAS. Please, someone could explain it to me, I would be so
glad.



These are merely subdivisions of total SASA.  Atoms are either polar 
(hydrophilic) or nonpolar (hydrophobic) and thus contribute to the total SASA 
based on their characteristics.


-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
* Only plain text messages are allowed!
* 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] question to mdp file

2012-07-03 Thread Justin A. Lemkul



On 7/3/12 8:25 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Ah okey. Yes sorry that was a typo...I ment nsteps.
So but is there a possibility to define a minimal step size so that the
minimization ends when the energy does not changes much any more?



Just set emtol to some unreasonably low value, and mdrun will end with this 
message (which is not really an error):


http://www.gromacs.org/Documentation/Errors#Stepsize_too_small.2c_or_no_change_in_energy._Converged_to_machine_precision.2c_but_not_to_the_requested_precision

You can set emstep to whatever size you like to reach this point.

-Justin




On 7/3/12 8:11 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
is it correct when I set the
nstep = -1  and
emtol = $number

that the minimization goes as long as the energy difference between the
previous step and this step is not lower as $number. And that there is
no
maximal stepsize?



No.  The value of emtol is the target for the maximum force to define
convergence, not the difference between previous and current steps, and
it's not
an energy term, it's force.  The maximum step size is always set in
emstep.
What's more, "nstep" is not a correct keyword, but "nsteps" is.

-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
* Only plain text messages are allowed!
* 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







--
============

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
* Only plain text messages are allowed!
* 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] Umbrella - Force constant

2012-07-03 Thread Justin A. Lemkul



On 7/3/12 8:41 AM, Steven Neumann wrote:

Dear Gmx Users,

Do you know or can you suggest some results based on the comparison of
the force constant in Umbrell Sampling? Any literature?



That would be lovely, but I've never seen such a thing.  One could probably 
write a book with all the test cases that would be required.  My gut tells me 
that you can't generalize too much in terms of pulling simulations - the 
approach depends on what is being pulled (small molecule, peptide, large 
protein), what the medium is (water, membrane, etc), and what the interacting 
partner is (protein surface, ion channel, binding pocket).



As far as I understand when you use the same staring coordinates (from
the same pulling simulation) for windows but you just change the force
constant (e.g. from 500 to 2000 kJ/mol nm2) you should increase number
of windows (for f=2000) as smaller force constant will cover wider
neigboruring distances - that makes sense.
I am curious whether the final result will be the same? I guess with
stronger force it will converge faster but more windows are required.
is it the only one difference?



Without a systematic comparison, it's hard to say, but in theory if one samples 
sufficiently and has good overlap between neighboring windows, the results 
should converge to the same answer.


If someone knows of some applicable literature that has done such comparisons, 
please post a reference.  I'd love to see it.  Most SMD and US methodology is 
written with hand-waving explanations as to what the authors did and why it 
worked, and I have a suspicion that most reviewers don't have a better idea so 
they can't refute such claims.


-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
* Only plain text messages are allowed!
* 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] Implicit solvent setup for large protein

2012-07-03 Thread Justin A. Lemkul



On 7/3/12 8:57 AM, Ramon Crehuet Simon wrote:

Dear all,
I know several questions about implict solvent have already been asked in this 
list. I think and hope the question I have has not been raised. Forgive me if I 
am wrong.
I have read that all-to-all kernels are the best option when doing implicit 
solvent. Otherwise one should use large cutoffs.
I wanto to simulate a tetramer which has more than 1200 residues. I'm afraid 
all-to-all interactions are much too expensive in this case. So I should use 
largecutoffs. My question is, how large? Or, to be more precise, what should I 
monitor to check the the cutoffs are large enough?
As an extra questions, if I don't use all-to-all kernels, is it still worth 
using the GPU to accelerate the calculation?


In my experience, using long cutoffs leads to unstable trajectories and very 
poor energy conservation.  The all-vs-all kernels are the only ones I use for 
implicit calculations.


-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
* Only plain text messages are allowed!
* 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] question to mdp file

2012-07-03 Thread Justin A. Lemkul



On 7/3/12 8:11 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
is it correct when I set the
nstep = -1  and
emtol = $number

that the minimization goes as long as the energy difference between the
previous step and this step is not lower as $number. And that there is no
maximal stepsize?



No.  The value of emtol is the target for the maximum force to define 
convergence, not the difference between previous and current steps, and it's not 
an energy term, it's force.  The maximum step size is always set in emstep. 
What's more, "nstep" is not a correct keyword, but "nsteps" is.


-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
* Only plain text messages are allowed!
* 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] mdrun no structural output

2012-07-03 Thread Justin A. Lemkul



On 7/3/12 5:40 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
I wanted to do a minimization with mdrun but the only output I get is:
3m71_minim.edr
3m71_minim.log
3m71_minim.trr
But no structure file like .pdb i.e.

There was no error in the step before where I prepared the input file with
grompp. My .mdp file looks like this:

define  = -DPOSRES
integrator  = steep
emtol   = 10
nsteps  = 5000
nstenergy   = 1
energygrps  = System
coulombtype = PME
rcoulomb= 0.9
rvdw= 0.9
rlist   = 0.9
fourierspacing  = 0.12
pme_order   = 4
ewald_rtol  = 1e-5
pbc = xyz



And also in this step there was no error.
The end of the 3m71_minim.log looks like this:

Step   Time Lambda
2647 2647.00.0

Energies (kJ/mol)
Bond  AngleProper Dih.  Improper Dih.  LJ-14
 5.35530e+022.26340e+031.17875e+041.31106e+024.73257e+03
  Coulomb-14LJ (SR)   Coulomb (SR)   Coul. recip. Position Rest.
 6.58188e+049.41710e+04   -7.10195e+05   -1.62296e+059.89235e+02
   Potential Pressure (bar)
-6.92062e+05   -3.32594e+03


And the stdout looks like this:

Step= 2646, Dmax= 3.6e-03 nm, Epot= -6.92039e+05 Fmax= 5.12625e+03, atom=
1022
Step= 2647, Dmax= 4.3e-03 nm, Epot= -6.92098e+05 Fmax= 3.72457e+03, atom=
1022
Step= 2647, Dmax= 4.3e-03 nm, Epot= -6.92062e+05 Fmax= 1.51933e+03, atom=
1022
Step= 2648, Dmax= 5.2e-03 nm, Epot= -6.92099e+05 Fmax= 4.25221e+03, atom=
1022
Step= 2648, Dmax= 5.2e-03 nm, Epot= -6.92041e+05 Fmax= 6.45781e+03, atom=
1022

The command for the mdrun was:

mpirun -n 2 $gromacsPath/mdrun_mpi -c $path/3m71_minim.pdb -compact
-deffnm $path/3m71_minim -s $path/3m71_minim_ion.tpr  -v
2>>$path/minLogErr 1>>$path/minLogOut

Can you please tell me whats wrong?



When EM is done, mdrun prints very clear messages indicating convergence (or 
lack thereof) with information about maximum force and potential energy.  If 
you're not seeing this information, mdrun isn't done or somehow got killed or 
hung up.


-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
* Only plain text messages are allowed!
* 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] bindng energy & H bond energy calculation

2012-07-03 Thread Justin A. Lemkul



On 7/3/12 2:24 AM, Ravi Raja Merugu wrote:

Dear all,
Can any one help e regarding
1.)  I want to find the binding energy of  ligand with the protein , and


Umbrella sampling or free energy calculations can determine binding energies.


2.)  Hydrogen Bond energy of the protein ligand interactions  (for
H-bond ).. and


This term will probably be related to electrostatic interactions, but I don't 
know of any way to directly extract a "hydrogen bonding energy" in a literal 
sense.  Most force fields don't include such a term.



3 )  if possible H bond distances along the time scale(say 10ns).
Is there any way


g_dist and/or g_hbond


4 ) to measure the contribution of hydrophobic / philic interactions
of the protein active site with ligands binding energy..


Free energy calculations, decoupling van der Waals and Coulombic terms 
separately, will tell you this information.



5) to calculate the distance between the residues(let say Tyr and Arg)
in my protein. (Actually I want to check a  key interaction between
them)


g_dist

-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
* Only plain text messages are allowed!
* 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] problems with editconf

2012-07-03 Thread Justin A. Lemkul



On 7/3/12 5:10 AM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi everybody,
I want to put my protein in the box with editconf but when I look at it it
is always at the border of the box and not at the center. I tried it with
those two commands:

editconf -f 3m71.gro -o 3m71_box.gro -center 4.59340   4.59470   5.17330
-c -bt dodecahedron -d 1.0 2>>logErr 1>>logOut

editconf -f 3m71.gro -o 3m71_box.gro -c -bt dodecahedron -d 1.0 2>>logErr
1>>logOut

With the first command it is at the upper border and with the second
command it is at the right border. The coordinates of the firs command are
from the 3m71.gro file.

Can you please tell me what is wrong?



Nothing.  When visualizing the unit cell, it will appear as a triclinic box with 
the protein in a seemingly random location.  You can re-wrap the unit cell with 
trjconv -pbc mol -ur compact (which requires a .tpr file) to see the 
dodecahedral unit cell with everything placed as you would expect.


-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
* Only plain text messages are allowed!
* 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] PMF trails off to infinity.

2012-07-02 Thread Justin A. Lemkul



On 7/2/12 4:51 PM, Laura Kingsley wrote:

Basically what I'm doing is pulling a ligand out of a protein towards a dummy
atom, which has a mass of 1 and no charge. I've attached the a portion .mdp
files for both the smd portion and the umbrella sampling. I know that the ligand
gets very close possibly crashing into the dummy atom. So from what you're
saying, I'm thinking  this might be the source of the problem.



Sounds like it to me.

-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
* Only plain text messages are allowed!
* 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] PMF trails off to infinity.

2012-07-02 Thread Justin A. Lemkul



On 7/2/12 4:30 PM, Laura Kingsley wrote:

Just for clarification, the PMF is read from right to left and the force profile
is read from left to right.



The dramatic change in the magnitude and sign of the force, coupled with the 
steady increase in PMF, indicates to me that some elements of your system are 
crashing into one another.  In the absence of an accompanying explanation of 
what you're doing (description of system, procedure with .mdp parameters, etc) 
that's the best I can offer.


-Justin


On 07/02/2012 04:27 PM, Kingsley, Laura J wrote:

Here is a link to both the PMF profile and the force profile:

http://s1064.photobucket.com/albums/u370/laurakingsley/?action=view¤t=pull_fig.jpg



-Original Message-
From: gmx-users-boun...@gromacs.org [mailto:gmx-users-boun...@gromacs.org] On
Behalf Of Justin A. Lemkul
Sent: Monday, July 02, 2012 3:56 PM
To: Discussion list for GROMACS users
Subject: Re: [gmx-users] PMF trails off to infinity.



On 7/2/12 3:53 PM, Laura Kingsley wrote:

Basically what I'm doing is pulling a ligand toward a dummy atom. I pull the
ligand until its COM is very close (a few angstroms) from the dummy atom. Could
this be causing this behavior? What information would help? I'm having a tough

Collisions between any elements of your system would explain the problem, if it
looks like what I'm picturing in my head :)


time getting the figures attached but may be I could email them directly to you?
Let me know.


Bullet point #4 here provides a suggestion:

http://www.gromacs.org/Support/Mailing_Lists#Mailing_List_Etiquette

-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
* Only plain text messages are allowed!
* 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] PMF trails off to infinity.

2012-07-02 Thread Justin A. Lemkul



On 7/2/12 3:53 PM, Laura Kingsley wrote:

Basically what I'm doing is pulling a ligand toward a dummy atom. I pull the
ligand until its COM is very close (a few angstroms) from the dummy atom. Could
this be causing this behavior? What information would help? I'm having a tough


Collisions between any elements of your system would explain the problem, if it 
looks like what I'm picturing in my head :)



time getting the figures attached but may be I could email them directly to you?
Let me know.



Bullet point #4 here provides a suggestion:

http://www.gromacs.org/Support/Mailing_Lists#Mailing_List_Etiquette

-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
* Only plain text messages are allowed!
* 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] PMF trails off to infinity.

2012-07-02 Thread Justin A. Lemkul



On 7/2/12 9:54 AM, Laura Kingsley wrote:

Hello,

I am using steered MD and umbrella sampling to generate a PMF profile for
pulling a small ligand 3 nm.  As I pull the ligand from 3A toward 0 A, the PMF
starts  a dramatic uphill climb. This does not agree with the force profile
which seems to level out. Any ideas about what might be going wrong here? I am
thinking this probably isn't correct, but I don't know where I've messed up.
Thanks! I can attach the graphs if necessary.



We'll probably need a better description of what you're doing, along with any 
applicable figures.  It sounds to me like your ligand is crashing into something 
in your system, but that's just a guess based on incomplete information.


-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
* Only plain text messages are allowed!
* 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] charmm36 parameters

2012-07-02 Thread Justin A. Lemkul
HBL 1   1.0080000.10A   0.235197261589  0.092048 ; Not 
in C27
CCL 6   12.011000.34A   0.356359487256  0.29288  ; Not 
in C27
CL  6   12.011000.90A   0.356359487256  0.29288
CTL16   12.011000.14A   0.405358916754  0.08368
CTL26   12.01100-0.18   A   0.358141284692  0.234304
CTL36   12.01100-0.27   A   0.363486677001  0.326352
CTL56   12.01100-0.35   A   0.367050271874  0.33472
CEL16   12.01100-0.15   A   0.372395664183  0.284512
CEL26   12.011000.000   A   0.370613866746  0.267776 ; 
partial charge def not found
OBL 8   15.999400   -0.52   A   0.302905564168  0.50208
OCL 8   15.999400   -0.76   A   0.302905564168  0.50208
O2L 8   15.999400   -0.78   A   0.302905564168  0.50208
OHL 8   15.999400   -0.66   A   0.315378146222  0.6363864
OSL 8   15.999400   -0.49   A   0.293996576986  0.4184 ; 
Different to C27
OSLP8   15.999400   -0.57   A   0.293996576986  0.4184 ; Not in 
C27
NH3L7   14.00700-0.30   A   0.329632525712  0.8368
NTL 7   14.00700-0.60   A   0.329632525712  0.8368
SL  16  32.06   1.33A   0.374177461619  1.96648
PL  15  30.974000   1.50A   0.38308644882.44764
; The following atom types are NOT part of the CHARMM distribution
; atomtypes for additional water models
#ifdef HEAVY_H
OWT38   9.951400-0.834  A   3.15058e-01 6.36386e-01 
; TIP3p O
HWT31   4.0320000.417   A   0.0 0.0 ; TIP3p H
OWT48   9.9514000.0 A   3.15365e-01 6.48520e-01 
; TIP4p O
HWT41   4.0320000.52A   0.0 0.0 ; TIP4p H
MWT40   0.00-1.04   A   0.0 0.0 ; TIP4p vsite
OWT58   9.9514000.0 A   3.12000e-01 6.69440e-01 
; TIP5p O
HWT51   4.0320000.241   A   0.0 0.0 ; TIP5p H
MWT50   0.00-0.241  A   0.0 0.0 ; TIP5p vsite
OW  8   9.951400-0.82   A   3.16557e-01 6.50194e-01 
; SPC 0
HW  1   4.0320000.41A   0.0 0.0; SPC H
#else
OWT38   15.999400   -0.834  A   3.15058e-01 6.36386e-01 
; TIP3p O
HWT31   1.0080000.417   A   0.0 0.0 ; TIP3p H
OWT48   15.999400   0.0 A   3.15365e-01 6.48520e-01 
; TIP4p O
HWT41   1.0080000.52A   0.0 0.0 ; TIP4p H
MWT40   0.00-1.04   A   0.0 0.0 ; TIP4p vsite
OWT58   15.999400   0.0 A   3.12000e-01 6.69440e-01 
; TIP5p O
HWT51   1.0080000.241   A   0.0 0.0 ; TIP5p H
MWT50   0.00-0.241  A   0.0 0.0 ; TIP5p vsite
OW  8   15.999400   -0.82   A   3.16557e-01 6.50194e-01 
; SPC O
HW  1   1.0080000.41A   0.0 0.0 ; SPC H
#endif
; special dummy-type particles
MNH30   0.000.00A   0.0 0.0
MNH20   0.000.00A   0.0 0.0
MCH30   0.000.00A   0.0 0.0
MCH3S   0   0.000.00A   0.0 0.0
; Ions and noble gases (useful for tutorials)
Cu2+29  63.546002.00A   2.08470e-01 4.76976e+00
Ar  18  39.948000.00A   3.41000e-01 2.74580e-02
; Added by DvdS 05/2005 copied from GROMACS force field.
SI  14  28.080000.00A   3.38550e-01 2.44704e+00




--
====

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
* Only plain text messages are allowed!
* 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] number of contact

2012-07-02 Thread Justin A. Lemkul



On 7/2/12 8:18 AM, Turgay Cakmak wrote:

Dear all,

Does the "number of contact calculated by g_mindist" mean that
interaction between two groups? If so, could you kindly check the
reliability of the case I have written below.



If the atoms are in contact within a given cutoff distance, then you can pretty 
safely assume they are interacting (provided the cutoff is not exorbitantly 
long, but 0.4 nm seems reasonable).


-Justin


I have 10 peptides (5 of them are X-peptide and 5 of them are
Y-peptide) in a box filled with water. After MD simulation (production
run), I saw all the peptides come together. Now, in order to
understand whether
hydrophobic forces are important for this aggregation or not, I
calculated number of contacts between 2 index groups (First group
consist of hydrophobic residues (FFAA) belong to X-peptide and
second group again consist of hydrophobic residues (FFAA) belong to
Y-peptide) using the following command;

g_mindist -f .xtc -s .tpr -n .ndx -on -d 0.4

Then, when I looked the number of contact versus time graph, there is
an increase in the number of contact over time. (I think this makes
sense in the aggregation case.)
Thanks in advance,

Turgay



--
====

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
* Only plain text messages are allowed!
* 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: COM Pulling

2012-07-01 Thread Justin A. Lemkul



On 7/1/12 12:02 AM, Raj wrote:

Dear Justin,

I've tried what you have suggested. I have used the pull code and gave 4
amino acid residues as a reference group. The location of the groups are one
below the ligand (in the protein core) and 3 were on the above ( towards the
active site gorge). When i used the code



This sounds like a poor approach.  If you have some residues that are internal 
and some that are external, mdrun won't do any productive pulling.  Define one 
group in a concerted area as a reference group and pull with respect to it to 
move the ligand.


-Justin


pull = umbrella
pull_geometry  = distance  ; simple distance increase
pull_start  = yes   ; define initial COM distance > 0
pull_ngroups= 1
pull_group0 = DRG
pull_group1 = reference groups
pull_rate1   = 0.003  ; 0.01 nm per ps = 10 nm per ns
pull_k1   = 1000  ; kJ mol^-1 nm^-2

the ligand migrated more and more towards the protein , not coming out of
the protein through the channel i defined by referring the amino acid groups

when I tried with pull_geometry = position, the pull_vec i gave as 0 0 1 and
pull_initial = 0.1 but the grompp ends up with an error saying pull_vec can
not be zero.

when i was trying to use pull_geometry = direction the system blowed up. but
the pull_geometry = distance worked well.

please give me a suggestion to resolve the issue


--
View this message in context: 
http://gromacs.5086.n6.nabble.com/COM-Pulling-tp4998944p4998989.html
Sent from the GROMACS Users Forum mailing list archive at Nabble.com.



--
========

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
* Only plain text messages are allowed!
* 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] Berger lipid

2012-06-30 Thread Justin A. Lemkul



On 6/30/12 11:03 AM, Shima Arasteh wrote:

Thanks Justin.
Yes, you are right. You wrote the tutorial for DPPC and I know that.

ok, I saw the temperatures upper than 271 K is proper for POPC. It's acceptable that 310 
mentioned in Peter's paper is correct. But I'm wondering how they chose "310" 
K? In his article, he explains that their study was done in mammalian cells.




310 K is physiological temperature for the human body.  Peter, please correct me 
if this was not the intent of your study.


-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
* Only plain text messages are allowed!
* 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] Berger lipid

2012-06-30 Thread Justin A. Lemkul



On 6/30/12 10:33 AM, Shima Arasteh wrote:

Dear Peter, Thanks for your link and the article. I'd like to know more about
your paper. You've mentioned in it that the temperature of POPC equilibrated,
is 310 K. As I saw in Justin's tutorial , 323 K is proper, however it was a


I'll let Peter address the question directed to him, but I want to refute 
something stated here.  One should not say that for membrane simulations "323 K 
is proper."  I explain in the tutorial why this particular temperature is used 
in the context of DPPC only.  I also provide a table of phase transition 
temperatures for several lipids to explain the reason why the elevated 
temperature was required in this case.


-Justin


different system simulated and also many parameters are not the same as your
system . Would you telling me about the reason of 310 K?

Thanks in advance


Sincerely, Shima


- Original Message - From: Peter C. Lai  To: Shima
Arasteh ; Discussion list for GROMACS users
 Cc: Sent: Friday, June 29, 2012 2:54 PM Subject: Re:
[gmx-users] Berger lipid

yes
http://www.frontiersin.org/Bioinformatics_and_Computational_Biology/10.3389/fgene.2012.00061/abstract

 The files are here: http://uab.hyperfine.info/~pcl/files/popc36/

On 2012-06-28 09:58:26PM -0700, Shima Arasteh wrote:

Yes, I remember now...you are right :) But I didn't know
the linked you sent me, was your own output! However  I wanted to know if
it is necessary to produce the .itp file on my own or not.

I still have this link, so will cite to you. It would be a good idea to see
its package in lipidbook too. Thanks Peter




Sincerely, Shima


 From: Peter Lai  To:
Discussion list for GROMACS users  Sent: Friday,
June 29, 2012 7:14 AM Subject: RE: [gmx-users] Berger lipid

Uh didn't we go through all of this like more than a month ago? I published
a paper using C36 POPC and even a linked to my popc.itp for it on this
list...

Of course Shima is welcome to pdb2gmx his own POPC, which I am fairly
certain will result in an identical file...

Lipidbook seems to only have C36 POPE. I guess maybe I will upload ours.
 From:
gmx-users-boun...@gromacs.org [gmx-users-boun...@gromacs.org] on behalf of
Justin A. Lemkul [jalem...@vt.edu] Sent: Thursday, June 28, 2012 7:56 PM
To: Discussion list for GROMACS users Subject: Re: [gmx-users] Berger
lipid

On 6/28/12 8:54 PM, Shima Arasteh wrote:

Yes, I know that as studied the Kalp15 tutorial. Sorry, the last question
:) :

DO I need to run pdb2gmx to get the top file of POPC in CHARMM36? Is it
ok?  Because I see that POPC.itp is also required for simulation of
protein in bilayer.



You need a topology of some sort.  It depends on what parameters you have
on hand.  If you do not have popc.itp from anywhere, then you need to
generate it somehow.  If it is present in the .rtp file for CHARMM36 that
you have, then you can run pdb2gmx on it.

-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 * Only plain text
messages are allowed! * 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 list
gmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users *
Only plain text messages are allowed! * 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 list
gmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users *
Only plain text messages are allowed! * 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




--


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

Re: [gmx-users] Re: .top file incoherent with the values in the ffbonded.itp file

2012-06-29 Thread Justin A. Lemkul



On 6/29/12 6:10 PM, sreeta.g wrote:

Hi Justin

I have changed all the force field parameters as below, (I have already
shown my ffbonded.itp file)
ffnonbonded.itp
; name  bond_typemasscharge   ptype  sigma  epsilon
  opls_966   CA612.01100  0.240   A3.5000e-01
3.3600e-01
  opls_967   CA612.01100  0.180   A3.5000e-01
3.3600e-01
  opls_968   HB1   11.00800   0.060   A2.5700e-01
2.1000e-01
  opls_969   HB2   11.00800   0.060   A2.5700e-01
2.1000e-01
  opls_970   HB3   11.00800   0.060   A2.5700e-01
2.1000e-01
  opls_971   HA1   11.00800   0.060   A2.5700e-01
2.1000e-01
  opls_972   HA2   11.00800   0.060   A2.5700e-01
2.1000e-01

atomname2type.n2t
Copls_136-0.120 12.011  4C 0.153   H 0.11H 0.110   H 0.110
Copls_966 0.240 12.011  4C 0.153   H 0.11H 0.110   O 0.143
Copls_135-0.180 12.011  4C 0.153   H 0.11H 0.110   H 0.110
Oopls_154-0.700 15.999  2C 0.143   H 0.097

aminoacid.rtp
[ Pva ]
   [ atoms ] ; see atomtypes.atp for explainations
 CBopls_136 -0.120 1
 HB1   opls_968 0.060 1
 HB2   opls_969 0.060 1
 CAopls_966 0.240 2
. (contd)

so I have maintained consistency in my atom types in the ffbonded.itp files.
When I am running Grompp, it gives 7988 errors and I have 9088 atoms in my
polymer chain. All the errors mentioned are No default Ryckaert-Bell. types,
which I explained to you in my first mail.

So I am still unsure as to what my problem is and how I should deal with it.



The only way to diagnose the problem is to check your self-consistency. 
Probably most of the errors are redundant.  Choose one of them (it tells you the 
line where the error occurs), write down the atom numbers involved and then look 
up which atom types (not names) to which these numbers pertain.  That is the 
missing interaction type that needs to be defined in ffbonded.itp.


I can't tell what's wrong based on the snippets of files shown.  You can, based 
on having a very careful look through all of your files.


-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
* Only plain text messages are allowed!
* 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] error with grompp

2012-06-29 Thread Justin A. Lemkul



On 6/29/12 5:38 PM, reising...@rostlab.informatik.tu-muenchen.de wrote:

Hi Justin,

yes I removed all the old resulting files and did everything again. So now
there is the topology and coordinate file with only NA and CL and not NA+
or CL-.
I also checked whether the molecules are listed in the same order as in
the .gro file and it is the case. So that is also correct.

What do you mean with:


What does your [molecules] directive specify?


My [molecules] part in the topology file looks like this:

[ molecules ]
; Compound#mols
Protein_chain_A 1
DUM 20088
SOL 13428
NA   29
CL   29



I see no reason this would not work.  However, I just noticed from your previous 
message:



; Include water topology
#include "amber03.ff/tip3p.itp"



; Include topology for ions
#include "amber03.ff/ions.itp"

#include "amber03.ff/spc.itp"
#include "amber03.ff/dum.itp"



You're using two different water models, so things are getting overridden there. 
 With AMBER03, you should be using TIP3P, not SPC.  The conflicting water 
models suggest you've made manual modifications to the topology.  Perhaps there 
is some error as a result.


-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
* Only plain text messages are allowed!
* 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: .top file incoherent with the values in the ffbonded.itp file

2012-06-29 Thread Justin A. Lemkul



On 6/29/12 4:44 PM, sreeta.g wrote:

Hi Justin
Thank you for your reply.
However, when I am using the grompp command, the topol.tpr file is not being
formed due to a fatal error. This fatal error is the cumulative of ' No
default Ryckaert-Bell.' types for many atoms in the polymer chain.


And as I suggested before, you should investigate which lines are causing 
complaints to diagnose why grompp is failing.



Also, regarding the comment on the atom types and atom names, I have tried
Google searching to find the difference, but in vain.


An atom name is an arbitrary label.  An atom type is a designator that supplies 
various nonbonded parameters (listed in the [atomtypes] directive of 
ffnonbonded.itp.  The files you posted before seem to define your bonded 
parameters based on the names present in the coordinate file.  This is an 
incorrect approach.  Refer to existing force fields and the manual for examples.


-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
* Only plain text messages are allowed!
* 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] .top file incoherent with the values in the ffbonded.itp file

2012-06-29 Thread Justin A. Lemkul



On 6/29/12 3:55 PM, sreeta.g wrote:

Hello
My A.gro file is:
  1PvB CB1   3.109   2.784   0.803  0.  0.  0.
 1PvBHB12   3.156   2.880   0.843  0.  0.  0.
 1PvBHB23   2.997   2.797   0.819  0.  0.  0.
 1PvBHB34   3.129   2.765   0.693  0.  0.  0.
 1PvB CA5   3.173   2.673   0.879  0.  0.  0.
... (contd.)

My ffbonded.itp file has:
[ bondtypes ]
;   ij   func   b0(bond length)  kb (force const)
CBHB1  10.11   284512
CBHB2  10.11   284512
CBHB3  10.11   284512.0
... (contd)

[ angletypes ]
;  i jk   func th0cth
   HOOH   CA1   105320
   CACB   CA1   109.45 482.3
   CBCA   CB1   109.45 482.3
   OHCA   CB1   107.8  460
... (contd)

[ dihedraltypes ]
;  ijkl   funcRB coefficients
OH  OH   CA   CB3 3.09.0 0.0   -12.0
0.0   0. ;
CA  CB   CA   CB3 5.75000   17.25000 0.0   -23.0
0.0   0. ;
. (contd)

However, my topol.top file looks like this:
[ bonds ]
;  aiaj functc0c1c2c3
 1 2 1
 1 3 1
 1 4 1
 1 5 1
 112 1
(contd)

[ angles ]
;  aiajak functc0c1c2
c3
 2 1 3 1
 2 1 4 1
 2 1 5 1
 2 112 1
 3 1 4 1
 3 1 5 1
. (contd)

[ dihedrals ]
;  aiajakal functc0c1c2
c3c4c5
 2 1 5 6 3
 2 1 5 7 3
 3 1 5 6 3
 3 1 5 7 3
 4 1 5 6 3
 4 1 5 7 3
... (contd)

I have two issues:
1) I have mentioned the bond lengths and bond constants in the [bondtypes]
and [angletypes ]directory of the ffbonded.itp file, however in the [bonds]
[ angles ]directory in the topol.top file the Ryckaert-Bellemans
coefficients are used, which are assigned as blanks. This is giving the
error No default Ryckaert-Bell. types (multiple times).


The "blanks" in the topology are totally normal.  When assembling the .tpr file, 
grompp will look up the appropriate values in ffbonded.itp that match with the 
atom types it finds.  If you are getting errors about missing parameters, you 
have failed to assign some interaction type.  The error message should contain 
the offending lines.  Note, too, that the parameters specified in ffbonded.itp 
are assigned by atom *type* and not atom *name* - perhaps this is the origin of 
your problem.



2) Also, the RB coefficients mentioned in the [dihedraltypes] of the
ffbonded.itp are not used in the [dihedrals] directory of the topol.top
file. This is again writing blanks for the RB coefficients.
What is the issue and why is my topl.top file inconsistent with the inputs
given in the ffbonded.itp file?


The blanks are normal; see above comments.

-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
* Only plain text messages are allowed!
* 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] pdb file of polymer

2012-06-29 Thread Justin A. Lemkul



On 6/29/12 2:21 PM, Parvez khan wrote:

Hi, I am trying to do polymer simulation with gromacs. I am new to
gromacs and trying to construct topology for a system of polymer
chains. My problem is that i am facing difficulties to creat pdb file
for polymer chain containing 1000 monomers. I have used PRODRG server
but it gives me a pdb and topology file up to maximum 41 monomers
chain. I am not understanding what wrong with PRODRG server. Is there


Nothing is wrong.  It's just that PRODRG limits how many atoms are allowed in 
the input structure.



any server or tool for generating db file.


There are lots of ways to generate a coordinate file.  Here are some 
suggestions:

http://www.gromacs.org/Documentation/File_Formats/Coordinate_File#Sources

-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
* Only plain text messages are allowed!
* 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] Structure optimization failure

2012-06-29 Thread Justin A. Lemkul



On 6/29/12 11:10 AM, massimo sandal wrote:


On 29 Jun 2012 15:09, "Justin A. Lemkul" mailto:jalem...@vt.edu>> wrote:
 >
 >
 >
 > On 6/29/12 9:04 AM, massimo sandal wrote:
 >>
 >>
 >> On 29 Jun 2012 14:08, "Justin A. Lemkul" mailto:jalem...@vt.edu> <mailto:jalem...@vt.edu <mailto:jalem...@vt.edu>>>
 >>
 >>  >
 >>  > The settings in my tutorial are for use with OPLS-AA and are thus not
 >> suitable for a simulation with CHARMM.  Cutoffs and other aspects will be
different.
 >>
 >> This is interesting. In general, where can you find optimal mdp settings for
 >> each force field?
 >>
 >
 > By reading the primary literature for how each force field was derived.  The
parameters are only guaranteed to work within the context of correct cutoffs,
electrostatic schemes, neighborlist update interval, water model used, etc. Some
older force fields have been demonstrated to work well with PME even if
originally derived with RF or some other method, but generally the other
settings must adhere to the original parameterization.  Some force fields can be
very sensitive to these incorrect settings.
 >
 >

Yes, I know, but I hoped someone could have put a table online or wrote a book
chapter somewhere about this.

Wouldn't it be nice to create a table of "standard settings" for each forcefield
in the gmx documentation (with lit references of course)?



Well, anyone is welcome to submit anything they feel would be useful... ;)

I have considered this in the past, but have rejected the thought every time it 
comes to mind.  It would be nice if there was a central repository of such 
parameters, but then likely a simulation becomes "do this, not that" with no 
required thought from the end user.  If you make simulations a black box, 
everyone expects them to automatically work without fail.


The other objection I have always had is the continual improvement of the field. 
 If Gromacs puts out an "official" or "standard" list of what to do, it is 
always subject to change and then becomes incumbent upon a Gromacs contributor 
to verify its accuracy frequently.  My personal fear is that this becomes an 
untenable task.  Then in the case of an error, someone's whole Ph.D. could go 
out the window...


In the end, I think it's always safe to ask the user to do a bit of legwork to 
understand the force field he or she wishes to use.  The knowledge gained from 
an hour of reading will save a lot of potential headaches.


I would reject any attempt to include such information in the manual (lest its 
settings become "official"), but if someone wants to put up a wiki page on the 
website with clear warnings and disclaimers, they are welcome to do so.


Just my $0.02.

-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
* Only plain text messages are allowed!
* 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] problem with git/4.6

2012-06-29 Thread Justin A. Lemkul



On 6/29/12 11:08 AM, Michael Brunsteiner wrote:




You should have it.  In CMakeLists.txt, PROJECT_VERSION should be set to
"4.6-dev" so you can check that.


ithat what i looks like ... i now get:

prompt> mdrun_d
[...]
:-)  VERSION 4.6-dev-20120629-9c6be1c  (-:
[...]

which gives me:
g_bar_d -f mdv*.xvg -b 100
[...]
total   0.000 -  1.000,   DG -9.00 +/-  0.15

with vanilla 4.5.5 i get:
g_bar_d -f mdv*.xvg -b 100
[...]
total   0.000 -  1.000,   DG -9.00 +/-  0.15

there's two points:

1) this value is fine, as, once combined with the coulomb part, it is very 
close to literature and
exptl data. however, in both cases  get this warning ...  Second Law of 
Thermodynamics etc...
and riduculously high s_B values at lambda=0.



This indeed sounds like a bug worth fixing since it is present in the 
release-4-6 branch.



2) the reason i saw different DeltaG values earlier was that i had included (by 
mistake)
the dhdl file from a single additional window at lambda=0 that was produced 
with exactly the same
input apart from the value of 0.3 (instead of 0.5) for sc_sigma ... as far as i 
understand things (not very
well perhaps) the resulting free energy difference should NOT depend on this 
parameter ... so neither should
the result of g_bar ... but maybe it does when using foreign lambdas as i did 
here ? if this is so
it might be a good idea making g_bar read tpr files to give a warning in such a 
case ...



Whether or not sc_sigma has an effect depends on the topology.  If you change 
the value, you can change the result of the calculation.  g_bar has no need to 
read .tpr files, so I doubt it would be very useful to include them as an 
additional input.



i guess for now i'll wait for a stable 4.6 - any ideas when this will be there?



No clue.  Maybe a few weeks, maybe not.  There are several very large changes 
being put in place.  But note that if you're getting the same error with the 
development version of the free energy code, you're likely to see the exact same 
thing in the 4.6 release.  The free energy code merge is considered done and 
ready for production.  If this is not the case, a bug report needs to be filed 
ASAP.  Please submit a report at redmine.gromacs.org with a minimal test case.


-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
* Only plain text messages are allowed!
* 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


  1   2   3   4   5   6   7   8   9   10   >