Previous Thread
Next Thread
Print Thread
Page 1 of 6 1 2 3 4 5 6
generating parameters for pyruvate
#23873 04/02/10 11:30 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Hi,

This is my first time attempting the protocol outlined on http://dogmans.umaryland.edu/~kenno/#CGenFF.

I am working on generating parameters for pyruvate. And I want to check to see if my approach is correct.

This is what I have done:
*Note: I wasn't able to follow everything in the tutorial and I grossly simplified many steps.

1. I created the rtf file by analogy.
2. I ran generate_tutorial.inp
3. I put the resultant gjf file into Gaussian 09.
4. I looked at the Hessian written in the internal coordinates and extracted the diagonal entries as my force constants.
5. I scaled the values into kcal/(mol Angstroms)
6. For the dihedral terms I divided by multiplicity squared, and I also used the minimized coordinates to get delta (by multiplying Xmin by multiplicity and subtracting 180 degrees).
7. I entered only the missing parameters into the rtf.
8. I then ran gaussian again with scrf=(pcm,solvent=water) to simulate the presence of water.
9. I used these Mulliken charges in my rtf.
10. I ran a minimization (without waters) and got a distorted structure.
11. I adjusted the charges so that the methyl group had .09 on each hydrogen and .21 on the carbon and reran the minimization (without waters).
12. There was less distortion but the structure is not perfect.


So at this point I have a few questions:
1. Will this approach work?
2. Can I just keep tweaking charges?
3. What should my goal be? Should I just try to match the minimized gaussian structure as closely as possible using the following criteria: "typically 0.03 for the bonds and 3 for the angles".

Thanks

Re: generating parameters for pyruvate
RD2015 #23888 04/05/10 06:53 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Originally Posted By: RD2015
I am working on generating parameters for pyruvate. And I want to check to see if my approach is correct.
Thank you for your interest. This looks like a worthwile endeavor. Just for your reference, here are the suggested atom types for pyruvate (it's important to start from good atom types).
Code:
ATOM C1     CG2O3
ATOM O11    OG2D2
ATOM O12    OG2D2
ATOM C2     CG2O5
ATOM O2     OG2D3
ATOM C3     CG331
ATOM H31    HGA3
ATOM H32    HGA3
ATOM H33    HGA3
Also, don't forget to define the two improper dihedrals!

Originally Posted By: RD2015

This is what I have done:
*Note: I wasn't able to follow everything in the tutorial and I grossly simplified many steps.
Can you please elaborate on which steps you were not able to follow and why? I'd like to improve the tutorial but we get precious little feedback.



Originally Posted By: RD2015
4. I looked at the Hessian written in the internal coordinates and extracted the diagonal entries as my force constants.
This generally doesn't work well because the diagonal elements of the Hessian in most cases cannot be trivially correlated with the parameters.

Originally Posted By: RD2015
6. For the dihedral terms I divided by multiplicity squared,
confused What did you divide? By what multiplicity? Why? I don't really understand what you're doing but it's unlikely to be correct. Please perform a potential energy scan around the C-C-C-O dihedral and fit a single 2-fold term with phase 180 (or 0) to it. This is the only correct way to fit this dihedral.

Originally Posted By: RD2015
and I also used the minimized coordinates to get delta (by multiplying Xmin by multiplicity and subtracting 180 degrees).
Danger, Will Robinson! How about reading the FAQ?

Originally Posted By: RD2015
8. I then ran gaussian again with scrf=(pcm,solvent=water) to simulate the presence of water.
9. I used these Mulliken charges in my rtf.
This is wrong on so many levels.
  • Using an implicit solvent to generate target data for parameterizing a force field that is predominantly meant to work with explicit solvent is inconsistent. It is like using a condensed phase force field to generate target data for parameterizing a DFT method...
  • Mulliken is already not such a great method for calculating charges (basis set dependence,...), and combining it with PCM takes the last bit of relevance out of it.
  • Straight Mulliken charges are generally not useful in force field design; force fields that use QM-based charges use an electrostatic fitting method such as RESP.
  • CHARMM is not one of those force fields. Again, read the fine FAQ.


Originally Posted By: RD2015
10. I ran a minimization (without waters) and got a distorted structure.
11. I adjusted the charges so that the methyl group had .09 on each hydrogen and .21 on the carbon and reran the minimization (without waters).
12. There was less distortion but the structure is not perfect.
Given everything you did wrong, it is hardly surprising you get a "distorted structure". Also, you probably forgot to define improper dihedrals on the carbonyl groups.

Originally Posted By: RD2015
So at this point I have a few questions:
1. Will this approach work?
2. Can I just keep tweaking charges?
3. What should my goal be? Should I just try to match -snip-
No, no and no. We didn't make up our parameterization procedure overnight on a napkin in the bar after drinking a bit to much. It is the result of decades of experience and extensive validation by many research groups. Please re-read the CGenFF paper and also look up some of the references (especially the MacKerell and the Jorgensen ones, if you don't believe me). Then, start all over again and try to follow the tutorial as closely as you can. I'll understand if you can't figure out how to use molvib on your own, and in the specific case of pyruvate, you may get away with it, but you cannot optimize the charges without water interactions, and there's no way around doing a potential energy scan for the central dihedral. You have access to Gaussian 09; all the calculations should be light enough to run on any desktop PC produced in the last 6 years.

Edit: oh yeah, make sure you use the downloadable version of the tutorial, else you won't have all the files.
http://dogmans.umaryland.edu/~kenno/cgenff/download.html

Re: generating parameters for pyruvate
Kenno #23894 04/05/10 10:37 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Thank you for your advice. I'll start over and get back to you with any new questions.

Re: generating parameters for pyruvate
RD2015 #23905 04/07/10 12:00 AM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Hi,

I restarted using the atom types you suggested. Everything works fine but I am now stuck at the "Target data" portion of the tutorial.

Here's what I have done so far:
1. Redid atom types as you suggested above.
2a. Changed zeed to 1 C1 1 C2 1 C3
2b. Ran generate tutorial twice to find out what parameters were missing
2c. The missing parameters were:

ANGLES
CG2O3 CG2O5 CG331 35.00 116.00 ! adapted from "ACO, acetone adm 11/08"

