[gmx-users] RDF calculation.
Hi, I am studying the dynamics of a double stranded DNA in presence of small molecules. I want to study the preferential binding of these molecules in the grooves and backbones of the DNA and also calculate the feasibility of these binding in terms of calculating the free energy change.My question is- 1. while calculating radial distribution functions say between DNA minor groove and choline, should I consider the center-off-mass atom of the minor groove by using the -COM option of g_rdf? 2. Do you think calculation of PMF from g(r) for a specific ion binding with the DNA can provide reliable free energy change or do I have to switch to umbrella sampling? PS: I only want to investigate the small molecule binding, not the stability of the DNA duplex as a whole. Thanks for your time in advance. Soumadwip Ghosh Research Fellow Indian Institute of Technology Bombay India -- 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] RDF calculation.
Well, since you are not interested in the duplex stability, I would not use -COM. See, during the simulation, the center-of-mass of such a small helix part can drastically move, leading to a questionable result. Check your system and make sure that you will not face this situation. About the PMF, that depends on the ions you are using and the forcefields. It is already known that some of them are poorly described, so check that out too. Hope it can help! Cheers! 2014-12-27 7:17 GMT-02:00 soumadwip ghosh : > Hi, > I am studying the dynamics of a double stranded DNA in presence of > small molecules. I want to study the preferential binding of these > molecules in the grooves and backbones of the DNA and also calculate the > feasibility of these binding in terms of calculating the free energy > change.My question is- > > 1. while calculating radial distribution functions say between DNA minor > groove and choline, should I consider the center-off-mass atom of the minor > groove by using the -COM option of g_rdf? > > 2. Do you think calculation of PMF from g(r) for a specific ion binding > with the DNA can provide reliable free energy change or do I have to switch > to umbrella sampling? > > PS: I only want to investigate the small molecule binding, not the > stability of the DNA duplex as a whole. > > Thanks for your time in advance. > > Soumadwip Ghosh > Research Fellow > Indian Institute of Technology Bombay > India > -- > 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. > -- Marcelo Depólo Polêto Biochemicist University of Viçosa - Brazil -- 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] length of beta sheet
On 12/25/14 1:43 AM, pratibha kapoor wrote:
Hi everyone,
I would like to calculate length of beta sheet. Is this any automatic tool
to serve the purpose?
To elaborate, I am aware that dssp can give me whether the residue is in
beta sheet ("E") or not at each time frame but I want to know exactly how
many hydrogen bonds are formed by each residue in E (and not just E if
minimum criterion of two hydrogen bonds is satisfied by DSSP). From this I
would like to calculate length of beta sheet.
Any help is highly appreciated.
You can extract this sort of information from the .xpm file that do_dssp writes
out. Just parse it with whatever scripting language you like.
-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.
Re: [gmx-users] Pulling
On 12/25/14 2:49 PM, asasa qsqs wrote: Dear Justin A. Lemkul, I want to use pulling code for calculate the PMF of a molecule at the across the DMPC, but i can't obtained a large pull across the membrane my pulling code was: ; Pull code pull= umbrella pull_geometry= position pull_dim= N N Y pull_start = yes; define initial COM distance > 0 pull_ngroups= 1 pull_group0= DMPC pull_group1= QBH pull_k1= 1000 ; kJ mol^-1 nm^-2 pull_vec1= 0 0 0.5 after "perl distances.pl" step i obtained a summary.dat file. but in this file pulling was very low, just about 0.3nm. My molecule name is QBH and is out of membrane at bulk of the water and i want to pull it until center of bilayer. what must i do? In principle, that should work, but you haven't said how long you conduct the pulling or shown what the distance time series is doing. Is QBH approaching the DMPC membrane over time, or just staying at some distance? -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.
Re: [gmx-users] How to choose the best pose of ligand from the trajectory file?
On 12/25/14 5:03 PM, Hassan Aaryapour wrote: Dear Gromacs Users How can I find and score protein-ligand binding poses and choose the best pose or frame from the trajectory file after doing a Molecular dynamics simulation of 40ns? "Best" is a qualitative assessment that should be driven by objective, quantitative metrics. Decide how you want/need to analyze the interactions of the ligand with the protein, and the pose that is "best" is probably the one that captures most/all of the principle features (contacts, hydrogen bonds, dihedral angles, etc). -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.
Re: [gmx-users] Umbrella sampling along a dihedral angle
On 12/26/14 6:11 PM, Nash, Anthony wrote: Dear Gromacs community, Using enforced rotation potentials I have generated a really smooth rotation of my ligand the catalytic binding domain of my protein. I then put together a simple perl script to calculate the dihedral angle of four particular atoms at each time step using g_angle. Taking a 6 degree interval (for a total of 60 windows to span the 360 degree rotation) I modified the topology with a [ dihedral_restraint ] with the respective dihedral angle for that particular configuation. I plan on then running each window for 10 ns to begin with (according to the PMF profile I can add more time). It is from this point where I am a little stuck. I have gone through a number of posts going back to 2011 with regards to using umbrella sampling along a dihedral angle reaction coordinate. The precise post was from Justin: "You'll have to build various configurations that correspond to different dihedral angles (which form the sampling windows), then restrain them. The energy attributed to the restraints is then stored in the .edr file. From these energies, you should be able to construct the energy curve over the sampling windows. There are examples of this in the literature, so I suspect you should be able to find some demonstrations of how it's applied." I was wondering whether gromacs now support this procedure natively? I noticed g_wham has an option -cycl for dihedrals but I am still not sure how this fits in with the traditional gromacs means of generating pull force output as I have only been using distance based rather than angular based umbrella sampling. g_wham is still only really designed to process the output from the pull code, not any generalized restraint. That would be a nice features, but AFAIK it is not in the pipeline. You can de-bias the results like I suggested above, but it's probably more work than it's worth. I would suggest just using Alan Grossfield's WHAM implementation. There, you just need to provide the raw numbers (in this case, the dihedral values) and the force constants used. It will then give you your PMF. We have done this recently; it is very easy. -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.
Re: [gmx-users] How to plot a droplet radius profile as a function of z-direction?
On 12/27/14 2:05 AM, Kester Wong wrote: Dear gmx-users, Q1) I would like to plot the radius profile of a water droplet on graphene, as a function of z-direction, as per below: http://journals.aps.org/prl/article/10.1103/PhysRevLett.109.176101/figures/2/large Which GROAMCS utility should I use, and how do I go about plotting the profile? It will require several steps. Not something I've ever done, but I would think g_traj for coordinate extraction, then post-processing to calculate the diameter or radius, would be the first step. Presumably you would have to do this for various slices (assuming you're looking at cross-sections in the x-y plane for varying values of z), so g_select to define atoms within that slice will be needed before even running g_traj. Then you construct the profile from this information. Q2) Another question that I have is, suppose I have equilibrated (NVT) a water droplet for 10ns, and would like to plot the H-bond density (along the z-direction), I type the following command: g_density_mpi -f traj_comp.xtc -n system.ndx -s topol.tpr -b 9000 -e 1 -sl 1000 -dens number -o number-density-OH.xvg I was only able to select: System Graphene (my substrate) Water How can I create a system.ndx such that it allows me to choose the OH (instead of water)? make_ndx or g_select create index files. If you're defining "hydrogen bond density" as "number of hydrogen bonds per slice," along with the above, you need g_select and then you need g_hbond. g_density is not the right program here. Q3) Regarding Q2, the g_density_mpi utility read the frames contained within the time frame of 9 - 10 ns, but does it average the frames to give a time-averaged density.xvg? Yes. But again, g_density is probably not what you want here; see the help description. -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.
Re: [gmx-users] RDF calculation.
On 12/27/14 4:17 AM, soumadwip ghosh wrote: Hi, I am studying the dynamics of a double stranded DNA in presence of small molecules. I want to study the preferential binding of these molecules in the grooves and backbones of the DNA and also calculate the feasibility of these binding in terms of calculating the free energy change.My question is- 1. while calculating radial distribution functions say between DNA minor groove and choline, should I consider the center-off-mass atom of the minor groove by using the -COM option of g_rdf? I'm not sure it's useful. Choline is a rather flexible molecule, so simply considering the COM of the minor groove vs. all atoms of choline is probably going to give you a very noisy profile. 2. Do you think calculation of PMF from g(r) for a specific ion binding with the DNA can provide reliable free energy change or do I have to switch to umbrella sampling? This doesn't generally work well for ions, but the quality depends on the force field. Depending on what the ion is, the exchange rate may be impossible to observe in an atomistic simulation. The exchange rates for most common ions can range from nanosecond (feasible with MD) to microsecond (not reasonable, even if the force field is very accurate). -Justin PS: I only want to investigate the small molecule binding, not the stability of the DNA duplex as a whole. Thanks for your time in advance. Soumadwip Ghosh Research Fellow Indian Institute of Technology Bombay India -- == 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.
Re: [gmx-users] Umbrella sampling along a dihedral angle
Thanks for the help Justin. I've found the site: http://membrane.urmc.rochester.edu/content/wham Cheers Anthony Dr Anthony Nash Department of Chemistry University College London 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: 27 December 2014 14:49 To: gmx-us...@gromacs.org Subject: Re: [gmx-users] Umbrella sampling along a dihedral angle On 12/26/14 6:11 PM, Nash, Anthony wrote: > Dear Gromacs community, > > Using enforced rotation potentials I have generated a really smooth rotation > of my ligand the catalytic binding domain of my protein. I then put together > a simple perl script to calculate the dihedral angle of four particular atoms > at each time step using g_angle. Taking a 6 degree interval (for a total of > 60 windows to span the 360 degree rotation) I modified the topology with a [ > dihedral_restraint ] with the respective dihedral angle for that particular > configuation. I plan on then running each window for 10 ns to begin with > (according to the PMF profile I can add more time). It is from this point > where I am a little stuck. > > I have gone through a number of posts going back to 2011 with regards to > using umbrella sampling along a dihedral angle reaction coordinate. The > precise post was from Justin: > > "You'll have to build various configurations that correspond to different > dihedral angles (which form the sampling windows), then restrain them. > > The energy attributed to the restraints is then stored in the .edr file. From > these energies, you should be able to construct the energy curve over the > sampling windows. There are examples of this in the literature, so I suspect > you should be able to find some demonstrations of how it's applied." > > I was wondering whether gromacs now support this procedure natively? I > noticed g_wham has an option -cycl for dihedrals but I am still not sure how > this fits in with the traditional gromacs means of generating pull force > output as I have only been using distance based rather than angular based > umbrella sampling. > g_wham is still only really designed to process the output from the pull code, not any generalized restraint. That would be a nice features, but AFAIK it is not in the pipeline. You can de-bias the results like I suggested above, but it's probably more work than it's worth. I would suggest just using Alan Grossfield's WHAM implementation. There, you just need to provide the raw numbers (in this case, the dihedral values) and the force constants used. It will then give you your PMF. We have done this recently; it is very easy. -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. -- 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] core dumped
Dear users when i simulated Aspirin on Human serum albumin (pdb: 2i30) in equilibration process (NVT), i encountered with below message and my simulation stopped. i dont know why? *pdfco@pdfco-X:~/Aspirin$ tail nvt.job-nex int0 Number of random exchanges to carry out eachexchange interval (N^3 is one suggestion). -nexzero or not specified gives neighbor replica exchange.-reseed int-1 Seed for replica exchange, -1 is generate a seed-[no]ionize bool no Do a simulation including the effect of an X-Raybombardment on your system* *Reading file nvt.tpr, VERSION 4.6.5 (single precision)Using 1 MPI thread[1]+ Segmentation fault (core dumped) mdrun -v -deffnm nvt &> nvt.job* my *NVT.mdp* content was: title = Protein NVT equilibration define = -DPOSRES ; position restrain the protein ; Run parameters integrator = md ; leap-frog integrator nsteps = 10 ; 0.002 * 10 = 200 ps dt = 0.002 ; 2 fs ; Output control nstxout = 1000 ; save coordinates every 0.2 ps nstvout = 1000 ; save velocities every 0.2 ps nstenergy = 1000 ; save energies every 0.2 ps nstlog = 1000 ; update log file every 0.2 ps ; Bond parameters continuation = no ; first dynamics run constraint_algorithm = lincs ; holonomic constraints constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained lincs_iter = 1 ; accuracy of LINCS lincs_order = 4 ; also related to accuracy ; Neighborsearching ns_type = grid ; search neighboring grid cells nstlist = 5 ; 10 fs rlist = 1.0 ; short-range neighborlist cutoff (in nm) rcoulomb = 1.0 ; short-range electrostatic cutoff (in nm) rvdw = 1.0 ; short-range van der Waals cutoff (in nm) ; Electrostatics coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics pme_order = 4 ; cubic interpolation fourierspacing = 0.16 ; grid spacing for FFT ; Temperature coupling is on tcoupl = v-rescale tc-grps = Protein_SAL SOL_CL ; two coupling groups - more accurate tau_t = 0.1 0.1 ; time constant, in ps ref_t = 310 310 ; reference temperature, one for each group, in K ; Pressure coupling is off pcoupl = no ; no pressure coupling in NVT ; Periodic boundary conditions pbc = XYZ ; Dispersion correction DispCorr = EnerPres ; account for cut-off vdW scheme ; Velocity generation gen_vel = yes ; assign velocities from Maxwell distribution gen_temp = 300 ; temperature for Maxwell distribution gen_seed = -1 ; generate a random seed i really appriciate for your hint regards -- 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] core dumped
On 12/27/14 11:42 AM, elham tazikeh wrote: Dear users when i simulated Aspirin on Human serum albumin (pdb: 2i30) in equilibration process (NVT), i encountered with below message and my simulation stopped. i dont know why? http://www.gromacs.org/Documentation/Terminology/Blowing_Up A seg fault leaves nothing to diagnose. Follow the advice in that link. -Justin *pdfco@pdfco-X:~/Aspirin$ tail nvt.job-nex int0 Number of random exchanges to carry out eachexchange interval (N^3 is one suggestion). -nexzero or not specified gives neighbor replica exchange.-reseed int-1 Seed for replica exchange, -1 is generate a seed-[no]ionize bool no Do a simulation including the effect of an X-Raybombardment on your system* *Reading file nvt.tpr, VERSION 4.6.5 (single precision)Using 1 MPI thread[1]+ Segmentation fault (core dumped) mdrun -v -deffnm nvt &> nvt.job* my *NVT.mdp* content was: title = Protein NVT equilibration define = -DPOSRES ; position restrain the protein ; Run parameters integrator = md ; leap-frog integrator nsteps = 10 ; 0.002 * 10 = 200 ps dt = 0.002 ; 2 fs ; Output control nstxout = 1000 ; save coordinates every 0.2 ps nstvout = 1000 ; save velocities every 0.2 ps nstenergy = 1000 ; save energies every 0.2 ps nstlog = 1000 ; update log file every 0.2 ps ; Bond parameters continuation = no ; first dynamics run constraint_algorithm = lincs ; holonomic constraints constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained lincs_iter = 1 ; accuracy of LINCS lincs_order = 4 ; also related to accuracy ; Neighborsearching ns_type = grid ; search neighboring grid cells nstlist = 5 ; 10 fs rlist = 1.0 ; short-range neighborlist cutoff (in nm) rcoulomb = 1.0 ; short-range electrostatic cutoff (in nm) rvdw = 1.0 ; short-range van der Waals cutoff (in nm) ; Electrostatics coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics pme_order = 4 ; cubic interpolation fourierspacing = 0.16 ; grid spacing for FFT ; Temperature coupling is on tcoupl = v-rescale tc-grps = Protein_SAL SOL_CL ; two coupling groups - more accurate tau_t = 0.1 0.1 ; time constant, in ps ref_t = 310 310 ; reference temperature, one for each group, in K ; Pressure coupling is off pcoupl = no ; no pressure coupling in NVT ; Periodic boundary conditions pbc = XYZ ; Dispersion correction DispCorr = EnerPres ; account for cut-off vdW scheme ; Velocity generation gen_vel = yes ; assign velocities from Maxwell distribution gen_temp = 300 ; temperature for Maxwell distribution gen_seed = -1 ; generate a random seed i really appriciate for your hint regards -- == 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.
Re: [gmx-users] Umbrella sampling along a dihedral angle
Hi Justin, I think I've shot myself in the foot and I'm trawling the gromacs manual with little success at the moment. I'm trying to define my [ dihedral_restraint ] over my ligand-enzyme complex. The topology of the ligand is in one .itp file and the topology of the enzyme is in a separate .itp file. Therefore, in my top level topology file (system.top) I have: #include "../enzyme.itp" #include "../ligand.itp" Typically if I want restraints to help equilibrate the system I would just define a posre_NAME.itp with the respective force constants beneath each #include directive. But given that each topology file starts with atom 1 I am a bit unsure how I can define a set of restraints over two topology files. Can this be done? I'm really hoping I won't need to clump all topologies into one big file. Thanks as ever, Anthony Dr Anthony Nash Department of Chemistry University College London 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: 27 December 2014 14:49 To: gmx-us...@gromacs.org Subject: Re: [gmx-users] Umbrella sampling along a dihedral angle On 12/26/14 6:11 PM, Nash, Anthony wrote: > Dear Gromacs community, > > Using enforced rotation potentials I have generated a really smooth rotation > of my ligand the catalytic binding domain of my protein. I then put together > a simple perl script to calculate the dihedral angle of four particular atoms > at each time step using g_angle. Taking a 6 degree interval (for a total of > 60 windows to span the 360 degree rotation) I modified the topology with a [ > dihedral_restraint ] with the respective dihedral angle for that particular > configuation. I plan on then running each window for 10 ns to begin with > (according to the PMF profile I can add more time). It is from this point > where I am a little stuck. > > I have gone through a number of posts going back to 2011 with regards to > using umbrella sampling along a dihedral angle reaction coordinate. The > precise post was from Justin: > > "You'll have to build various configurations that correspond to different > dihedral angles (which form the sampling windows), then restrain them. > > The energy attributed to the restraints is then stored in the .edr file. From > these energies, you should be able to construct the energy curve over the > sampling windows. There are examples of this in the literature, so I suspect > you should be able to find some demonstrations of how it's applied." > > I was wondering whether gromacs now support this procedure natively? I > noticed g_wham has an option -cycl for dihedrals but I am still not sure how > this fits in with the traditional gromacs means of generating pull force > output as I have only been using distance based rather than angular based > umbrella sampling. > g_wham is still only really designed to process the output from the pull code, not any generalized restraint. That would be a nice features, but AFAIK it is not in the pipeline. You can de-bias the results like I suggested above, but it's probably more work than it's worth. I would suggest just using Alan Grossfield's WHAM implementation. There, you just need to provide the raw numbers (in this case, the dihedral values) and the force constants used. It will then give you your PMF. We have done this recently; it is very easy. -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. -- 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] Error during entropy calculation
Hi all, I am trying to calculate entropy using g_anaeig (Gromacs version 5.0.1), after using g_nmeig to calculate eigenvectors and eigenvalues (after computing hessian matrix): g_nmeig_mpi -f hessian.mtx -s run.tpr -nom g_anaeig_mpi -v eigenvec.trr -entropy But I am getting the following error: Fatal error: Can not calculate entropies from mass-weighted eigenvalues, redo the analysis without mass-weighting. -- Thanks and Regards, Bipin Singh -- 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] charge distribution
Dear gromacs users i d like to know about charge distribution in .itp file that how can fix the charge of atoms for instance, i simulated Aspirin in Human serum albumin, and i want to know about their charges in an article charge distribution was as below, but in my .itp file (aquired by PRODRG server) i had 11 atoms instead of 13 atoms please help me how can fix partial charges? thanks for any hint Atom numberAtom TypeNameCharge groupPartial ChargeMass -- 1OMOL1−0.635 15.99942CCK10.27 12.0113OMOM1−0.635 15.99944CCJ20 12.0115CR1CI 20 13.0196CR1CH20 13.0197CR1CG30 13.0198CR1CF30 13.0199CCE30 12.01110OAOD40 15.999411CCB40.38 12.01112OOC4−0.38 15.999413CH3CA50 15.035 -- 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] core dumped
Dear Justin thanks for your advise only my remained question is, if my problem was for my topology file (as i recieved this error a few times: your topology file has non-zero charge = -0.61) by editing my topology file (.itp file), can i going on my simulation? please tell me, what is the best way for distribute the charges on a molecule (Aspirin)? thanks for your help -- 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.
