hello sorry for writing again and again as from the previous parameters, now I had removed the -ve sign from the pull_coord1_vec =0 0 -1 and pull_coord1_rate = -0.002 as follow
; Pull code pull = umbrella pull_geometry = direction pull_dim = N N Y pull_coord1_vec =0 0 1 pull_start = yes pull_ngroups =2 pull-ncoords =1 pull_coord1_groups = 1 2 pull_group1_name = NPROT_ref pull_group2_name = LIG pull_coord1_rate = 0.002 ; 0.002 nm per ps = 2 nm per ns pull_coord1_k = 1000 ; kJ mol^-1 nm^-2 pull_nstxout = 50 ; every 1 ps pull_nstfout = 50 ; every 1 ps Now I think simulation and the results looking correct and no strange behavior was observed. So can you kindly comments on the setting which I have mentioned above and their correctness. (there may be chances that I may have understood the parameters wrongly), and from this simulation I want to do the umbrella sampling for free energy calculation. with regards Mahender Singh > Hi > Christopher Neale and Justin > > I used now following pull code (what i understood from the previous emails) > to pull the drug across the membrane > ; Pull code > pull = umbrella > pull_geometry = direction > pull_dim = N N Y > pull_coord1_vec =0 0 -1 > pull_start = yes > pull_ngroups =2 > pull-ncoords =1 > pull_coord1_groups = 1 2 > pull_group1_name = NPROT_ref > pull_group2_name = LIG > pull_coord1_rate = -0.002 ; 0.002 nm per ps = 2 nm per ns > pull_coord1_k = 1000 ; kJ mol^-1 nm^-2 > pull_nstxout = 50 ; every 1 ps > pull_nstfout = 50 ; every 1 ps > > there was no error and the molecule was able to cross the membrane from one > side to the other side in a single SMD. > But my concern here is the pullf.xvg output, because in the graph all the > force is in negative sign. but if I consider the magnitude of the force, it's > giving me the expected graph as there will be increase in the force when drug > will enter in the membrane. Can you kindly comment of the statement. > > with regards > Mahender Singh > > > > Thank you dear Justin and Christopher Neale for your help. > > > > Solute is at -z relative to the bilayer and wants to pull toward +z from > > the bulk water on the -Z to the +z bulk water through the membrane. > > I observed that the + and the - sign in the pull_coord1_rate is associated > > with increase or decrease in the COM distance between the pull group and > > the reference group, respectively. > > I am trying to split simulation into two parts i.e. each leaflet as > > separate SMD simulation. > > When solute is moving from the -z to the center of the bilayer, it will > > stop here as pull_coord1_rate = -0.0020 (decrease in the COM, pull force > > will be - in sign, am I right that sign in force values is negative because > > decrease in the distance between COM of two group with time is happening), > > still in the last frame, COM of the molecule will not moves to the +z side > > and I am not able to do the second part of the simulation by setting > > pull_start=yes , as if I will change the pull_coord1_rate = 0.0020 it will > > pull it back in the -z direction (pull force will be in positive sign) > > (pull_geometry= distance). > > > > @Christopher Neale > > I will try the following setting and will let you know the results. > > > > @ Justin > > When I am pulling the molecule from the center of bilayer to the bulk water > > and the same molecule from the bulk water to the center of bilayer (by > > changing the sign in the pull_coord1_rate = -/+0.0020 ). pull-force v/s > > time graphs are looking like following (first graph is the pulling of > > molecule from water to the center of bilayer and second graph is from the > > center of bilayer to the bulk water). My query was the - and the + sign of > > the pull force, which I think is due to the increase and decrease in the > > COM distance between the pull group and the reference group, am I right? > > > > > > > > with regards > > Mahender singh > > > > > > > Dear Justin: > > > > > > I think you are correct. It may not even really be that much of a hassle > > > though one may need to be very careful in system setup for umbrellas near > > > the center of the bilayer when using pull_start=yes rather than > > > pull_coord1_init (i.e., one may need to make sure that the solute COM is > > > on the side of the bilayer that one specifies with pull_coord1_vec, > > > though having not used it myself perhaps this is not even an issue?). > > > Separately, I presume that pull_coord1_init=0 should be identical with > > > pull_coord1_vec = 0 0 1 or pull_coord1_vec = 0 0 -1, though if I was > > > using the new code that would be something I would check at the outset. > > > > > > Dear Mahender: > > > > > > Based on Justin's suggestion and your initial mdp file, you should try > > > something like the following. I'm still using gromacs 4.6.7 so I can not > > > verify the accuracy of your other parameters (or even be entirely sure > > > about these ones, but they are worth a try). > > > > > > ;; Presuming that you start with the solute at +z relative to the bilayer > > > and want to pull toward -z > > > pull = umbrella > > > pull_geometry = direction > > > pull_dim = N N Y > > > pull_coord1_vec = 0 0 1 > > > pull_coord1_rate = -0.002 > > > pull_start = yes > > > > > > ;; I presume that > > > ;; pull_coord1_vec = 0 0 1 and pull_coord1_rate = -0.002 > > > ;; is identical to > > > ;; pull_coord1_vec = 0 0 -1 and pull_coord1_rate = 0.002 > > > > > > ;; I am not sure if pull_dim = N N Y is necessary, but I don't see how it > > > could hurt. > > > ;; Presumably the distance and the vector on which the pull_coord1_rate > > > acts are already going to be entirely along z due to pull_coord1_vec , > > > but I'd test that too if not using pull_dim = N N Y > > > > > > Please report back whether this works or whether you get strange > > > behaviour when you hit the bilayer center. > > > > > > ________________________________________ > > > From: [email protected] > > > <[email protected]> on behalf of Justin > > > Lemkul <[email protected]> > > > Sent: 22 April 2015 08:12 > > > To: [email protected] > > > Subject: Re: [gmx-users] Steered MD simulation of drug molecules across > > > the plasma membrane in gromacs 5.0.4 > > > > > > On 4/21/15 6:54 AM, Christopher Neale wrote: > > > > Dear Mahender: > > > > > > > > It's a real pity that the pull_geometry = position has been removed. So > > > > now it's impossible to do umbrella sampling of a solute across a lipid > > > > bilayer properly? Anyway, I see a note about this removal here ( > > > > http://redmine.gromacs.org/issues/1346 ), though no reason is given. > > > > I'd suggest that you revert to a more functional version of gromacs ;) > > > > > > > > Hopefully I'm missing something and somebody else can point you to a > > > > better solution. > > > > > > > > > > I think the behavior can be recovered with "direction" geometry, but that > > > requires each half of the symmetric sampling windows to be set up > > > separately, I > > > think, rather than just a continuous vector like "position" used to > > > provide. > > > That's a shame if true. It's a workaround, but it's much less convenient. > > > > > > -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 > > > > > > [email protected] | (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 [email protected]. > > > > ------------------------------ > > Message: 4 > Date: Thu, 23 Apr 2015 01:18:01 -0600 > From: Alex <[email protected]> > To: "[email protected]" > <[email protected]>, <[email protected]> > Subject: Re: [gmx-users] Problems in the modelling of interaction > between peptide and copper > Message-ID: <[email protected]> > Content-Type: text/plain; charset=us-ascii > > Hi, > > Your model is not reasonable. First of all, let's forget about > simulations for a moment, because GMX is a wildly inappropriate choice. > Even if the energy was not zero (which in your case is probably due to > some LJ parameters set up improperly, or some really good GMX code > giving you a hint that this is a bad idea), the results would describe a > protein next to a collection of LJ entities put together in a crystal > lattice. Here is why: > > When exposed to air or water, copper surface > would undergo quick oxidation. That aside, let us say we have a magic > copper surface, which is not oxidated. This is even more problematic > for MD-type methods, because metal surfaces don't really have > conventional van der Waals due to presence of Fermi gas in > metallic crystals, as opposed to the type of polarization in covalent > dielectrics or even semiconductors. > > The only method somewhat reasonable for metals in MD I am aware of is > the Modified Embedded Atom Method (MEAM), which, I believe, is implemented > in LAMMPS. At the same time, it won't handle the protein. > > Alternatively, you can use DFTMD in CP2K and actually place > the short region of your protein close to the surface of copper and > see what happens. This in fact appears to be a half-decent idea and > I'd try it. Depending on the exchange correlation functions and the > basis you choose, you could probably do about a 100 atoms for a few > picoseconds within reasonable time with parallelization over 40-50 > decent cores. One thing to keeo in mind is that > the last time a checked (a few months back), CP2K had no electron > energy sampling aside from the Gamma point (VASP can do that, but it > is expensive). If that doesn't bother you, I would recommend going the DFT > route. You could even try to > build the energy surface for your protein-metal interaction and then > try to use it to parameterize a lower level method such as MD. > > Hope this helps. > > Alex > > > MT> Dear all, > > MT> I got a problem when modelling the interaction of a peptide and > MT> copper. I am using the force field of opls-aa. The charge of > MT> copper atoms are defined as zero and the non-bonded parameters are > MT> obtained from the CU ions, which can be > MT> found in opls-aa. > > MT> I applied smd and found that the peptide is approaching the > MT> copper surface and stopped when it is very close to the surface, > MT> due to some abnormal interaction. The whole process seems quite > MT> reasonable, however, when I rerun the > MT> simulation to obtain the interaction energy, for the vdw > MT> interaction between peptide and copper is zero, absolute zero, in the > whole process. > > MT> Can anyone give me some suggestions? Is my model reasonable? I am > MT> not so confident with modelling strategy of using the LJ > MT> parameters for copper bulk material. Also the zero interaction energy is > another big issue. > > MT> Cheers, > MT> Tng > > > > > > ------------------------------ > > Message: 5 > Date: Thu, 23 Apr 2015 01:18:01 -0600 > From: Alex <[email protected]> > To: "[email protected]" > <[email protected]>, <[email protected]> > Subject: Re: [gmx-users] Problems in the modelling of interaction > between peptide and copper > Message-ID: <[email protected]> > Content-Type: text/plain; charset=us-ascii > > Hi, > > Your model is not reasonable. First of all, let's forget about > simulations for a moment, because GMX is a wildly inappropriate choice. > Even if the energy was not zero (which in your case is probably due to > some LJ parameters set up improperly, or some really good GMX code > giving you a hint that this is a bad idea), the results would describe a > protein next to a collection of LJ entities put together in a crystal > lattice. Here is why: > > When exposed to air or water, copper surface > would undergo quick oxidation. That aside, let us say we have a magic > copper surface, which is not oxidated. This is even more problematic > for MD-type methods, because metal surfaces don't really have > conventional van der Waals due to presence of Fermi gas in > metallic crystals, as opposed to the type of polarization in covalent > dielectrics or even semiconductors. > > The only method somewhat reasonable for metals in MD I am aware of is > the Modified Embedded Atom Method (MEAM), which, I believe, is implemented > in LAMMPS. At the same time, it won't handle the protein. > > Alternatively, you can use DFTMD in CP2K and actually place > the short region of your protein close to the surface of copper and > see what happens. This in fact appears to be a half-decent idea and > I'd try it. Depending on the exchange correlation functions and the > basis you choose, you could probably do about a 100 atoms for a few > picoseconds within reasonable time with parallelization over 40-50 > decent cores. One thing to keeo in mind is that > the last time a checked (a few months back), CP2K had no electron > energy sampling aside from the Gamma point (VASP can do that, but it > is expensive). If that doesn't bother you, I would recommend going the DFT > route. You could even try to > build the energy surface for your protein-metal interaction and then > try to use it to parameterize a lower level method such as MD. > > Hope this helps. > > Alex > > > MT> Dear all, > > MT> I got a problem when modelling the interaction of a peptide and > MT> copper. I am using the force field of opls-aa. The charge of > MT> copper atoms are defined as zero and the non-bonded parameters are > MT> obtained from the CU ions, which can be > MT> found in opls-aa. > > MT> I applied smd and found that the peptide is approaching the > MT> copper surface and stopped when it is very close to the surface, > MT> due to some abnormal interaction. The whole process seems quite > MT> reasonable, however, when I rerun the > MT> simulation to obtain the interaction energy, for the vdw > MT> interaction between peptide and copper is zero, absolute zero, in the > whole process. > > MT> Can anyone give me some suggestions? Is my model reasonable? I am > MT> not so confident with modelling strategy of using the LJ > MT> parameters for copper bulk material. Also the zero interaction energy is > another big issue. > > MT> Cheers, > MT> Tng > > > > > > ------------------------------ > > -- > 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 [email protected]. > > End of gromacs.org_gmx-users Digest, Vol 132, Issue 106 > ******************************************************* -- 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 [email protected].