DIHEDRALS
CG2O3 CG2O5 CG331 HGA3 0.1000 3 0.00 !adapted from "BTON, butanone
OG2D2 CG2O3 CG2O5 CG331 0.3000 2 180.00 !adapted from "BIPHENYL ANALOGS unmodified,peml"

3. Ran gaussian on pyr_opt_freq_mp2.gjf

4. Assigned .09 charge to each hydrogen. And rescaled Mulliken charges:
i) x = (Sum of Mulliken charges for hydrogens) - (# of hydrogens * .09)
ii) x/(remaining # of atoms)=y
iii) (mulliken charge on remaining atom) +y= Charge inputed into pyruvate.str file

Now I have a few questions:

1. In the tutorial you prescribe using Merz-Kollman charges. I can't find them in the output gaussian log file. Is it okay to use Mulliken charges?

2. I am really confused about what to do next. Am I supposed to generate different monohydrate pyruvates and then minimize in gaussian to obtain an energy? (If so how many do I need? 3?)

3. Assuming that's the right idea, how do I generate these files that are analogous to "prld_water_hf_H1_0.gjf"?

4. I was also looking at the next script "water_constr_tutorial.inp" and I noticed that it requires a file called @residue_mp2.crd. I realize this comes from "pyr_opt_freq_mp2" and is somewhere under "Input Orientation" in the log file. What is the standard way to generate the crd from the log file?

Thanks again

EDIT: I forgot to ask about the improper terms.

I am not sure I understand when the improper terms need to be included. Are the two impropers I need to include: C2 attached to C1 C3 O2, and C1 attached to O11, O12, C2?

Last edited by RD2015; 04/07/10 01:00 AM. Reason: forgot impropers
Re: generating parameters for pyruvate
RD2015 #23913 04/07/10 05:13 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Originally Posted By: RD2015
-snip-
2c. The missing parameters were:
Code:
ANGLES
CG2O3  CG2O5  CG331    35.00    116.00               ! adapted from "ACO, acetone adm 11/08"

DIHEDRALS
CG2O3  CG2O5  CG331  HGA3       0.1000  3     0.00   !adapted from "BTON, butanone
OG2D2  CG2O3  CG2O5  CG331      0.3000  2   180.00   !adapted from "BIPHENYL ANALOGS unmodified,peml"
These look like good initial guesses to me.

Originally Posted By: RD2015
1. In the tutorial you prescribe using Merz-Kollman charges. I can't find them in the output gaussian log file. Is it okay to use Mulliken charges?
Nope, still not OK. Read what I wrote about Mulliken charges in my previous post. We haven't changed our mind overnight.

I don't know about Gaussian 09, but Gaussian 03 produces an output section titled "Electrostatic Properties Using The MP2 Density" which contains (among other things) the Merz-Kollman charges, under "Charges from ESP fit" (if you used the input from the tutorial). Be careful, it does this 2 times, once for the input geometry and a second time for the minimized geometry. It's the second one you want to look at!

Also note that we use the Merz-Kollman charges only as an initial guess. If you're not going to generate water interactions and optimize the charges based on those, then forget about QM-based charges altogether and just use charges by analogy to existing compounds (but then again, the analogy will not be very good). Again, see this FAQ.

Originally Posted By: RD2015
4. Assigned .09 charge to each hydrogen. And rescaled Mulliken charges:
i) x = (Sum of Mulliken charges for hydrogens) - (# of hydrogens * .09)
ii) x/(remaining # of atoms)=y
iii) (mulliken charge on remaining atom) +y= Charge inputed into pyruvate.str file
Holy contrivance, batman! I have to repeat that you can't just make up your own scheme for assigning charges. Also, what is wrong with the + and - operators? Scaling charges is not a good idea because you're also scaling the dipole moment and all polar interactions the molecule would make in condensed phase. Here's what you should do instead (using the Merz-Kollman charges):
  • Sum the 4 charges on the CH3 group together.
  • Subtract (3 x .09 = .27) to get the charge on the carbon.
  • Assign .09 to each hydrogen.
In other words: the charge on the whole CH3 group should remain the same. Easy, isn't it?

Originally Posted By: RD2015
2. I am really confused about what to do next. Am I supposed to generate different monohydrate pyruvates and then minimize in gaussian to obtain an energy? (If so how many do I need? 3?)
The more water interactions you have, the less your system is underdefined and the easier it gets to generate relevant charges. Look at figure 2 in the CGenFF paper: we put waters wherever we saw possibility for a hydrogen-bond-like interaction. However, for pyruvate, it's slightly more complicated because you have carbonyl groups. For carbonyl groups, we usually generate two water interactions per oxygen: one where the oxygen donates a hydrogen bond in the lone pair direction and one where the oxygen donates a hydrogen bond in the direction of the C=O bond (ie. linear). Out of the two lone pairs each oxygen has, you should choose the one that suffers the least contamination by other parts of the molecule. For example, on the C2=O carbonyl, you should choose the lone pair that points away from the carboxylate. In total, you would have about 8 water interactions (2 per oxygen + 2 for the methyl hydrogens, assuming that the molecule is planar so that the 3rd methyl hydrogen is equivalent to one of the others for symmetry reasons).

Originally Posted By: RD2015
3. Assuming that's the right idea, how do I generate these files that are analogous to "prld_water_hf_H1_0.gjf"?
Since you're the first one asking this question on the forums and it's not a totally trivial thing to do, I'm willing to generate the input files for you so that other people can use them as an example. So if you can't figure it out on yourself, just use the forum's "File Manager" feature to attach your toppar file and your gaussian output to your next post, and I'll take care of it when I find a little hole in my schedule.

Originally Posted By: RD2015
4. I was also looking at the next script "water_constr_tutorial.inp" and I noticed that it requires a file called @residue_mp2.crd. I realize this comes from "pyr_opt_freq_mp2" and is somewhere under "Input Orientation" in the log file. What is the standard way to generate the crd from the log file?
The standard way is to take an existing crd file for the same molecule and copy-pasting the cartesian coordinates from the gaussian log file into it. Beware: the fields in the crd format are fixed-width. I have a script that does this but I may have forgotten to include it in the "downloadable" tutorial.

