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).
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.
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
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.
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.
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.
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.