[gmx-users] Question on pressure coupling
Dear Gromacs users,
Inverted hexagonal phases are usually simulated (and should be) using a
triclinic box with either semi-isotropic or anisotropic pressure couplings.
I am not sure what a semi-isotropic and anisotropic means in a "Triclinic"
box, exactly?
The .log files look the same for the Triclinic box and rectangular box.
Here is an example:
Semi-isotropic (.mdp file):
compressibility = 4.5e-5 4.6e-5
ref_p = 1.5 1.6
The .log file looks like:
compressibility (3x3):
compressibility[0]={ 4.5e-05, 0.0e+00, 0.0e+00}
compressibility[1]={ 0.0e+00, 4.5e-05, 0.0e+00}
compressibility[2]={ 0.0e+00, 0.0e+00, 4.6e-05}
ref-p (3x3):
ref-p[0]={ 1.5e+00, 0.0e+00, 0.0e+00}
ref-p[1]={ 0.0e+00, 1.5e+00, 0.0e+00}
ref-p[2]={ 0.0e+00, 0.0e+00, 1.6e+00}
and for Anisotropic (.mdp file)
compressibility = 4.5e-5 4.6e-5 4.7e-5 4.8e-5 4.9e-5 4.1e-5
ref_p = 1.5 1.6 1.7 1.8 1.9 1.1
The .log file looks like:
compressibility (3x3):
compressibility[0]={ 4.5e-05, 4.8e-05, 4.9e-05}
compressibility[1]={ 4.8e-05, 4.6e-05, 4.1e-05}
compressibility[2]={ 4.9e-05, 4.1e-05, 4.7e-05}
ref-p (3x3):
ref-p[0]={ 1.5e+00, 1.8e+00, 1.9e+00}
ref-p[1]={ 1.8e+00, 1.6e+00, 1.1e+00}
ref-p[2]={ 1.9e+00, 1.1e+00, 1.7e+00}
How the scaling is done differently for Triclinic and rectangular boxes?
(The question is also valid for isotropic coupling, but for semi and
aniso is more confusing)
Or, is there any difference in the implementation of pressure coupling for
different PBC types at all?
Thanks for your comments in advance,
Cheers,
Mohsen
--
*Rewards work better than punishment ...*
--
Gromacs Users mailing list
* Please search the archive at
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a
mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Fluctuations in HII phase
Thanks Justin! It sounds an interesting approach. I will give it a try. Cheers, Mohsen On Fri, May 4, 2018 at 8:15 AM, Justin Lemkul wrote: > > > On 5/3/18 10:36 PM, Mohsen Ramezanpour wrote: > >> Hi Gromacs users, >> >> I have a cylinder composed of lipids with waters inside (HII phase). >> This lipid tube fluctuates along the Z-axis. >> That being said, if we slice the cylinder in the Z direction, the center >> of >> mass (COM) for each slice is moving around in the XY plane. You can >> imagine >> it as a 1D undulations in the Z direction. >> >> I want to quantitively show that some of these cylinders are more dynamics >> than others in nature. >> >> Snapshots from the system are one of the solutions but can not be used for >> quantitative comparisons between several systems with a different level of >> undulations. It just shows these systems are bent. >> >> The average size of the system in the Z direction and its fluctuations >> might be another naive quantity to be considered. If the Z size is fixed >> and the system is behaving like a wave, this would not give what I want >> either. >> >> So, I would like to know your opinion on what would be the best way to >> show >> such undulations in a quantitative way? >> > > I would plot the x-y COM position as a function of some selection of atoms > along z, and calculate the fluctuation in COM position about the mean. If > the structure is indeed more dynamic, that slice should fluctuate more than > another system that is more rigid. Think of it as the RMSF of the COM in > the plane. > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Assistant Professor > Virginia Tech Department of Biochemistry > > 303 Engel Hall > 340 West Campus Dr. > Blacksburg, VA 24061 > > jalem...@vt.edu | (540) 231-3129 > http://www.thelemkullab.com > > == > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Fluctuations in HII phase
Hi Gromacs users, I have a cylinder composed of lipids with waters inside (HII phase). This lipid tube fluctuates along the Z-axis. That being said, if we slice the cylinder in the Z direction, the center of mass (COM) for each slice is moving around in the XY plane. You can imagine it as a 1D undulations in the Z direction. I want to quantitively show that some of these cylinders are more dynamics than others in nature. Snapshots from the system are one of the solutions but can not be used for quantitative comparisons between several systems with a different level of undulations. It just shows these systems are bent. The average size of the system in the Z direction and its fluctuations might be another naive quantity to be considered. If the Z size is fixed and the system is behaving like a wave, this would not give what I want either. So, I would like to know your opinion on what would be the best way to show such undulations in a quantitative way? Thanks in advance for your comments. Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] gmx densmap
Hi Gromacs users, I have a cylinder composed of lipids with waters and salts inside (HII phase) This cylinder is aligned in the Z direction. I would like to calculate the density of each group (lipid P atoms, water, salt) in the XY plane and show them all in one figure with different colors. gmx densmap can give the 2D map for each group separately. The output is an .xpm file. Using this module for each group, I will get several .xpm files. To combine them, I am using the gmx xpm2ps module and use -combine mult option. Unfortunately, it does not give what I want. I was wondering if you might have a solution for that? Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] gmx traj output
Thanks, Justin. Cheers, Mohsen On Tue, Apr 24, 2018 at 7:09 PM, Justin Lemkul wrote: > > > On 4/23/18 11:08 PM, Mohsen Ramezanpour wrote: > >> Hi Gromacs users, >> >> I am using gmx traj to get the box size (X, Y, Z). >> I use the following command: >> >> echo 0 | gmx traj -f trajectory.xtc -s md.tpr -ob box.xvg -xvg none >> >> The output has 77 columns. >> 1st is the time step, >> >> 2nd, 3rd, and 4th seem to be the X, Y, and Z size of the box. >> I was wondering what are the 5th, 6th, and 7th columns giving, exactly? >> > > The output is the half-matrix of the box tensor. The first three elements > are the diagonal, the remaining three are the off-diagonal elements. > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Assistant Professor > Virginia Tech Department of Biochemistry > > 303 Engel Hall > 340 West Campus Dr. > Blacksburg, VA 24061 > > jalem...@vt.edu | (540) 231-3129 > http://www.thelemkullab.com > > == > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] gmx traj output
Hi Gromacs users, I am using gmx traj to get the box size (X, Y, Z). I use the following command: echo 0 | gmx traj -f trajectory.xtc -s md.tpr -ob box.xvg -xvg none The output has 77 columns. 1st is the time step, 2nd, 3rd, and 4th seem to be the X, Y, and Z size of the box. I was wondering what are the 5th, 6th, and 7th columns giving, exactly? Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] PBC fix for visualization
Hi Dallas and other Gromacs users, I used -pbc whole and -ur compact in the first step "System" index group And then, used the output file for -pbc cluster. Choosing the "System" index for clustering gave the best result I got. (Although there are still few lipids which are not completely in the places they should be). Anyway, the outputs are reasonable. However, I have ca. 100 simulation systems each of ca. 500 ns. (each system composed of ca. 4000 molecules) Doing this clustering on the whole system takes a considerable time. I thought I might be able to do the clustering on a group of atoms (an index group composed of P in lipids and O in water molecules). This, however, did not work and gmx trjconv complained about: "Molecule 1 marked for clustering but not atom 1 in it - check your index!" I was wondering what would be your solution to this? Cheers, Mohsen On Mon, May 22, 2017 at 11:45 AM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Dallas, > > Thanks for your reply. > I did try -pbc cluster for waters. It could fix it somehow but not > completely. > After that, I had to use -pbc center to fix it. Still, I do not get what I > want. > Unfortunately, some waters and lipids are appearing from the other side of > the box. > > Cheers, > Mohsen > > > On Sun, May 21, 2017 at 8:12 PM, Dallas Warren > wrote: > >> I have found the cluster option of -pbc to work well for putting >> aggregates back together correctly. Some times you do need an index >> file and appropriate groups to assist with it getting it right. >> >> gmx trjconv -pbc cluster >> Catch ya, >> >> Dr. Dallas Warren >> Drug Delivery, Disposition and Dynamics >> Monash Institute of Pharmaceutical Sciences, Monash University >> 381 Royal Parade, Parkville VIC 3052 >> dallas.war...@monash.edu >> - >> When the only tool you own is a hammer, every problem begins to resemble >> a nail. >> >> >> On 16 May 2017 at 07:50, Mohsen Ramezanpour >> wrote: >> > Dear Gromacs users, >> > >> > I have an HII phase made of one inverted cylinder (and waters inside) >> in a >> > triclinic box with 90, 90, 60 angles. After running the simulation, this >> > cylinder become bent like a curve. I.e. is not a perfect cylinder >> anymore. >> > As a result, some water molecules and lipids pass the box sides and >> enter >> > from the other side of box because of PBC. >> > Now, I want to make the cylinder again but I am not sure how to do so. >> > >> > The best I could do was to use "-pbc mol -ur compact" options in >> trjconv. >> > However, here are still some lipids and molecules which are not part of >> the >> > cylinder. >> > >> > Any idea how could I fix the effect of these PBC in visualization? >> > >> > Thanks >> > Mohsen >> > >> > -- >> > *Rewards work better than punishment ...* >> > -- >> > Gromacs Users mailing list >> > >> > * Please search the archive at http://www.gromacs.org/Support >> /Mailing_Lists/GMX-Users_List before posting! >> > >> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >> > >> > * For (un)subscribe requests visit >> > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >> send a mail to gmx-users-requ...@gromacs.org. >> -- >> Gromacs Users mailing list >> >> * Please search the archive at http://www.gromacs.org/Support >> /Mailing_Lists/GMX-Users_List before posting! >> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >> >> * For (un)subscribe requests visit >> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >> send a mail to gmx-users-requ...@gromacs.org. >> > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] gangle for a plane and z axis
Hi Gromacs users, I would like to calculate the angle distribution between a plane and z axis using gangle. In a 5-member ring, there are three atoms of O1-C-O2 which form an angle. So, I will choose these three atom names to define the plane. I am not sure how the g_angle makes the vector normal for this plane, though. I know that it first make two vectors, and then does a cross product of these two initial vectors to make the normal vector to the plane. But, what is the direction of the vector normal to this plane? Looking at the calc_vec function (lines 652 to 690) code from: https://github.com/gromacs/gromacs/blob/master/src/gromacs/trajectoryanalysis/modules/angle.cpp It is using functions cprod and rvec_sub. I found cprod definition here: http://manual.gromacs.org/documentation/2018-beta3/doxygen/html-lib/group__module__simd.xhtml#ga54124af0ff118d3f171b0eef07105c76 cprod (v1,v2,xout) means C = v1 * v2 (so the Right Hand Rule can be used to get the direction) I could not find rvec_sub function, though. However, I think: rvec_sub(x[1], x[0], v1) gives the v1 as v1= x[1]-x[0]. i.e. a vector starting from x[0] and ending on x[1]. Are these correct? Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] inter- and intra-molecular hbond analysis
Hi Gromacs users, I would like to know your opinion on this. Thanks in advance for your comments Cheers, Mohsen On Fri, Jan 19, 2018 at 5:16 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Dear Gromacs users, > > Given a mixed bilayer composed of two lipids (type A and B), I would like > to calculate: > > 1) "Intramolecular" h-bonds for a lipid type to see if part of lipid (NH3+ > of PE) is interacting with some other part of itself (PO4- of the same > lipid). > > 2) "intermolecular" hbonds between A-A, A-B, and B-B. > > If I understood correctly, gmx hbond does not distinguish between the > inter- and intra-molecular hbonds. So, choosing A and A lipids will give > both the intramolecular and intermolecular hbonds between lipids type A in > the system. > > One solution could be doing this analysis on each individual lipid of type > A first, and then make an average over all the lipid type A in the system. > This will give the intramolecular hbonds (lets call it H-intra) > > Next, calculating the hbonds between lipid A and A as gmx hbond does by > defualt. This will give both intermolecular and intramolecular > hbonds between A lipids in the system (Lets call it H-Both). > Thus, (H-both) - (H-intra) = (H-intermolecular hbond between lipids of the > same type) > > I think there should be an easier way to do this (especially H-intra) and > I do not know. > I was wondering if anyone has any suggestion? > > This is the same problem I am facing when I want to calculate the min-dist > between two groups, i.e. the intermolecular and intramolecular min-distance > between two groups of atoms. > > Thanks in advance for your comments, > Cheers, > Mohsen > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] inter- and intra-molecular hbond analysis
Dear Gromacs users, Given a mixed bilayer composed of two lipids (type A and B), I would like to calculate: 1) "Intramolecular" h-bonds for a lipid type to see if part of lipid (NH3+ of PE) is interacting with some other part of itself (PO4- of the same lipid). 2) "intermolecular" hbonds between A-A, A-B, and B-B. If I understood correctly, gmx hbond does not distinguish between the inter- and intra-molecular hbonds. So, choosing A and A lipids will give both the intramolecular and intermolecular hbonds between lipids type A in the system. One solution could be doing this analysis on each individual lipid of type A first, and then make an average over all the lipid type A in the system. This will give the intramolecular hbonds (lets call it H-intra) Next, calculating the hbonds between lipid A and A as gmx hbond does by defualt. This will give both intermolecular and intramolecular hbonds between A lipids in the system (Lets call it H-Both). Thus, (H-both) - (H-intra) = (H-intermolecular hbond between lipids of the same type) I think there should be an easier way to do this (especially H-intra) and I do not know. I was wondering if anyone has any suggestion? This is the same problem I am facing when I want to calculate the min-dist between two groups, i.e. the intermolecular and intramolecular min-distance between two groups of atoms. Thanks in advance for your comments, Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] tilt angle for POPC
Thanks, I agree. I will do as you recommended. Cheers, Mohsen On Fri, Jan 19, 2018 at 9:23 AM, Justin Lemkul wrote: > > > On 1/19/18 11:16 AM, Mohsen Ramezanpour wrote: > >> On Fri, Jan 19, 2018 at 5:55 AM, Justin Lemkul wrote: >> >> >>> On 1/18/18 1:28 PM, Mohsen Ramezanpour wrote: >>> >>> Hi Justin, >>>> >>>> On Thu, Jan 18, 2018 at 5:42 AM, Justin Lemkul wrote: >>>> >>>> >>>> On 1/18/18 12:24 AM, Mohsen Ramezanpour wrote: >>>>> >>>>> Dear Gromacs users, >>>>> >>>>>> I am interested in calculation of tilt angle for the POPC headgroup >>>>>> (angle >>>>>> distribution between the P-N vector and Z axis). >>>>>> I am not sure if my approach is correct as my angle distribution does >>>>>> not >>>>>> seem reasonable. >>>>>> >>>>>> Given a bilayer with 200 lipids (100 lipids in each leaflet with >>>>>> resides 1-100 and 101-200 for upper and lower leaflets, respectively) >>>>>> simulated for 200 ns: >>>>>> >>>>>> gmx make_ndx -f SYSTEM.gro -o index.ndx >>>>>> keep 0 >>>>>> r 1-100 >>>>>> name 1 upperleaflet >>>>>> 1 & a P >>>>>> 1 & a N >>>>>> 2 | 3 >>>>>> name 4 vector >>>>>> q >>>>>> >>>>>> Next, I use: >>>>>> gmx gangle -f md.xtc -s md.tpr -n index.ndx -g1 vector -g2 z >>>>>>-b >>>>>> 10 -group1 -oh histogram.xvg -binw 0.01 >>>>>> >>>>>> and choose index group 4 and then Ctrl+D. >>>>>> >>>>>> Please let me know your opinion. I think I am doing something wrong, >>>>>> especially with the construction of P-N vector. >>>>>> >>>>>> What results do you get from this approach? >>>>>> >>>>> A normal distribution centered around 105 with a high pick standing out >>>>> >>>> around 90. >>>> I expected it to be around ~73-77 >>>> >>>> Upon what are you basing your expectation? Previous literature? Some >>> other >>> calculation or visual inspection? >>> >> It is mentioned upon the last paragraph of left column at page (3777) in >> this article: >> Saches et al. Biophys J >> <https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1304278/#>. 2004 Jun; >> 86(6): >> 3772–3782 >> > > They're using a different (outdated and worse in general) force field than > you. You can't expect the same results. You need to know if there's > precedent for how *your* chosen force field performs in this regard. > > >> >>> What happens if you try to analyze only a single lipid? >>>> >>>>> I get an approximate normal distribution centered around 100 degrees. >>>>> >>>> >>>> I have two systems, one is salt free and the other one has 0.15 M NaCL >>>> salt. >>>> The shape of histograms is the same. >>>> >>>> I had to separate the leaflets, right? >>>> >>>> If the membrane is homogeneous, I don't see any reason why. >>> >> Choosing z axis to calculate the angle distribution for, P->N vector has >> opposite directions (in average) between two leaflets. >> for one leaflet one makes (Just as an example) a normal distribution >> centered around 73, and the lower leaflet will make a normal distribution >> centered around (180-73). >> However, if we consider both, we will get a normal distribution centered >> around ~90. >> This is why I think we must separate those, even for a homogenous bilayer. >> Is this reasonable? >> > > OK, makes sense. > > -Justin > > > >>> -Justin >>> >>> >>> So, any idea why this is happening? >>> >>>> -Justin >>>> >>>>> -- >>>>> == >>>>> >>>>> Justin A. Lemkul, Ph.D. >>>>> Assistant Professor >>>>> Virginia Tech Department of Biochemistry >>>>> >>>>> 303 Engel Hall >>>>> 340 West Campus Dr. >>>>> Blacksburg, VA
Re: [gmx-users] tilt angle for POPC
On Fri, Jan 19, 2018 at 5:55 AM, Justin Lemkul wrote: > > > On 1/18/18 1:28 PM, Mohsen Ramezanpour wrote: > >> Hi Justin, >> >> On Thu, Jan 18, 2018 at 5:42 AM, Justin Lemkul wrote: >> >> >>> On 1/18/18 12:24 AM, Mohsen Ramezanpour wrote: >>> >>> Dear Gromacs users, >>>> >>>> I am interested in calculation of tilt angle for the POPC headgroup >>>> (angle >>>> distribution between the P-N vector and Z axis). >>>> I am not sure if my approach is correct as my angle distribution does >>>> not >>>> seem reasonable. >>>> >>>> Given a bilayer with 200 lipids (100 lipids in each leaflet with >>>> resides 1-100 and 101-200 for upper and lower leaflets, respectively) >>>> simulated for 200 ns: >>>> >>>> gmx make_ndx -f SYSTEM.gro -o index.ndx >>>> keep 0 >>>> r 1-100 >>>> name 1 upperleaflet >>>> 1 & a P >>>> 1 & a N >>>> 2 | 3 >>>> name 4 vector >>>> q >>>> >>>> Next, I use: >>>> gmx gangle -f md.xtc -s md.tpr -n index.ndx -g1 vector -g2 z >>>> -b >>>>10 -group1 -oh histogram.xvg -binw 0.01 >>>> >>>> and choose index group 4 and then Ctrl+D. >>>> >>>> Please let me know your opinion. I think I am doing something wrong, >>>> especially with the construction of P-N vector. >>>> >>>> What results do you get from this approach? >>> >>> A normal distribution centered around 105 with a high pick standing out >> around 90. >> I expected it to be around ~73-77 >> > > Upon what are you basing your expectation? Previous literature? Some other > calculation or visual inspection? It is mentioned upon the last paragraph of left column at page (3777) in this article: Saches et al. Biophys J <https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1304278/#>. 2004 Jun; 86(6): 3772–3782 > > >> What happens if you try to analyze only a single lipid? >>> >>> I get an approximate normal distribution centered around 100 degrees. >> >> >> I have two systems, one is salt free and the other one has 0.15 M NaCL >> salt. >> The shape of histograms is the same. >> >> I had to separate the leaflets, right? >> > > If the membrane is homogeneous, I don't see any reason why. Choosing z axis to calculate the angle distribution for, P->N vector has opposite directions (in average) between two leaflets. for one leaflet one makes (Just as an example) a normal distribution centered around 73, and the lower leaflet will make a normal distribution centered around (180-73). However, if we consider both, we will get a normal distribution centered around ~90. This is why I think we must separate those, even for a homogenous bilayer. Is this reasonable? > > > -Justin > > > So, any idea why this is happening? >> >> -Justin >>> >>> -- >>> == >>> >>> Justin A. Lemkul, Ph.D. >>> Assistant Professor >>> Virginia Tech Department of Biochemistry >>> >>> 303 Engel Hall >>> 340 West Campus Dr. >>> Blacksburg, VA 24061 >>> >>> jalem...@vt.edu | (540) 231-3129 >>> http://www.biochem.vt.edu/people/faculty/JustinLemkul.html >>> >>> == >>> >>> -- >>> Gromacs Users mailing list >>> >>> * Please search the archive at http://www.gromacs.org/Support >>> /Mailing_Lists/GMX-Users_List before posting! >>> >>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >>> >>> * For (un)subscribe requests visit >>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >>> send a mail to gmx-users-requ...@gromacs.org. >>> >>> >> >> > -- > == > > Justin A. Lemkul, Ph.D. > Assistant Professor > Virginia Tech Department of Biochemistry > > 303 Engel Hall > 340 West Campus Dr. > Blacksburg, VA 24061 > > jalem...@vt.edu | (540) 231-3129 > http://www.biochem.vt.edu/people/faculty/JustinLemkul.html > > == > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] tilt angle for POPC
I did more tests and here are the conclusions for the pick mentioned earlier. The pick appears when the -binw gets smaller and smaller. For the default value of gangle (which is 1), one cannot detect such a pick and you will get a normal distribution. Even for a -binw of 0.1 everything looks fine. However, for -binw of lower than 0.1, e.g. 0.5 or so, the picks appear. Well, is this a problem with POPC simulation, POPC topology, POPC parameters, the algorithm behind gangle, or it is what is really happening and we do not detect it due to the NOT enough Small value of -binw in the histogram? I vote for the last option. Cheers, On Thu, Jan 18, 2018 at 12:07 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Looking at index file made in the above mentioned way, as way as the > following way: > > gmx make_ndx -f SYSTEM.gro -o index.ndx > keep 0 > r 1-100 > name 1 upperleaflet > 1 & a P N > name 2 vector > q > > gave the same results. > Interestingly, even if I make the index file like: > gmx make_ndx -f SYSTEM.gro -o index.ndx > keep 0 > r 1-100 > name 1 upperleaflet > 1 & a N P > name 2 vector > q > > > The index file (i.e. the order and numbers in the index file) for group 2 > is exactly the same, i.e. will give the N->P vector, not P->N. > So, I think this is why there was a shift to more than 90 degrees. > doing (180-alpha) gave the P->N angle distribution. > > STILL, I do not know why there is a pick standing out around 90 degrees. > It is very sharp narrow pick. > Other tilt angle distributions I have seem do not show such a pick. > > Thanks in advance for your comments and help > > > On Thu, Jan 18, 2018 at 11:29 AM, Mohsen Ramezanpour < > ramezanpour.moh...@gmail.com> wrote: > >> And, I used Charmm36 FF for simulations in Gromacs v. 2016.3 >> >> Cheers, >> Mohsen >> >> On Thu, Jan 18, 2018 at 11:28 AM, Mohsen Ramezanpour < >> ramezanpour.moh...@gmail.com> wrote: >> >>> Hi Justin, >>> >>> On Thu, Jan 18, 2018 at 5:42 AM, Justin Lemkul wrote: >>> >>>> >>>> >>>> On 1/18/18 12:24 AM, Mohsen Ramezanpour wrote: >>>> >>>>> Dear Gromacs users, >>>>> >>>>> I am interested in calculation of tilt angle for the POPC headgroup >>>>> (angle >>>>> distribution between the P-N vector and Z axis). >>>>> I am not sure if my approach is correct as my angle distribution does >>>>> not >>>>> seem reasonable. >>>>> >>>>> Given a bilayer with 200 lipids (100 lipids in each leaflet with >>>>> resides 1-100 and 101-200 for upper and lower leaflets, respectively) >>>>> simulated for 200 ns: >>>>> >>>>> gmx make_ndx -f SYSTEM.gro -o index.ndx >>>>> keep 0 >>>>> r 1-100 >>>>> name 1 upperleaflet >>>>> 1 & a P >>>>> 1 & a N >>>>> 2 | 3 >>>>> name 4 vector >>>>> q >>>>> >>>>> Next, I use: >>>>> gmx gangle -f md.xtc -s md.tpr -n index.ndx -g1 vector -g2 >>>>> z -b >>>>> 10 -group1 -oh histogram.xvg -binw 0.01 >>>>> >>>>> and choose index group 4 and then Ctrl+D. >>>>> >>>>> Please let me know your opinion. I think I am doing something wrong, >>>>> especially with the construction of P-N vector. >>>>> >>>> >>>> What results do you get from this approach? >>>> >>> A normal distribution centered around 105 with a high pick standing out >>> around 90. >>> I expected it to be around ~73-77 >>> >>> >>>> >>>> What happens if you try to analyze only a single lipid? >>>> >>> I get an approximate normal distribution centered around 100 degrees. >>> >>> >>> I have two systems, one is salt free and the other one has 0.15 M NaCL >>> salt. >>> The shape of histograms is the same. >>> >>> I had to separate the leaflets, right? >>> So, any idea why this is happening? >>> >>>> >>>> -Justin >>>> >>>> -- >>>> == >>>> >>>> Justin A. Lemkul, Ph.D. >>>> Assistant Professor >>>> Virginia Tech Department of Biochemistry >>>> >>>> 303 Engel Hall &
Re: [gmx-users] tilt angle for POPC
Looking at index file made in the above mentioned way, as way as the following way: gmx make_ndx -f SYSTEM.gro -o index.ndx keep 0 r 1-100 name 1 upperleaflet 1 & a P N name 2 vector q gave the same results. Interestingly, even if I make the index file like: gmx make_ndx -f SYSTEM.gro -o index.ndx keep 0 r 1-100 name 1 upperleaflet 1 & a N P name 2 vector q The index file (i.e. the order and numbers in the index file) for group 2 is exactly the same, i.e. will give the N->P vector, not P->N. So, I think this is why there was a shift to more than 90 degrees. doing (180-alpha) gave the P->N angle distribution. STILL, I do not know why there is a pick standing out around 90 degrees. It is very sharp narrow pick. Other tilt angle distributions I have seem do not show such a pick. Thanks in advance for your comments and help On Thu, Jan 18, 2018 at 11:29 AM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > And, I used Charmm36 FF for simulations in Gromacs v. 2016.3 > > Cheers, > Mohsen > > On Thu, Jan 18, 2018 at 11:28 AM, Mohsen Ramezanpour < > ramezanpour.moh...@gmail.com> wrote: > >> Hi Justin, >> >> On Thu, Jan 18, 2018 at 5:42 AM, Justin Lemkul wrote: >> >>> >>> >>> On 1/18/18 12:24 AM, Mohsen Ramezanpour wrote: >>> >>>> Dear Gromacs users, >>>> >>>> I am interested in calculation of tilt angle for the POPC headgroup >>>> (angle >>>> distribution between the P-N vector and Z axis). >>>> I am not sure if my approach is correct as my angle distribution does >>>> not >>>> seem reasonable. >>>> >>>> Given a bilayer with 200 lipids (100 lipids in each leaflet with >>>> resides 1-100 and 101-200 for upper and lower leaflets, respectively) >>>> simulated for 200 ns: >>>> >>>> gmx make_ndx -f SYSTEM.gro -o index.ndx >>>> keep 0 >>>> r 1-100 >>>> name 1 upperleaflet >>>> 1 & a P >>>> 1 & a N >>>> 2 | 3 >>>> name 4 vector >>>> q >>>> >>>> Next, I use: >>>> gmx gangle -f md.xtc -s md.tpr -n index.ndx -g1 vector -g2 z >>>> -b >>>> 10 -group1 -oh histogram.xvg -binw 0.01 >>>> >>>> and choose index group 4 and then Ctrl+D. >>>> >>>> Please let me know your opinion. I think I am doing something wrong, >>>> especially with the construction of P-N vector. >>>> >>> >>> What results do you get from this approach? >>> >> A normal distribution centered around 105 with a high pick standing out >> around 90. >> I expected it to be around ~73-77 >> >> >>> >>> What happens if you try to analyze only a single lipid? >>> >> I get an approximate normal distribution centered around 100 degrees. >> >> >> I have two systems, one is salt free and the other one has 0.15 M NaCL >> salt. >> The shape of histograms is the same. >> >> I had to separate the leaflets, right? >> So, any idea why this is happening? >> >>> >>> -Justin >>> >>> -- >>> == >>> >>> Justin A. Lemkul, Ph.D. >>> Assistant Professor >>> Virginia Tech Department of Biochemistry >>> >>> 303 Engel Hall >>> 340 West Campus Dr. >>> Blacksburg, VA 24061 >>> >>> jalem...@vt.edu | (540) 231-3129 >>> http://www.biochem.vt.edu/people/faculty/JustinLemkul.html >>> >>> == >>> >>> -- >>> Gromacs Users mailing list >>> >>> * Please search the archive at http://www.gromacs.org/Support >>> /Mailing_Lists/GMX-Users_List before posting! >>> >>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >>> >>> * For (un)subscribe requests visit >>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >>> send a mail to gmx-users-requ...@gromacs.org. >>> >> >> >> >> -- >> *Rewards work better than punishment ...* >> > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] tilt angle for POPC
And, I used Charmm36 FF for simulations in Gromacs v. 2016.3 Cheers, Mohsen On Thu, Jan 18, 2018 at 11:28 AM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Justin, > > On Thu, Jan 18, 2018 at 5:42 AM, Justin Lemkul wrote: > >> >> >> On 1/18/18 12:24 AM, Mohsen Ramezanpour wrote: >> >>> Dear Gromacs users, >>> >>> I am interested in calculation of tilt angle for the POPC headgroup >>> (angle >>> distribution between the P-N vector and Z axis). >>> I am not sure if my approach is correct as my angle distribution does not >>> seem reasonable. >>> >>> Given a bilayer with 200 lipids (100 lipids in each leaflet with >>> resides 1-100 and 101-200 for upper and lower leaflets, respectively) >>> simulated for 200 ns: >>> >>> gmx make_ndx -f SYSTEM.gro -o index.ndx >>> keep 0 >>> r 1-100 >>> name 1 upperleaflet >>> 1 & a P >>> 1 & a N >>> 2 | 3 >>> name 4 vector >>> q >>> >>> Next, I use: >>> gmx gangle -f md.xtc -s md.tpr -n index.ndx -g1 vector -g2 z >>> -b >>> 10 -group1 -oh histogram.xvg -binw 0.01 >>> >>> and choose index group 4 and then Ctrl+D. >>> >>> Please let me know your opinion. I think I am doing something wrong, >>> especially with the construction of P-N vector. >>> >> >> What results do you get from this approach? >> > A normal distribution centered around 105 with a high pick standing out > around 90. > I expected it to be around ~73-77 > > >> >> What happens if you try to analyze only a single lipid? >> > I get an approximate normal distribution centered around 100 degrees. > > > I have two systems, one is salt free and the other one has 0.15 M NaCL > salt. > The shape of histograms is the same. > > I had to separate the leaflets, right? > So, any idea why this is happening? > >> >> -Justin >> >> -- >> == >> >> Justin A. Lemkul, Ph.D. >> Assistant Professor >> Virginia Tech Department of Biochemistry >> >> 303 Engel Hall >> 340 West Campus Dr. >> Blacksburg, VA 24061 >> >> jalem...@vt.edu | (540) 231-3129 >> http://www.biochem.vt.edu/people/faculty/JustinLemkul.html >> >> == >> >> -- >> Gromacs Users mailing list >> >> * Please search the archive at http://www.gromacs.org/Support >> /Mailing_Lists/GMX-Users_List before posting! >> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >> >> * For (un)subscribe requests visit >> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >> send a mail to gmx-users-requ...@gromacs.org. >> > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] tilt angle for POPC
Hi Justin, On Thu, Jan 18, 2018 at 5:42 AM, Justin Lemkul wrote: > > > On 1/18/18 12:24 AM, Mohsen Ramezanpour wrote: > >> Dear Gromacs users, >> >> I am interested in calculation of tilt angle for the POPC headgroup (angle >> distribution between the P-N vector and Z axis). >> I am not sure if my approach is correct as my angle distribution does not >> seem reasonable. >> >> Given a bilayer with 200 lipids (100 lipids in each leaflet with >> resides 1-100 and 101-200 for upper and lower leaflets, respectively) >> simulated for 200 ns: >> >> gmx make_ndx -f SYSTEM.gro -o index.ndx >> keep 0 >> r 1-100 >> name 1 upperleaflet >> 1 & a P >> 1 & a N >> 2 | 3 >> name 4 vector >> q >> >> Next, I use: >> gmx gangle -f md.xtc -s md.tpr -n index.ndx -g1 vector -g2 z >> -b >> 10 -group1 -oh histogram.xvg -binw 0.01 >> >> and choose index group 4 and then Ctrl+D. >> >> Please let me know your opinion. I think I am doing something wrong, >> especially with the construction of P-N vector. >> > > What results do you get from this approach? > A normal distribution centered around 105 with a high pick standing out around 90. I expected it to be around ~73-77 > > What happens if you try to analyze only a single lipid? > I get an approximate normal distribution centered around 100 degrees. I have two systems, one is salt free and the other one has 0.15 M NaCL salt. The shape of histograms is the same. I had to separate the leaflets, right? So, any idea why this is happening? > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Assistant Professor > Virginia Tech Department of Biochemistry > > 303 Engel Hall > 340 West Campus Dr. > Blacksburg, VA 24061 > > jalem...@vt.edu | (540) 231-3129 > http://www.biochem.vt.edu/people/faculty/JustinLemkul.html > > == > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] tilt angle for POPC
Dear Gromacs users, I am interested in calculation of tilt angle for the POPC headgroup (angle distribution between the P-N vector and Z axis). I am not sure if my approach is correct as my angle distribution does not seem reasonable. Given a bilayer with 200 lipids (100 lipids in each leaflet with resides 1-100 and 101-200 for upper and lower leaflets, respectively) simulated for 200 ns: gmx make_ndx -f SYSTEM.gro -o index.ndx keep 0 r 1-100 name 1 upperleaflet 1 & a P 1 & a N 2 | 3 name 4 vector q Next, I use: gmx gangle -f md.xtc -s md.tpr -n index.ndx -g1 vector -g2 z -b 10 -group1 -oh histogram.xvg -binw 0.01 and choose index group 4 and then Ctrl+D. Please let me know your opinion. I think I am doing something wrong, especially with the construction of P-N vector. Thanks in advance, Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Electron density using g_density
Hi All, I did a short test on this. Following is the results. I made two electron.dat files, one like: 1 N = 7 and the second one like: 2 N = 7 N = 7 I also made an index.ndx with these groups: r lipid1 & a N r lipid2 & a N Next, I calculated the electron density for these two groups using two electron.dat files, separately. The results are exactly the same. Cheers On Wed, Jan 17, 2018 at 9:56 AM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Gromacs users, > > My system is composed of three lipid molecules. > When I made the topology file for new lipids in the system, I chose some > atom names by my choice. > Now, when I want to make an "electron.dat" file which includes the atomic > number for each atom name in the system, I have more than one line for some > atom names, e.g. for N, C, and H atoms. > In an ideal case, the atom names should be different in order to the > interpret the gmx density output without any confusion. > However, I think I could/should keep just one line for each atom name > because there is no difference in their atomic number between the lipids. > Based on the code, it seems that when a group is selected from > an index.ndx file, the partial charges for that specific atom (i.e. atom > name N from residue lipid1) is taken into account and added to the atomic > number defined in the "electron.dat". This way, it seems both having one > line and three lines for an atom name (e.g. N) will give the same results, > as gmx density can recognize this based on index group. > > I would like to know your opinion on this. > Thanks in advance for your comments > > Cheers, > Mohsen > > > On Tue, Jan 16, 2018 at 10:16 AM, Mohsen Ramezanpour < > ramezanpour.moh...@gmail.com> wrote: > >> Thanks, Justin. >> >> I will as you recommended. >> The last approach was described in a book chapter for bilayer simulation. >> I had difficulty to justify that approach. >> >> Cheers, >> Mohsen >> >> >> >> On Tue, Jan 16, 2018 at 9:42 AM, Justin Lemkul wrote: >> >>> >>> >>> On 1/16/18 11:37 AM, Mohsen Ramezanpour wrote: >>> >>>> Dear Gromacs users, >>>> >>>> In the man page of g_density for calculation of electron density, it is >>>> mentioned that: >>>> >>>> "The number of electrons for each atom is modified by its atomic partial >>>> charge." >>>> >>>> https://linux.die.net/man/1/g_density >>>> >>>> In some publications, it seems the authors just took the atomic number >>>> as >>>> the number of electrons for each atom name. >>>> For instance, for POPC: they took N=7, P=15, O=8, C=6, and H=1 >>>> >>>> In some others, they did not consider the partial charge. >>>> For instance, for POPC: they took N = 6, P = 16 and all the O atoms = 8, >>>> C=6, and H=1 >>>> >>>> So, I was wondering which approach do you think is better? >>>> Will the difference between these approaches matter for making >>>> conclusions? >>>> >>> >>> The code accounts for the partial charges (see calc_electron_density() >>> in gmx_density.cpp). You should specify the number of electrons on the >>> elemental form of the atom. One cannot ascribe the entire +1 charge to N or >>> -1 to P, anyway, so that's a false approach. >>> >>> -Justin >>> >>> -- >>> == >>> >>> Justin A. Lemkul, Ph.D. >>> Assistant Professor >>> Virginia Tech Department of Biochemistry >>> >>> 303 Engel Hall >>> 340 West Campus Dr. >>> Blacksburg, VA 24061 >>> >>> jalem...@vt.edu | (540) 231-3129 >>> http://www.biochem.vt.edu/people/faculty/JustinLemkul.html >>> >>> == >>> >>> -- >>> Gromacs Users mailing list >>> >>> * Please search the archive at http://www.gromacs.org/Support >>> /Mailing_Lists/GMX-Users_List before posting! >>> >>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >>> >>> * For (un)subscribe requests visit >>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >>> send a mail to gmx-users-requ...@gromacs.org. >>> >> >> >> >> -- >> *Rewards work better than punishment ...* >> > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Electron density using g_density
Hi Gromacs users, My system is composed of three lipid molecules. When I made the topology file for new lipids in the system, I chose some atom names by my choice. Now, when I want to make an "electron.dat" file which includes the atomic number for each atom name in the system, I have more than one line for some atom names, e.g. for N, C, and H atoms. In an ideal case, the atom names should be different in order to the interpret the gmx density output without any confusion. However, I think I could/should keep just one line for each atom name because there is no difference in their atomic number between the lipids. Based on the code, it seems that when a group is selected from an index.ndx file, the partial charges for that specific atom (i.e. atom name N from residue lipid1) is taken into account and added to the atomic number defined in the "electron.dat". This way, it seems both having one line and three lines for an atom name (e.g. N) will give the same results, as gmx density can recognize this based on index group. I would like to know your opinion on this. Thanks in advance for your comments Cheers, Mohsen On Tue, Jan 16, 2018 at 10:16 AM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Thanks, Justin. > > I will as you recommended. > The last approach was described in a book chapter for bilayer simulation. > I had difficulty to justify that approach. > > Cheers, > Mohsen > > > > On Tue, Jan 16, 2018 at 9:42 AM, Justin Lemkul wrote: > >> >> >> On 1/16/18 11:37 AM, Mohsen Ramezanpour wrote: >> >>> Dear Gromacs users, >>> >>> In the man page of g_density for calculation of electron density, it is >>> mentioned that: >>> >>> "The number of electrons for each atom is modified by its atomic partial >>> charge." >>> >>> https://linux.die.net/man/1/g_density >>> >>> In some publications, it seems the authors just took the atomic number as >>> the number of electrons for each atom name. >>> For instance, for POPC: they took N=7, P=15, O=8, C=6, and H=1 >>> >>> In some others, they did not consider the partial charge. >>> For instance, for POPC: they took N = 6, P = 16 and all the O atoms = 8, >>> C=6, and H=1 >>> >>> So, I was wondering which approach do you think is better? >>> Will the difference between these approaches matter for making >>> conclusions? >>> >> >> The code accounts for the partial charges (see calc_electron_density() in >> gmx_density.cpp). You should specify the number of electrons on the >> elemental form of the atom. One cannot ascribe the entire +1 charge to N or >> -1 to P, anyway, so that's a false approach. >> >> -Justin >> >> -- >> == >> >> Justin A. Lemkul, Ph.D. >> Assistant Professor >> Virginia Tech Department of Biochemistry >> >> 303 Engel Hall >> 340 West Campus Dr. >> Blacksburg, VA 24061 >> >> jalem...@vt.edu | (540) 231-3129 >> http://www.biochem.vt.edu/people/faculty/JustinLemkul.html >> >> == >> >> -- >> Gromacs Users mailing list >> >> * Please search the archive at http://www.gromacs.org/Support >> /Mailing_Lists/GMX-Users_List before posting! >> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >> >> * For (un)subscribe requests visit >> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >> send a mail to gmx-users-requ...@gromacs.org. >> > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Electron density using g_density
Thanks, Justin. I will as you recommended. The last approach was described in a book chapter for bilayer simulation. I had difficulty to justify that approach. Cheers, Mohsen On Tue, Jan 16, 2018 at 9:42 AM, Justin Lemkul wrote: > > > On 1/16/18 11:37 AM, Mohsen Ramezanpour wrote: > >> Dear Gromacs users, >> >> In the man page of g_density for calculation of electron density, it is >> mentioned that: >> >> "The number of electrons for each atom is modified by its atomic partial >> charge." >> >> https://linux.die.net/man/1/g_density >> >> In some publications, it seems the authors just took the atomic number as >> the number of electrons for each atom name. >> For instance, for POPC: they took N=7, P=15, O=8, C=6, and H=1 >> >> In some others, they did not consider the partial charge. >> For instance, for POPC: they took N = 6, P = 16 and all the O atoms = 8, >> C=6, and H=1 >> >> So, I was wondering which approach do you think is better? >> Will the difference between these approaches matter for making >> conclusions? >> > > The code accounts for the partial charges (see calc_electron_density() in > gmx_density.cpp). You should specify the number of electrons on the > elemental form of the atom. One cannot ascribe the entire +1 charge to N or > -1 to P, anyway, so that's a false approach. > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Assistant Professor > Virginia Tech Department of Biochemistry > > 303 Engel Hall > 340 West Campus Dr. > Blacksburg, VA 24061 > > jalem...@vt.edu | (540) 231-3129 > http://www.biochem.vt.edu/people/faculty/JustinLemkul.html > > == > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Electron density using g_density
Dear Gromacs users, In the man page of g_density for calculation of electron density, it is mentioned that: "The number of electrons for each atom is modified by its atomic partial charge." https://linux.die.net/man/1/g_density In some publications, it seems the authors just took the atomic number as the number of electrons for each atom name. For instance, for POPC: they took N=7, P=15, O=8, C=6, and H=1 In some others, they did not consider the partial charge. For instance, for POPC: they took N = 6, P = 16 and all the O atoms = 8, C=6, and H=1 So, I was wondering which approach do you think is better? Will the difference between these approaches matter for making conclusions? Thanks in advance for your comments and replies Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] g_order for all-atom simulations for double bonds
Dear Gromacs users, I have a bilayer composed of lipids with two unsaturated bonds in each tail. I did my simulations in Gromacs v. 2016.3 in Charmm36FF. BUT, I am doing the analysis in Gromacs v. 2016.2. I want to calculate the order parameter for this lipid. Here is what I do: 1) I calculate the S_cd for all the carbon atoms in each tail 2) I calculate the S_cd for two carbon atoms in the first double bond (with -unsat option) 3) I calculate the S_cd for two carbon atoms in the second double bond (with -unsat option) 4) I replace the values from steps 2 and 3 in the output files from step 1. This was the well-accepted approach for calculation of unsaturated tails S_cd using g_order. (please correct me if I am wrong) Today, I came across a recent publication (I would recommend everyone to look at) by Thomas Piggot et al. JCTC, 2017 which challenges this approach. Briefly speaking, the standard g_order (even version 2016) seems to give incorrect S_cd for double bonds. The modified version of the g_order in this work seems to provide to give much better values. However, it is recommended to be used when a UA force field has been used for the simulation. I was wondering what tools do you think I should use for this case, and why? 1) NMRlipids script for all_atom S_cd calculation 2) installing and using this modified version 3) if these modifications have been already implemented in the Gromacs 2016.3 or maybe newer versions if any. Thanks in advance for your comments Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Vacuum in the hydrophobic part of bilayers
Hi Gromacs users, In the hydrophobic part of the bilayers, there are some spaces (especially in the mid-layer) which is not completely filled with atoms (so, vacuum). I was wondering if there is any tool which is capable of measuring this (i.e. the volume of this space) in bilayers? Thanks in advance for your reply, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] lipid bilayer simulation-simulation setup
Hi Gormacs users, I was wondering if you have any opinion regarding these questions? Thanks in advance for your comments Cheers, Mohsen On Fri, Dec 29, 2017 at 5:55 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Everyone, > > Lipid bilayers are (almost all the times) simulated in a rectangular or > cubic boxes where all the angles are 90 between the PBC box sides. They are > usually equilibrated with semi-isotropic Berendsen P-coupling whereas the > production runs use semi-isotropic Parrinello-Rahman p-coupling (this part > make sense). I have two questions: > > 1) What is the problem if we use a triclinic box with 90 90 120 angles for > bilayer simulation: > If we choose the z-direction as the normal vector to the bilayer. and a > semi-isotropic p-coupling will be applied to the XY plane? > > > 2) There are some publications who used the Berendsen coupling for > production run as well. They do not usually explain why they prefer this > barostat for the production run. It is known that Berendsen barostat is not > producing the target ensemble as good as Parrinello-Rahman. > It might be, however, because of using anisotropic p-coupling where both > the box sides lengths and angles are allowed to change. > > Any idea? > Thanks in advance for your input on this questions. > > Cheers, > Mohsen > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] lipid bilayer simulation-simulation setup
Hi Everyone, Lipid bilayers are (almost all the times) simulated in a rectangular or cubic boxes where all the angles are 90 between the PBC box sides. They are usually equilibrated with semi-isotropic Berendsen P-coupling whereas the production runs use semi-isotropic Parrinello-Rahman p-coupling (this part make sense). I have two questions: 1) What is the problem if we use a triclinic box with 90 90 120 angles for bilayer simulation: If we choose the z-direction as the normal vector to the bilayer. and a semi-isotropic p-coupling will be applied to the XY plane? 2) There are some publications who used the Berendsen coupling for production run as well. They do not usually explain why they prefer this barostat for the production run. It is known that Berendsen barostat is not producing the target ensemble as good as Parrinello-Rahman. It might be, however, because of using anisotropic p-coupling where both the box sides lengths and angles are allowed to change. Any idea? Thanks in advance for your input on this questions. Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] lipid molecular director axis definition
Thanks, Justin. I think so. Sure, I will. On Fri, Nov 10, 2017 at 5:56 AM, Justin Lemkul wrote: > > > On 11/9/17 6:21 PM, Mohsen Ramezanpour wrote: > >> Hi Everyone, >> >> I have a system formed of lipids (DOPE). I would like to know the average >> angle of the "molecular director axis" with z-direction. Using g_angle it >> should be easy. However, my question is regarding the molecular director >> axis. >> >> Is there any definition for how to determine a "lipid molecular axis"? >> >> If yes, is this definition applicable to all lipids, e.g. PC, PE, PS >> lipids, disregards the angle the headgroups makes with the lipid/water >> interface? >> >> >> I looked at the literature but I could not find any definition. >> For instance, such an angle was calculated in the following article. >> >> https://www.ncbi.nlm.nih.gov/pubmed/17598103 >> >> However, I could not find how the molecular director axis was defined. >> > > It appears to be the principal axis of the lipid molecule, but probably > the authors of that study can tell you how it was done, either via GROMACS > tools or their own code. I'm not sure how GROMACS tools would calculate > that quantity. > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Assistant Professor > Virginia Tech Department of Biochemistry > > 303 Engel Hall > 340 West Campus Dr. > Blacksburg, VA 24061 > > jalem...@vt.edu | (540) 231-3129 > http://www.biochem.vt.edu/people/faculty/JustinLemkul.html > > == > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] lipid molecular director axis definition
Hi Everyone, I have a system formed of lipids (DOPE). I would like to know the average angle of the "molecular director axis" with z-direction. Using g_angle it should be easy. However, my question is regarding the molecular director axis. Is there any definition for how to determine a "lipid molecular axis"? If yes, is this definition applicable to all lipids, e.g. PC, PE, PS lipids, disregards the angle the headgroups makes with the lipid/water interface? I looked at the literature but I could not find any definition. For instance, such an angle was calculated in the following article. https://www.ncbi.nlm.nih.gov/pubmed/17598103 However, I could not find how the molecular director axis was defined. Thanks in advance for your reply Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Order Parameter for HII phase
Hi Everyone, Like to know your thoughts on this. I appreciate any comment on this in advance. Cheers, Mohsen On Wed, Oct 18, 2017 at 7:26 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Justin, Tom, Antonio, > > I read the articles recommended in the previous emails. > To get the order parameter, for any system and by any technique (e.g. > g_angle), we need to calculate the C-D vector angle with respect to the > "local normal vector" at the place of that lipid. > > For the bilayer, it is the Z direction which is the normal vector to the > bilayer surface. Easy. > > For the HII, however, this "local normal vector" is challenging. > The cylinders in HII (which is long enough) are not perfect cylinders and > are bent (in a shape like a Tilda "~" ). Even when the cylinders are not > bent (i.e.they are like "__", still a "local normal vector" need to be > calculated in this calculation. > > > Here are four solutions I am thinking about: > > 1) Using the "radial" option in g_order, but I am not sure if this will > give what I want. > > 2) Slicing the HII cylinders in Z direction (i.e. the cylinder axis), > calculating a center of mass for those lipids in each slice, and then using > the g_angle and its "-g2 sphnorm". In this way, the COM of that slice will > be taken as the center of the sphere (an approximation). > > 3) Using g_angle but using "-g1 plane": in this case, the coordinates of > three atoms in that lipid at each frame will be used to make the local > normal vector. > > 4) Maybe a combination of g_sgangle and g_angle/g_order. The > g_sgangle seems to give the local normal vector for each time frame. This > can be subsequently used in g_order or g_angle for angle calculation. > > One of the main problems in methods 3 and 4 is that: the way the three > atoms on the lipid are selected, the normal vector will be different, > especially if the lipid is tilted. > > Is there any other way which I could calculate the C-D angle and order > parameter in HII phase? > > Please let me know your opinion. Any help is appreciated in advance. > > Cheers, > Mohsen > > On Thu, Feb 23, 2017 at 4:55 PM, Mohsen Ramezanpour < > ramezanpour.moh...@gmail.com> wrote: > >> Thanks Antonio. >> I will give them a try :-) >> >> Cheers >> Mohsen >> >> On Thu, Feb 23, 2017 at 4:15 PM, Antonio Baptista >> wrote: >> >>> Hi Mohsen, >>> >>> I suggest you have a look at the QRB 1977 review by Seelig, that Tom >>> already mentioned. Like for the planar case, they discuss how spectra of >>> cylindrical lipid phases are related to the order parameter tensor (but I >>> never looked into the details for that case). Anyway, I don't know if you >>> want to compare your simulations with experimental data or are just looking >>> for a convenient structural parameter. But, whatever the case is, reading >>> the literature and thinking about the geometry of your system (Duliez's >>> papers are a good train for that) should give you some hints for a relevant >>> parameter for your system. >>> >>> Once you have selected a parameter, you can compute it from stratch >>> using other GROMACS tools besides g_order. In particular, you can use >>> g_gangle to get several sorts of angles and then process them yourself. As >>> an example, you can use this approach to compute SCD using the equation >>> already mentioned by Tom, and then check if it agrees with the SCD given by >>> g_order (I once did that, just to be sure). >>> >>> Good luck! :) >>> >>> Best, >>> Antonio >>> >>> >>> >>> On Thu, 23 Feb 2017, Justin Lemkul wrote: >>> >>> >>>> >>>> On 2/23/17 1:25 PM, Mohsen Ramezanpour wrote: >>>> >>>>> And I agree with Piggot. The paper by Chau is about on option in >>>>> g_order. >>>>> >>>>> >>>> Yes, and if memory serves (been a while since I used the program, so >>>> perhaps I am confusing something), this is what -o provides you so I >>>> thought it was relevant. -od gives you the deuterium order parameters. >>>> The documentation for gmx order is indeed very sparse. >>>> >>>> Anyway, a general question: >>>>> Can we expect to find a published article for each/some module(s) in >>>>> Gromacs? >>>>> >>>> >>>> No. Many of the tools are just
Re: [gmx-users] Bug?! Lost particles while sorting
Hi Mark,
Thanks for replying on this email.
I forgot to reply again and mention that it was not a bug.
I could fix it simply by following some instructions in:
https://docs.computecanada.ca/wiki/GROMACS
It is working perfectly fine now (GPU with version 2016.3)
Here are two things I changed in my script when submitted the jobs on the
cluster:
1) replacing "srun" with "mpiexec"
2) including the following line and removing the -ntomp 1 in mdrun command
export OMP_NUM_THREADS="${SLURM_CPUS_PER_TASK:-1}"
Cheers,
Mohsen
On Fri, Oct 13, 2017 at 4:02 AM, Mark Abraham
wrote:
> Hi,
>
> Unfortunately that will just get lost in my inbox, and anyway I will have
> to go and make the redmine issue. We need that (rather than email) because
> multiple people might find the issue and need to be able to search for the
> answer. Also, somebody else fixing the issue won't be able to contact you
> easily if we need to follow up to get further information, or check that
> the issue is resolved.
>
> To file a request, please follow the link on https://redmine.gromacs.org/
> after
> you register an account. All quite easy, I hope! Do let us know if some
> aspect is hard.
>
> Mark
>
> On Thu, Aug 31, 2017 at 3:04 PM Mohsen Ramezanpour <
> ramezanpour.moh...@gmail.com> wrote:
>
> > Hi Mark,
> >
> > Thanks for your reply.
> > I tried but I was not sure how to do so.
> > May I send the corresponding files to your personal email, instead?
> >
> > Cheers,
> > Mohsen
> >
> >
> > On Wed, Aug 30, 2017 at 6:25 PM, Mark Abraham
> > wrote:
> >
> > > Hi,
> > >
> > > That particular output indicates that the code is not working as
> > intended.
> > > Please open an issue on the GROMACS redmine and attach the tpr, the log
> > > file from your run, and and instructions on how to reproduce.
> > >
> > > Thanks,
> > >
> > > Mark
> > >
> > > On Wed, 30 Aug 2017 19:34 Mohsen Ramezanpour <
> > ramezanpour.moh...@gmail.com
> > > >
> > > wrote:
> > >
> > > > Dear Gromacs users,
> > > >
> > > > I am running simulations using Gromacs version 2016.3 while using
> GPU.
> > > >
> > > > I get an error in my simulations as follows:
> > > >
> > > > Program: gmx mdrun, version 2016.3
> > > > Source file: src/gromacs/mdlib/nbnxn_grid.cpp (line 404)
> > > > MPI rank:15 (out of 16)
> > > >
> > > > Software inconsistency error:
> > > > Lost particles while sorting
> > > >
> > > > For more information and tips for troubleshooting, please check the
> > > GROMACS
> > > > website at http://www.gromacs.org/Documentation/Errors
> > > >
> > > >
> > > > Googling the error I found that there were relevant bugs before
> > > (discussed
> > > > by Dr. Berk Hess):
> > > >
> > > > https://redmine.gromacs.org/issues/1379
> > > >
> > > > https://redmine.gromacs.org/issues/1153
> > > >
> > > > However, they should have been fixed in the newer versions as
> mentioned
> > > in
> > > > the above links.
> > > >
> > > > Can someone please let me know what the problem could be in this
> case?
> > > >
> > > > Many thanks in advance,
> > > >
> > > > Cheers,
> > > > Mohsen
> > > > --
> > > > *Rewards work better than punishment ...*
> > > > --
> > > > Gromacs Users mailing list
> > > >
> > > > * Please search the archive at
> > > > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > > > posting!
> > > >
> > > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> > > >
> > > > * For (un)subscribe requests visit
> > > > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users
> or
> > > > send a mail to gmx-users-requ...@gromacs.org.
> > > >
> > > --
> > > Gromacs Users mailing list
> > >
> > > * Please search the archive at http://www.gromacs.org/
> > > Support/Mailing_Lists/GMX-Users_List before posting!
> > >
> > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> > >
> > > * For (un)subscribe requests visit
> > > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > &g
Re: [gmx-users] Order Parameter for HII phase
Hi Justin, Tom, Antonio, I read the articles recommended in the previous emails. To get the order parameter, for any system and by any technique (e.g. g_angle), we need to calculate the C-D vector angle with respect to the "local normal vector" at the place of that lipid. For the bilayer, it is the Z direction which is the normal vector to the bilayer surface. Easy. For the HII, however, this "local normal vector" is challenging. The cylinders in HII (which is long enough) are not perfect cylinders and are bent (in a shape like a Tilda "~" ). Even when the cylinders are not bent (i.e.they are like "__", still a "local normal vector" need to be calculated in this calculation. Here are four solutions I am thinking about: 1) Using the "radial" option in g_order, but I am not sure if this will give what I want. 2) Slicing the HII cylinders in Z direction (i.e. the cylinder axis), calculating a center of mass for those lipids in each slice, and then using the g_angle and its "-g2 sphnorm". In this way, the COM of that slice will be taken as the center of the sphere (an approximation). 3) Using g_angle but using "-g1 plane": in this case, the coordinates of three atoms in that lipid at each frame will be used to make the local normal vector. 4) Maybe a combination of g_sgangle and g_angle/g_order. The g_sgangle seems to give the local normal vector for each time frame. This can be subsequently used in g_order or g_angle for angle calculation. One of the main problems in methods 3 and 4 is that: the way the three atoms on the lipid are selected, the normal vector will be different, especially if the lipid is tilted. Is there any other way which I could calculate the C-D angle and order parameter in HII phase? Please let me know your opinion. Any help is appreciated in advance. Cheers, Mohsen On Thu, Feb 23, 2017 at 4:55 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Thanks Antonio. > I will give them a try :-) > > Cheers > Mohsen > > On Thu, Feb 23, 2017 at 4:15 PM, Antonio Baptista > wrote: > >> Hi Mohsen, >> >> I suggest you have a look at the QRB 1977 review by Seelig, that Tom >> already mentioned. Like for the planar case, they discuss how spectra of >> cylindrical lipid phases are related to the order parameter tensor (but I >> never looked into the details for that case). Anyway, I don't know if you >> want to compare your simulations with experimental data or are just looking >> for a convenient structural parameter. But, whatever the case is, reading >> the literature and thinking about the geometry of your system (Duliez's >> papers are a good train for that) should give you some hints for a relevant >> parameter for your system. >> >> Once you have selected a parameter, you can compute it from stratch using >> other GROMACS tools besides g_order. In particular, you can use g_gangle to >> get several sorts of angles and then process them yourself. As an example, >> you can use this approach to compute SCD using the equation already >> mentioned by Tom, and then check if it agrees with the SCD given by g_order >> (I once did that, just to be sure). >> >> Good luck! :) >> >> Best, >> Antonio >> >> >> >> On Thu, 23 Feb 2017, Justin Lemkul wrote: >> >> >>> >>> On 2/23/17 1:25 PM, Mohsen Ramezanpour wrote: >>> >>>> And I agree with Piggot. The paper by Chau is about on option in >>>> g_order. >>>> >>>> >>> Yes, and if memory serves (been a while since I used the program, so >>> perhaps I am confusing something), this is what -o provides you so I >>> thought it was relevant. -od gives you the deuterium order parameters. >>> The documentation for gmx order is indeed very sparse. >>> >>> Anyway, a general question: >>>> Can we expect to find a published article for each/some module(s) in >>>> Gromacs? >>>> >>> >>> No. Many of the tools are just implementations of simple algorithms >>> used widely in simulations. Some do have specific references and those are >>> generally printed to the screen output in bold blocks of text. >>> >>> With regards to deuterium order parameters, there is a well accepted and >>> ubiquitously used equation, which Tom posted. >>> >>> I mean how/where can we figure out the underlying algorithm and >>>> comprehensive description of each analysis tool? >>>> >>>> >>> That's why GROMACS is open source :) Efforts are certainly always made >>> to pro
Re: [gmx-users] p-coupling for inverted hexagonal phase
I did this test already, actually. The main problem I am facing right now is that the semi-isotropic test systems get equilibrated (judging on box size for now) after ca. 50 ns or less. BUT for the anisotropic coupling, the box sizes 1) fluctuate a lot and 2) does not seem to be equilibrated even after 300 ns. This makes me less confident on the results, especially with the point you mentioned on the effect of anisotropic after a long simulation time. Cheers, Mohsen On Thu, Sep 28, 2017 at 5:51 AM, Justin Lemkul wrote: > > > On 9/28/17 12:31 AM, Mohsen Ramezanpour wrote: > >> Thanks Justin for your comment, >> >> I think applying the semi-isotropic will put too much symmetry on the >> system while it might not be the case and system should be free to change >> if it is energetically favourable. >> The semi-isotropic, in fact, is scaling both the X and Y sides of the >> simulation box with the same factor, if I understood correctly. So, it >> will >> keep the original geometry of the system the same through the simulation >> time. >> Think of situations where the water channels prefer to be elongated in >> either X or Y direction. This will prevent this, right? >> Also, if there is any phase transition possible (e.g. from HII to >> Lamellar), semi-isotropic will probably limit this transition more than >> what an anisotropic coupling will do. >> >> Please let me know your thoughts on this. >> > > It seems like you already have a rationale for what to do. In the absence > of literature precedent, you should probably test both to confirm your > suspicions and validate your approach. Just be aware of potential artifacts > due to anisotropic pressure coupling. They can be quite extreme over the > course of a long simulation. > > -Justin > > > Cheers, >> Mohsen >> >> >> On Wed, Sep 27, 2017 at 6:26 PM, Justin Lemkul wrote: >> >> >>> On 9/27/17 10:57 AM, Mohsen Ramezanpour wrote: >>> >>> Hi Everyone, >>>> >>>> Please let me know your thoughts on this. >>>> >>>> >>>> I don't see why these would be anisotropic. The lateral dimension is >>> symmetric, so it seems that semiisotropic would be most appropriate. >>> Anisotropic pressure coupling can lead to significant distortions of the >>> box, and is most applicable in crystalline systems. >>> >>> -Justin >>> >>> Thanks in advance, >>> >>>> Cheers, >>>> Mohsen >>>> >>>> On Fri, Sep 22, 2017 at 11:13 AM, Mohsen Ramezanpour < >>>> ramezanpour.moh...@gmail.com> wrote: >>>> >>>> Dear Gromacs users, >>>> >>>>> I am doing a simulation on inverted hexagonal (HII) phase composed of >>>>> just >>>>> lipids and water. >>>>> >>>>> (please have a look at HII phase here: >>>>> https://openi.nlm.nih.gov/detailedresult.php?img= >>>>> PMC2695813_1757-5036-2-3-3&req=4) >>>>> >>>>> My question is regarding the pressure coupling for such systems. I am >>>>> using the anisotropic p-coupling at the moment: >>>>> >>>>> pcoupl = Parrinello-Rahman; replaced >>>>> pcoupltype = anisotropic >>>>> tau_p = 5.0 >>>>> compressibility = 4.5e-5 4.5e-5 4.5e-5 0 0 0 >>>>> ref_p = 1.0 1.0 1.0 0 0 0 >>>>> >>>>> I think this is the correct way to treat this type of system. >>>>> >>>>> Some might argue that a semi-isotropic p-coupling is the correct >>>>> treatment for that. i.e. >>>>> >>>>> pcoupl = Parrinello-Rahman; replaced >>>>> pcoupltype = anisotropic >>>>> tau_p = 5.0 >>>>> compressibility = 4.5e-5 4.5e-5 >>>>> ref_p = 1.0 1.0 >>>>> >>>>> where the first is for Z (the cylindrical axis of HII phase) and the >>>>> second parameters control and scale the box size in both X and Y >>>>> directions >>>>> equally. >>>>> >>>>> I just like to know your opinions on this. Which one do you think is >>>>> the >>>>> better way to treat such systems? >>>>> >>>>> Thanks in advance for your comments. >>>>> Mohsen &
Re: [gmx-users] p-coupling for inverted hexagonal phase
Thanks Justin for your comment, I think applying the semi-isotropic will put too much symmetry on the system while it might not be the case and system should be free to change if it is energetically favourable. The semi-isotropic, in fact, is scaling both the X and Y sides of the simulation box with the same factor, if I understood correctly. So, it will keep the original geometry of the system the same through the simulation time. Think of situations where the water channels prefer to be elongated in either X or Y direction. This will prevent this, right? Also, if there is any phase transition possible (e.g. from HII to Lamellar), semi-isotropic will probably limit this transition more than what an anisotropic coupling will do. Please let me know your thoughts on this. Cheers, Mohsen On Wed, Sep 27, 2017 at 6:26 PM, Justin Lemkul wrote: > > > On 9/27/17 10:57 AM, Mohsen Ramezanpour wrote: > >> Hi Everyone, >> >> Please let me know your thoughts on this. >> >> > I don't see why these would be anisotropic. The lateral dimension is > symmetric, so it seems that semiisotropic would be most appropriate. > Anisotropic pressure coupling can lead to significant distortions of the > box, and is most applicable in crystalline systems. > > -Justin > > Thanks in advance, >> Cheers, >> Mohsen >> >> On Fri, Sep 22, 2017 at 11:13 AM, Mohsen Ramezanpour < >> ramezanpour.moh...@gmail.com> wrote: >> >> Dear Gromacs users, >>> >>> I am doing a simulation on inverted hexagonal (HII) phase composed of >>> just >>> lipids and water. >>> >>> (please have a look at HII phase here: >>> https://openi.nlm.nih.gov/detailedresult.php?img= >>> PMC2695813_1757-5036-2-3-3&req=4) >>> >>> My question is regarding the pressure coupling for such systems. I am >>> using the anisotropic p-coupling at the moment: >>> >>> pcoupl = Parrinello-Rahman; replaced >>> pcoupltype = anisotropic >>> tau_p = 5.0 >>> compressibility = 4.5e-5 4.5e-5 4.5e-5 0 0 0 >>> ref_p = 1.0 1.0 1.0 0 0 0 >>> >>> I think this is the correct way to treat this type of system. >>> >>> Some might argue that a semi-isotropic p-coupling is the correct >>> treatment for that. i.e. >>> >>> pcoupl = Parrinello-Rahman; replaced >>> pcoupltype = anisotropic >>> tau_p = 5.0 >>> compressibility = 4.5e-5 4.5e-5 >>> ref_p = 1.0 1.0 >>> >>> where the first is for Z (the cylindrical axis of HII phase) and the >>> second parameters control and scale the box size in both X and Y >>> directions >>> equally. >>> >>> I just like to know your opinions on this. Which one do you think is the >>> better way to treat such systems? >>> >>> Thanks in advance for your comments. >>> Mohsen >>> >>> >>> >>> -- >>> *Rewards work better than punishment ...* >>> >>> >> >> >> > -- > == > > Justin A. Lemkul, Ph.D. > Assistant Professor > Virginia Tech Department of Biochemistry > > 303 Engel Hall > 340 West Campus Dr. > Blacksburg, VA 24061 > > jalem...@vt.edu | (540) 231-3129 > http://www.biochem.vt.edu/people/faculty/JustinLemkul.html > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] p-coupling for inverted hexagonal phase
Hi Everyone, Please let me know your thoughts on this. Thanks in advance, Cheers, Mohsen On Fri, Sep 22, 2017 at 11:13 AM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Dear Gromacs users, > > I am doing a simulation on inverted hexagonal (HII) phase composed of just > lipids and water. > > (please have a look at HII phase here: > https://openi.nlm.nih.gov/detailedresult.php?img= > PMC2695813_1757-5036-2-3-3&req=4) > > My question is regarding the pressure coupling for such systems. I am > using the anisotropic p-coupling at the moment: > > pcoupl = Parrinello-Rahman; replaced > pcoupltype = anisotropic > tau_p = 5.0 > compressibility = 4.5e-5 4.5e-5 4.5e-5 0 0 0 > ref_p = 1.0 1.0 1.0 0 0 0 > > I think this is the correct way to treat this type of system. > > Some might argue that a semi-isotropic p-coupling is the correct > treatment for that. i.e. > > pcoupl = Parrinello-Rahman; replaced > pcoupltype = anisotropic > tau_p = 5.0 > compressibility = 4.5e-5 4.5e-5 > ref_p = 1.0 1.0 > > where the first is for Z (the cylindrical axis of HII phase) and the > second parameters control and scale the box size in both X and Y directions > equally. > > I just like to know your opinions on this. Which one do you think is the > better way to treat such systems? > > Thanks in advance for your comments. > Mohsen > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] p-coupling for inverted hexagonal phase
Dear Gromacs users, I am doing a simulation on inverted hexagonal (HII) phase composed of just lipids and water. (please have a look at HII phase here: https://openi.nlm.nih.gov/detailedresult.php?img=PMC2695813_1757-5036-2-3-3&req=4 ) My question is regarding the pressure coupling for such systems. I am using the anisotropic p-coupling at the moment: pcoupl = Parrinello-Rahman; replaced pcoupltype = anisotropic tau_p = 5.0 compressibility = 4.5e-5 4.5e-5 4.5e-5 0 0 0 ref_p = 1.0 1.0 1.0 0 0 0 I think this is the correct way to treat this type of system. Some might argue that a semi-isotropic p-coupling is the correct treatment for that. i.e. pcoupl = Parrinello-Rahman; replaced pcoupltype = anisotropic tau_p = 5.0 compressibility = 4.5e-5 4.5e-5 ref_p = 1.0 1.0 where the first is for Z (the cylindrical axis of HII phase) and the second parameters control and scale the box size in both X and Y directions equally. I just like to know your opinions on this. Which one do you think is the better way to treat such systems? Thanks in advance for your comments. Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Bug?! Lost particles while sorting
Hi Mark, Thanks for your reply. I tried but I was not sure how to do so. May I send the corresponding files to your personal email, instead? Cheers, Mohsen On Wed, Aug 30, 2017 at 6:25 PM, Mark Abraham wrote: > Hi, > > That particular output indicates that the code is not working as intended. > Please open an issue on the GROMACS redmine and attach the tpr, the log > file from your run, and and instructions on how to reproduce. > > Thanks, > > Mark > > On Wed, 30 Aug 2017 19:34 Mohsen Ramezanpour > > wrote: > > > Dear Gromacs users, > > > > I am running simulations using Gromacs version 2016.3 while using GPU. > > > > I get an error in my simulations as follows: > > > > Program: gmx mdrun, version 2016.3 > > Source file: src/gromacs/mdlib/nbnxn_grid.cpp (line 404) > > MPI rank:15 (out of 16) > > > > Software inconsistency error: > > Lost particles while sorting > > > > For more information and tips for troubleshooting, please check the > GROMACS > > website at http://www.gromacs.org/Documentation/Errors > > > > > > Googling the error I found that there were relevant bugs before > (discussed > > by Dr. Berk Hess): > > > > https://redmine.gromacs.org/issues/1379 > > > > https://redmine.gromacs.org/issues/1153 > > > > However, they should have been fixed in the newer versions as mentioned > in > > the above links. > > > > Can someone please let me know what the problem could be in this case? > > > > Many thanks in advance, > > > > Cheers, > > Mohsen > > -- > > *Rewards work better than punishment ...* > > -- > > Gromacs Users mailing list > > > > * Please search the archive at > > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before > > posting! > > > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > > > * For (un)subscribe requests visit > > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > > send a mail to gmx-users-requ...@gromacs.org. > > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/ > Support/Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Bug?! Lost particles while sorting
Dear Gromacs users, I am running simulations using Gromacs version 2016.3 while using GPU. I get an error in my simulations as follows: Program: gmx mdrun, version 2016.3 Source file: src/gromacs/mdlib/nbnxn_grid.cpp (line 404) MPI rank:15 (out of 16) Software inconsistency error: Lost particles while sorting For more information and tips for troubleshooting, please check the GROMACS website at http://www.gromacs.org/Documentation/Errors Googling the error I found that there were relevant bugs before (discussed by Dr. Berk Hess): https://redmine.gromacs.org/issues/1379 https://redmine.gromacs.org/issues/1153 However, they should have been fixed in the newer versions as mentioned in the above links. Can someone please let me know what the problem could be in this case? Many thanks in advance, Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] reading topology by Gromacs
Dear Gromacs users, I want to add the modified parameters for several proper dihedrals. So, I want to add them to the section [ dihedraltypes ] in ffbonded.itp. There are two dihedrals already for the combination I want: C HC HC NC2 2 0.00 0.00 . . . NC2XXC 2 0.00 376.56 If I have the following dihedral (just for the sake of example): C HC HC NC2 2 angle force_constant Where shall I add this line exactly? Since the angles can be either 0 or 180 in Charmm36, the parameters for: C HC HC NC2 should be the same with: NC2 HC HC C In other words, does Gromacs take the first or the last match it finds in ffbonded.itp for the assigned atom types in desire dihedral? I had a look at .tpr file using gmx dump, and I just can see the the type number, e.g. 574, but I do not know where can I see the parameters for these types (for gromacs version of Charmm36FF). Thanks in advance for your comments Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Atom types comparison between CGenFF and Charm36FF
Thanks, Justin. Sure. On Wed, Jun 28, 2017 at 10:21 AM, Justin Lemkul wrote: > > > On 6/28/17 9:10 AM, Mohsen Ramezanpour wrote: > >> Thanks Justin for your comment. >> >> I have a bit of difficulty for the finding the analog parts to these two >> parts. >> >> Is there any database (preferably with shapes) for the molecules available >> in charmmFF? >> >> > The CHARMM topology files have full residue names and chemical structure > drawings. Within the toppar directory (download the tarball from our > site), just grep -r "RESI" * > > If you have additional questions about this or CHARMM parametrization, > contact me off list as this is not really a GROMACS issue. > > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Atom types comparison between CGenFF and Charm36FF
Thanks Justin for your comment. I have a bit of difficulty for the finding the analog parts to these two parts. Is there any database (preferably with shapes) for the molecules available in charmmFF? Cheers, Mohsen On Tue, Jun 27, 2017 at 7:51 PM, Justin Lemkul wrote: > > > On 6/27/17 3:27 PM, Mohsen Ramezanpour wrote: > >> Dear Gromacs Users, >> >> I am trying to parameterize a molecule in Charmm36FF. >> >> As part of this molecule, there is a "neutral trimethylamine nitrogen" and >> its protonated case. i.e. C-N-C(2) and C-N(H)-C(2), respectively. >> >> To have some idea about the appropriate atom types, I first used the >> GAAMP server. >> >> Here are the assigned atom types by GAAMP (based on CGenFF force field): >> >> >> *For Neutral case:* >> NG301 for N, >> CG3AM0 for two methyl group connected to the N, i.e for bolded Carbons in >> C-N-*C(2).* >> HGA3 for each Hydrogen atom in CG3AM0 methyl groups, >> CG321 for the carbon atom connecting the NC(2) to the rest of molecule. >> i.e. for bolded Carbon in *C*-N-C(2), >> and HGA2 for the Hydrogen atoms connected to CG321. >> ### >> >> *For Protonated case:* >> NG3P1 for N, >> CG334 for two methyl group connected to the N, i.e for bolded Carbons in >> C-N(H)-*C(2).* >> HGA3 for each Hydrogen atom in CG334 methyl groups, >> CG324 for the carbon atom connecting the NC(2) to the rest of molecule. >> i.e. for bolded Carbon in *C*-N(H)-C(2), >> and HGA2 for the Hydrogen atoms connected to CG321. >> >> >> Now, I want to find the best equivalent atom type in Charmm36FF for each >> case. >> >> For the protonated case, I think it is easier as it is something between >> PE >> and PC lipid headgroups. For the neutral case, however, it is difficult to >> find similar atom types in Charmm36. >> >> My approach was to check the LJ and bonded parameters for the assigned >> atom >> types (which are in CGenFF) and find the atom types in Charmm36 with the >> exact same values. Although possible for some cases, but there are many >> problems with this approach: >> >> 1) the exact values are not found in Charmm36FF. >> 2) If it is found, there are a lot of parameters missing in the Charmm36FF >> force field. Not all the bonds, angles, and dihedrals are defined in >> Charmm36FF. >> >> Thanks in advance for any comment or suggestion. >> >> > The "G" in CGenFF is for "general," which means the parameters are not > necessarily the same as the parent CHARMM force field, and are subsequently > a compromise between being highly optimized (e.g. CHARMM36) and being > broadly applicable (general) such that the types can be used across > different molecules. > > If you want to parametrize something for CHARMM, e.g. to be merged into a > larger molecule, you're wasting your time with trying to generate a CGenFF > topology and try to find exact matches in CHARMM36. By definition, you > won't. If your goal is a CHARMM topology, then start from existing CHARMM > atom types, import charges by analogy, and do a full parametrization > procedure (well described in numerous places in the literature). > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Atom types comparison between CGenFF and Charm36FF
Dear Gromacs Users, I am trying to parameterize a molecule in Charmm36FF. As part of this molecule, there is a "neutral trimethylamine nitrogen" and its protonated case. i.e. C-N-C(2) and C-N(H)-C(2), respectively. To have some idea about the appropriate atom types, I first used the GAAMP server. Here are the assigned atom types by GAAMP (based on CGenFF force field): *For Neutral case:* NG301 for N, CG3AM0 for two methyl group connected to the N, i.e for bolded Carbons in C-N-*C(2).* HGA3 for each Hydrogen atom in CG3AM0 methyl groups, CG321 for the carbon atom connecting the NC(2) to the rest of molecule. i.e. for bolded Carbon in *C*-N-C(2), and HGA2 for the Hydrogen atoms connected to CG321. ### *For Protonated case:* NG3P1 for N, CG334 for two methyl group connected to the N, i.e for bolded Carbons in C-N(H)-*C(2).* HGA3 for each Hydrogen atom in CG334 methyl groups, CG324 for the carbon atom connecting the NC(2) to the rest of molecule. i.e. for bolded Carbon in *C*-N(H)-C(2), and HGA2 for the Hydrogen atoms connected to CG321. Now, I want to find the best equivalent atom type in Charmm36FF for each case. For the protonated case, I think it is easier as it is something between PE and PC lipid headgroups. For the neutral case, however, it is difficult to find similar atom types in Charmm36. My approach was to check the LJ and bonded parameters for the assigned atom types (which are in CGenFF) and find the atom types in Charmm36 with the exact same values. Although possible for some cases, but there are many problems with this approach: 1) the exact values are not found in Charmm36FF. 2) If it is found, there are a lot of parameters missing in the Charmm36FF force field. Not all the bonds, angles, and dihedrals are defined in Charmm36FF. Thanks in advance for any comment or suggestion. Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Gromacs on GPU
Thanks Stephane, Sure, but how if I am in the middle of simulation and must switch to GPU? I have 50 systems each around 100 ns, but I need 500 ns of each. :-) What would you suggest in this case? Cheers, Mohsen On Fri, Jun 16, 2017 at 12:20 PM, Téletchéa Stéphane < stephane.teletc...@univ-nantes.fr> wrote: > Le 16/06/2017 à 20:07, Mohsen Ramezanpour a écrit : > >> Thanks Justin. >> >> So, can we say that simulation on CPU and GPU (as far as we use the same >> version of Gromacs) are compatible? >> >> If yes, is that okay to continue a simulation which was done with CPU (say >> till 100 ns) to 500 ns (using GPU)? >> Or I should start from t=0 with GPU? >> >> It is important for my case as the allocations on supercomputers change >> from CPU to GPU. So, not sure if I should start all again or it is fine to >> continue. >> >> Cheers >> > > Dear Moshen, > > What I'm doing daily is using my "small" workstation with CPU then > continue once equilibration and primary production has worked well. > > Go for GPU for a boost in performance (big) once the primary steps are ok > :-) > > Best, > > Stéphane > > -- > Assistant Professor in BioInformatics, UFIP, UMR 6286 CNRS, Team Protein > Design In Silico > UFR Sciences et Techniques, 2, rue de la Houssinière, Bât. 25, 44322 > Nantes cedex 03, France > Tél : +33 251 125 636 / Fax : +33 251 125 632 > http://www.ufip.univ-nantes.fr/ - http://www.steletch.org > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Gromacs on GPU
Thanks Justin. So, can we say that simulation on CPU and GPU (as far as we use the same version of Gromacs) are compatible? If yes, is that okay to continue a simulation which was done with CPU (say till 100 ns) to 500 ns (using GPU)? Or I should start from t=0 with GPU? It is important for my case as the allocations on supercomputers change from CPU to GPU. So, not sure if I should start all again or it is fine to continue. Cheers On Fri, Jun 16, 2017 at 12:00 PM, Justin Lemkul wrote: > > > On 6/16/17 1:55 PM, Mohsen Ramezanpour wrote: > >> Dear Gromacs users, >> >> I would like to use *Gromacs on GPU on supercomputers*. On my desktop, I >> am >> running on GPU by default. I use the same .mdp parameters and usual >> Gromacs >> commands which I use for CPUs. Everything looks fine to the best of my >> knowledge. >> >> I was reading a presentation by Dr. Lindahl on using Gromacs on GPU which >> makes me worried a bit. >> >> http://on-demand.gputechconf.com/gtc/2013/webinar/gromacs-ke >> pler-gpus-gtc-express-webinar.pdf >> >> >> I want to use the charmm36 FF in Gromacs version 5-1-4 for my project. >> I would prefer to keep the parameters as the recommended parameters for >> using Charmm36 in Gromacs package. >> >> 1) I was wondering if there is any special change required to be made in >> .mdp file? >> >> > No. > > 2) Is that necessary to do all the steps (EM, NVT, ...) with GPU or I can >> simply equilibrate my system on my Desktop (say CPU) and then use GPU for >> just the production run? >> >> > You do not have to do everything on the GPU. You can do any or all parts > of your protocol (though I've never tried EM on a GPU). > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Gromacs on GPU
Dear Gromacs users, I would like to use *Gromacs on GPU on supercomputers*. On my desktop, I am running on GPU by default. I use the same .mdp parameters and usual Gromacs commands which I use for CPUs. Everything looks fine to the best of my knowledge. I was reading a presentation by Dr. Lindahl on using Gromacs on GPU which makes me worried a bit. http://on-demand.gputechconf.com/gtc/2013/webinar/gromacs-kepler-gpus-gtc-express-webinar.pdf I want to use the charmm36 FF in Gromacs version 5-1-4 for my project. I would prefer to keep the parameters as the recommended parameters for using Charmm36 in Gromacs package. 1) I was wondering if there is any special change required to be made in .mdp file? 2) Is that necessary to do all the steps (EM, NVT, ...) with GPU or I can simply equilibrate my system on my Desktop (say CPU) and then use GPU for just the production run? Thanks in advance for your replies, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] bad contacts on lipid tails
Hi Justin, I could solve the problem. Recognizing the lipids which had overlap and subsequent removal of one of them in each pair could solve it. Now, it can be run with a time step of 2 without any problem. Cheers, Mohsen On Wed, Jun 14, 2017 at 3:26 PM, Justin Lemkul wrote: > > > On 6/14/17 5:23 PM, Mohsen Ramezanpour wrote: > >> I have done equilibration for almost 45 ns or even 60 ns for one test >> system. >> >> I expected with this energy and this equilibration time, it works fine. >> But >> it fails. >> >> I will construct the systems as you suggested to see what I get. >> >> > With such a long equilibration, there's no reason it should fail with dt = > 2 fs if the clashes have been eliminated, because I would expect they > should be gone after a few ps. I'd suspect something else in the topology > or .mdp file but you'll have to do a more thorough diagnosis (e.g. > http://www.gromacs.org/Documentation/Terminology/Blowing_Up# > Diagnosing_an_Unstable_System) > > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] bad contacts on lipid tails
I have done equilibration for almost 45 ns or even 60 ns for one test system. I expected with this energy and this equilibration time, it works fine. But it fails. I will construct the systems as you suggested to see what I get. Thanks On Wed, Jun 14, 2017 at 3:20 PM, Justin Lemkul wrote: > > > On 6/14/17 5:18 PM, Mohsen Ramezanpour wrote: > >> Hi Justin, >> >> Thanks for your reply, >> >> Here are the outputs for EM step: >> >> Steepest Descents converged to Fmax < 100 in 5198 steps >> Potential Energy = -2.6418572e+05 >> Maximum force = 9.5666359e+01 on atom 54770 >> Norm of force = 5.7127428e+00 >> >> It seems okay to me. >> >> I chose the setting you mentioned. >> >> Interestingly, it also works for MD with dt=1 fs. Not for dt=2 fs, though. >> Based on visualization, some lipid chains are crossed which I think this >> cause the problem. Not sure how can I fix them. >> >> > Surprising. With such a low force, I wouldn't expect there to be any > overlap or detrimental close contact. An initial equilibration with dt = 1 > fs is a possibility, or a more robust method for constructing the system, > whatever that may be. > > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] bad contacts on lipid tails
Hi Justin, Thanks for your reply, Here are the outputs for EM step: Steepest Descents converged to Fmax < 100 in 5198 steps Potential Energy = -2.6418572e+05 Maximum force = 9.5666359e+01 on atom 54770 Norm of force = 5.7127428e+00 It seems okay to me. I chose the setting you mentioned. Interestingly, it also works for MD with dt=1 fs. Not for dt=2 fs, though. Based on visualization, some lipid chains are crossed which I think this cause the problem. Not sure how can I fix them. Cheers, Mohsen On Wed, Jun 14, 2017 at 3:02 PM, Justin Lemkul wrote: > > > On 6/14/17 1:05 PM, Mohsen Ramezanpour wrote: > >> Dear Gromacs users, >> >> I am doing simulations in HII phase in absence of any sugar or alkane to >> fill the gaps (in other words, HII in water) >> In few of my systems (which I made), there are lipids which their tails >> have bad contacts. Based on visualization, the CH2-CH2 bond from one tail >> in lipid 1 cross the CH2-CH2 bond in the tail from an adjacent lipid. They >> kinda stuck there. >> >> Interestingly, these systems are running properly during EM step, as well >> as MD with a time step of 1 fs but fail with 2 fs, even after 45 ns of >> simulation with 1 fs. The lipid-lipid contacts cause LINCS errors for >> angles on those bonds and involved hydrogens, and eventually, the system >> will blow up. >> >> This is what I understand after many tests to figure out the problem for >> this failure (segmentation fault (core dumped)). The topology of the >> molecule is correct and there is not water molecule which causes the >> problem. I also did simulations on one lipid both in the vacuum and in >> water. Everything is fine. The system setup is the only thing which I am >> suspicious about. >> >> Is there any trick to: >> >> 1) recognize these bad contacts? (I deleted few lipids but still, there >> are >> other lipids which cause the problem) >> >> > Nope, aside from your eyeballs and inspection of maximum forces following > EM that might point to a problem. > > 2) remove these bad contacts? (applying stronger bond and angle force >> constant, changing the LINCS parameters, or something like these?) >> >> > Changing the force field is unwise. Changing LINCS parameters can make > the constraint algorithm more permissive and can overcome problems, but > that can mask problems. > > If your system is crashing despite "successful" EM, that means you > probably only found a metastable state. One thing that we have recently > found is that the default settings for the GROMACS steepest descent > minimizer can lead to premature termination but a decrease of initial step > size to e.g. 0.002 (instead of 0.01, which is the default) does a much > better job, especially when we compare it to the CHARMM minimizer (which is > pretty bulletproof and does a great job of resolving nasty clashes that > GROMACS doesn't handle well with default settings). > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] bad contacts on lipid tails
Dear Gromacs users, I am doing simulations in HII phase in absence of any sugar or alkane to fill the gaps (in other words, HII in water) In few of my systems (which I made), there are lipids which their tails have bad contacts. Based on visualization, the CH2-CH2 bond from one tail in lipid 1 cross the CH2-CH2 bond in the tail from an adjacent lipid. They kinda stuck there. Interestingly, these systems are running properly during EM step, as well as MD with a time step of 1 fs but fail with 2 fs, even after 45 ns of simulation with 1 fs. The lipid-lipid contacts cause LINCS errors for angles on those bonds and involved hydrogens, and eventually, the system will blow up. This is what I understand after many tests to figure out the problem for this failure (segmentation fault (core dumped)). The topology of the molecule is correct and there is not water molecule which causes the problem. I also did simulations on one lipid both in the vacuum and in water. Everything is fine. The system setup is the only thing which I am suspicious about. Is there any trick to: 1) recognize these bad contacts? (I deleted few lipids but still, there are other lipids which cause the problem) 2) remove these bad contacts? (applying stronger bond and angle force constant, changing the LINCS parameters, or something like these?) I highly appreciate your help in advance, Cheers Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] auto correlation function and time
Hi guys, I have a quick question on calculating the autocorrelation time: I have the file for a dihedral angle vs time (named angle.xvg). Now, I need to calculate the autocorrelation time for this angle. There are two problems, though: 1) If I use the angle.xvg with g_angle itself (i.e. using the -oc option), the graph for the autocorrelation function is different than if I use g_analyze -ac option. Here are the commands I used: gmx angle -f md.xtc -n index.ndx -ov test.xvg -oc test-corr1.xvg -type dihedral gmx analyze -f test.xvg-ac test-corr2.xvg 2) given the output, either test-corr1.xvg or test-corr2.xvg, how can I calculate the "correlation time" for this property? please consider that this is an angle. Googling, I found different things. One told to integrate the autocorrelation function, and somewhere else, I found that I should first fit the autocorrelation function with an exponential function and then take the decay time as the correlation time. Not sure what is the correct way to do it. Thanks in advance for your reply and comments. Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] long equilibration for HII phase
Dear Gromacs users, I am using Charmm36 in Gromacs for doing MD on HII phase. I use the recommended optimized parameters by Charmm-GUI for using Charmm36 in Gromacs, although I made the HII phase systems. all the systems are running without any problem with a time step of 1 fs. The systems seem to be equilibrated (they have been run for almost 45 ns with such time step). However, still, when I used a dt of 2 fs, it crashes and some water molecules jump away. What would you suggest for this situation? Shall I use 1fs for the rest or there is another trick? Cheers, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] RMSF per residue
Hi Justin, Here is what I understand. I would like to know if I understood correctly. Given a selection, which has all the heavy atoms (5 atoms) in the side chain of one specific residue (say Lysine), I want to have only "one value" for the average fluctuation of this selection. Thus, I will use "-res" option in gmx_rmsf command. If I understood correctly, the RMSF for each atom in this selection is calculated first. So, we will have the rmsf(i) values, which i=1 to 5. Next, a mass-weighted average will be taken over all of these atoms. Thus, the final value I get as output is the following: RMSF per selection = (Sum (m(i)*rmsf(i)))/(Sum m(i)) where m(i) is the mass for atom i. Is this correct? If yes, why the RMSF per selection should be a mass-weighted? In other words, whey the usual average is not used, e.g. (Sum (rmsf(i))/5) where 5 is the total number of atoms in this selection. Thanks Mohsen On Thu, May 25, 2017 at 12:54 PM, Justin Lemkul wrote: > > > On 5/24/17 4:12 PM, Mohsen Ramezanpour wrote: > >> Dear Gromacs users, >> >> I have a question on RMSF calculation: >> >> I am not sure how -res works in RMSF and I could not find any useful >> explanation for it. >> >> It will give the average fluctuations per residue, fine. but how exactly? >> >> 1) it calculates the RMSF for each atom in that residue first. Then, it >> averages over all of these atoms and reports it as the average for the >> residue. >> >> 2) it takes the center of mass or specific atom of the residue and >> calculates the RMSF for that. >> >> Or something else? >> >> >> I am interested in the average fluctuation for sidechain of one residue >> (excluding H atoms from this calculation). I want to have only one value >> as >> the average fluctuation for this group. >> >> I choose the residue-sidechain-H group (which I made in an index file) >> when >> prompted in the following command: >> >> gmx rmsf -f md.xtc -fit -s md.tpr -n index.ndx -o >> residues-sidechain-H.xvg -res >> >> >> an alternative would be to use the following command and make an average >> over all atoms manually: >> >> gmx rmsf -f md.xtc -fit -s md.tpr -n index.ndx -o >> residues-sidechain-H.xvg >> >> >> The results are NOT the same as expected. I think I should use the second >> command and make the average manually. What do you think? >> >> > The total RMSF of the selection is computed, then a mass-weighted average > over each of the atoms is computed and that subsequent value assigned for > output. See the "average_residues" function in > src/gromacs/gmxana/gmx_rmsf.cpp > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] RMSF per residue
Hi Justin, I have a naive question on tool development for gromacs. If I change part of the code, or if I write a new analysis tool in c++ and add it to this directory (after compiling with g++), shall I still recompile the whole Gromacs again? Cheers, Mohsen On Thu, May 25, 2017 at 1:43 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Thanks Justin! > > On Thu, May 25, 2017 at 12:54 PM, Justin Lemkul wrote: > >> >> >> On 5/24/17 4:12 PM, Mohsen Ramezanpour wrote: >> >>> Dear Gromacs users, >>> >>> I have a question on RMSF calculation: >>> >>> I am not sure how -res works in RMSF and I could not find any useful >>> explanation for it. >>> >>> It will give the average fluctuations per residue, fine. but how exactly? >>> >>> 1) it calculates the RMSF for each atom in that residue first. Then, it >>> averages over all of these atoms and reports it as the average for the >>> residue. >>> >>> 2) it takes the center of mass or specific atom of the residue and >>> calculates the RMSF for that. >>> >>> Or something else? >>> >>> >>> I am interested in the average fluctuation for sidechain of one residue >>> (excluding H atoms from this calculation). I want to have only one value >>> as >>> the average fluctuation for this group. >>> >>> I choose the residue-sidechain-H group (which I made in an index file) >>> when >>> prompted in the following command: >>> >>> gmx rmsf -f md.xtc -fit -s md.tpr -n index.ndx -o >>> residues-sidechain-H.xvg -res >>> >>> >>> an alternative would be to use the following command and make an average >>> over all atoms manually: >>> >>> gmx rmsf -f md.xtc -fit -s md.tpr -n index.ndx -o >>> residues-sidechain-H.xvg >>> >>> >>> The results are NOT the same as expected. I think I should use the second >>> command and make the average manually. What do you think? >>> >>> >> The total RMSF of the selection is computed, then a mass-weighted average >> over each of the atoms is computed and that subsequent value assigned for >> output. See the "average_residues" function in >> src/gromacs/gmxana/gmx_rmsf.cpp >> >> -Justin >> >> -- >> == >> >> Justin A. Lemkul, Ph.D. >> Ruth L. Kirschstein NRSA Postdoctoral Fellow >> >> Department of Pharmaceutical Sciences >> School of Pharmacy >> Health Sciences Facility II, Room 629 >> University of Maryland, Baltimore >> 20 Penn St. >> Baltimore, MD 21201 >> >> jalem...@outerbanks.umaryland.edu | (410) 706-7441 >> http://mackerell.umaryland.edu/~jalemkul >> >> == >> -- >> Gromacs Users mailing list >> >> * Please search the archive at http://www.gromacs.org/Support >> /Mailing_Lists/GMX-Users_List before posting! >> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >> >> * For (un)subscribe requests visit >> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >> send a mail to gmx-users-requ...@gromacs.org. >> > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] RMSF per residue
Thanks Justin! On Thu, May 25, 2017 at 12:54 PM, Justin Lemkul wrote: > > > On 5/24/17 4:12 PM, Mohsen Ramezanpour wrote: > >> Dear Gromacs users, >> >> I have a question on RMSF calculation: >> >> I am not sure how -res works in RMSF and I could not find any useful >> explanation for it. >> >> It will give the average fluctuations per residue, fine. but how exactly? >> >> 1) it calculates the RMSF for each atom in that residue first. Then, it >> averages over all of these atoms and reports it as the average for the >> residue. >> >> 2) it takes the center of mass or specific atom of the residue and >> calculates the RMSF for that. >> >> Or something else? >> >> >> I am interested in the average fluctuation for sidechain of one residue >> (excluding H atoms from this calculation). I want to have only one value >> as >> the average fluctuation for this group. >> >> I choose the residue-sidechain-H group (which I made in an index file) >> when >> prompted in the following command: >> >> gmx rmsf -f md.xtc -fit -s md.tpr -n index.ndx -o >> residues-sidechain-H.xvg -res >> >> >> an alternative would be to use the following command and make an average >> over all atoms manually: >> >> gmx rmsf -f md.xtc -fit -s md.tpr -n index.ndx -o >> residues-sidechain-H.xvg >> >> >> The results are NOT the same as expected. I think I should use the second >> command and make the average manually. What do you think? >> >> > The total RMSF of the selection is computed, then a mass-weighted average > over each of the atoms is computed and that subsequent value assigned for > output. See the "average_residues" function in > src/gromacs/gmxana/gmx_rmsf.cpp > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] RMSF per residue
Dear Gromacs users, I have a question on RMSF calculation: I am not sure how -res works in RMSF and I could not find any useful explanation for it. It will give the average fluctuations per residue, fine. but how exactly? 1) it calculates the RMSF for each atom in that residue first. Then, it averages over all of these atoms and reports it as the average for the residue. 2) it takes the center of mass or specific atom of the residue and calculates the RMSF for that. Or something else? I am interested in the average fluctuation for sidechain of one residue (excluding H atoms from this calculation). I want to have only one value as the average fluctuation for this group. I choose the residue-sidechain-H group (which I made in an index file) when prompted in the following command: gmx rmsf -f md.xtc -fit -s md.tpr -n index.ndx -o residues-sidechain-H.xvg -res an alternative would be to use the following command and make an average over all atoms manually: gmx rmsf -f md.xtc -fit -s md.tpr -n index.ndx -o residues-sidechain-H.xvg The results are NOT the same as expected. I think I should use the second command and make the average manually. What do you think? Cheers Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] PBC fix for visualization
Hi Dallas, Thanks for your reply. I did try -pbc cluster for waters. It could fix it somehow but not completely. After that, I had to use -pbc center to fix it. Still, I do not get what I want. Unfortunately, some waters and lipids are appearing from the other side of the box. Cheers, Mohsen On Sun, May 21, 2017 at 8:12 PM, Dallas Warren wrote: > I have found the cluster option of -pbc to work well for putting > aggregates back together correctly. Some times you do need an index > file and appropriate groups to assist with it getting it right. > > gmx trjconv -pbc cluster > Catch ya, > > Dr. Dallas Warren > Drug Delivery, Disposition and Dynamics > Monash Institute of Pharmaceutical Sciences, Monash University > 381 Royal Parade, Parkville VIC 3052 > dallas.war...@monash.edu > - > When the only tool you own is a hammer, every problem begins to resemble a > nail. > > > On 16 May 2017 at 07:50, Mohsen Ramezanpour > wrote: > > Dear Gromacs users, > > > > I have an HII phase made of one inverted cylinder (and waters inside) in > a > > triclinic box with 90, 90, 60 angles. After running the simulation, this > > cylinder become bent like a curve. I.e. is not a perfect cylinder > anymore. > > As a result, some water molecules and lipids pass the box sides and enter > > from the other side of box because of PBC. > > Now, I want to make the cylinder again but I am not sure how to do so. > > > > The best I could do was to use "-pbc mol -ur compact" options in > trjconv. > > However, here are still some lipids and molecules which are not part of > the > > cylinder. > > > > Any idea how could I fix the effect of these PBC in visualization? > > > > Thanks > > Mohsen > > > > -- > > *Rewards work better than punishment ...* > > -- > > Gromacs Users mailing list > > > > * Please search the archive at http://www.gromacs.org/ > Support/Mailing_Lists/GMX-Users_List before posting! > > > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > > > * For (un)subscribe requests visit > > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/ > Support/Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] PBC fix for visualization
Dear Gromacs users, I have an HII phase made of one inverted cylinder (and waters inside) in a triclinic box with 90, 90, 60 angles. After running the simulation, this cylinder become bent like a curve. I.e. is not a perfect cylinder anymore. As a result, some water molecules and lipids pass the box sides and enter from the other side of box because of PBC. Now, I want to make the cylinder again but I am not sure how to do so. The best I could do was to use "-pbc mol -ur compact" options in trjconv. However, here are still some lipids and molecules which are not part of the cylinder. Any idea how could I fix the effect of these PBC in visualization? Thanks Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] GAAMP parameterization
Thanks Justin. I used CGenFF as initial parameters and the output file from GAAMP was not good. By "not good" I mean that: the QM PS was not matched with the fitted parameters. In fact, the relative depth of local minima and even the location of minima was not the same as the QM one. When I use GAFF as initial parameters, everything is fine and all the dihedrals have been fitted to the QM PS. So, I thought maybe I could use these charges and dihedral parameters instead. Cheers, Mohsen On Tue, May 2, 2017 at 9:50 AM, Justin Lemkul wrote: > > > On 5/2/17 11:46 AM, Mohsen Ramezanpour wrote: > >> Hi David, >> >> My bad. >> GAAMP server is a server (developed by Dr. Roux at Chicago university) >> which does an automated QM calculations for charge and dihedral fitting >> for >> both Amber and Charmm36 ff. >> Just to start, it takes the initial parameters from either GAFF or CGenFF >> and tries to make them optimized. >> Here is the link for others who might be interested: >> >> http://gaamp.lcrc.anl.gov/ >> >> I read the article but I still do not know the answer for the questions I >> asked. >> >> > GAFF parameters are to be used with AMBER, CGenFF with CHARMM. So if you > used GAFF to try to parametrize something ultimately for use with CHARMM, I > wouldn't use it. Do it properly with CGenFF. Note that the CGenFF server > is still operable and you should get a result almost immediately. That > doesn't mean it won't still need some work, but all the target data from > GAAMP (e.g. the QM) is still valid for refining the topology. > > -Justin > > > Cheers, >> Mohsen >> >> >> On Mon, May 1, 2017 at 11:06 PM, David van der Spoel < >> sp...@xray.bmc.uu.se> >> wrote: >> >> On 02/05/17 06:06, Mohsen Ramezanpour wrote: >>> >>> Hi Everyone, >>>> >>>> I would like to know your opinion on these questions. >>>> >>>> I highly appreciate your comments in advance. >>>> >>>> You could at least at the address of the server and tell in one line >>> what >>> it is supposed to do. >>> >>> You should never mix force fields unless you really really know what you >>> are doing. Not even change water models - most people trying that don't >>> really know what they are doing. >>> >>> >>>> Cheers, >>>> Mohsen >>>> >>>> On Thu, Apr 27, 2017 at 10:23 PM, Mohsen Ramezanpour < >>>> ramezanpour.moh...@gmail.com> wrote: >>>> >>>> Dear Gromacs users, >>>> >>>>> >>>>> Using GAAMP server for parameterization, if the GAFF for >>>>> "initial parameters" is used instead of "CGenFF": >>>>> >>>>> 1) Can the results still be used for using simulation in "Charmm36" >>>>> force >>>>> field? Assuming that all the other values are chosen as default ones in >>>>> the >>>>> GAAMP webpage. >>>>> >>>>> 2) Will the result be the same (or approximately the same) with when >>>>> "CGenFF" is used as initial parameters? >>>>> >>>>> 3) Will the units for angles, dihedrals, distances, ... be in CGenFF or >>>>> GAFF? I was wondering if I can still use the cgen_charmm2gmx.py script >>>>> for >>>>> conversions to use it in Gromacs package? >>>>> >>>>> As I understood from the paper, GAAMP can give output for >>>>> parameterization >>>>> in both Charmm and Amber force fields. >>>>> >>>>> This could be easily tested by running the same job two times. >>>>> Unfortunately, GAAMP server has some problems at the moment so it will >>>>> take >>>>> too much time to check this. >>>>> >>>>> I appreciate your comments in advance. >>>>> >>>>> Best, >>>>> Mohsen >>>>> >>>>> >>>>> -- >>>>> *Rewards work better than punishment ...* >>>>> >>>>> >>>>> >>>> >>>> >>>> >>> -- >>> David van der Spoel, Ph.D., Professor of Biology >>> Head of Department, Cell & Molecular Biology, Uppsala University. >>> Box 596, SE-75124 Uppsala, Sweden. Phone: +46184714205. >>> http://www.icm.uu.se >>> -- >>>
Re: [gmx-users] GAAMP parameterization
Hi David, My bad. GAAMP server is a server (developed by Dr. Roux at Chicago university) which does an automated QM calculations for charge and dihedral fitting for both Amber and Charmm36 ff. Just to start, it takes the initial parameters from either GAFF or CGenFF and tries to make them optimized. Here is the link for others who might be interested: http://gaamp.lcrc.anl.gov/ I read the article but I still do not know the answer for the questions I asked. Cheers, Mohsen On Mon, May 1, 2017 at 11:06 PM, David van der Spoel wrote: > On 02/05/17 06:06, Mohsen Ramezanpour wrote: > >> Hi Everyone, >> >> I would like to know your opinion on these questions. >> >> I highly appreciate your comments in advance. >> > You could at least at the address of the server and tell in one line what > it is supposed to do. > > You should never mix force fields unless you really really know what you > are doing. Not even change water models - most people trying that don't > really know what they are doing. > >> >> Cheers, >> Mohsen >> >> On Thu, Apr 27, 2017 at 10:23 PM, Mohsen Ramezanpour < >> ramezanpour.moh...@gmail.com> wrote: >> >> Dear Gromacs users, >>> >>> Using GAAMP server for parameterization, if the GAFF for >>> "initial parameters" is used instead of "CGenFF": >>> >>> 1) Can the results still be used for using simulation in "Charmm36" force >>> field? Assuming that all the other values are chosen as default ones in >>> the >>> GAAMP webpage. >>> >>> 2) Will the result be the same (or approximately the same) with when >>> "CGenFF" is used as initial parameters? >>> >>> 3) Will the units for angles, dihedrals, distances, ... be in CGenFF or >>> GAFF? I was wondering if I can still use the cgen_charmm2gmx.py script >>> for >>> conversions to use it in Gromacs package? >>> >>> As I understood from the paper, GAAMP can give output for >>> parameterization >>> in both Charmm and Amber force fields. >>> >>> This could be easily tested by running the same job two times. >>> Unfortunately, GAAMP server has some problems at the moment so it will >>> take >>> too much time to check this. >>> >>> I appreciate your comments in advance. >>> >>> Best, >>> Mohsen >>> >>> >>> -- >>> *Rewards work better than punishment ...* >>> >>> >> >> >> > > -- > David van der Spoel, Ph.D., Professor of Biology > Head of Department, Cell & Molecular Biology, Uppsala University. > Box 596, SE-75124 Uppsala, Sweden. Phone: +46184714205. > http://www.icm.uu.se > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] GAAMP parameterization
Hi Everyone, I would like to know your opinion on these questions. I highly appreciate your comments in advance. Cheers, Mohsen On Thu, Apr 27, 2017 at 10:23 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Dear Gromacs users, > > Using GAAMP server for parameterization, if the GAFF for > "initial parameters" is used instead of "CGenFF": > > 1) Can the results still be used for using simulation in "Charmm36" force > field? Assuming that all the other values are chosen as default ones in the > GAAMP webpage. > > 2) Will the result be the same (or approximately the same) with when > "CGenFF" is used as initial parameters? > > 3) Will the units for angles, dihedrals, distances, ... be in CGenFF or > GAFF? I was wondering if I can still use the cgen_charmm2gmx.py script for > conversions to use it in Gromacs package? > > As I understood from the paper, GAAMP can give output for parameterization > in both Charmm and Amber force fields. > > This could be easily tested by running the same job two times. > Unfortunately, GAAMP server has some problems at the moment so it will take > too much time to check this. > > I appreciate your comments in advance. > > Best, > Mohsen > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] GAAMP parameterization
Dear Gromacs users, Using GAAMP server for parameterization, if the GAFF for "initial parameters" is used instead of "CGenFF": 1) Can the results still be used for using simulation in "Charmm36" force field? Assuming that all the other values are chosen as default ones in the GAAMP webpage. 2) Will the result be the same (or approximately the same) with when "CGenFF" is used as initial parameters? 3) Will the units for angles, dihedrals, distances, ... be in CGenFF or GAFF? I was wondering if I can still use the cgen_charmm2gmx.py script for conversions to use it in Gromacs package? As I understood from the paper, GAAMP can give output for parameterization in both Charmm and Amber force fields. This could be easily tested by running the same job two times. Unfortunately, GAAMP server has some problems at the moment so it will take too much time to check this. I appreciate your comments in advance. Best, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] gmx hbond
Thanks Erik, It helped. Cheers On Sat, Apr 15, 2017 at 9:30 AM, Erik Marklund wrote: > Dear Moshen, > > I doubt the difference in versions will cause any problems in this case. > > The second column is ill-described. It contains the number of pairs within > 0.3 nm that don’t fulfil the angle criterion. > > Kind regards, > Erik > __ > Erik Marklund, PhD, Marie Skłodowska Curie INCA Fellow > Department of Chemistry – BMC, Uppsala University > +46 (0)18 471 4539 > erik.markl...@kemi.uu.se<mailto:erik.markl...@kemi.uu.se> > > On 14 Apr 2017, at 03:20, Mohsen Ramezanpour <mailto:ramezanpour.moh...@gmail.com>> wrote: > > Hi Gromacs users, > > I have a question regarding the output file from g_hbond: > > First:I have done my simulations with version 4.6.7 and doing the analysis > with version 2016.2 > Will this cause any hidden problem in analysis? the analysis is working > fine but I am asking about the correctness of results because of this > change. > > > Second and the main question: > > When I use this command: > gmx hbond -f md.xtc -n index.ndx -s file.tpr -num hbond.xvg -g test.log -r > 0.30 > > I get this: > > ... > @ s0 legend "Hydrogen bonds" > @ s1 legend "Pairs within 0.3 nm" >400 0 > * 40.1 1 0* > 40.2 1 3 > 40.3 2 2 > 40.4 0 2 > 40.5 1 1 > * 40.6 2 0* > ... > the first column is time frame > the second is number of h-bonds > the third is "Pairs within 0.3 nm" > > > by -r 0.30 I wanted to change the cut-off for h-bond definition. > i.e. report a h-bond if and only if the D...A distance is less or equal to > 0.30 ns, AND, the angle of H..D...A is less or equal than 35 (this angle is > default in 2016.2 version). > > Regarding the highlighted lines: > How it is possible that the there is not any pair in this distance but we > still have hbond? > If they are not in this cut-off, thy should not be recognized as hbond > either. > > I assume the pair means the D...A pair because I did NOT use -noda option. > > Also, it does not give any test.log file, and other problems :-) > > Please let me know your opinion about this results. > > Thanks in advance for your reply > Cheers > Mohsen > > > -- > *Rewards work better than punishment ...* > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/ > Support/Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org<mailto:gmx-users-request@ > gromacs.org>. > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/ > Support/Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] gmx hbond
Hi Gromacs users, I have a question regarding the output file from g_hbond: First:I have done my simulations with version 4.6.7 and doing the analysis with version 2016.2 Will this cause any hidden problem in analysis? the analysis is working fine but I am asking about the correctness of results because of this change. Second and the main question: When I use this command: gmx hbond -f md.xtc -n index.ndx -s file.tpr -num hbond.xvg -g test.log -r 0.30 I get this: ... @ s0 legend "Hydrogen bonds" @ s1 legend "Pairs within 0.3 nm" 400 0 * 40.1 1 0* 40.2 1 3 40.3 2 2 40.4 0 2 40.5 1 1 * 40.6 2 0* ... the first column is time frame the second is number of h-bonds the third is "Pairs within 0.3 nm" by -r 0.30 I wanted to change the cut-off for h-bond definition. i.e. report a h-bond if and only if the D...A distance is less or equal to 0.30 ns, AND, the angle of H..D...A is less or equal than 35 (this angle is default in 2016.2 version). Regarding the highlighted lines: How it is possible that the there is not any pair in this distance but we still have hbond? If they are not in this cut-off, thy should not be recognized as hbond either. I assume the pair means the D...A pair because I did NOT use -noda option. Also, it does not give any test.log file, and other problems :-) Please let me know your opinion about this results. Thanks in advance for your reply Cheers Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Problem with combining the .edr files
Thanks Mark :-) Problem solved. On Wed, Mar 15, 2017 at 1:10 PM, Mark Abraham wrote: > Hi, > > Judging from your start times, you're concatenating apples with oranges :-) > In this case, probably NVT equilibration with NPT production, or similar. > Arguably, concatenating the trr files in that case is also conceptually > wrong, even if the file format is not robust enough to hint that you're > trying to do something that will leave you with a file whose contents have > an interpretation that is somewhere between unclear and wrong. edr is > actually smart enough to perhaps stop you doing something that isn't > useful. > > Mark > > On Wed, Mar 15, 2017 at 7:26 PM Mohsen Ramezanpour < > ramezanpour.moh...@gmail.com> wrote: > > > Dear Gromacs users, > > > > I did a simulation in parts using -noappend command, so I have a few .trr > > and .edr files which I wish to merge. > > > > I have difficulty with that and I could not find any related post on that > > except this one: > > http://comments.gmane.org/gmane.science.biology.gromacs.user/68806 > > > > Using Gromacs 2016.2 installed on my PC (without MPI), and the following > > command: > > gmx trjcat -f *.trr -o fixed.trr > > I will have a reasonable .trr file. However, when I use: > > > > gmx eneconv -f *.edr -o fixed.edr > > > > I get the followings: > > > > Opened NJOB_1.edr as single precision energy file > > Last energy frame read 0 time0.000 > > Opened NJOB_2.edr as single precision energy file > > Energy files don't match, different number of energies: > > NJOB_1.edr: 55 > > NJOB_2.edr: 55 > > > > Continue conversion using only the first 55 terms (n/y)? > > (you should be sure that the energy terms match) > > y > > Reading energy frame 0 time0.000 > > Opened NJOB_3.edr as single precision energy file > > Energy files don't match, different number of energies: > > NJOB_2.edr: 55 > > NJOB_3.edr: 58 > > > > Continue conversion using only the first 55 terms (n/y)? > > (you should be sure that the energy terms match) > > y > > Reading energy frame 0 time0.000 > > Opened NJOB_4.edr as single precision energy file > > Energy files don't match, different number of energies: > > NJOB_3.edr: 58 > > NJOB_4.edr: 57 > > > > Continue conversion using only the first 55 terms (n/y)? > > (you should be sure that the energy terms match) > > y > > Reading energy frame 0 time 4.000 > > Opened NJOB_5.edr as single precision energy file > > Energy files don't match, different number of energies: > > NJOB_4.edr: 57 > > NJOB_5.edr: 60 > > > > Continue conversion using only the first 55 terms (n/y)? > > (you should be sure that the energy terms match) > > y > > Reading energy frame 0 time 8.000 > > Opened system_MD.5.edr.edr as single precision energy file > > Reading energy frame 0 time 12.000 > > Opened system_MD.5.part0003.edr.edr as single precision energy file > > Reading energy frame 0 time 32.000 > > Opened system_MD.5.part0004.edr.edr as single precision energy file > > Reading energy frame 0 time 57.000 > > Opened system_MD.5.part0005.edr.edr as single precision energy file > > Reading energy frame 0 time 82.000 > > > > Summary of files and start times used: > > > > FileStart time > > - > >NJOB_1.edr0.000 > >NJOB_2.edr0.000 > >NJOB_3.edr0.000 > >NJOB_4.edr4.000 > >NJOB_5.edr8.000 > > system_MD.5.edr.edr 12.000 > > system_MD.5.part0003.edr.edr 32.000 > > system_MD.5.part0004.edr.edr 57.000 > > system_MD.5.part0005.edr.edr 82.000 > > > > Opened NJOB_1.edr as single precision energy file > > Reading energy frame 0 time0.000 > > Continue writing frames from t=0, step=0 > > Last energy frame read 0 time0.000 > > Last step written from NJOB_1.edr: t 0, step 0 > > > > Opened NJOB_2.edr as single precision energy file > > Reading energy frame 1 time2.500 > > Last step written from NJOB_2.edr: t 0, step 0 > > > > Opened NJOB_3.edr as single precision energy file > > Reading energy frame 1 time 20.000 > > Continue writing frames from t=20, step=1 > > Last energ
[gmx-users] Problem with combining the .edr files
Dear Gromacs users, I did a simulation in parts using -noappend command, so I have a few .trr and .edr files which I wish to merge. I have difficulty with that and I could not find any related post on that except this one: http://comments.gmane.org/gmane.science.biology.gromacs.user/68806 Using Gromacs 2016.2 installed on my PC (without MPI), and the following command: gmx trjcat -f *.trr -o fixed.trr I will have a reasonable .trr file. However, when I use: gmx eneconv -f *.edr -o fixed.edr I get the followings: Opened NJOB_1.edr as single precision energy file Last energy frame read 0 time0.000 Opened NJOB_2.edr as single precision energy file Energy files don't match, different number of energies: NJOB_1.edr: 55 NJOB_2.edr: 55 Continue conversion using only the first 55 terms (n/y)? (you should be sure that the energy terms match) y Reading energy frame 0 time0.000 Opened NJOB_3.edr as single precision energy file Energy files don't match, different number of energies: NJOB_2.edr: 55 NJOB_3.edr: 58 Continue conversion using only the first 55 terms (n/y)? (you should be sure that the energy terms match) y Reading energy frame 0 time0.000 Opened NJOB_4.edr as single precision energy file Energy files don't match, different number of energies: NJOB_3.edr: 58 NJOB_4.edr: 57 Continue conversion using only the first 55 terms (n/y)? (you should be sure that the energy terms match) y Reading energy frame 0 time 4.000 Opened NJOB_5.edr as single precision energy file Energy files don't match, different number of energies: NJOB_4.edr: 57 NJOB_5.edr: 60 Continue conversion using only the first 55 terms (n/y)? (you should be sure that the energy terms match) y Reading energy frame 0 time 8.000 Opened system_MD.5.edr.edr as single precision energy file Reading energy frame 0 time 12.000 Opened system_MD.5.part0003.edr.edr as single precision energy file Reading energy frame 0 time 32.000 Opened system_MD.5.part0004.edr.edr as single precision energy file Reading energy frame 0 time 57.000 Opened system_MD.5.part0005.edr.edr as single precision energy file Reading energy frame 0 time 82.000 Summary of files and start times used: FileStart time - NJOB_1.edr0.000 NJOB_2.edr0.000 NJOB_3.edr0.000 NJOB_4.edr4.000 NJOB_5.edr8.000 system_MD.5.edr.edr 12.000 system_MD.5.part0003.edr.edr 32.000 system_MD.5.part0004.edr.edr 57.000 system_MD.5.part0005.edr.edr 82.000 Opened NJOB_1.edr as single precision energy file Reading energy frame 0 time0.000 Continue writing frames from t=0, step=0 Last energy frame read 0 time0.000 Last step written from NJOB_1.edr: t 0, step 0 Opened NJOB_2.edr as single precision energy file Reading energy frame 1 time2.500 Last step written from NJOB_2.edr: t 0, step 0 Opened NJOB_3.edr as single precision energy file Reading energy frame 1 time 20.000 Continue writing frames from t=20, step=1 Last energy frame read 2000 time 4.000 Writing frame time 4 Last step written from NJOB_3.edr: t 4, step 2000 Opened NJOB_4.edr as single precision energy file Reading energy frame 1 time 40020.000 Continue writing frames from t=40020, step=2001 Last energy frame read 2000 time 8.000 Writing frame time 8 Last step written from NJOB_4.edr: t 8, step 4000 Opened NJOB_5.edr as single precision energy file Reading energy frame 1 time 80020.000 Continue writing frames from t=80020, step=4001 Last energy frame read 2000 time 12.000 Writing frame time 12 Last step written from NJOB_5.edr: t 12, step 6000 Opened system_MD.5.edr.edr as single precision energy file Reading energy frame 1 time 120020.000 Continue writing frames from t=120020, step=6001 Last energy frame read 1 time 32.000 riting frame time 32 Last step written from system_MD.5.edr.edr: t 32, step 16000 Opened system_MD.5.part0003.edr.edr as single precision energy file Reading energy frame 1 time 320020.000 Continue writing frames from t=320020, step=16001 Last energy frame read 12500 time 57.000 riting frame time 56 Last step written from system_MD.5.part0003.edr.edr: t 57, step 28500 Opened system_MD.5.part0004.edr.edr as single precision energy file Reading energy frame 1 time 570020.000 Continue writing frames from t=570020, step=28501 Last energy frame read 12500 time 82.000 riting frame time 82 Last step written from system_MD.5.part0004.edr.edr: t 82, step 41000 Opened system_MD.5.part0005.edr.edr as single precision energy file Reading energy frame 1 time 820020.00
Re: [gmx-users] Gromac4.6.7 installation problem
Thanks Mark. I tried gcc 4.9 and used it as follows: cmake .. -DGMX_BUILD_OWN_FFTW=ON -DREGRESSIONTEST_DOWNLOAD=ON -DGMX_GPU=ON -DCMAKE_INSTALL_PREFIX=/home/mohsen/programs-mohsen -DCMAKE_CXX_LINK_FLAGS="-Wl,-rpath,/usr/lib/gcc/x86_64-linux-gnu/4.9 -L/usr/lib/gcc/x86_64-linux-gnu/4.9" I tried simple make for both -DGMX_GPU=ON and -DGMX_GPU=off for -DGMX_GPU=ON it failed and the same problem happened, although I used gcc 4.9 When I made -DGMX_GPU=off, it said: -- Configuring done -- Generating done CMake Warning: Manually-specified variables *were not used* by the project: *CMAKE_CXX_LINK_FLAGS* again, some of the previous errors happened. for instance: CMake Warning (dev) at cmake/gmxSetBuildInformation.cmake:156 (set): Policy CMP0053 is not set: Simplify variable reference and escape sequence evaluation. Run "cmake --help-policy CMP0053" for policy details. Use the cmake_policy command to set the policy and suppress this warning. For input: '@OUTPUT_CPU_FEATURES@' the old evaluation rules produce: 'aes apic avx avx2 clfsh cmov cx8 cx16 f16c fma htt lahf_lm mmx msr nonstop_tsc pcid pclmuldq pdcm pdpe1gb popcnt pse rdrnd rdtscp sse2 sse3 sse4.1 sse4.2 ssse3 tdt x2apic' but the new evaluation rules produce: '@OUTPUT_CPU_FEATURES@' Using the old result for compatibility since the policy is not set. The, I did a make check, and I got this. [ 1%] Built target fftwBuild [ 64%] Built target gmx [ 83%] Built target md [ 94%] Built target gmxana [ 94%] Built target editconf [ 98%] Built target gmxpreprocess [100%] Built target mdrun [100%] Built target pdb2gmx [100%] Built target gmxcheck [100%] Built target grompp Scanning dependencies of target gmxtests [100%] Built target gmxtests Scanning dependencies of target check Test project /home/mohsen/Downloads/gromacs-4.6.7/build Start 1: regressiontests/simple 1/6 Test #1: regressiontests/simple ... Passed1.17 sec Start 2: regressiontests/complex 2/6 Test #2: regressiontests/complex .. Passed2.84 sec Start 3: regressiontests/kernel 3/6 Test #3: regressiontests/kernel ... Passed 26.26 sec Start 4: regressiontests/freeenergy 4/6 Test #4: regressiontests/freeenergy ...***Failed4.83 sec Mdrun cannot use the requested (or automatic) number of cores, retrying with 8. Mdrun cannot use the requested (or automatic) number of cores, retrying with 8. Mdrun cannot use the requested (or automatic) number of cores, retrying with 8. FAILED. Check checkpot.out (18 errors) files in relative Mdrun cannot use the requested (or automatic) number of cores, retrying with 8. 1 out of 9 freeenergy tests FAILED Start 5: regressiontests/pdb2gmx 5/6 Test #5: regressiontests/pdb2gmx .. Passed 11.53 sec Start 6: regressiontests/rotation 6/6 Test #6: regressiontests/rotation . Passed1.87 sec 83% tests passed, 1 tests failed out of 6 Total Test time (real) = 48.50 sec The following tests FAILED: 4 - regressiontests/freeenergy (Failed) Errors while running CTest CMakeFiles/check.dir/build.make:57: recipe for target 'CMakeFiles/check' failed make[3]: *** [CMakeFiles/check] Error 8 CMakeFiles/Makefile2:131: recipe for target 'CMakeFiles/check.dir/all' failed make[2]: *** [CMakeFiles/check.dir/all] Error 2 CMakeFiles/Makefile2:138: recipe for target 'CMakeFiles/check.dir/rule' failed make[1]: *** [CMakeFiles/check.dir/rule] Error 2 Makefile:221: recipe for target 'check' failed make: *** [check] Error 2 Is this critical or I can continue installing? Thanks in advance for your reply and help Cheers Mohsen -- Build files have been written to: /home/mohsen/Downloads/gromacs-4.6.7/build The problem still exist. Cheers Mohsen On Fri, Mar 3, 2017 at 12:19 AM, Mark Abraham wrote: > Hi, > > Using parallel make obscures problems when they exist, so run again with > simple make. Probably cuda does not support your compiler, so you will run > into problems later. I suggest gcc 4.9 > > Mark > > On Thu, 2 Mar 2017 06:58 Mohsen Ramezanpour > wrote: > > > Dear Gromacs users, > > > > I am trying to install gromacs 4.6.7 on Ubuntu 16.04. > > Unfortunately, I get an Error which I do not know how to solve it. > > > > After > > mkdir build > > cd build > > > > *I use:* > > cmake .. -DGMX_BUILD_OWN_FFTW=ON -DREGRESSIONTEST_DOWNLOAD=ON > -DGMX_GPU=ON > > -DCMAKE_INSTALL_PREFIX=/home/mohsen/programs > > -DCMAKE_CXX_LINK_FLAGS="-Wl,-rpath,/usr/lib/gcc/x86_64-linux-gnu/5 > > -L/usr/lib/gcc/x86_64-linux-gnu/5" > > > > *I get this:* > > -- The C compiler identification is GNU 5.4.0 > > -- Check for working C compiler: /usr/bin/cc > > -- Check for working C compiler: /usr/bin/cc -- works > > -- Detecting
[gmx-users] Gromac4.6.7 installation problem
Dear Gromacs users,
I am trying to install gromacs 4.6.7 on Ubuntu 16.04.
Unfortunately, I get an Error which I do not know how to solve it.
After
mkdir build
cd build
*I use:*
cmake .. -DGMX_BUILD_OWN_FFTW=ON -DREGRESSIONTEST_DOWNLOAD=ON -DGMX_GPU=ON
-DCMAKE_INSTALL_PREFIX=/home/mohsen/programs
-DCMAKE_CXX_LINK_FLAGS="-Wl,-rpath,/usr/lib/gcc/x86_64-linux-gnu/5
-L/usr/lib/gcc/x86_64-linux-gnu/5"
*I get this:*
-- The C compiler identification is GNU 5.4.0
-- Check for working C compiler: /usr/bin/cc
-- Check for working C compiler: /usr/bin/cc -- works
-- Detecting C compiler ABI info
-- Detecting C compiler ABI info - done
-- Detecting C compile features
-- Detecting C compile features - done
-- Looking for NVIDIA GPUs present in the system
-- Number of NVIDIA GPUs detected: 1
-- Looking for pthread.h
-- Looking for pthread.h - found
-- Looking for pthread_create
-- Looking for pthread_create - not found
-- Looking for pthread_create in pthreads
-- Looking for pthread_create in pthreads - not found
-- Looking for pthread_create in pthread
-- Looking for pthread_create in pthread - found
-- Found Threads: TRUE
-- Found CUDA: /usr (found suitable version "7.5", minimum required is
"3.2")
-- The CXX compiler identification is GNU 5.4.0
-- Check for working CXX compiler: /usr/bin/c++
-- Check for working CXX compiler: /usr/bin/c++ -- works
-- Detecting CXX compiler ABI info
-- Detecting CXX compiler ABI info - done
-- Detecting CXX compile features
-- Detecting CXX compile features - done
-- Checking for GCC x86 inline asm
-- Checking for GCC x86 inline asm - supported
-- Detecting best acceleration for this CPU
CMake Warning (dev) at cmake/gmxDetectAcceleration.cmake:65 (set):
Policy CMP0053 is not set: Simplify variable reference and escape sequence
evaluation. Run "cmake --help-policy CMP0053" for policy details. Use
the
cmake_policy command to set the policy and suppress this warning.
For input:
'@GCC_INLINE_ASM_DEFINE@ -I${CMAKE_SOURCE_DIR}/include
-DGMX_CPUID_STANDALONE'
the old evaluation rules produce:
'-DGMX_X86_GCC_INLINE_ASM
-I/home/mohsen/Downloads/gromacs-4.6.7/include -DGMX_CPUID_STANDALONE'
but the new evaluation rules produce:
'@GCC_INLINE_ASM_DEFINE@ -I/home/mohsen/Downloads/gromacs-4.6.7/include
-DGMX_CPUID_STANDALONE'
Using the old result for compatibility since the policy is not set.
Call Stack (most recent call first):
cmake/gmxDetectAcceleration.cmake:98 (gmx_suggest_x86_acceleration)
CMakeLists.txt:184 (gmx_detect_acceleration)
This warning is for project developers. Use -Wno-dev to suppress it.
*and similar blocks like this again and again.*
*it ends up with:*.
.
.
-- [download 96% complete]
-- [download 97% complete]
-- [download 98% complete]
-- [download 99% complete]
-- [download 100% complete]
-- Configuring done
-- Generating done
-- Build files have been written to:
/home/mohsen/Downloads/gromacs-4.6.7/build
*So, I use:*
*make -j 16*
*and I get lots of lines which end with:*
.
.
.
Making install in neon
Making install in reodft
Making install in api
/usr/bin/install -c -m 644
/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/src/fftwBuild/api/fftw3.h
/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/src/fftwBuild/api/fftw3.f
/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/src/fftwBuild/api/fftw3l.f03
/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/src/fftwBuild/api/fftw3q.f03
'/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/include'
/usr/bin/install -c -m 644 fftw3.f03
'/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/include'
Making install in libbench2
Making install in .
/bin/bash ./libtool --mode=install /usr/bin/install -c libfftw3f.la
'/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/lib'
/usr/bin/install -c -m 644 fftw3f.pc
'/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/lib/pkgconfig'
libtool: install: /usr/bin/install -c .libs/libfftw3f.lai
/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/lib/
libfftw3f.la
libtool: install: /usr/bin/install -c .libs/libfftw3f.a
/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/lib/libfftw3f.a
libtool: install: chmod 644
/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/lib/libfftw3f.a
libtool: install: ranlib
/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/lib/libfftw3f.a
libtool: finish:
PATH="/home/mohsen/bin:/home/mohsen/.local/bin:/usr/local/sbin:/usr/local/bin:/usr/sbin:/usr/bin:/sbin:/bin:/usr/games:/usr/local/games:/snap/bin:/sbin"
ldconfig -n
/home/mohsen/Downloads/gromacs-4.6.7/build/src/contrib/fftw/fftwBuild-prefix/lib
--
Libraries have been installed in:
/home/m
Re: [gmx-users] Order Parameter for HII phase
Thanks Antonio. I will give them a try :-) Cheers Mohsen On Thu, Feb 23, 2017 at 4:15 PM, Antonio Baptista wrote: > Hi Mohsen, > > I suggest you have a look at the QRB 1977 review by Seelig, that Tom > already mentioned. Like for the planar case, they discuss how spectra of > cylindrical lipid phases are related to the order parameter tensor (but I > never looked into the details for that case). Anyway, I don't know if you > want to compare your simulations with experimental data or are just looking > for a convenient structural parameter. But, whatever the case is, reading > the literature and thinking about the geometry of your system (Duliez's > papers are a good train for that) should give you some hints for a relevant > parameter for your system. > > Once you have selected a parameter, you can compute it from stratch using > other GROMACS tools besides g_order. In particular, you can use g_gangle to > get several sorts of angles and then process them yourself. As an example, > you can use this approach to compute SCD using the equation already > mentioned by Tom, and then check if it agrees with the SCD given by g_order > (I once did that, just to be sure). > > Good luck! :) > > Best, > Antonio > > > > On Thu, 23 Feb 2017, Justin Lemkul wrote: > > >> >> On 2/23/17 1:25 PM, Mohsen Ramezanpour wrote: >> >>> And I agree with Piggot. The paper by Chau is about on option in g_order. >>> >>> >> Yes, and if memory serves (been a while since I used the program, so >> perhaps I am confusing something), this is what -o provides you so I >> thought it was relevant. -od gives you the deuterium order parameters. >> The documentation for gmx order is indeed very sparse. >> >> Anyway, a general question: >>> Can we expect to find a published article for each/some module(s) in >>> Gromacs? >>> >> >> No. Many of the tools are just implementations of simple algorithms used >> widely in simulations. Some do have specific references and those are >> generally printed to the screen output in bold blocks of text. >> >> With regards to deuterium order parameters, there is a well accepted and >> ubiquitously used equation, which Tom posted. >> >> I mean how/where can we figure out the underlying algorithm and >>> comprehensive description of each analysis tool? >>> >>> >> That's why GROMACS is open source :) Efforts are certainly always made >> to provide references when possible. >> >> -Justin >> >> Cheers >>> Mohsen >>> >>> >>> On Thu, Feb 23, 2017 at 10:06 AM, Mohsen Ramezanpour < >>> ramezanpour.moh...@gmail.com> wrote: >>> >>> Hi Justin, Piggot, >>>> >>>> Thanks for your replies. >>>> I agree with that. The problem is that the situation is straightforward >>>> for bilayers as bilayers are usually in specific angles regarding the >>>> magnetic field (e.g. they are perpendicular, as it is the case in >>>> Vermeer study.) >>>> For example the normal to the bilayer is z and we do not need be worry >>>> about anything else. >>>> However, in the case of HII phase, there are cylinders which are not >>>> perfectly cylinder. >>>> Even if they are perfect cylinders, I think I should use "membrane >>>> radial >>>> axis" to mimic the local normal vector to the lipid surface. >>>> >>>> Cheers >>>> Mohsen >>>> >>>> >>>> On Thu, Feb 23, 2017 at 8:08 AM, Piggot T. >>>> wrote: >>>> >>>> Hi, >>>>> >>>>> I don't believe this Chau paper is for deuterium order parameters (as I >>>>> think you are most likely referring to) but is related to some of the >>>>> other >>>>> order parameter options in gmx order. The Vermeer paper you mention >>>>> gives a >>>>> good overview of deuterium order parameters and you can find more >>>>> details >>>>> in the papers referenced therein (e.g. the Douliez et al. and the >>>>> Seelig et >>>>> al. papers in particular). >>>>> >>>>> For united-atom systems (as gmx order assumes), I believe the >>>>> calculation >>>>> is performed in gmx order using the "SCD = (2/3) Sxx + (1/3) Syy" >>>>> equation >>>>> (which is nicely derived in the 1998 Douliez paper: >>>&
Re: [gmx-users] Order Parameter for HII phase
And I agree with Piggot. The paper by Chau is about on option in g_order. Anyway, a general question: Can we expect to find a published article for each/some module(s) in Gromacs? I mean how/where can we figure out the underlying algorithm and comprehensive description of each analysis tool? Cheers Mohsen On Thu, Feb 23, 2017 at 10:06 AM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Justin, Piggot, > > Thanks for your replies. > I agree with that. The problem is that the situation is straightforward > for bilayers as bilayers are usually in specific angles regarding the > magnetic field (e.g. they are perpendicular, as it is the case in > Vermeer study.) > For example the normal to the bilayer is z and we do not need be worry > about anything else. > However, in the case of HII phase, there are cylinders which are not > perfectly cylinder. > Even if they are perfect cylinders, I think I should use "membrane radial > axis" to mimic the local normal vector to the lipid surface. > > Cheers > Mohsen > > > On Thu, Feb 23, 2017 at 8:08 AM, Piggot T. wrote: > >> Hi, >> >> I don't believe this Chau paper is for deuterium order parameters (as I >> think you are most likely referring to) but is related to some of the other >> order parameter options in gmx order. The Vermeer paper you mention gives a >> good overview of deuterium order parameters and you can find more details >> in the papers referenced therein (e.g. the Douliez et al. and the Seelig et >> al. papers in particular). >> >> For united-atom systems (as gmx order assumes), I believe the calculation >> is performed in gmx order using the "SCD = (2/3) Sxx + (1/3) Syy" equation >> (which is nicely derived in the 1998 Douliez paper: >> http://aip.scitation.org/doi/pdf/10.1063/1.476823) . The calculation >> requires an appropriate definition of the molecular frame of the system, >> where the z-axis is the bilayer normal. >> >> Note that if you are calculating order parameters for unsaturated carbons >> (e.g. as in carbon-carbon double bonds), gmx order doesn't currently do >> this correctly as it requires a different calculation with a different >> molecular frame. >> >> Cheers >> >> Tom >> >> >> From: gromacs.org_gmx-users-boun...@maillist.sys.kth.se [ >> gromacs.org_gmx-users-boun...@maillist.sys.kth.se] on behalf of Justin >> Lemkul [jalem...@vt.edu] >> Sent: 23 February 2017 12:16 >> To: gmx-us...@gromacs.org >> Subject: Re: [gmx-users] Order Parameter for HII phase >> >> On 2/22/17 9:51 PM, Mohsen Ramezanpour wrote: >> > Dear Gromacs users, >> > >> > Unfortunately, I did not get any reply on this post. >> > >> > I was wondering if anyone is aware of the original paper on g_order so >> that >> > I can understand the underlying algorithm of the options? >> > I am not sure how to use it for HII phase order parameters. >> > >> > I found this article but I am not sure if these are the ones, as they do >> > not describe g_order. >> > >> > "A new order parameter for tetrahedral configurations" by Chau and >> Hardwick >> > (1998) >> >> This is the reference given at the end of the gmx order help description. >> >> -Justin >> >> > and >> > "Acyl chain order parameter profiles in phospholipid bilayers: >> computation >> > from molecular dynamics simulations and comparison with 2 H NMR >> experiments" >> > by Vermeer et al. >> > >> > Thanks in advance for your comments >> > Mohsen >> > >> > >> > >> > On Sun, Feb 19, 2017 at 10:18 AM, Mohsen Ramezanpour < >> > ramezanpour.moh...@gmail.com> wrote: >> > >> >> Hi Everyone, >> >> >> >> I was wondering if someone has any suggestion on the following post, >> >> please? >> >> >> >> Cheers >> >> Mohsen >> >> >> >> -- Forwarded message -- >> >> From: Mohsen Ramezanpour >> >> Date: Thu, Feb 16, 2017 at 3:13 PM >> >> Subject: Order Parameter for HII phase >> >> To: Discussion list for GROMACS users >> >> >> >> >> >> Dear Gromacs users, >> >> >> >> I am interested in calculation of order parameter of lipids in a HII >> phase. >> >> >> >> I am not sure about the correct command to use. Which one o
Re: [gmx-users] Order Parameter for HII phase
Hi Justin, Piggot, Thanks for your replies. I agree with that. The problem is that the situation is straightforward for bilayers as bilayers are usually in specific angles regarding the magnetic field (e.g. they are perpendicular, as it is the case in Vermeer study.) For example the normal to the bilayer is z and we do not need be worry about anything else. However, in the case of HII phase, there are cylinders which are not perfectly cylinder. Even if they are perfect cylinders, I think I should use "membrane radial axis" to mimic the local normal vector to the lipid surface. Cheers Mohsen On Thu, Feb 23, 2017 at 8:08 AM, Piggot T. wrote: > Hi, > > I don't believe this Chau paper is for deuterium order parameters (as I > think you are most likely referring to) but is related to some of the other > order parameter options in gmx order. The Vermeer paper you mention gives a > good overview of deuterium order parameters and you can find more details > in the papers referenced therein (e.g. the Douliez et al. and the Seelig et > al. papers in particular). > > For united-atom systems (as gmx order assumes), I believe the calculation > is performed in gmx order using the "SCD = (2/3) Sxx + (1/3) Syy" equation > (which is nicely derived in the 1998 Douliez paper: > http://aip.scitation.org/doi/pdf/10.1063/1.476823) . The calculation > requires an appropriate definition of the molecular frame of the system, > where the z-axis is the bilayer normal. > > Note that if you are calculating order parameters for unsaturated carbons > (e.g. as in carbon-carbon double bonds), gmx order doesn't currently do > this correctly as it requires a different calculation with a different > molecular frame. > > Cheers > > Tom > > > From: gromacs.org_gmx-users-boun...@maillist.sys.kth.se [ > gromacs.org_gmx-users-boun...@maillist.sys.kth.se] on behalf of Justin > Lemkul [jalem...@vt.edu] > Sent: 23 February 2017 12:16 > To: gmx-us...@gromacs.org > Subject: Re: [gmx-users] Order Parameter for HII phase > > On 2/22/17 9:51 PM, Mohsen Ramezanpour wrote: > > Dear Gromacs users, > > > > Unfortunately, I did not get any reply on this post. > > > > I was wondering if anyone is aware of the original paper on g_order so > that > > I can understand the underlying algorithm of the options? > > I am not sure how to use it for HII phase order parameters. > > > > I found this article but I am not sure if these are the ones, as they do > > not describe g_order. > > > > "A new order parameter for tetrahedral configurations" by Chau and > Hardwick > > (1998) > > This is the reference given at the end of the gmx order help description. > > -Justin > > > and > > "Acyl chain order parameter profiles in phospholipid bilayers: > computation > > from molecular dynamics simulations and comparison with 2 H NMR > experiments" > > by Vermeer et al. > > > > Thanks in advance for your comments > > Mohsen > > > > > > > > On Sun, Feb 19, 2017 at 10:18 AM, Mohsen Ramezanpour < > > ramezanpour.moh...@gmail.com> wrote: > > > >> Hi Everyone, > >> > >> I was wondering if someone has any suggestion on the following post, > >> please? > >> > >> Cheers > >> Mohsen > >> > >> -- Forwarded message -- > >> From: Mohsen Ramezanpour > >> Date: Thu, Feb 16, 2017 at 3:13 PM > >> Subject: Order Parameter for HII phase > >> To: Discussion list for GROMACS users > >> > >> > >> Dear Gromacs users, > >> > >> I am interested in calculation of order parameter of lipids in a HII > phase. > >> > >> I am not sure about the correct command to use. Which one of the > following > >> should I use to be able to compare my data with deutirated NMR > experiments? > >> > >> 1) gmx order -f md.xtc -n index.ndx -s md.tpr -o order.xvg -od > >> deuter.xvg -d z -szonly > >> > >> or > >> > >> 2) gmx order -f md.xtc -n index.ndx -s md.tpr -o order.xvg -od > >> deuter.xvg -d z > >> (In this case, there are four columns of data (x, y, z for each carbon > >> atom). > >> The column 3 seems to be the same with when we include -szonly) > >> > >> or anything else? > >> > >> Corresponding available experimental data, I assume the second command > >> should be correct in my case . However, I am not sure x or y or what to > use > >> for
Re: [gmx-users] Order Parameter for HII phase
Dear Gromacs users, Unfortunately, I did not get any reply on this post. I was wondering if anyone is aware of the original paper on g_order so that I can understand the underlying algorithm of the options? I am not sure how to use it for HII phase order parameters. I found this article but I am not sure if these are the ones, as they do not describe g_order. "A new order parameter for tetrahedral configurations" by Chau and Hardwick (1998) and "Acyl chain order parameter profiles in phospholipid bilayers: computation from molecular dynamics simulations and comparison with 2 H NMR experiments" by Vermeer et al. (2007) Thanks in advance for your comments Mohsen On Sun, Feb 19, 2017 at 10:18 AM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Everyone, > > I was wondering if someone has any suggestion on the following post, > please? > > Cheers > Mohsen > > -- Forwarded message -- > From: Mohsen Ramezanpour > Date: Thu, Feb 16, 2017 at 3:13 PM > Subject: Order Parameter for HII phase > To: Discussion list for GROMACS users > > > Dear Gromacs users, > > I am interested in calculation of order parameter of lipids in a HII phase. > > I am not sure about the correct command to use. Which one of the following > should I use to be able to compare my data with deutirated NMR experiments? > > 1) gmx order -f md.xtc -n index.ndx -s md.tpr -o order.xvg -od > deuter.xvg -d z -szonly > > or > > 2) gmx order -f md.xtc -n index.ndx -s md.tpr -o order.xvg -od > deuter.xvg -d z > (In this case, there are four columns of data (x, y, z for each carbon > atom). > The column 3 seems to be the same with when we include -szonly) > > or anything else? > > Corresponding available experimental data, I assume the second command > should be correct in my case . However, I am not sure x or y or what to use > for judging. > > Cylinders in HII phase are parallel to the Z axis in my system. > > Many thanks in advance for your comments > Mohsen > -- > *Rewards work better than punishment ...* > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Fwd: Order Parameter for HII phase
Hi Everyone, I was wondering if someone has any suggestion on the following post, please? Cheers Mohsen -- Forwarded message -- From: Mohsen Ramezanpour Date: Thu, Feb 16, 2017 at 3:13 PM Subject: Order Parameter for HII phase To: Discussion list for GROMACS users Dear Gromacs users, I am interested in calculation of order parameter of lipids in a HII phase. I am not sure about the correct command to use. Which one of the following should I use to be able to compare my data with deutirated NMR experiments? 1) gmx order -f md.xtc -n index.ndx -s md.tpr -o order.xvg -od deuter.xvg -d z -szonly or 2) gmx order -f md.xtc -n index.ndx -s md.tpr -o order.xvg -od deuter.xvg -d z (In this case, there are four columns of data (x, y, z for each carbon atom). The column 3 seems to be the same with when we include -szonly) or anything else? Corresponding available experimental data, I assume the second command should be correct in my case . However, I am not sure x or y or what to use for judging. Cylinders in HII phase are parallel to the Z axis in my system. Many thanks in advance for your comments Mohsen -- *Rewards work better than punishment ...* -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Order Parameter for HII phase
Dear Gromacs users, I am interested in calculation of order parameter of lipids in a HII phase. I am not sure about the correct command to use. Which one of the following should I use to be able to compare my data with deutirated NMR experiments? 1) gmx order -f md.xtc -n index.ndx -s md.tpr -o order.xvg -od deuter.xvg -d z -szonly or 2) gmx order -f md.xtc -n index.ndx -s md.tpr -o order.xvg -od deuter.xvg -d z (In this case, there are four columns of data (x, y, z for each carbon atom). The column 3 seems to be the same with when we include -szonly) or anything else? Corresponding available experimental data, I assume the second command should be correct in my case . However, I am not sure x or y or what to use for judging. Cylinders in HII phase are parallel to the Z axis in my system. Many thanks in advance for your comments Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Question regarding equilibration of the system
Hi Ali, There are several parameters to check. It is, however, highly dependent to the system. Please describe your system, is it a protein-ligand system, a lipid bilayer, ... ? For lipid bilayers, the box size, and consequently, the are per lipid is a good one to check for. Cheers On Thu, Feb 9, 2017 at 4:12 PM, wrote: > Dear Gromacs users, > > Is there any criteria I can guarantee that system has reached its > equilibrium state? > > Eagerly looking forward to hear back from you. > > Sincerely, > > Ali > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/ > Support/Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] dihedral MM scan
Regarding the CGenFF article by Vanommeslaeghe et al. (2009): It is mentioned that: "... are correlated, as are their parameters. In such cases, the individual dihedrals are still scanned at the QM level with all other degrees of freedom, including other dihedrals that have to be fit allowed to relax. However, when scanning the corresponding MM surfaces for a given dihedral *it is necessary* that the remaining dihedrals associated with parameters that are being fit *be restrained in the QM optimized structures from the QM PES*. The use of these additional restraints assures that the energy contributions from those dihedrals are taken into account when performing the fitting" (Vanommeslaeghe et al. (2009)) When I did the MM scanning, I ignored the other associated dihedrals (some of them were highly correlated) which are in their local minima in QM scanning. However, I did an EM step with a tolerance of 50 and I assumed that those associated dihedrals will go into their local minima which should be similar with QM case. If I understood correctly from this article, I should use exactly the QM optimized structure for each specific angle, and apply restraints to all of them (except the one that I am trying to scan) during zero-step MD energy calculation. That being said, taking the structure from that and restraining, I do not need the EM step by MM force field any more. Is this right? Thanks in advance for your opinion Cheers Mohsen On Wed, Jan 25, 2017 at 10:14 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Justin, > > I did a dihedral population analysis in both bilayer and HII phase. > Interestingly, the population analysis is consistent with the QM profile > from GAAMP. (i.e. I see more population in deeper local minima in QM > profile, and the population is 0 at QM barriers) > There is only one exception, though. and it is when the Nitrogen atom is > protonated (there is a (CH2)-(CH2)-(N)-2(CH3) group at the top part of > lipid head group), > If I do a dihedral population analysis in vacuum, it matches the QM > profile pretty reasonable. > However, if I use this cationic lipid in a system (either bilayer or HII > phase), the population profile changes significantly (I mean the shape > change completely). > This one is not consistent with QM profile. I think this might be because > of strong interaction with PS lipids which are negatively charged. > > So, I have made a mistake in my MM scan. I will check it again to see what > is the problem. Based on these analysis I expect to get quite good MM PES > profiles. > Cheers > Mohsen > > On Tue, Jan 24, 2017 at 1:48 PM, Mohsen Ramezanpour < > ramezanpour.moh...@gmail.com> wrote: > >> Hi Justin, >> >> Sorry, I am just doing some dihedral population analysis to see which >> part is wrong. >> I agree. Most likely I am making some mistakes which is hidden to my >> eyes. I will prepare and share those with you with the result of population >> analysis when prepared. >> >> Cheers >> Mohsen >> >> On Thu, Jan 19, 2017 at 7:48 PM, Justin Lemkul wrote: >> >>> >>> >>> On 1/19/17 9:40 PM, Mohsen Ramezanpour wrote: >>> >>>> Thanks Justin, >>>> >>>> I did the scan for the whole range from -180 to 180 in 10 intervals. >>>> I used emtol = 50, and restrain force of 10 to have a better >>>> relaxation >>>> and restraining, respectively. >>>> Unfortunately, my MM profiles do not look similar to GAAMP QM or MM >>>> profiles, although I have used the suggested parameters. >>>> >>>> To check the method I am using and the parameters in the .mdp files, I >>>> did >>>> the MM scanning on a rotatable bond (CH2-CH2-CH2-CH3) similar to one >>>> can be >>>> found in a lipid tail. >>>> This dihedral was detected as a soft dihedral by GAAMP and there was a >>>> good >>>> agreement between QM and MM profiles. >>>> Here again, my MM profile from scanning does not seem reasonable. >>>> Most importantly, there is a very large barrier at angle 0 (ca. 300 >>>> times >>>> of energies in other angles). Besides, there are sudden increase or >>>> decrease in energy just by 10 degree of rotation. >>>> >>>> Any suggestion is appreciated in advance. >>>> >>> >>> There's really nothing I can suggest without actually seeing the files >>> (plots of the energy, the QM inputs and your topology). Something is >>> inconsistent but it's not clear what. >>> >>> >>> -Justin >>> >&g
Re: [gmx-users] dihedral MM scan
Hi Justin, I did a dihedral population analysis in both bilayer and HII phase. Interestingly, the population analysis is consistent with the QM profile from GAAMP. (i.e. I see more population in deeper local minima in QM profile, and the population is 0 at QM barriers) There is only one exception, though. and it is when the Nitrogen atom is protonated (there is a (CH2)-(CH2)-(N)-2(CH3) group at the top part of lipid head group), If I do a dihedral population analysis in vacuum, it matches the QM profile pretty reasonable. However, if I use this cationic lipid in a system (either bilayer or HII phase), the population profile changes significantly (I mean the shape change completely). This one is not consistent with QM profile. I think this might be because of strong interaction with PS lipids which are negatively charged. So, I have made a mistake in my MM scan. I will check it again to see what is the problem. Based on these analysis I expect to get quite good MM PES profiles. Cheers Mohsen On Tue, Jan 24, 2017 at 1:48 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Justin, > > Sorry, I am just doing some dihedral population analysis to see which part > is wrong. > I agree. Most likely I am making some mistakes which is hidden to my eyes. > I will prepare and share those with you with the result of population > analysis when prepared. > > Cheers > Mohsen > > On Thu, Jan 19, 2017 at 7:48 PM, Justin Lemkul wrote: > >> >> >> On 1/19/17 9:40 PM, Mohsen Ramezanpour wrote: >> >>> Thanks Justin, >>> >>> I did the scan for the whole range from -180 to 180 in 10 intervals. >>> I used emtol = 50, and restrain force of 10 to have a better >>> relaxation >>> and restraining, respectively. >>> Unfortunately, my MM profiles do not look similar to GAAMP QM or MM >>> profiles, although I have used the suggested parameters. >>> >>> To check the method I am using and the parameters in the .mdp files, I >>> did >>> the MM scanning on a rotatable bond (CH2-CH2-CH2-CH3) similar to one can >>> be >>> found in a lipid tail. >>> This dihedral was detected as a soft dihedral by GAAMP and there was a >>> good >>> agreement between QM and MM profiles. >>> Here again, my MM profile from scanning does not seem reasonable. >>> Most importantly, there is a very large barrier at angle 0 (ca. 300 times >>> of energies in other angles). Besides, there are sudden increase or >>> decrease in energy just by 10 degree of rotation. >>> >>> Any suggestion is appreciated in advance. >>> >> >> There's really nothing I can suggest without actually seeing the files >> (plots of the energy, the QM inputs and your topology). Something is >> inconsistent but it's not clear what. >> >> >> -Justin >> >> -- >> == >> >> Justin A. Lemkul, Ph.D. >> Ruth L. Kirschstein NRSA Postdoctoral Fellow >> >> Department of Pharmaceutical Sciences >> School of Pharmacy >> Health Sciences Facility II, Room 629 >> University of Maryland, Baltimore >> 20 Penn St. >> Baltimore, MD 21201 >> >> jalem...@outerbanks.umaryland.edu | (410) 706-7441 >> http://mackerell.umaryland.edu/~jalemkul >> >> == >> -- >> Gromacs Users mailing list >> >> * Please search the archive at http://www.gromacs.org/Support >> /Mailing_Lists/GMX-Users_List before posting! >> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >> >> * For (un)subscribe requests visit >> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >> send a mail to gmx-users-requ...@gromacs.org. >> > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Dihedral Population Analysis -Index file
Thanks Mark. It works now. I did the same for other dihedrals and everything looks fine. Cheers Mohsen On Wed, Jan 25, 2017 at 3:10 AM, Mark Abraham wrote: > Hi, > > Permutations are tricky to describe, but the code looks like it works the > way > http://manual.gromacs.org/documentation/5.1/onlinehelp/ > selections.html#permuting-selections-permute > describes. You want 4 2 1 3 > > Mark > > On Tue, 24 Jan 2017 21:46 Mohsen Ramezanpour > > wrote: > > > Thanks Mark. That is great. > > > > I tried that and I made an Index file of what I want. > > > > "DIH" resname DOPC and atomname O C12 C23 C > > > > > > However, looking at index files, the order of atom numbers is not right > and > > give me a wrong dihedral for analysis with g_angle. > > > > 8101 8116 8117 8126 > > > > > > *I want it to be like:* > > *8117 8116 8126 8101* (based on atom numbers in the .gro file) > > because the order matters for the dihedrals. > > > > So, I used the following one but it did not work neither: > > "DIH" resname DOPC and atomname O C12 C23 C permute 1 2 3 4 > > To keep them in the order they are (I assumed 1, 2, 3, 4 represents the > > position of atom names in command), and the out put is again: > > > > 8101 8116 8117 8126 which is wrong again. > > > > If I assume that the numbers are the positions for atom numbers (not the > > atomnames), then 1 would be for 8101, 2 for 8116, 3 for 8117, and 4 for > > 8126 > > "DIH" resname DOPC and atomname O C12 C23 C permute 3 2 4 1 > > > > and the output will be: > > 8126 8116 8101 8117 > > > > Which is not what I want. > > Could you please clarify how to make the order as I want? > > I read but could not figure it out. > > > > Thanks. > > > > > > On Tue, Jan 24, 2017 at 3:56 AM, Mark Abraham > > wrote: > > > > > Hi, > > > > > > On Tue, Jan 24, 2017 at 12:35 AM Mohsen Ramezanpour < > > > ramezanpour.moh...@gmail.com> wrote: > > > > > > > Dear Gromacs users, > > > > > > > > I have a trajectory file which I wish to do a dihedral analysis for > > > > specific dihedral. > > > > Based on what I understood from mailing list and manual, I made an > > > > index.ndx manually for the dihedral of interest in the molecule. > > > > > > > > Usually, when we make an index file, say for atom P in the system, > > > make_ndx > > > > makes an index file which includes the numbers for all the P atom in > > the > > > > system. These numbers are the index numbers in the .gro file for the > > > whole > > > > system. > > > > > > > > In my index file, I made it based on the atom numbers in the topology > > > file > > > > for the lipid. > > > > > > > > > > That forces you to use a trajectory file that contains only a single > > lipid. > > > With suitable index groups (for which I recommend gmx select), there's > no > > > need to do that. > > > > > > > > > > For doing analysis, I first used trjconv and saved only the .xtc file > > for > > > > only the lipid of interest, and then used g_angle to calculate the > > > dihedral > > > > population. I think what I get is only for ONE lipid over time NOT an > > > > average over time and lipids in the system. > > > > > > > > > > Perforce. You only gave gmx angle one lipid and used an index group > with > > > one dihedral. Whether you got a time series, an average or a > distribution > > > depends on what you did with gmx angle. > > > > > > > > > > The obvious solution would be to include all the dihedrals for all > the > > > > lipids in it. > > > > BUT, this is not easy to do manually. I tried mk_angndx but it is not > > > > giving me what I want. > > > > Is there any way to do this? > > > > > > > > > > Sure, see > > > http://manual.gromacs.org/documentation/5.1/onlinehelp/selections.html > > > > > > Something along the lines of > > > > > > resname POPC and name C2 C3 C4 C5 > > > > > > will produce a group that has that dihedral from every POPC. > > > > > > Mark > > > > > > > > > > Thanks in advance for any comment and suggestion > > > > > >
Re: [gmx-users] dihedral MM scan
Hi Justin, Sorry, I am just doing some dihedral population analysis to see which part is wrong. I agree. Most likely I am making some mistakes which is hidden to my eyes. I will prepare and share those with you with the result of population analysis when prepared. Cheers Mohsen On Thu, Jan 19, 2017 at 7:48 PM, Justin Lemkul wrote: > > > On 1/19/17 9:40 PM, Mohsen Ramezanpour wrote: > >> Thanks Justin, >> >> I did the scan for the whole range from -180 to 180 in 10 intervals. >> I used emtol = 50, and restrain force of 10 to have a better >> relaxation >> and restraining, respectively. >> Unfortunately, my MM profiles do not look similar to GAAMP QM or MM >> profiles, although I have used the suggested parameters. >> >> To check the method I am using and the parameters in the .mdp files, I did >> the MM scanning on a rotatable bond (CH2-CH2-CH2-CH3) similar to one can >> be >> found in a lipid tail. >> This dihedral was detected as a soft dihedral by GAAMP and there was a >> good >> agreement between QM and MM profiles. >> Here again, my MM profile from scanning does not seem reasonable. >> Most importantly, there is a very large barrier at angle 0 (ca. 300 times >> of energies in other angles). Besides, there are sudden increase or >> decrease in energy just by 10 degree of rotation. >> >> Any suggestion is appreciated in advance. >> > > There's really nothing I can suggest without actually seeing the files > (plots of the energy, the QM inputs and your topology). Something is > inconsistent but it's not clear what. > > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Dihedral Population Analysis -Index file
Thanks Mark. That is great. I tried that and I made an Index file of what I want. "DIH" resname DOPC and atomname O C12 C23 C However, looking at index files, the order of atom numbers is not right and give me a wrong dihedral for analysis with g_angle. 8101 8116 8117 8126 *I want it to be like:* *8117 8116 8126 8101* (based on atom numbers in the .gro file) because the order matters for the dihedrals. So, I used the following one but it did not work neither: "DIH" resname DOPC and atomname O C12 C23 C permute 1 2 3 4 To keep them in the order they are (I assumed 1, 2, 3, 4 represents the position of atom names in command), and the out put is again: 8101 8116 8117 8126 which is wrong again. If I assume that the numbers are the positions for atom numbers (not the atomnames), then 1 would be for 8101, 2 for 8116, 3 for 8117, and 4 for 8126 "DIH" resname DOPC and atomname O C12 C23 C permute 3 2 4 1 and the output will be: 8126 8116 8101 8117 Which is not what I want. Could you please clarify how to make the order as I want? I read but could not figure it out. Thanks. On Tue, Jan 24, 2017 at 3:56 AM, Mark Abraham wrote: > Hi, > > On Tue, Jan 24, 2017 at 12:35 AM Mohsen Ramezanpour < > ramezanpour.moh...@gmail.com> wrote: > > > Dear Gromacs users, > > > > I have a trajectory file which I wish to do a dihedral analysis for > > specific dihedral. > > Based on what I understood from mailing list and manual, I made an > > index.ndx manually for the dihedral of interest in the molecule. > > > > Usually, when we make an index file, say for atom P in the system, > make_ndx > > makes an index file which includes the numbers for all the P atom in the > > system. These numbers are the index numbers in the .gro file for the > whole > > system. > > > > In my index file, I made it based on the atom numbers in the topology > file > > for the lipid. > > > > That forces you to use a trajectory file that contains only a single lipid. > With suitable index groups (for which I recommend gmx select), there's no > need to do that. > > > > For doing analysis, I first used trjconv and saved only the .xtc file for > > only the lipid of interest, and then used g_angle to calculate the > dihedral > > population. I think what I get is only for ONE lipid over time NOT an > > average over time and lipids in the system. > > > > Perforce. You only gave gmx angle one lipid and used an index group with > one dihedral. Whether you got a time series, an average or a distribution > depends on what you did with gmx angle. > > > > The obvious solution would be to include all the dihedrals for all the > > lipids in it. > > BUT, this is not easy to do manually. I tried mk_angndx but it is not > > giving me what I want. > > Is there any way to do this? > > > > Sure, see > http://manual.gromacs.org/documentation/5.1/onlinehelp/selections.html > > Something along the lines of > > resname POPC and name C2 C3 C4 C5 > > will produce a group that has that dihedral from every POPC. > > Mark > > > > Thanks in advance for any comment and suggestion > > > > Cheers > > Mohsen > > > > -- > > *Rewards work better than punishment ...* > > -- > > Gromacs Users mailing list > > > > * Please search the archive at > > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before > > posting! > > > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > > > * For (un)subscribe requests visit > > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > > send a mail to gmx-users-requ...@gromacs.org. > > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/ > Support/Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Dihedral Population Analysis -Index file
Dear Gromacs users, I have a trajectory file which I wish to do a dihedral analysis for specific dihedral. Based on what I understood from mailing list and manual, I made an index.ndx manually for the dihedral of interest in the molecule. Usually, when we make an index file, say for atom P in the system, make_ndx makes an index file which includes the numbers for all the P atom in the system. These numbers are the index numbers in the .gro file for the whole system. In my index file, I made it based on the atom numbers in the topology file for the lipid. For doing analysis, I first used trjconv and saved only the .xtc file for only the lipid of interest, and then used g_angle to calculate the dihedral population. I think what I get is only for ONE lipid over time NOT an average over time and lipids in the system. The obvious solution would be to include all the dihedrals for all the lipids in it. BUT, this is not easy to do manually. I tried mk_angndx but it is not giving me what I want. Is there any way to do this? Thanks in advance for any comment and suggestion Cheers Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] dihedral MM scan
Thanks Justin, I did the scan for the whole range from -180 to 180 in 10 intervals. I used emtol = 50, and restrain force of 10 to have a better relaxation and restraining, respectively. Unfortunately, my MM profiles do not look similar to GAAMP QM or MM profiles, although I have used the suggested parameters. To check the method I am using and the parameters in the .mdp files, I did the MM scanning on a rotatable bond (CH2-CH2-CH2-CH3) similar to one can be found in a lipid tail. This dihedral was detected as a soft dihedral by GAAMP and there was a good agreement between QM and MM profiles. Here again, my MM profile from scanning does not seem reasonable. Most importantly, there is a very large barrier at angle 0 (ca. 300 times of energies in other angles). Besides, there are sudden increase or decrease in energy just by 10 degree of rotation. Any suggestion is appreciated in advance. Cheers Mohsen On Wed, Jan 18, 2017 at 7:55 PM, Justin Lemkul wrote: > > > On 1/18/17 9:51 PM, Mohsen Ramezanpour wrote: > >> Dear Gromacs users, >> >> I want to do a MM scan on a dihedral and compare with the available QM >> profile using Gromacs 5.1.3 >> >> I used the following (updated compared to my last email on parameters) >> em.mdp and md.mdp files for relaxation and zero-step md, respectively. >> >> *em.mdp:* >> define = -DPOSRESDIH >> integrator = steep >> emtol = 100.0 >> nsteps = 5 >> vdwtype = Cut-off >> vdw-modifier= none >> coulombtype = Cut-off >> constraints= h-bonds >> constraint_algorithm= LINCS >> nstlist= 0 >> ns-type = simple >> pbc = no >> cutoff-scheme = group >> rcoulomb = 0 >> rvdw = 0 >> rvdw_switch = 0 >> rlist= 0 >> >> Which, by -DPOSRES I am applying position restraint on specific dihedral >> of >> interest. >> The following lines were included at the end of molecule.itp >> >> #ifdef POSRES >> >> [ dihedral_restraints ] >> >> 118 8 9 1 120 0 1 >> >> #endif >> >> >> >> *as for md.mdp:* >> >> define = >> integrator = md >> dt= 0.001 >> nsteps = 0 >> nstlog = 0 >> nstxout = 0 >> nstvout = 0 >> nstfout = 0 >> nstcalcenergy = 0 >> nstenergy = 0 >> nstxout-compressed = 10 >> compressed-x-precision = 1000 >> coulombtype = Cut-off >> vdwtype = Cut-off >> vdw-modifier = none >> constraints = h-bonds >> constraint_algorithm= LINCS >> nstcomm = 100 >> comm_mode = linear >> comm_grps = system >> refcoord_scaling = com >> nstlist = 0 >> ns-type= simple >> pbc= no >> cutoff-scheme = group >> rcoulomb= 0 >> rvdw = 0 >> rvdw_switch= 0 >> rlist = 0 >> continuation= yes >> >> It seems to work fine, at least because grompp do not complain and it >> gives >> me the desired files. >> >> There is one problem: >> Looking at the Potential Energy which zero-step md for angle = 120 gives >> me >> a value of 384.682 KJ/mol (ca. 91 Kcal/mol) which is almost 20 times the >> value from QM profile (4.5 Kcal/mol) by GAAMP server (as I understood, the >> units for dihedral profiles by GAAMP are in Kcal/mol). >> >> I used the parameters GAAMP gave me (partial charges and dihedral >> parameters). >> >> Briefly, the same molecule which I uploaded to GAAMP was used for MM scan >> in vaccume. No box and infinite cutoff was used. I did an EM step by above >> em.mdp and let the system to relax. Then, I removed the restraints and did >> the zero-step md part and got this energy. >> >> It is worth mentioning that the molecule is positively charges. BUT, as I >> understood, this is what GAAMP also uses for the QM scan. >> >> What could be the problem? >> >> I highly appreciate your comm
[gmx-users] dihedral MM scan
Dear Gromacs users, I want to do a MM scan on a dihedral and compare with the available QM profile using Gromacs 5.1.3 I used the following (updated compared to my last email on parameters) em.mdp and md.mdp files for relaxation and zero-step md, respectively. *em.mdp:* define = -DPOSRESDIH integrator = steep emtol = 100.0 nsteps = 5 vdwtype = Cut-off vdw-modifier= none coulombtype = Cut-off constraints= h-bonds constraint_algorithm= LINCS nstlist= 0 ns-type = simple pbc = no cutoff-scheme = group rcoulomb = 0 rvdw = 0 rvdw_switch = 0 rlist= 0 Which, by -DPOSRES I am applying position restraint on specific dihedral of interest. The following lines were included at the end of molecule.itp #ifdef POSRES [ dihedral_restraints ] 118 8 9 1 120 0 1 #endif *as for md.mdp:* define = integrator = md dt= 0.001 nsteps = 0 nstlog = 0 nstxout = 0 nstvout = 0 nstfout = 0 nstcalcenergy = 0 nstenergy = 0 nstxout-compressed = 10 compressed-x-precision = 1000 coulombtype = Cut-off vdwtype = Cut-off vdw-modifier = none constraints = h-bonds constraint_algorithm= LINCS nstcomm = 100 comm_mode = linear comm_grps = system refcoord_scaling = com nstlist = 0 ns-type= simple pbc= no cutoff-scheme = group rcoulomb= 0 rvdw = 0 rvdw_switch= 0 rlist = 0 continuation= yes It seems to work fine, at least because grompp do not complain and it gives me the desired files. There is one problem: Looking at the Potential Energy which zero-step md for angle = 120 gives me a value of 384.682 KJ/mol (ca. 91 Kcal/mol) which is almost 20 times the value from QM profile (4.5 Kcal/mol) by GAAMP server (as I understood, the units for dihedral profiles by GAAMP are in Kcal/mol). I used the parameters GAAMP gave me (partial charges and dihedral parameters). Briefly, the same molecule which I uploaded to GAAMP was used for MM scan in vaccume. No box and infinite cutoff was used. I did an EM step by above em.mdp and let the system to relax. Then, I removed the restraints and did the zero-step md part and got this energy. It is worth mentioning that the molecule is positively charges. BUT, as I understood, this is what GAAMP also uses for the QM scan. What could be the problem? I highly appreciate your comments as expert in this. Cheers Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] MM dihedral scanning
On Wed, Jan 18, 2017 at 10:35 AM, Justin Lemkul wrote: > > > On 1/17/17 10:37 PM, Mohsen Ramezanpour wrote: > >> Hi Gromacs users, >> >> I want to follow recommended procedure (by Justin) for MM scanning of >> rotatable dihedrals and do it manually by Gromacs. >> >> I have made the following em.mdp (for restraining desired rotatable bond >> and relaxing the rest of molecule) and md-zero.mdp for the single-point >> energy calculation. >> I was wondering if these parameters seems reasonable. >> >> >> *em.mdp:* >> >> define = -DPOSRES >> > > What are you restraining here? The specific dihedral of interest for scan. explained more at the end. > > > integrator = steep >> emtol= 1000.0 >> nsteps = 5000 >> vdwtype= Cut-off >> vdw-modifier = Force-switch >> > > There is no force switching because the cutoffs are infinite. You'll > probably get an error. > Thanks. I belive I should use "None" then. Based on manual: "None: Use an unmodified Van der Waals potential. With the *group* scheme this means *no exact cut-off* is used, energies and forces are calculated for all pairs in the neighborlist. " And I have used the cutoff-scheme = group belove. > > coulombtype = pme >> ; >> constraints = h-bonds >> constraint_algorithm = LINCS >> ; to do a MM scan >> nstlist = 0 >> ns-type = simple >> pbc = no >> cutoff-scheme = group >> rcoulomb = 0 >> rvdw= 0 >> rvdw_switch = 0 >> rlist = 0 >> >> >> >> >> *md-zero.mdp:* >> >> define = >> integrator = md >> dt = 0.001 >> nsteps = 0 >> nstlog = 0 >> nstxout = 0 >> nstvout = 0 >> nstfout = 0 >> nstcalcenergy = 0 >> nstenergy = 0 >> ; >> nstxout-compressed= 10 >> compressed-x-precision = 1000 >> coulombtype = pme >> vdwtype = Cut-off >> vdw-modifier= Force-switch >> > > See above. > > constraints = h-bonds >> constraint_algorithm= LINCS >> nstcomm = 100 >> comm_mode = linear >> comm_grps = system >> refcoord_scaling= com >> ; >> nstlist = 0 >> ns-type = simple >> pbc = no >> cutoff-scheme = group >> rcoulomb = 0 >> rvdw= 0 >> rvdw_switch = 0 >> rlist = 0 >> >> *I will use this command to calculate it:* >> >> mdrun -s input.tpr -rerun configuration.pdb -nsteps 0 >> > > You don't need -nsteps. Sure. Regarding the restraints (I forward my last email here): If I understood correctly, NONE of the following options are necessary (not present) in Gromacs version 5.1.3: dihre = yes dihre_fc = 100 ; kJ/(mol rad^2) dihre_tau = 0.0 nstdihreout = 250 And there is not any similar parameter in this version too. However, the kfac in the following section should be the "real force" for restraining. So, for dihedral restraining during em part, including following lines in the topol.top (or at the end of molecule.itp) should be okay with above em.mdp parameters. #ifdef POSRES [ dihedral_restraints ] ; ai aj ak al type label phi dphi kfac A B C D 11 120 0 1 #endif Which means that I want to apply a force of 1 kj/mol/rad^2 to a specific dihedral of ABCD to an angle of 120 for the purpose of scanning. I chose dphi as "0" and I did not include "power" (because of version 5) Cheers > > > -Justin > > > >> Thanks in advance for your comments, >> Mohsen >> >> On Fri, Jan 13, 2017 at 2:04 PM, Justin Lemkul wrote: >> >> >>> >>> On 1/13/17 3:26 PM, Mohsen Ramezanpour wrote: >>> >>> Thanks. interesting! >>>> So, for now, lets focus on one of those dihedrals which has another >>>> equivalent. >>>> If I do PES on a rotatable bond in one ethyl, the other equivalent >>>> rotatable bond in other ethyl will be free to move. If I understood >>>&
Re: [gmx-users] MM dihedral scanning
If I understood correctly, NONE of the following options are necessary (not present) in Gromacs version 5.1.3: dihre = yes dihre_fc = 100 ; kJ/(mol rad^2) dihre_tau = 0.0 nstdihreout = 250 And there is not any similar parameter in this version too. However, the kfac in the following section should be the "real force" for restraining. So, for dihedral restraining during em part, including following lines in the topol.top (or at the end of molecule.itp) should be okay with above em.mdp parameters. #ifdef POSRES [ dihedral_restraints ] ; ai aj ak al type label phi dphi kfac A B C D 11 120 0 1 #endif Which means that I want to apply a force of 1 kj/mol/rad^2 to a specific dihedral of ABCD to an angle of 120 for the purpose of scanning. I chose dphi as "0" and I did not include "power" (because of version 5) Cheers On Tue, Jan 17, 2017 at 8:39 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > And oh, none of the Kenno's web pages for tutorial and software is working > at the moment. > > On Tue, Jan 17, 2017 at 8:37 PM, Mohsen Ramezanpour < > ramezanpour.moh...@gmail.com> wrote: > >> Hi Gromacs users, >> >> I want to follow recommended procedure (by Justin) for MM scanning of >> rotatable dihedrals and do it manually by Gromacs. >> >> I have made the following em.mdp (for restraining desired rotatable bond >> and relaxing the rest of molecule) and md-zero.mdp for the single-point >> energy calculation. >> I was wondering if these parameters seems reasonable. >> >> >> *em.mdp:* >> >> define = -DPOSRES >> integrator = steep >> emtol= 1000.0 >> nsteps = 5000 >> vdwtype= Cut-off >> vdw-modifier = Force-switch >> coulombtype = pme >> ; >> constraints = h-bonds >> constraint_algorithm = LINCS >> ; to do a MM scan >> nstlist = 0 >> ns-type = simple >> pbc = no >> cutoff-scheme = group >> rcoulomb = 0 >> rvdw= 0 >> rvdw_switch = 0 >> rlist = 0 >> >> >> >> >> *md-zero.mdp:* >> >> define = >> integrator = md >> dt = 0.001 >> nsteps = 0 >> nstlog = 0 >> nstxout = 0 >> nstvout = 0 >> nstfout = 0 >> nstcalcenergy = 0 >> nstenergy = 0 >> ; >> nstxout-compressed= 10 >> compressed-x-precision = 1000 >> coulombtype = pme >> vdwtype = Cut-off >> vdw-modifier= Force-switch >> constraints = h-bonds >> constraint_algorithm= LINCS >> nstcomm = 100 >> comm_mode = linear >> comm_grps = system >> refcoord_scaling= com >> ; >> nstlist = 0 >> ns-type = simple >> pbc = no >> cutoff-scheme = group >> rcoulomb = 0 >> rvdw= 0 >> rvdw_switch = 0 >> rlist = 0 >> >> *I will use this command to calculate it:* >> >> mdrun -s input.tpr -rerun configuration.pdb -nsteps 0 >> >> Thanks in advance for your comments, >> Mohsen >> >> On Fri, Jan 13, 2017 at 2:04 PM, Justin Lemkul wrote: >> >>> >>> >>> On 1/13/17 3:26 PM, Mohsen Ramezanpour wrote: >>> >>>> Thanks. interesting! >>>> So, for now, lets focus on one of those dihedrals which has another >>>> equivalent. >>>> If I do PES on a rotatable bond in one ethyl, the other equivalent >>>> rotatable bond in other ethyl will be free to move. If I understood >>>> correctly, you meant it does not matter and I can do PES with one of the >>>> two. >>>> >>>> Looking at QM and MM profiles by GAAMP, I figured out something >>>> interesting: >>>> If I shift the MM graph for 120 in phi direction, MM and QM will be >>>> similar >>>> in their general behaviour and it is much much better that before. >>>> It works for shift of 240 too, but 120 results in better agreement with >>>> QM >>>>
Re: [gmx-users] MM dihedral scanning
And oh, none of the Kenno's web pages for tutorial and software is working at the moment. On Tue, Jan 17, 2017 at 8:37 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Gromacs users, > > I want to follow recommended procedure (by Justin) for MM scanning of > rotatable dihedrals and do it manually by Gromacs. > > I have made the following em.mdp (for restraining desired rotatable bond > and relaxing the rest of molecule) and md-zero.mdp for the single-point > energy calculation. > I was wondering if these parameters seems reasonable. > > > *em.mdp:* > > define = -DPOSRES > integrator = steep > emtol= 1000.0 > nsteps = 5000 > vdwtype= Cut-off > vdw-modifier = Force-switch > coulombtype = pme > ; > constraints = h-bonds > constraint_algorithm = LINCS > ; to do a MM scan > nstlist = 0 > ns-type = simple > pbc = no > cutoff-scheme = group > rcoulomb = 0 > rvdw= 0 > rvdw_switch = 0 > rlist = 0 > > > > > *md-zero.mdp:* > > define = > integrator = md > dt = 0.001 > nsteps = 0 > nstlog = 0 > nstxout = 0 > nstvout = 0 > nstfout = 0 > nstcalcenergy = 0 > nstenergy = 0 > ; > nstxout-compressed= 10 > compressed-x-precision = 1000 > coulombtype = pme > vdwtype = Cut-off > vdw-modifier= Force-switch > constraints = h-bonds > constraint_algorithm= LINCS > nstcomm = 100 > comm_mode = linear > comm_grps = system > refcoord_scaling= com > ; > nstlist = 0 > ns-type = simple > pbc = no > cutoff-scheme = group > rcoulomb = 0 > rvdw= 0 > rvdw_switch = 0 > rlist = 0 > > *I will use this command to calculate it:* > > mdrun -s input.tpr -rerun configuration.pdb -nsteps 0 > > Thanks in advance for your comments, > Mohsen > > On Fri, Jan 13, 2017 at 2:04 PM, Justin Lemkul wrote: > >> >> >> On 1/13/17 3:26 PM, Mohsen Ramezanpour wrote: >> >>> Thanks. interesting! >>> So, for now, lets focus on one of those dihedrals which has another >>> equivalent. >>> If I do PES on a rotatable bond in one ethyl, the other equivalent >>> rotatable bond in other ethyl will be free to move. If I understood >>> correctly, you meant it does not matter and I can do PES with one of the >>> two. >>> >>> Looking at QM and MM profiles by GAAMP, I figured out something >>> interesting: >>> If I shift the MM graph for 120 in phi direction, MM and QM will be >>> similar >>> in their general behaviour and it is much much better that before. >>> It works for shift of 240 too, but 120 results in better agreement with >>> QM >>> PES. >>> >> >> This sounds like the phase and/or multiplicity are off. Perhaps a bad >> initial guess of the parameters threw off the fitting. >> >> So, there might be two possible reasons for this: >>> 1) parameters are not good and this is how it is >>> 2) GAAMP has started rotation from different starting points (with a >>> difference of 120) >>> 3) GAAMP made a mistake for taking one dihedral for QM and took the >>> equivalent one for MM PES. (i.e. taking one of oxygen atoms in QM and >>> other >>> Oxygen in MM for the same rotatable bond) >>> >>> I should read the deatils of how GAAMP did these exactly to rule out the >>> options 2 and 3 >>> >> >> That's the beauty of GAAMP - it gives you full input and output for all >> the QM and MM processes. It should be quite apparent what the programs >> did. Perhaps it did not deal elegantly with the symmetry of the molecule. >> Its programs are quite robust but no black box is perfect (automated >> parametrization is a long-standing goal but in truth there's a lot of work >> left to be done). >> >> -Justin >> >> >> Cheers >>> >>> >>> >>> On Fri, Jan 13, 2017 at 5:45 AM, Justin Lemkul wrote: >>> >
Re: [gmx-users] MM dihedral scanning
Hi Gromacs users, I want to follow recommended procedure (by Justin) for MM scanning of rotatable dihedrals and do it manually by Gromacs. I have made the following em.mdp (for restraining desired rotatable bond and relaxing the rest of molecule) and md-zero.mdp for the single-point energy calculation. I was wondering if these parameters seems reasonable. *em.mdp:* define = -DPOSRES integrator = steep emtol= 1000.0 nsteps = 5000 vdwtype= Cut-off vdw-modifier = Force-switch coulombtype = pme ; constraints = h-bonds constraint_algorithm = LINCS ; to do a MM scan nstlist = 0 ns-type = simple pbc = no cutoff-scheme = group rcoulomb = 0 rvdw= 0 rvdw_switch = 0 rlist = 0 *md-zero.mdp:* define = integrator = md dt = 0.001 nsteps = 0 nstlog = 0 nstxout = 0 nstvout = 0 nstfout = 0 nstcalcenergy = 0 nstenergy = 0 ; nstxout-compressed= 10 compressed-x-precision = 1000 coulombtype = pme vdwtype = Cut-off vdw-modifier= Force-switch constraints = h-bonds constraint_algorithm= LINCS nstcomm = 100 comm_mode = linear comm_grps = system refcoord_scaling= com ; nstlist = 0 ns-type = simple pbc = no cutoff-scheme = group rcoulomb = 0 rvdw= 0 rvdw_switch = 0 rlist = 0 *I will use this command to calculate it:* mdrun -s input.tpr -rerun configuration.pdb -nsteps 0 Thanks in advance for your comments, Mohsen On Fri, Jan 13, 2017 at 2:04 PM, Justin Lemkul wrote: > > > On 1/13/17 3:26 PM, Mohsen Ramezanpour wrote: > >> Thanks. interesting! >> So, for now, lets focus on one of those dihedrals which has another >> equivalent. >> If I do PES on a rotatable bond in one ethyl, the other equivalent >> rotatable bond in other ethyl will be free to move. If I understood >> correctly, you meant it does not matter and I can do PES with one of the >> two. >> >> Looking at QM and MM profiles by GAAMP, I figured out something >> interesting: >> If I shift the MM graph for 120 in phi direction, MM and QM will be >> similar >> in their general behaviour and it is much much better that before. >> It works for shift of 240 too, but 120 results in better agreement with QM >> PES. >> > > This sounds like the phase and/or multiplicity are off. Perhaps a bad > initial guess of the parameters threw off the fitting. > > So, there might be two possible reasons for this: >> 1) parameters are not good and this is how it is >> 2) GAAMP has started rotation from different starting points (with a >> difference of 120) >> 3) GAAMP made a mistake for taking one dihedral for QM and took the >> equivalent one for MM PES. (i.e. taking one of oxygen atoms in QM and >> other >> Oxygen in MM for the same rotatable bond) >> >> I should read the deatils of how GAAMP did these exactly to rule out the >> options 2 and 3 >> > > That's the beauty of GAAMP - it gives you full input and output for all > the QM and MM processes. It should be quite apparent what the programs > did. Perhaps it did not deal elegantly with the symmetry of the molecule. > Its programs are quite robust but no black box is perfect (automated > parametrization is a long-standing goal but in truth there's a lot of work > left to be done). > > -Justin > > > Cheers >> >> >> >> On Fri, Jan 13, 2017 at 5:45 AM, Justin Lemkul wrote: >> >> >>> >>> On 1/12/17 11:20 PM, Mohsen Ramezanpour wrote: >>> >>> On Thu, Jan 12, 2017 at 8:36 PM, Justin Lemkul wrote: >>>> >>>> >>>> >>>>> On 1/12/17 10:21 PM, Mohsen Ramezanpour wrote: >>>>> >>>>> On Thu, Jan 12, 2017 at 6:31 PM, Justin Lemkul >>>>> wrote: >>>>> >>>>>> >>>>>> >>>>>> >>>>>> On 1/12/17 7:47 PM, Mohsen Ramezanpour wrote: >>>>>>> >>>>>>> On Thu, Jan 12, 2017 at 4:24 PM, Justin Lemkul >>>>>>> wrote: >>>>>>> >>>>>>> >>>>>>>> >>&g
Re: [gmx-users] MM dihedral scanning
Thanks. interesting! So, for now, lets focus on one of those dihedrals which has another equivalent. If I do PES on a rotatable bond in one ethyl, the other equivalent rotatable bond in other ethyl will be free to move. If I understood correctly, you meant it does not matter and I can do PES with one of the two. Looking at QM and MM profiles by GAAMP, I figured out something interesting: If I shift the MM graph for 120 in phi direction, MM and QM will be similar in their general behaviour and it is much much better that before. It works for shift of 240 too, but 120 results in better agreement with QM PES. So, there might be two possible reasons for this: 1) parameters are not good and this is how it is 2) GAAMP has started rotation from different starting points (with a difference of 120) 3) GAAMP made a mistake for taking one dihedral for QM and took the equivalent one for MM PES. (i.e. taking one of oxygen atoms in QM and other Oxygen in MM for the same rotatable bond) I should read the deatils of how GAAMP did these exactly to rule out the options 2 and 3 Cheers On Fri, Jan 13, 2017 at 5:45 AM, Justin Lemkul wrote: > > > On 1/12/17 11:20 PM, Mohsen Ramezanpour wrote: > >> On Thu, Jan 12, 2017 at 8:36 PM, Justin Lemkul wrote: >> >> >>> >>> On 1/12/17 10:21 PM, Mohsen Ramezanpour wrote: >>> >>> On Thu, Jan 12, 2017 at 6:31 PM, Justin Lemkul wrote: >>>> >>>> >>>> >>>>> On 1/12/17 7:47 PM, Mohsen Ramezanpour wrote: >>>>> >>>>> On Thu, Jan 12, 2017 at 4:24 PM, Justin Lemkul >>>>> wrote: >>>>> >>>>>> >>>>>> >>>>>> >>>>>> On 1/12/17 3:18 PM, Mohsen Ramezanpour wrote: >>>>>>> >>>>>>> Dear Gromacs users, >>>>>>> >>>>>>> >>>>>>>> For parameterization of a molecule in Charmm36, I have got the QM >>>>>>>> scanning >>>>>>>> and partial charges from GAMMP server. However, the fitted >>>>>>>> parameters >>>>>>>> are >>>>>>>> not good enough. >>>>>>>> >>>>>>>> >>>>>>>> That's very surprising. What's wrong with what GAAMP gave you? >>>>>>>> >>>>>>> >>>>>>> The dihedral has two local minima in both QM and fitted ones both >>>>>>> from >>>>>>> >>>>>>> GAAMP. The angle for the minima are okay but the corresponding depths >>>>>> are >>>>>> not. >>>>>> In fact, the depth for the first local minimum is larger than second >>>>>> one >>>>>> in >>>>>> QM, while the situation is reverse in MM profile with fitted >>>>>> parameters. >>>>>> This makes the dihedral to be more (statistically) in wrong angle (in >>>>>> local >>>>>> minimum which is not the most favourable one). >>>>>> GAAMP, unfortunately, did not work well with my case (some critical >>>>>> partial >>>>>> charges and critical dihedrals). >>>>>> >>>>>> >>>>>> I decided to do the MM scanning and try to get better parameters for >>>>>> the >>>>>> >>>>>> >>>>>>> dihedral. >>>>>>> >>>>>>>> >>>>>>>> Unfortunately, I do not have any experience with this part, and I >>>>>>>> could >>>>>>>> not >>>>>>>> find any tutorial for how to do this in Gromacs. >>>>>>>> I was wondering if you are aware of any tutorial which could help me >>>>>>>> to >>>>>>>> overcome this challenging step. >>>>>>>> >>>>>>>> >>>>>>>> Tutorial (CGenFF theory is the same as CHARMM, by design): >>>>>>>> >>>>>>>> http://mackerell.umaryland.edu/~kenno/cgenff/download.php#tutor >>>>>>> >>>>>>> Fitting program and other resources: >>>>>>> http://mackerell.umaryland.edu/~kenno/lsfitpar/ >>>>>>> http://mackerell.umaryland.edu/ff_dev.shtml >>>>>>> >>>>>>> Obviously, these are all CHARMM-centric approaches and frankly the
Re: [gmx-users] MM dihedral scanning
On Thu, Jan 12, 2017 at 8:36 PM, Justin Lemkul wrote: > > > On 1/12/17 10:21 PM, Mohsen Ramezanpour wrote: > >> On Thu, Jan 12, 2017 at 6:31 PM, Justin Lemkul wrote: >> >> >>> >>> On 1/12/17 7:47 PM, Mohsen Ramezanpour wrote: >>> >>> On Thu, Jan 12, 2017 at 4:24 PM, Justin Lemkul wrote: >>>> >>>> >>>> >>>>> On 1/12/17 3:18 PM, Mohsen Ramezanpour wrote: >>>>> >>>>> Dear Gromacs users, >>>>> >>>>>> >>>>>> For parameterization of a molecule in Charmm36, I have got the QM >>>>>> scanning >>>>>> and partial charges from GAMMP server. However, the fitted parameters >>>>>> are >>>>>> not good enough. >>>>>> >>>>>> >>>>>> That's very surprising. What's wrong with what GAAMP gave you? >>>>> >>>>> The dihedral has two local minima in both QM and fitted ones both from >>>>> >>>> GAAMP. The angle for the minima are okay but the corresponding depths >>>> are >>>> not. >>>> In fact, the depth for the first local minimum is larger than second one >>>> in >>>> QM, while the situation is reverse in MM profile with fitted parameters. >>>> This makes the dihedral to be more (statistically) in wrong angle (in >>>> local >>>> minimum which is not the most favourable one). >>>> GAAMP, unfortunately, did not work well with my case (some critical >>>> partial >>>> charges and critical dihedrals). >>>> >>>> >>>> I decided to do the MM scanning and try to get better parameters for the >>>> >>>>> >>>>> dihedral. >>>>>> >>>>>> Unfortunately, I do not have any experience with this part, and I >>>>>> could >>>>>> not >>>>>> find any tutorial for how to do this in Gromacs. >>>>>> I was wondering if you are aware of any tutorial which could help me >>>>>> to >>>>>> overcome this challenging step. >>>>>> >>>>>> >>>>>> Tutorial (CGenFF theory is the same as CHARMM, by design): >>>>>> >>>>> http://mackerell.umaryland.edu/~kenno/cgenff/download.php#tutor >>>>> >>>>> Fitting program and other resources: >>>>> http://mackerell.umaryland.edu/~kenno/lsfitpar/ >>>>> http://mackerell.umaryland.edu/ff_dev.shtml >>>>> >>>>> Obviously, these are all CHARMM-centric approaches and frankly the >>>>> modules >>>>> within CHARMM make parametrization rather straightforward (not "easy," >>>>> mind >>>>> you, but straightforward). Since I began working with CHARMM, it has >>>>> become indispensable in my daily routine. >>>>> >>>>> If you want to do things in GROMACS, the main issue is that you will >>>>> have >>>>> to do MM scans in a more manual fashion, by restraining the target >>>>> dihedrals (very strongly) in a series of configurations (typically at >>>>> intervals of 15 degrees over a full rotation) while allowing the rest >>>>> of >>>>> the molecule to relax to match the QM. >>>>> >>>>> >>>> If I got it right, I must do EM on the molecule while only the desired >>>> dihedral is fixed in a specific angle. Which aspect should match with >>>> QM? >>>> If you mean structurally, what is the criteria for matching (RMSD?.)? >>>> >>>> >>>> "Match" in this case means "treat equivalently," therefore only one >>> constraint (restraint in the MM) should be applied while allowing the >>> rest >>> of the molecule to relax freely. You do have a difficult case because >>> each >>> of your molecules is symmetric; this means the same dihedral term is >>> affecting multiple torsions. >>> >>> So, probably I should 4 dihedrals and try to optimize all at once?!( >> because all 4 dihedrals of O-C-CH2-CH3 seems equivalent to me). Am I >> right? >> >> > No, there are 2 O-C-CH2-CH3 dihedrals, not 4. > How come? :-) there are two ethyls and two oxygens in the ring. Each ethyl forms two dihedrals (one with each oxygen). Thus, there
Re: [gmx-users] MM dihedral scanning
I was thinking of a potential solution for this. Please let me know your opinion as expert. Sorry, if it seems stupid :-) I have the QM scan profile already. 1-4 interactions for this dihedral must be included in charmm36. My idea: I was thinking to find the dihedral parameters *purely mathematically* (by using the dihedral potential function), so that the fitted mathematical curve will match as much as possible to the QM one. This should be doable by including more terms in dihedral fitting and playing with force constants, multiplicity, and angles.* BUT*, this time, I will exclude the 1-4 interactions *in this specific dihedral* when I am doing my simulations using charmm36. I think this approach might work. However, I am not sure if that is okay to exclude one specific dihedral in all the molecules of system (in all lipids of a bilayer) while other 1-4 are included (as it is the case for charmm36 ff). Please let me know your opinion. Cheers Mohsen On Thu, Jan 12, 2017 at 8:21 PM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > On Thu, Jan 12, 2017 at 6:31 PM, Justin Lemkul wrote: > >> >> >> On 1/12/17 7:47 PM, Mohsen Ramezanpour wrote: >> >>> On Thu, Jan 12, 2017 at 4:24 PM, Justin Lemkul wrote: >>> >>> >>>> >>>> On 1/12/17 3:18 PM, Mohsen Ramezanpour wrote: >>>> >>>> Dear Gromacs users, >>>>> >>>>> For parameterization of a molecule in Charmm36, I have got the QM >>>>> scanning >>>>> and partial charges from GAMMP server. However, the fitted parameters >>>>> are >>>>> not good enough. >>>>> >>>>> >>>> That's very surprising. What's wrong with what GAAMP gave you? >>>> >>>> The dihedral has two local minima in both QM and fitted ones both from >>> GAAMP. The angle for the minima are okay but the corresponding depths are >>> not. >>> In fact, the depth for the first local minimum is larger than second one >>> in >>> QM, while the situation is reverse in MM profile with fitted parameters. >>> This makes the dihedral to be more (statistically) in wrong angle (in >>> local >>> minimum which is not the most favourable one). >>> GAAMP, unfortunately, did not work well with my case (some critical >>> partial >>> charges and critical dihedrals). >>> >>> >>> I decided to do the MM scanning and try to get better parameters for the >>>> >>>>> dihedral. >>>>> >>>>> Unfortunately, I do not have any experience with this part, and I could >>>>> not >>>>> find any tutorial for how to do this in Gromacs. >>>>> I was wondering if you are aware of any tutorial which could help me to >>>>> overcome this challenging step. >>>>> >>>>> >>>>> Tutorial (CGenFF theory is the same as CHARMM, by design): >>>> http://mackerell.umaryland.edu/~kenno/cgenff/download.php#tutor >>>> >>>> Fitting program and other resources: >>>> http://mackerell.umaryland.edu/~kenno/lsfitpar/ >>>> http://mackerell.umaryland.edu/ff_dev.shtml >>>> >>>> Obviously, these are all CHARMM-centric approaches and frankly the >>>> modules >>>> within CHARMM make parametrization rather straightforward (not "easy," >>>> mind >>>> you, but straightforward). Since I began working with CHARMM, it has >>>> become indispensable in my daily routine. >>>> >>>> If you want to do things in GROMACS, the main issue is that you will >>>> have >>>> to do MM scans in a more manual fashion, by restraining the target >>>> dihedrals (very strongly) in a series of configurations (typically at >>>> intervals of 15 degrees over a full rotation) while allowing the rest of >>>> the molecule to relax to match the QM. >>>> >>> >>> If I got it right, I must do EM on the molecule while only the desired >>> dihedral is fixed in a specific angle. Which aspect should match with QM? >>> If you mean structurally, what is the criteria for matching (RMSD?.)? >>> >>> >> "Match" in this case means "treat equivalently," therefore only one >> constraint (restraint in the MM) should be applied while allowing the rest >> of the molecule to relax freely. You do have a difficult case because each >> of your molecules is symmetric; this means the same dihedral term is >> affecting
Re: [gmx-users] MM dihedral scanning
On Thu, Jan 12, 2017 at 6:31 PM, Justin Lemkul wrote: > > > On 1/12/17 7:47 PM, Mohsen Ramezanpour wrote: > >> On Thu, Jan 12, 2017 at 4:24 PM, Justin Lemkul wrote: >> >> >>> >>> On 1/12/17 3:18 PM, Mohsen Ramezanpour wrote: >>> >>> Dear Gromacs users, >>>> >>>> For parameterization of a molecule in Charmm36, I have got the QM >>>> scanning >>>> and partial charges from GAMMP server. However, the fitted parameters >>>> are >>>> not good enough. >>>> >>>> >>> That's very surprising. What's wrong with what GAAMP gave you? >>> >>> The dihedral has two local minima in both QM and fitted ones both from >> GAAMP. The angle for the minima are okay but the corresponding depths are >> not. >> In fact, the depth for the first local minimum is larger than second one >> in >> QM, while the situation is reverse in MM profile with fitted parameters. >> This makes the dihedral to be more (statistically) in wrong angle (in >> local >> minimum which is not the most favourable one). >> GAAMP, unfortunately, did not work well with my case (some critical >> partial >> charges and critical dihedrals). >> >> >> I decided to do the MM scanning and try to get better parameters for the >>> >>>> dihedral. >>>> >>>> Unfortunately, I do not have any experience with this part, and I could >>>> not >>>> find any tutorial for how to do this in Gromacs. >>>> I was wondering if you are aware of any tutorial which could help me to >>>> overcome this challenging step. >>>> >>>> >>>> Tutorial (CGenFF theory is the same as CHARMM, by design): >>> http://mackerell.umaryland.edu/~kenno/cgenff/download.php#tutor >>> >>> Fitting program and other resources: >>> http://mackerell.umaryland.edu/~kenno/lsfitpar/ >>> http://mackerell.umaryland.edu/ff_dev.shtml >>> >>> Obviously, these are all CHARMM-centric approaches and frankly the >>> modules >>> within CHARMM make parametrization rather straightforward (not "easy," >>> mind >>> you, but straightforward). Since I began working with CHARMM, it has >>> become indispensable in my daily routine. >>> >>> If you want to do things in GROMACS, the main issue is that you will have >>> to do MM scans in a more manual fashion, by restraining the target >>> dihedrals (very strongly) in a series of configurations (typically at >>> intervals of 15 degrees over a full rotation) while allowing the rest of >>> the molecule to relax to match the QM. >>> >> >> If I got it right, I must do EM on the molecule while only the desired >> dihedral is fixed in a specific angle. Which aspect should match with QM? >> If you mean structurally, what is the criteria for matching (RMSD?.)? >> >> > "Match" in this case means "treat equivalently," therefore only one > constraint (restraint in the MM) should be applied while allowing the rest > of the molecule to relax freely. You do have a difficult case because each > of your molecules is symmetric; this means the same dihedral term is > affecting multiple torsions. > So, probably I should 4 dihedrals and try to optimize all at once?!( because all 4 dihedrals of O-C-CH2-CH3 seems equivalent to me). Am I right? > > >> >> Deactivate the restraint, obtain the potential energy of the molecule >> via >> >>> mdrun -rerun and plot as a function of the dihedral. >>> >> >> This should be a zero step EM, right? The molecule should not be allowed >> to >> change its conformation. >> >> > No, a zero-step MD. EM actually changes the coordinates before step zero. > > A bit of shell scripting and careful topology modification and this can >> >>> be done. >>> >>> >> Two more questions on this part: >> 1) I am using the QM scanning data and partial charges from GAAMP. When I >> do this MM scanning, do I need to exclude any 1-4 interactions or I can >> behave this dihedral as other dihedrals? >> > > 1-4 interactions are always at full strength in CHARMM. > > 2) this dihedral is part of a lipid. >> Do I need to do these on only "one Lipid in vacuum" or >> OR >> on all lipids in "a bilayer in vacuum" or in a "bilayer in solvent" >> >> I think it should be "one Lipid in vacuum". >> >
Re: [gmx-users] MM dihedral scanning
On Thu, Jan 12, 2017 at 4:24 PM, Justin Lemkul wrote: > > > On 1/12/17 3:18 PM, Mohsen Ramezanpour wrote: > >> Dear Gromacs users, >> >> For parameterization of a molecule in Charmm36, I have got the QM scanning >> and partial charges from GAMMP server. However, the fitted parameters are >> not good enough. >> > > That's very surprising. What's wrong with what GAAMP gave you? > The dihedral has two local minima in both QM and fitted ones both from GAAMP. The angle for the minima are okay but the corresponding depths are not. In fact, the depth for the first local minimum is larger than second one in QM, while the situation is reverse in MM profile with fitted parameters. This makes the dihedral to be more (statistically) in wrong angle (in local minimum which is not the most favourable one). GAAMP, unfortunately, did not work well with my case (some critical partial charges and critical dihedrals). > I decided to do the MM scanning and try to get better parameters for the >> dihedral. >> >> Unfortunately, I do not have any experience with this part, and I could >> not >> find any tutorial for how to do this in Gromacs. >> I was wondering if you are aware of any tutorial which could help me to >> overcome this challenging step. >> >> > Tutorial (CGenFF theory is the same as CHARMM, by design): > http://mackerell.umaryland.edu/~kenno/cgenff/download.php#tutor > > Fitting program and other resources: > http://mackerell.umaryland.edu/~kenno/lsfitpar/ > http://mackerell.umaryland.edu/ff_dev.shtml > > Obviously, these are all CHARMM-centric approaches and frankly the modules > within CHARMM make parametrization rather straightforward (not "easy," mind > you, but straightforward). Since I began working with CHARMM, it has > become indispensable in my daily routine. > > If you want to do things in GROMACS, the main issue is that you will have > to do MM scans in a more manual fashion, by restraining the target > dihedrals (very strongly) in a series of configurations (typically at > intervals of 15 degrees over a full rotation) while allowing the rest of > the molecule to relax to match the QM. If I got it right, I must do EM on the molecule while only the desired dihedral is fixed in a specific angle. Which aspect should match with QM? If you mean structurally, what is the criteria for matching (RMSD?.)? Deactivate the restraint, obtain the potential energy of the molecule via > mdrun -rerun and plot as a function of the dihedral. This should be a zero step EM, right? The molecule should not be allowed to change its conformation. A bit of shell scripting and careful topology modification and this can > be done. > Two more questions on this part: 1) I am using the QM scanning data and partial charges from GAAMP. When I do this MM scanning, do I need to exclude any 1-4 interactions or I can behave this dihedral as other dihedrals? 2) this dihedral is part of a lipid. Do I need to do these on only "one Lipid in vacuum" or OR on all lipids in "a bilayer in vacuum" or in a "bilayer in solvent" I think it should be "one Lipid in vacuum". > > -Justin > > There are three parts of molecule which I am interested in its dihedral >> parameters that I am stuck in it for a while: >> 1) (*2,2-Diethyl-1,3-dioxolane*) >> http://www.chemspider.com/Chemical-Structure.217102.html? >> rid=378be046-1c14-4bce-ac46-6591776f7e08 >> >> dihedrals of O-C-CH2-CH3 >> >> 2) (DIMETHYLPROPYLAMINE) >> http://www.chemspider.com/Chemical-Structure.55178.html?rid= >> 16e49844-d25d-48c8-b1b5-89fb2f9372bd >> >> dihedral of CH3-N-CH2-CH2 >> >> Many thanks in advance for any suggestion. >> >> Cheers >> Mohsen >> >> > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] MM dihedral scanning
Dear Gromacs users, For parameterization of a molecule in Charmm36, I have got the QM scanning and partial charges from GAMMP server. However, the fitted parameters are not good enough. I decided to do the MM scanning and try to get better parameters for the dihedral. Unfortunately, I do not have any experience with this part, and I could not find any tutorial for how to do this in Gromacs. I was wondering if you are aware of any tutorial which could help me to overcome this challenging step. There are three parts of molecule which I am interested in its dihedral parameters that I am stuck in it for a while: 1) (*2,2-Diethyl-1,3-dioxolane*) http://www.chemspider.com/Chemical-Structure.217102.html?rid=378be046-1c14-4bce-ac46-6591776f7e08 dihedrals of O-C-CH2-CH3 2) (DIMETHYLPROPYLAMINE) http://www.chemspider.com/Chemical-Structure.55178.html?rid=16e49844-d25d-48c8-b1b5-89fb2f9372bd dihedral of CH3-N-CH2-CH2 Many thanks in advance for any suggestion. Cheers Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] rlist, rcolumb, and rvdw
I found the discussion and the links inside very useful. Thanks! I agree with was discussed above. I should not mess around with cut-offs :-) Regarding the main question on this post (Just for curiosity): If one wish to use 0.8-0.1 (mess around :-) ), shall we change the rlist and rcolumb to 0.1 as well? OR, as far as rlist is bigger than rvdw and rcolumb, we are good to to (except taking more computing time) Cheers On Wed, Jan 11, 2017 at 7:56 AM, Justin Lemkul wrote: > > > On 1/11/17 9:42 AM, Thomas Piggot wrote: > >> Well, IMHO, for a lipids only system I'd would actually change those to >> use >> rvdw-switch = 0.8. Not only is this consistent with the original CHARMM36 >> lipids >> publication but for DPPC/POPC (at least) produces better membrane >> properties in >> GROMACS. >> >> See https://github.com/NMRLipids/NmrLipidsCholXray/issues/4 and my >> comments/results in the discussion for more details. >> >> > Interesting discussion. Indeed, the 1.0/1.2 cutoff settings are preferred > for compatibility with everything in the CHARMM36 force field. In our > hands, when we tested DPPC, 0.8/1.2 and 1.0/1.2 produced results in > reasonable agreement (e.g. the CHARMM-GUI paper referenced several times in > that discussion) though compressibility increased with the shorter > switching region. We tested 0.8/1.0 because everyone complains about the > "long" cutoff required to use CHARMM :) If only everything were like > AMBER, which uses 0.8... > > All of this should become moot as we work to refine the force field with > the recently implemented dispersion PME. No promises on when that will be > done, of course. > > -Justin > > > Cheers >> >> Tom >> >> On 11/01/17 12:51, Justin Lemkul wrote: >> >>> >>> >>> On 1/11/17 2:57 AM, Mark Abraham wrote: >>> >>>> Hi, >>>> >>>> Those ancient comments pertain only to the deprecated "group" cutoff >>>> scheme. You should look at the extensive documentation of both schemes >>>> in >>>> the current reference manuals. Particularly vdw cutoffs are baked into >>>> the >>>> forcefield and should not be varied without extensive testing. See >>>> http://www.gromacs.org/Documentation/Terminology/Force_Fields/CHARMM >>>> >>>> >>> To emphasize, this is especially true for CHARMM lipid parameters, as is >>> the >>> case here. >>> >>> Repeat after me: "I will not mess around with cutoffs for lipids." :) >>> >>> -Justin >>> >>> Mark >>>> >>>> On Wed, 11 Jan 2017 04:51 Mohsen Ramezanpour < >>>> ramezanpour.moh...@gmail.com> >>>> wrote: >>>> >>>> Dear gromacs users, >>>>> >>>>> Reading through mailing list I found a nice discussion on the relation >>>>> between rlist, rcolumb, and rvdw: >>>>> >>>>> https://www.mail-archive.com/gmx-users@gromacs.org/msg08387.html >>>>> >>>>> >>>>> https://mailman-1.sys.kth.se/pipermail/gromacs.org_gmx-users >>>>> /2011-February/058507.html >>>>> >>>>> >>>>> If I understood correctly, >>>>> >>>>> rcoulomb < rlist and rvdw < rlist is the most accurate way, >>>>> >>>>> and >>>>> rlist=rcoulomb=rvdw is the commonly used way for these parameters. >>>>> >>>>> I was wondering what would be the case if I want to change the >>>>> cut-offs? >>>>> >>>>> For instance, to use charmm36 in Groamcs, it is recommended to use >>>>> 1.0-1.2 >>>>> >>>>> as the cut-off for LJ, and rlist=rcoulomb=rvdw >>>>> >>>>> However, this setting might not result in good behaviour for some >>>>> lipid bilayers. >>>>> >>>>> Therefore, I want to check if other cut-offs works better in gromacs. >>>>> >>>>> If I change "1.2 nm" to "A nm" for rvdw, do I need to change rcolumb >>>>> and >>>>> rlist to A as well? i.e. rlist=rcolumb=rvdw=A? >>>>> >>>>> I am using following parameters already: >>>>> >>>>> cutoff-scheme = Verlet >>>>> nstlist = 20 >>>>> rlist = *1.2* >>>>> coulombty
Re: [gmx-users] Fwd: LJ cut-offs
Thanks Justin, Good to know. I used version 5.0 as Lee et al. 2015. I agree with you in general about keeping cut-offs the same. My systems are pure PS lipid membranes. The study by Klauda in 2016 shows that the higher order in PS membranes might be because of PS-PS or PS-K+ interactions. However, this was not done in Gromacs. This order was intensified when I used Gromacs (version 5.0, though). I will do with 5.1.* to see how it improves (hope it works). If PS-PS and PS-K+ interactions are not well defined, which results in higher order parameters, what would be a good solution for it? 1) reparameterization of PS lipids 2) reparameterization of PS-cationic ions interactions 3) possibly changing the cut-offs, say, 8-10 A. What do you think for this specific case (PS lipids)? Thanks On Wed, Jan 11, 2017 at 11:19 AM, Justin Lemkul wrote: > > > On 1/11/17 12:46 PM, Mohsen Ramezanpour wrote: > >> Hi Guys, >> >> Thanks for your comments. My question was exactly what Dawid clarified. >> Sure, I will read those as you suggested. >> >> Dawid, regarding this: >> "You need to keep in mind however that each force field was optimized for >> given set of cut-offs and >> they virtually become part of that FF (just like FF formula and >> parameters), so it is strongly recommended >> not to change them" >> >> I agree. BUT, using charmm36 in gromacs or other software requires to >> change some critical parameters like cut-offs to get a better agreements >> with experiments as well as with simulations done with charmm/NAMD. >> This is why. >> >> > As of version GROMACS 5.1, after which a small bug in force-switching was > fixed (see the link that Tom posted), this should not be true. The > settings listed in the GROMACS website for CHARMM36 should correctly > reproduce what people do in NAMD or CHARMM using the same cutoffs. We did > a lot of work to make sure this was the case. Cutoffs are a function of > the force field, not the software. > > For a pure membrane, there is a valid comparison to be made between > 0.8/1.2 and 1.0/1.2 switching as Tom has noted, but for a membrane protein > 1.0/1.2 is what you should use (and from data I have collected and seen, > the differences between the two switching ranges are small, in any case). > > -Justin > > > Cheers >> Mohsen >> >> On Wed, Jan 11, 2017 at 5:54 AM, Justin Lemkul wrote: >> >> >>> >>> On 1/10/17 10:40 PM, Mohsen Ramezanpour wrote: >>> >>> Dear gromacs users, >>>> >>>> Please let me know your opinion on the following question: >>>> Thanks in advance for your comments >>>> >>>> >>>> dx.doi.org/10.1002/jcc.21287 >>> >>> While that paper (obviously) focuses on CHARMM, there is a vast amount of >>> general MD information described in it, as well as cited literature. The >>> section on nonbonded interactions probably addresses most of what you >>> want >>> to know. >>> >>> -Justin >>> >>> >>> -- Forwarded message -- >>>> From: Mohsen Ramezanpour >>>> Date: Thu, Jan 5, 2017 at 5:20 PM >>>> Subject: LJ cut-offs >>>> To: Discussion list for GROMACS users >>>> >>>> >>>> Dear Gromacs users, >>>> >>>> Every force field has been parametrized with a specific LJ cut-off which >>>> must be the same for simulations using that force field. >>>> However, I was wondering if there is any reason why people usually take >>>> even numbers (e.g. 8-10, 8-12, 10-12 all with a difference of 2) for LJ >>>> cut-offs in force field development? >>>> Is there any rule that prohibit the use of 9-10, 9-11, or 11-12 for LJ >>>> cut-off? >>>> >>>> Thanks >>>> >>>> Cheers >>>> Mohsen >>>> >>>> >>>> >>>> >>>> >>>> >>>> -- >>> == >>> >>> Justin A. Lemkul, Ph.D. >>> Ruth L. Kirschstein NRSA Postdoctoral Fellow >>> >>> Department of Pharmaceutical Sciences >>> School of Pharmacy >>> Health Sciences Facility II, Room 629 >>> University of Maryland, Baltimore >>> 20 Penn St. >>> Baltimore, MD 21201 >>> >>> jalem...@outerbanks.umaryland.edu | (410) 706-7441 >>> http://mackerell.umaryland.edu/~jalemkul >>> >
Re: [gmx-users] Fwd: LJ cut-offs
I found these two links from other post of mine which might be useful for this discussion too: People have used 11-12 switching range as well. So, if it works, there is no restriction to not choose such numbers with 0.1 nm difference in cutoffs. http://pubs.acs.org/doi/abs/10.1021/jp101759q There is a great discussion which is relevant and better to read. https://github.com/NMRLipids/NmrLipidsCholXray/issues/4 Cheers Mohsen On Wed, Jan 11, 2017 at 10:46 AM, Mohsen Ramezanpour < ramezanpour.moh...@gmail.com> wrote: > Hi Guys, > > Thanks for your comments. My question was exactly what Dawid clarified. > Sure, I will read those as you suggested. > > Dawid, regarding this: > "You need to keep in mind however that each force field was optimized for > given set of cut-offs and > they virtually become part of that FF (just like FF formula and > parameters), so it is strongly recommended > not to change them" > > I agree. BUT, using charmm36 in gromacs or other software requires to > change some critical parameters like cut-offs to get a better agreements > with experiments as well as with simulations done with charmm/NAMD. > This is why. > > Cheers > Mohsen > > On Wed, Jan 11, 2017 at 5:54 AM, Justin Lemkul wrote: > >> >> >> On 1/10/17 10:40 PM, Mohsen Ramezanpour wrote: >> >>> Dear gromacs users, >>> >>> Please let me know your opinion on the following question: >>> Thanks in advance for your comments >>> >>> >> dx.doi.org/10.1002/jcc.21287 >> >> While that paper (obviously) focuses on CHARMM, there is a vast amount of >> general MD information described in it, as well as cited literature. The >> section on nonbonded interactions probably addresses most of what you want >> to know. >> >> -Justin >> >> >>> -- Forwarded message -- >>> From: Mohsen Ramezanpour >>> Date: Thu, Jan 5, 2017 at 5:20 PM >>> Subject: LJ cut-offs >>> To: Discussion list for GROMACS users >>> >>> >>> Dear Gromacs users, >>> >>> Every force field has been parametrized with a specific LJ cut-off which >>> must be the same for simulations using that force field. >>> However, I was wondering if there is any reason why people usually take >>> even numbers (e.g. 8-10, 8-12, 10-12 all with a difference of 2) for LJ >>> cut-offs in force field development? >>> Is there any rule that prohibit the use of 9-10, 9-11, or 11-12 for LJ >>> cut-off? >>> >>> Thanks >>> >>> Cheers >>> Mohsen >>> >>> >>> >>> >>> >>> >> -- >> == >> >> Justin A. Lemkul, Ph.D. >> Ruth L. Kirschstein NRSA Postdoctoral Fellow >> >> Department of Pharmaceutical Sciences >> School of Pharmacy >> Health Sciences Facility II, Room 629 >> University of Maryland, Baltimore >> 20 Penn St. >> Baltimore, MD 21201 >> >> jalem...@outerbanks.umaryland.edu | (410) 706-7441 >> http://mackerell.umaryland.edu/~jalemkul >> >> == >> >> -- >> Gromacs Users mailing list >> >> * Please search the archive at http://www.gromacs.org/Support >> /Mailing_Lists/GMX-Users_List before posting! >> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists >> >> * For (un)subscribe requests visit >> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or >> send a mail to gmx-users-requ...@gromacs.org. >> > > > > -- > *Rewards work better than punishment ...* > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] Fwd: LJ cut-offs
Hi Guys, Thanks for your comments. My question was exactly what Dawid clarified. Sure, I will read those as you suggested. Dawid, regarding this: "You need to keep in mind however that each force field was optimized for given set of cut-offs and they virtually become part of that FF (just like FF formula and parameters), so it is strongly recommended not to change them" I agree. BUT, using charmm36 in gromacs or other software requires to change some critical parameters like cut-offs to get a better agreements with experiments as well as with simulations done with charmm/NAMD. This is why. Cheers Mohsen On Wed, Jan 11, 2017 at 5:54 AM, Justin Lemkul wrote: > > > On 1/10/17 10:40 PM, Mohsen Ramezanpour wrote: > >> Dear gromacs users, >> >> Please let me know your opinion on the following question: >> Thanks in advance for your comments >> >> > dx.doi.org/10.1002/jcc.21287 > > While that paper (obviously) focuses on CHARMM, there is a vast amount of > general MD information described in it, as well as cited literature. The > section on nonbonded interactions probably addresses most of what you want > to know. > > -Justin > > >> -- Forwarded message -- >> From: Mohsen Ramezanpour >> Date: Thu, Jan 5, 2017 at 5:20 PM >> Subject: LJ cut-offs >> To: Discussion list for GROMACS users >> >> >> Dear Gromacs users, >> >> Every force field has been parametrized with a specific LJ cut-off which >> must be the same for simulations using that force field. >> However, I was wondering if there is any reason why people usually take >> even numbers (e.g. 8-10, 8-12, 10-12 all with a difference of 2) for LJ >> cut-offs in force field development? >> Is there any rule that prohibit the use of 9-10, 9-11, or 11-12 for LJ >> cut-off? >> >> Thanks >> >> Cheers >> Mohsen >> >> >> >> >> >> > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] rlist, rcolumb, and rvdw
Dear gromacs users, Reading through mailing list I found a nice discussion on the relation between rlist, rcolumb, and rvdw: https://www.mail-archive.com/gmx-users@gromacs.org/msg08387.html https://mailman-1.sys.kth.se/pipermail/gromacs.org_gmx-users/2011-February/058507.html If I understood correctly, rcoulomb < rlist and rvdw < rlist is the most accurate way, and rlist=rcoulomb=rvdw is the commonly used way for these parameters. I was wondering what would be the case if I want to change the cut-offs? For instance, to use charmm36 in Groamcs, it is recommended to use 1.0-1.2 as the cut-off for LJ, and rlist=rcoulomb=rvdw However, this setting might not result in good behaviour for some lipid bilayers. Therefore, I want to check if other cut-offs works better in gromacs. If I change "1.2 nm" to "A nm" for rvdw, do I need to change rcolumb and rlist to A as well? i.e. rlist=rcolumb=rvdw=A? I am using following parameters already: cutoff-scheme = Verlet nstlist = 20 rlist = *1.2* coulombtype = *pme* rcoulomb = *1.2* vdwtype = Cut-off vdw-modifier= *Force-switch* rvdw_switch = *1.0* rvdw= *1.2* Thanks in advance for your comments Cheers Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] Fwd: LJ cut-offs
Dear gromacs users, Please let me know your opinion on the following question: Thanks in advance for your comments -- Forwarded message -- From: Mohsen Ramezanpour Date: Thu, Jan 5, 2017 at 5:20 PM Subject: LJ cut-offs To: Discussion list for GROMACS users Dear Gromacs users, Every force field has been parametrized with a specific LJ cut-off which must be the same for simulations using that force field. However, I was wondering if there is any reason why people usually take even numbers (e.g. 8-10, 8-12, 10-12 all with a difference of 2) for LJ cut-offs in force field development? Is there any rule that prohibit the use of 9-10, 9-11, or 11-12 for LJ cut-off? Thanks Cheers Mohsen -- *Rewards work better than punishment ...* -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] LJ cut-offs
Dear Gromacs users, Every force field has been parametrized with a specific LJ cut-off which must be the same for simulations using that force field. However, I was wondering if there is any reason why people usually take even numbers (e.g. 8-10, 8-12, 10-12 all with a difference of 2) for LJ cut-offs in force field development? Is there any rule that prohibit the use of 9-10, 9-11, or 11-12 for LJ cut-off? Thanks Cheers Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] adding new parameters to force field
Thanks Mark, :-) Maybe my question was not clear. I have a molecule and I am trying to parametrize it. I have assigned atom types to each atom and now trying to find the approperiate bonded and nonbonded parameters for it. Based on these atom types, there are parameters for it in FF already. However, the parameters (especially, the dihedrals) do not result in reasonable agreement with exp data. So, I had to optimize the parameters for a few specific dihedrals (keeping the same atom types). Now, I want to add new parameters into the force field. This is why I was curious to see if grompp will read the new inserted parameters or the original ones. In gromos force field, each bond, angle, dihedral has a code (e.g. gb_1, ga_1, and gd_1), thus, I can specifically define which one to use for my molecule. However, in Charmm36, there is not any similar code. So, I was not sure which parameter grompp will uses. One solution was to add new optimized parameters in corresponding line in molecule.itp file. (there are a lot of lines and I rather to included lines in ff than molecule.itp file) An alternative was to define new atom types with similar properties but new dihedral parameters compared to original ones. I checked it yesterday and interestingly, grompp gives warnings regarding duplicates, which one it takes, and if it overrides original parameters or not. Cheers Mohsen On Wed, Dec 21, 2016 at 11:06 PM, Mark Abraham wrote: > Hi, > > Sometimes. It depends exactly which functions and function types, and > whether it's part of a moleculetype. If you don't want something, don't > include it :-) > > Mark > > On Thu, 22 Dec 2016 10:07 Mohsen Ramezanpour > > wrote: > > > Dear Gromacs users, > > > > I have a new file with both bonded and nonbonded parameters for some atom > > types in it. I want to use these new parameters for simulation and ignore > > the parameters in force field ( if there is any parameter already exist). > > > > 1) If I add the lines from each section of new file at the end of > > corresponding section at ffbonded.itp and ffnonbonded.itp files, do I > still > > need to comment out all the present parameters or they will be all > > overwritten by new defined parameters? > > > > 2) How if I just #include new file after inclusion of normal force field? > > Will it overwrite previous parameters? > > > > Reading gromacs forum discussions and manual, I still was not sure of > that > > (I am using gromacs versions 4 and 5). > > > > > > *From manual 5.8.3.* > > *Adding atom types* > > > > As of GROMACS version 3.1.3, atom types can be added in an extra [ > > atomtypes ] section > > after the the inclusion of the normal force field. After the definition > of > > the *new atom type(s)*, ad- > > ditional non-bonded and pair parameters can be defined. In pre-3.1.3 > > versions of GROMACS, the > > new atom types needed to be added in the [ atomtypes ] section of the > force > > field files, be- > > cause all non-bonded parameters above the last [ atomtypes ] section > would > > be overwritten > > using the standard combination rules. > > > > > > Best, > > Mohsen > > > > > > > > -- > > *Rewards work better than punishment ...* > > -- > > Gromacs Users mailing list > > > > * Please search the archive at > > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before > > posting! > > > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > > > * For (un)subscribe requests visit > > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > > send a mail to gmx-users-requ...@gromacs.org. > > > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/ > Support/Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
[gmx-users] adding new parameters to force field
Dear Gromacs users, I have a new file with both bonded and nonbonded parameters for some atom types in it. I want to use these new parameters for simulation and ignore the parameters in force field ( if there is any parameter already exist). 1) If I add the lines from each section of new file at the end of corresponding section at ffbonded.itp and ffnonbonded.itp files, do I still need to comment out all the present parameters or they will be all overwritten by new defined parameters? 2) How if I just #include new file after inclusion of normal force field? Will it overwrite previous parameters? Reading gromacs forum discussions and manual, I still was not sure of that (I am using gromacs versions 4 and 5). *From manual 5.8.3.* *Adding atom types* As of GROMACS version 3.1.3, atom types can be added in an extra [ atomtypes ] section after the the inclusion of the normal force field. After the definition of the *new atom type(s)*, ad- ditional non-bonded and pair parameters can be defined. In pre-3.1.3 versions of GROMACS, the new atom types needed to be added in the [ atomtypes ] section of the force field files, be- cause all non-bonded parameters above the last [ atomtypes ] section would be overwritten using the standard combination rules. Best, Mohsen -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
Re: [gmx-users] using implicit water for atomistic simulation of lipid bilayer
Thanks Justin. Agreed. I wanted to make sure about that. Cheers On Sat, Dec 10, 2016 at 12:48 PM, Justin Lemkul wrote: > > > On 12/9/16 5:58 PM, Mohsen Ramezanpour wrote: > >> Dear Gromacs Users, >> >> I was wondering if there is any example of using Implicit water model for >> lipid bilayers in Gromacs? >> >> I have read on using implicit solvents in Gromacs and found interesting >> discussions: >> >> https://groups.google.com/forum/#!topic/archive-gmx-users/0HIFVop390I >> >> https://mailman-1.sys.kth.se/pipermail/gromacs.org_gmx-users >> /2011-December/066830.html >> >> >> As I understood, using implicit solvent assumes an infinite solvent around >> the solute. >> >> Which seems good for peptide/protein folding. I did not see any example >> for >> bilayers. >> >> Assuming that implicit solvent is well parameterized to work with lipid >> membranes, the one problem seems to be the PBC option in mdp file. >> >> PBC = no is appropriate for using implicit solvent while PBC = xyz is the >> usual case for bilayers. >> >> Is there any way to split the PBC treatment for water and bilayer? Like >> what we do in T-coupling for different groups. >> >> Also, I want to mix this implicit solvent with atomistic force field, like >> Gromos FF, not MARTINI. >> >> > I don't think GROMACS (or any other program) can handle this. Perhaps you > can get CHARMM to do it via the GBMW module, but even that (I'm pretty > sure) requires an implicit bilayer (e.g. layers of high and low dielectric > constant only, no mixture of explicit + implicit). > > -Justin > > -- > == > > Justin A. Lemkul, Ph.D. > Ruth L. Kirschstein NRSA Postdoctoral Fellow > > Department of Pharmaceutical Sciences > School of Pharmacy > Health Sciences Facility II, Room 629 > University of Maryland, Baltimore > 20 Penn St. > Baltimore, MD 21201 > > jalem...@outerbanks.umaryland.edu | (410) 706-7441 > http://mackerell.umaryland.edu/~jalemkul > > == > -- > Gromacs Users mailing list > > * Please search the archive at http://www.gromacs.org/Support > /Mailing_Lists/GMX-Users_List before posting! > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists > > * For (un)subscribe requests visit > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or > send a mail to gmx-users-requ...@gromacs.org. > -- *Rewards work better than punishment ...* -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