Note that because your molecule has a net charge, you should not apply the 1.16 scaling factors to your interaction energies, and also, you cannot compare the QM and CGenFF dipole moments, as the dipole moment is ill-defined for a charged entity. We do have a procedure to align the QM and MM geometry to the origin and to each other in a way that yields a reasonable comparison for most molecules, but this is not documented in the paper nor in the tutorial. I'm busy right now, but if you're interested, I'll try to post the CHARMM script one of the coming days.

Originally Posted By: RD2015
I am not sure I understand when the improper terms need to be included. Are the two impropers I need to include: C2 attached to C1 C3 O2, and C1 attached to O11, O12, C2?
That's correct.

Re: generating parameters for pyruvate
Kenno #23940 04/08/10 08:05 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Here are my files. Thanks for your help.

Attached Files
params-RD2015.tar.gz (71.48 KB, 374 downloads)
Re: generating parameters for pyruvate
RD2015 #24042 04/16/10 09:55 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Sorry for my late reply, it slipped a bit off my radar. Attached is your toppar stream file, after the following corrections:
  • completed your IC table
  • impropers must be defined both in the topology and parameter section
  • we're going to abolish wildcard impropers in a future CGenFF release so I changed the wildcard atoms into explicit ones
  • when editing machine-readable ASCII files, please use a font in which you can differentiate between O and 0... corrected several of these mistakes.
Also attached are the Gaussian input files for the water interactions. Note that there's some redundancy: you'll have to make a choice between O2_0 and O2_180 based on which one suffers the least contamination. This may become obvious once you have the interaction distances and energies. Same for O12_0 and O12_180 . Get back to me when you have a list of HF water interaction distances and energies, then I'll show you the next step (reproducing them with CHARMM).

Oh yeah, there's also a file named pyru_center.gjf . You can run it if you want to have dipole moments that are (somewhat) comparable with the CHARMM ones.

Attached Files
pyru_water_gauss.tgz (2.07 KB, 434 downloads)
Re: generating parameters for pyruvate
Kenno #24050 04/19/10 03:35 AM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Thanks again for the help.

I had a problem running g09 on pyru_water_hf_H31.gjf. The run terminates via linke 9999 error. I tried increasing resources and also put in the calcall option. And it still crashes.

I've attached the files.

Attached Files
pyru_water.tar.gz (40.56 KB, 415 downloads)
Re: generating parameters for pyruvate
RD2015 #24055 04/19/10 05:44 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
What happened is so obvious I'm not going to tell you. Carefully look at your output. Shotgun debugging almost never works!

Re: generating parameters for pyruvate
Kenno #24069 04/21/10 05:02 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
This doesn't bode well for your ability to bring this parameter optimization to a good end, but here's a hint. Your output contains the following:
Code:
 Optimization stopped.
    -- Number of steps exceeded,  NStep=  20
This is indeed the most common cause of Gaussian's "Error termination request processed by link 9999." Your geometry optimization ran out of cycles. Why did it run out of cycles? Look at the geometry. Or, in your case, since you're optimizing only one degree of freedom, you could just look at the value of that degree of freedom in your Gaussian output. What do you see? Why do you think this would happen?

Last edited by Kenno; 04/21/10 07:09 PM. Reason: spelling error
Re: generating parameters for pyruvate
Kenno #24080 04/22/10 03:45 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
It looks like the water failed to find a stable distance from H31. After each optimization cycle it kept flying further and further away.

Perhaps roh should have had a smaller initial value?

Edited:

I had something very strange occur. I played around with roh and set it to 2, 1.5, 1 and had no convergence. Then I set it to .8 and it converged to an optimized distance of 3.6287.

To check the value then I put 3.6287 value back into the initial script and it did not converge. Was the convergence I got for .8 some sort of artifact?

I also had another thought, maybe I should make the step size smaller and increase the number of steps?

Last edited by RD2015; 04/22/10 04:27 PM.
Re: generating parameters for pyruvate
RD2015 #24082 04/22/10 04:41 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Pyruvate is negatively charged and as we're only optimizing the distance, the water is always pointing with its oxygen atom towards the negative charge. This yields an unfavorable charge-dipole interaction. The interaction is repulsive, and as this is an entirely physical phenomenon, there's nothing to be done about it.

How do we use this information to optimize our charges? Either we just use the qualitative information ("our force field should also have a repulsive interaction here"), or we fix the water molecule at a certain distance both in the QM and the MM and we try to reproduce the (positive) interaction energy. The latter approach gives us more quantitative information, but has the important disadvantage that we introduce an arbitrary distance into the procedure. Furthermore, equilibrium distances in gas phase and in high-dielectric condensed phase are known to differ; see Figure 8 and the accompanying discussion in the CGenFF paper. This is not well-studied for repulsive interactions (it's not a trivial thing to do) so we have no idea whether we should apply and offset, and if yes, by how much.

In the end, here's what I'd recommend you to do:
  • Calculate all the QM water interaction energies
  • Check whether the two other interactions with the methyl group are repulsive as well. The optimization often seems to converge because the repulsive force becomes smaller than the convergence criterion at a certain distance, but you'll still get a positive interaction energy.
    • If all 3 interactions with the methyl group are repulsive, we might need to fix the distance. But talk to me first if this is the case.
    • If at least one of them is attractive (there are favorable dipole-dipole interactions that may overcome the unfavorable charge-dipole interaction for some interaction geometries), we have enough information to carry on with our charge optimization without having to bother about quantifying the repulsive interactions.

Re: generating parameters for pyruvate
Kenno #24084 04/22/10 06:34 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I think all three might be repulsive. The H31 and H33 script crash outright. The H32 script converges, but I calculated a SCF Done: E(RHF) = -415.972756709 for "far away" (roh=50) and I calculated a SCF Done: E(RHF) = -415.968272866 for the converged value (roh = 3.20987). So I think it may be repulsive as well.

I'll continue to run the other scripts in the meantime.

Edited:

By the way, I am not sure I need to model pyruvate extremely accurately for my experiments. Am I at the point where I have a "good enough" set of parameters so that I can start minimizing and heating my protein/ligand complex (at least for now)?

Last edited by RD2015; 04/22/10 07:00 PM.
Re: generating parameters for pyruvate
RD2015 #24091 04/22/10 07:42 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Then we'll just forget about the methyl group. On second thought, you have only 4 degrees of freedom in your charge optimization (9 atoms - 3 aliphatic hydrogens - 1 because the carboxylate oxygens are equivalent - 1 because of the total charge constraint), so as long as you have 4 or more acceptable (attractive) water interactions, your system won't be underdefined.

Next step is to calculate the interaction energies as (energy of the complex) - (energy of the model compound) - (energy of tip3 water = -76.010530 Hartree) and then convert to kcal/mol. When you're done with that, I'll show you how to do the CHARMM water interactions.

Originally Posted By: RD2015
By the way, I am not sure I need to model pyruvate extremely accurately for my experiments. Am I at the point where I have a "good enough" set of parameters so that I can start minimizing and heating my protein/ligand complex (at least for now)?
What do you mean? You don't have any charges yet! Also, minimizing and heating usually only takes a couple of hours.

Re: generating parameters for pyruvate
Kenno #24106 04/22/10 09:15 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I put the values in the .ods file inside the archive.

Attached Files
archive.tar.gz (623.94 KB, 701 downloads)
Re: generating parameters for pyruvate
RD2015 #24107 04/22/10 11:15 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
That looks pretty good. Here's the decisions I made:
  • Between O2_0 and O2_180, I chose O2_180 because O2_0 is somewhat hindered by the electrostatic effects of the carboxylate group (as apparent from the interaction distance and energy) and because O2_180 is closer to the methyl group, a region where we have limited coverage.
  • Between O12_0 and O12_180, I chose O12_0 even though O12_180 has a more favorable interaction energy. The reasons for doing this are:
    • We already have a water molecule at the carboxylate side: O11_180 . I don't expect O12_180 to contribute a lot more information.
    • O12_0 happens to be in the region between the carboxylate carbonyl and the ketone carbonyl; it might be interesting to have some coverage there and the distance from the ketone carbonyl is still long enough to not have too much contamination.

I prepared a CHARMM input file analogous to the water_constr.inp in the tutorial. You can use this to optimize your charges.

Tip: in cases where it's impossible to get two water interactions of the same atom right, the strongest water interaction has priority. Unless the difference between the two QM interaction energies is small (we sometimes use "smaller than 2 kcal/mol" as a rule of thumb), then it's better to go for a compromise between the two water interactions.

Edit: Now may be a good time to start your MP2 dihedral potential energy scan. See the tutorial for examples; don't "fix" any coordinates.

Attached Files
pyru_water_constr.inp.txt (22.24 KB, 757 downloads)
pyru_mp2.crd (1.34 KB, 718 downloads)
Last edited by Kenno; 04/23/10 04:06 PM. Reason: Oops... forgot attachments
Re: generating parameters for pyruvate
Kenno #24110 04/23/10 05:52 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I've attached the results of the run.

Is this the right idea?

i) The energy difference in 1 and 2 is less than 2kcal so I can keep the charge on O2?
ii) Increase the charge on O11 to get ene diff below 2kcal.
iii) keep the same charge on the methyl hydrogens?

Attached Files
archive.tar.gz (22.52 KB, 656 downloads)
Last edited by RD2015; 04/23/10 06:14 PM.
Re: generating parameters for pyruvate
RD2015 #24113 04/23/10 07:53 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
You have to do some things by yourself. This is largely a matter of experimentation: adjust your charges a bit, see how the water interactions are affected, adjust them again, repeat ad nauseum... As the CHARMM run takes mere seconds, you can do several of these cycles per minute if you organize yourself a bit.

Deviation of 2 kcal/mol are not acceptable as such. In general, we aim to get all the interaction energies within 0.2 kcal/mol, as mentioned in the CGenFF paper and demonstrated in the tutorial. But there are a few exceptions. If you have two water interactions for the same HB acceptor, then in some cases it is impossible to get them all right, because any adjustment of the charges will affect both of them in the same direction. If that's the case, then one should either go for the strongest interaction, or make a compromise between the two interactions. To decide whether to make a compromise or go for the strongest, you look at the difference between the two QM interaction energies, and see if it's lower or higer than 2 kcal/mol; see my previous post. But there are cases in which you can get both of them right - I'd recommend you play with the balance of charge between the 3 carbons and see what the effect is.

Originally Posted By: RD2015
iii) keep the same charge on the methyl hydrogens?
Yes, they should be kept at 0.09. And the two carboxylate oxygens should always be given the same charge because they're chemically equivalent. And the total charge of the molecule should be kept at -1 at all time. This can be most easily accomplished by shifting charge around. For instance, you could try to subtract 0.02 charge units from each of the carboxylate oxygens (making them more negative) and move the resulting 2x0.02=0.04 charge units onto the carboxylate carbon (making it more positive).

As this leaves you only 4 degrees of freedom, this is a really easy assignment to begin with.

The initial guess water interactions look not too bad to begin with; in fact, the ones on O2 and O12 look pretty good. But do play around with the charges a bit, you might get them better. Keep an eye on the dipole moments while doing so; we want to overestimate them a bit (although they're less important than the interaction energies, especially for charged systems). And try to keep the charges physical; it's sometimes tempting to run off in an unphysical direction. If you're not sure about what's physical, try to not deviate more than 0.1 or 0.2 charge units from your Merz-Kollman initial guess.

Last edited by Kenno; 04/23/10 07:58 PM. Reason: Typos...
Re: generating parameters for pyruvate
Kenno #24167 05/03/10 08:04 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
It's been a while... any luck with either the charge optimization or the dihedral potential energy scan?

If you can't find a way to improve the charges, well, in about 2% of cases, the Merz-Kollman charges are as good as it gets, so you might just be "lucky" wink

Edit: however, I would recommend rounding the charges to two digits after the comma. This is to preserve compatibility with some softwares that do not have a higher precision.

Last edited by Kenno; 05/03/10 08:08 PM. Reason: added note about rounding.
Re: generating parameters for pyruvate
Kenno #24259 05/11/10 02:56 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Hi. Sorry I've been swamped with exams and a presentation (I just finished today). I'll post my results soon.

Last edited by RD2015; 05/11/10 02:57 PM.
Re: generating parameters for pyruvate
RD2015 #24266 05/11/10 09:04 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I've been working on the optimization today and this the best I've been able to get.

For example, I am not sure how to decrease the energy difference for O11...HOH, 180. DEGREES without raising for O11...HOH LINEAR. If I increase charge on O11 I'll get a better 3) energy difference but a worse difference for 4).

What do you think?

Also note: .bak files are the initial values.

Attached Files
chargeoptimize.tar.gz (23.08 KB, 1151 downloads)
Last edited by RD2015; 05/11/10 09:06 PM.
Re: generating parameters for pyruvate
RD2015 #24284 05/12/10 06:41 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
I'm in a time crunch myself now, so I can only have a look at it next week. In the meanwhile, you may want to start the Gaussian dihedral potential energy scan; you don't need optimized charges for the QM part.

What you're saying about improving one interaction making another one worse is pretty common; it's very well possible that you just can't get it better than that.

Re: generating parameters for pyruvate
Kenno #24289 05/12/10 08:05 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Ok. I'll try to get started on it.

But I also have a more general question. The protein I am trying to simulate is LDH. There is a NADH and a pyruvate.

Is it correct to use the following parameter sets: par_all27_prot_na.prm, top_all27_prot_na.rtf, and toppar_all27_na_nad_ppi.str?

Re: generating parameters for pyruvate
RD2015 #24291 05/12/10 08:28 PM
Joined: Sep 2003
Posts: 4,783
Likes: 2
Forum Member
Online Content
Forum Member
Joined: Sep 2003
Posts: 4,783
Likes: 2
Yes, if all residues in LDH, NADH and pyruvate are contained in the files you mention; after reading these files CHARMM would complain if the entities in the sequence and coordinate files are not recognized.


Lennart Nilsson
Karolinska Institutet
Stockholm, Sweden
Re: generating parameters for pyruvate
lennart #24298 05/13/10 06:17 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
We've just had a discussion about a similar problem here. mooctopus needed to read in the CHARMM protein force field, carbohydrate force field and CGenFF; it sounds like you need the CHARMM protein force field, nucleic acid force field and CGenFF. The way to do this is analogous.

Re: generating parameters for pyruvate
Kenno #24327 05/19/10 12:56 AM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I ran potential energy scans of c1-c2-c3-h32 and o11-c1-c2-c3 in gaussian. From -180 to +180 in increments of 1 degree.

Do I need to do this for the two impropers too?

I am not sure I understand the charmm script for the PES. What line in the code specifies the dihedral angle and increment? Is this even the correct question; or is the dihedral angle and increment being pulled from an external initial pdb file?

I'd like to attach my log files but they are too big, what should I do?

***
I tried running the program "pes" using a modified code extract_pyruvate.

The attached output files are now included.

Attached Files
input.zip (900 Bytes, 534 downloads)
c1-c2-h32.tar.gz (36.37 KB, 615 downloads)
o11-c1-c2-c3.tar.gz (5.37 KB, 527 downloads)
Last edited by RD2015; 05/19/10 04:37 PM.
Re: generating parameters for pyruvate
RD2015 #24328 05/19/10 05:53 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
About the charge optimization: that looks like an improvement.

About the scans:
  • The methyl group is not all that critical so you probably would have got away without scanning it, but it surely doesn't hurt.
  • Scanning in increments of 1° is a bit overkill, but if your computers can handle it in a reasonable amount of time, it doesn't hurt.
  • Given the symmetry, you could have scanned the C(sp²)-C(sp²) dihedral from -180 to 0; your energy profile is probably symmetric. Yet again, this doesn't hurt.
  • No, you don't need to scan the impropers. We only scan impropers when they're floppy (typically for nitrogens that are not involved in a double bond but are planar because of resonance). The rigid ones (like in your molecule), we ordinarily fit using molvib, but the parameters by analogy are probably good enough in your particular case. If you want to optimize them nevertheless, and if you don't want to try molvib, then scanning them would be the easiest way to go.
  • There's a python script in the tutorial that extracts all scan points from the Gaussian output and saves them as pdb files. It also makes an s_file.txt , each line of which is the full path to the pdb file of the corresponding scan point. The CHARMM script reads this s_file line by line, and for each line, it reads the geometry from the pdb file, constrains the dihedral(s) being scanned, minimizes all other degrees of freedom, and writes the energy (not including the constraint energy).
  • Note that for simple cases like this, we recommend fitting the dihedrals manually rather than using the automated fitting procedure, which only becomes interesting for more complicated cases.
  • Apart from the tutorial, the paper and the previous discussions on the public CHARMM forums, here's a mirror of a private discussion thread that might be helpful.


Edit: running the python script "pes" through a modified version of "extract_tutorial" or "extract_pes" indeed is the way to go.

Last edited by Kenno; 05/19/10 07:23 PM.
Re: generating parameters for pyruvate
Kenno #24330 05/19/10 08:41 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I am not sure I understand the purpose of loop_mead. Do I even need it?

Re: generating parameters for pyruvate
RD2015 #24331 05/19/10 10:12 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
The loop is for doing measurements on selected dihedrals. This is needed for the automatic fitting procedure but indeed not for the manual fitting you're going to do. However, do not mess with that script; it ain't broken yet, and given the robustness of the CHARMM "programming language", you may easily break it. The script is written in a way that all you ever need to change is the list of variables and "open" commands at the top. For instance, to bypass the measurements, simply set "meas" to 0 . That's why the "if @meas eq 0 goto loop_pes" line is there.

Generally spoken, running the tutorial and looking at the outputs it creates is a good way to get more insight in what the scripts do.

Re: generating parameters for pyruvate
Kenno #24334 05/19/10 11:31 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I've been ramping up the kpsi on o11-c1-c2-c3 and I am getting a somewhat better fit.

My files are attached.

Attached Files
scan.tar.gz (424.79 KB, 655 downloads)
graph.jpg (28.9 KB, 815 downloads)
Last edited by RD2015; 05/20/10 03:29 PM. Reason: added jpg
Re: generating parameters for pyruvate
RD2015 #24363 05/21/10 01:50 AM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
OK, we're getting somewhere. Here are a few tips:
  • For the OCCC, the shape could be better. I would try adding a 4-fold contribution to improve it. However, be careful not to add contributions with an odd phase. For background information, see the CGenFF FAQ and the private discussion thread I refered to earlier.
  • For the CCCH:
    • You see that the left part of your graph has a constant offset, which makes it more difficult to fit. One way to solve this is to split it in two: one scan for the CCCH (in which only CCCH is constrained) and another one for the OCCC, both of them aligned to their own minimum. The other (more challenging) solution is to first get the shape of the OCCC right (if you take your QM minimum in the OCCC scan and look at the MM energy at that same dihedral angle, there you have your offset...)
    • You currently only have a 3-fold contribution for the CCCH. Keep it that way.

Re: generating parameters for pyruvate
Kenno #24364 05/21/10 03:21 AM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I removed the constraint for the c-c-c-h like you said. Fit looks a little better?

Changing multiplicity for the other dihedral produced strange results. I started getting peaks inside of the wells.

***
By the way, I need to set up explicit solvation and periodic boundary conditions. Offhand can you recommend any good resources?

Thanks again.

Attached Files
scan.tar.gz (634.26 KB, 536 downloads)
scan.jpg (27.38 KB, 704 downloads)
Last edited by RD2015; 05/21/10 04:22 PM.
Re: generating parameters for pyruvate
RD2015 #24366 05/21/10 05:17 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Originally Posted By: RD2015
Changing multiplicity for the other dihedral produced strange results. I started getting peaks inside of the wells.
Yes. We're trying to reproduce the QM. The QM has peaks inside the wells. One thing you could do is visualize the QM scan and see whether any of the C(sp²) centers became slightly nonplanar around the energy minimum. That is a possible source of the peaks, and in that case, we don't necessarily want to reproduce them.

At any rate, I didn't say changing multiplicity - I said adding a multiplicity, so that you have a 2-fold and a 4-fold term for the same dihedral. There are tons of resources (most importantly, the tutorial) that show how this works. But check the planarity first.

Yesterday, I forgot to look at the scale of the graph. I now discovered it's from 0 to 1 kcal/mol; in other words, the discrepancies around the minimum are are small in terms of energy. As far as I am concerned, the fit is acceptable. But I still think you can get the OCCC even better by adding a 4-fold (and/or possibly a 6-fold, but let's not go there for now).

Originally Posted By: RD2015
By the way, I need to set up explicit solvation and periodic boundary conditions. Offhand can you recommend any good resources?
The script archive comes to mind...

Re: generating parameters for pyruvate
Kenno #24372 05/21/10 10:50 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I visualized the qm pdb files. You were right there is a loss of sp2 planarity of about 3 degrees.

Are these parameters good to use now?

***

I have been looking at the scripts from the archive as well as schlick's book. I was just wondering if there is anything else to read for background.

Re: generating parameters for pyruvate
RD2015 #24378 05/24/10 05:17 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Can you gzip the gaussian output file of the OCCC scan and post it?

Re: generating parameters for pyruvate
Kenno #24380 05/24/10 05:22 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
It's just over the limit for the forum. So i posted on an external site.

http://ifile.it/nsk86br

Re: generating parameters for pyruvate
RD2015 #24382 05/24/10 05:59 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
While the tiny deviation from planarity might be part of the cause of the "peak inside the well", we don't want to try to reproduce this deviation, because doing so would likely compromise the vibrational spectrum and the transferability. Nevertheless, I'd prefer the shape of the well to be reproduced a bit better. You can use the parameters as they are if you think they're good enough for your purpose, but I'd recommend trying with the additional term(s). It's not a lot of work.

Re: generating parameters for pyruvate
Kenno #24389 05/24/10 07:55 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Thanks.

Are you saying I should add in an additional term like:
OG2D2 CG2O3 CG2O5 CG331 (?) 4 0.00

And see if it helps the fit?

***
I added the additional term. Picture is attached.

Attached Files
newscan.jpg (25.59 KB, 694 downloads)
Last edited by RD2015; 05/24/10 08:06 PM. Reason: added the term
Re: generating parameters for pyruvate
RD2015 #24401 05/25/10 04:59 AM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Looks like a dramatic improvement, doesn't it? smile

However, something else is bothering me now. Although the picture is very small, it looks to me like the symmetry is broken. What are the final parameters for that dihedral?

Re: generating parameters for pyruvate
Kenno #24406 05/25/10 03:16 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
CG2O3 CG2O5 CG331 HGA3 0.1200 3 0.00
OG2D2 CG2O3 CG2O5 CG331 0.6450 2 180.00
OG2D2 CG2O3 CG2O5 CG331 0.0800 4 0.00

Re: generating parameters for pyruvate
RD2015 #24411 05/25/10 09:23 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
OK, I took a very close look at the whole case, and discovered that the line
Code:
IMPR C1   O11  O12  C2
is missing from your topology. Unfortunately, this means that you'll have to redo your dihedral fit. For starters, try 0.44 for your 2-fold and 0.06 for your 3-fold. Also, if you're still suffering asymmetry problems, a possible solution would be to not constrain CCCH while fitting OCCC.

Edit: beware: it is not impossible that CCCH needs to be re-fit as well.

Last edited by Kenno; 05/25/10 09:28 PM.
Re: generating parameters for pyruvate
Kenno #24412 05/25/10 09:40 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I thought that improper was accounted for by

CG2O3 OG2D2 OG2D2 CG2O5 96.0000 0 0.00
.

Complete file is attached

Attached Files
pyruvate-rd2015.tar.gz (844 Bytes, 390 downloads)
Re: generating parameters for pyruvate
RD2015 #24414 05/25/10 10:29 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Originally Posted By: RD2015
I thought that improper was accounted for by
Code:
CG2O3  OG2D2  OG2D2  CG2O5     96.0000  0     0.00
Every improper needs to be defined both in the topology and in the parameter section of the stream file, as discussed elaborately earlier in this thread. If you don't believe me, just try it.

Re: generating parameters for pyruvate
Kenno #24416 05/25/10 11:04 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Sorry I didn't even realize it was missing there.

The scan & rtf are attached.

Attached Files
pyruvate-rd2015.tar.gz (847 Bytes, 364 downloads)
scan.jpg (98.48 KB, 663 downloads)
Last edited by RD2015; 05/26/10 12:16 AM. Reason: rotated pic
Re: generating parameters for pyruvate
RD2015 #24418 05/26/10 01:39 AM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
You see, the major asymmetry around the peaks is gone now. There is still a minor asymmetry around the minima, which most probably is caused by the small offset in the location of the minima in the CCCH scan, and which can be resolved by not constraining CCCH while scanning OCCC. I suggest you try that now. Once that problem is resolved, you might want to tweak your amplitudes a bit, although it may not be necessary.

Re: generating parameters for pyruvate
Kenno #24419 05/26/10 01:56 AM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Ok I removed the other constraint. Picture is attached.

Attached Files
scan2.jpg (98.05 KB, 641 downloads)
Last edited by RD2015; 05/26/10 01:57 AM.
Re: generating parameters for pyruvate
RD2015 #24420 05/26/10 02:20 AM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Nice! Next step in the procedure is to take the input file for the water interactions, change "stage" to 5, run it, and see if the water interaction energies are still as good. If not, the charges may be tweaked.

Re: generating parameters for pyruvate
Kenno #24426 05/26/10 05:06 AM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
stage 5 result is attached.

Attached Files
stage5.tar.gz (23.23 KB, 519 downloads)
Re: generating parameters for pyruvate
RD2015 #24427 05/26/10 04:14 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Are there significant differences in the water interaction energies, compared to the constrained case? If yes, you might want to consider tweaking the charges. However, if you decide to do so, you'd have to tweak the dihedrals again after tweaking the charges. Then, you have to calculate the charges again and see if tweaking the dihedrals made any difference. Usually, you won't see big changes at that point (meaning that the force field is self-consistent).

Re: generating parameters for pyruvate
Kenno #24428 05/26/10 04:20 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
I'm not sure. I'll put the data inside of the post.

constrained:
Code:
PYRU WATER INTERACTION
######## USING MP2 GEOMETRY
QM DIPOLE (DEBYE)
MP2/6-31G*: X= 4.6182 Y= 1.9920 Z= 0.6033 TOT= 5.0655
EMPIRICAL DIPOLE: 6.08792 1.98516 0.682021 6.43962
1) O2...HOH, 180. DEGREES
A.I. -10.21 1.94
EMP. -9.51408 -10.7018 1.18769 1.7
ENE DIFF: 0.69592 DIST DIFF: -0.24
2) O2...HOH LINEAR, HOH: 90. DEGREES
A.I. -8.86 1.97
EMP. -9.51903 -10.7673 1.24829 1.68
ENE DIFF: -0.65903 DIST DIFF: -0.29
3) O11...HOH, 180. DEGREES
A.I. -13.45 1.87
EMP. -11.8093 -13.3237 1.51444 1.72
ENE DIFF: 1.6407 DIST DIFF: -0.15
4) O11...HOH LINEAR, HOH: 90. DEGREES
A.I. -12.97 1.84
EMP. -14.0462 -16.0888 2.04265 1.68
ENE DIFF: -1.0762 DIST DIFF: -0.16
5) O12...HOH, 0. DEGREES
A.I. -11.03 2.01
EMP. -9.8125 -11.1797 1.36724 1.77
ENE DIFF: 1.2175 DIST DIFF: -0.24
6) O12...HOH LINEAR, HOH: 90. DEGREES
A.I. -12.65 1.86
EMP. -14.0059 -16.1521 2.14621 1.68
ENE DIFF: -1.3559 DIST DIFF: -0.18
7) C3H..OHH 90.0 DEGREES
A.I. 0 4.50
EMP. 3.25087 3.3397 -8.882918E-02 4.49
ENE DIFF: 3.25087 DIST DIFF: -1E-02
8) C3H..OHH 90.0 DEGREES
A.I. 0 4.50
EMP. 3.07274 3.42676 -0.354018 2.89
ENE DIFF: 3.07274 DIST DIFF: -1.61
9) C3H..OHH 90.0 DEGREES
A.I. 0 4.50
EMP. 2.17095 2.23453 -6.357937E-02 4.49
ENE DIFF: 2.17095 DIST DIFF: -1E-02
AVEDIFF, RMSDIFF, AVERAGE ABSOLUTE ERROR
0.995283 1.62949 1.6822


unconstrained
Code:
PYRU WATER INTERACTION
QM DIPOLE (DEBYE)
MP2/6-31G*: X= 4.6182 Y= 1.9920 Z= 0.6033 TOT= 5.0655
EMPIRICAL DIPOLE: 6.51472 2.04064 0.629351 6.85579
1) O2...HOH, 180. DEGREES
A.I. -10.21 1.94
EMP. -9.10737 -10.4329 1.32555 1.72
ENE DIFF: 1.10263 DIST DIFF: -0.22
2) O2...HOH LINEAR, HOH: 90. DEGREES
A.I. -8.86 1.97
EMP. -9.37746 -10.6219 1.24447 1.68
ENE DIFF: -0.51746 DIST DIFF: -0.29
3) O11...HOH, 180. DEGREES
A.I. -13.45 1.87
EMP. -11.8425 -13.3543 1.51179 1.72
ENE DIFF: 1.6075 DIST DIFF: -0.15
4) O11...HOH LINEAR, HOH: 90. DEGREES
A.I. -12.97 1.84
EMP. -14.1356 -16.3738 2.23822 1.67
ENE DIFF: -1.1656 DIST DIFF: -0.17
5) O12...HOH, 0. DEGREES
A.I. -11.03 2.01
EMP. -10.4721 -11.9395 1.46745 1.74
ENE DIFF: 0.5579 DIST DIFF: -0.27
6) O12...HOH LINEAR, HOH: 90. DEGREES
A.I. -12.65 1.86
EMP. -14.1965 -16.3492 2.15268 1.68
ENE DIFF: -1.5465 DIST DIFF: -0.18
7) C3H..OHH 90.0 DEGREES
A.I. 0 4.50
EMP. 2.16734 2.23126 -6.391633E-02 4.49
ENE DIFF: 2.16734 DIST DIFF: -1E-02
8) C3H..OHH 90.0 DEGREES
A.I. 0 4.50
EMP. 2.75731 3.08851 -0.331209 2.82
ENE DIFF: 2.75731 DIST DIFF: -1.68
9) C3H..OHH 90.0 DEGREES
A.I. 0 4.50
EMP. 3.08129 3.16696 -8.56656E-02 4.49
ENE DIFF: 3.08129 DIST DIFF: -1E-02
AVEDIFF, RMSDIFF, AVERAGE ABSOLUTE ERROR
0.893823 1.58932 1.6115


So unconstrained is actually a little better?

Re: generating parameters for pyruvate
RD2015 #24429 05/26/10 07:08 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Don't get blinded by the averages on the last line, they are the reason why (at least some) automatic fitting programs perform so badly.
  • The interaction energies with the CH3 group got better, but they're totally irrelevant because they're repulsive.
  • Interaction 1) got substantially worse; I'd try to fix that up.
  • Interaction 5) got substantially better, at the cost of interaction 6), which got a little bit worse. Unfortunately, out of these two interactions, interaction 6) is more important (lower AI energy). Yet, improving interaction 6) should not lead to a worsening of interaction 3). The carboxylate group is a bit messy and I'm not sure you can improve it, but you can try. I'm not a psychic.

When you're done with that, yes, it's perfectly possible that the final unconstrained values are better than the constrained ones.

Re: generating parameters for pyruvate
Kenno #24495 06/07/10 03:48 PM
Joined: Jan 2010
Posts: 85
R
RD2015 Offline OP
Forum Member
OP Offline
Forum Member
R
Joined: Jan 2010
Posts: 85
Sorry I haven't posted my results in a while. I was on break.

I adjusted the charges and replotted the scan.

Code:
PYRU WATER INTERACTION
QM DIPOLE (DEBYE)
MP2/6-31G*: X= 4.6182 Y= 1.9920 Z= 0.6033 TOT= 5.0655
EMPIRICAL DIPOLE: 6.17438 2.37572 0.68391 6.65092
1) O2...HOH, 180. DEGREES
A.I. -10.21 1.94
EMP. -9.98738 -11.5515 1.56408 1.7
ENE DIFF: 0.22262 DIST DIFF: -0.24
2) O2...HOH LINEAR, HOH: 90. DEGREES
A.I. -8.86 1.97
EMP. -10.2175 -11.673 1.45541 1.66
ENE DIFF: -1.3575 DIST DIFF: -0.31
3) O11...HOH, 180. DEGREES
A.I. -13.45 1.87
EMP. -11.6132 -13.1247 1.51148 1.72
ENE DIFF: 1.8368 DIST DIFF: -0.15
4) O11...HOH LINEAR, HOH: 90. DEGREES
A.I. -12.97 1.84
EMP. -13.9382 -16.0165 2.07834 1.68
ENE DIFF: -0.9682 DIST DIFF: -0.16
5) O12...HOH, 0. DEGREES
A.I. -11.03 2.01
EMP. -9.98958 -11.3142 1.32461 1.75
ENE DIFF: 1.04042 DIST DIFF: -0.26
6) O12...HOH LINEAR, HOH: 90. DEGREES
A.I. -12.65 1.86
EMP. -14.0225 -16.1759 2.1534 1.68
ENE DIFF: -1.3725 DIST DIFF: -0.18
7) C3H..OHH 90.0 DEGREES
A.I. 0 4.50
EMP. 2.22844 2.29245 -6.401592E-02 4.49
ENE DIFF: 2.22844 DIST DIFF: -1E-02
8) C3H..OHH 90.0 DEGREES
A.I. 0 4.50
EMP. 2.67073 2.99919 -0.328457 2.8
ENE DIFF: 2.67073 DIST DIFF: -1.7
9) C3H..OHH 90.0 DEGREES
A.I. 0 4.50
EMP. 3.03498 3.12064 -8.56574E-02 4.49
ENE DIFF: 3.03498 DIST DIFF: -1E-02
AVEDIFF, RMSDIFF, AVERAGE ABSOLUTE ERROR
0.815088 1.64993 1.63691

Attached Files CCCH.pngOCCC.png
Re: generating parameters for pyruvate
RD2015 #24498 06/07/10 06:16 PM
Joined: Dec 2005
Posts: 1,535
Forum Member
Offline
Forum Member
Joined: Dec 2005
Posts: 1,535
Hmm, the charges are not convincingly improved. This seems to be a bit of a tough case; the ranking of interactions 1-6 is significantly different in the QM and the MM. Maybe we should go for a compromise after all, like getting the average difference of 1 and 2 as close to 0 as possible, and same for the average difference of 3,4,5, and 6.

The scans are quite nice.

Re: generating parameters for pyruvate
Kenno #34581 11/21/14 02:34 AM
Joined: Nov 2014
Posts: 1
Y
yxb Offline
Forum Member
Offline
Forum Member
Y
Joined: Nov 2014
Posts: 1
Hi Kenno,I'm generating parameters for my molecular following the "CGenFF tutorial" (file:///E:/documents/ug-tutorial/Charmm/cgenff-tutorial/cgenff_tutorial_old/supp/tutorial.html) and your guidance.

I have finished the part "Target Part".When I do the next "Initial charge optimization using the QM geometry",I'm confused with how to get the file water_constr.inp. Is there any generator to get it? Or should I write it manually? If it's the latter case, how can I prepare "the IC table for surface of interest"? Is there any principles to write this part?

Sorry for so many questions,maybe some are silly. frown I'm new for CHARMM and I'm really trapped here. Looking forward for your reply.

Re: generating parameters for pyruvate
RD2015 #37644 08/29/19 09:48 AM
Joined: Aug 2019
Posts: 2
T
Forum Member
Offline
Forum Member
T
Joined: Aug 2019
Posts: 2
Can you please explain how the dummy atoms are chosen to be placed in the script ?

In this folder :pyru_water_gauss.tgz

Last edited by thapa; 08/29/19 09:51 AM.
Page 1 of 6 1 2 3 4 5 6

Moderated by  alex, lennart, rmv 

Link Copied to Clipboard
Powered by UBB.threads™ PHP Forum Software 7.7.4
(Release build 20200307)
Responsive Width:

PHP: 5.6.33-0+deb8u1 Page Time: 0.034s Queries: 140 (0.017s) Memory: 1.4694 MB (Peak: 2.0035 MB) Data Comp: Off Server Time: 2020-08-07 18:30:00 UTC
Valid HTML 5 and Valid CSS