The directory contains examples which can be used
to test features of GCEM and provide examples of
the necessary program input files and settings.
The examples may also be useful as primers for
setting up own calculations.

For the sake of simplicity, the system is always
named by the single character "d".

Each GCEM calculation consist of at least two program runs.
The preprocessing run generates sidechain rotamers,
calculates molecular mechanics energy terms, writes
a restar file for the postprocessing run and
sets up the necessary input for the continuum
electrostatics calculations that are done by
independent MEAD programs (my_xyz_solver).
This approach allows a very efficient and simple
parallelization without communication overhead
(see the PERL scripts rst-gcem-xyz.pl for an example).
The postprocessing run reads the restart file and
the results of the continuum electrostatics calculations
and creates the necessary input for GMCT.
The postprocessing can also be run several times
alternating with recomputation of the electrostatic
energies to eliminate energetically unfavorable conformers.
Thereby, the continuum electrostatics calculations
are refined resulting in a decrease of the error introduced
by the inflation of the low dielectric molecule interior
by the high-energy conformers.


Hen-egg white lysozyme (HEWL)
=============================

Hen-egg white lysozyme (HEWL is often used as
model system for protonation state calculations.
See for example:

Bashford, D.; and Karplus, M.; Biochemistry,
1990, 29, 10219-10225

Webb, H.; Tynan-Connolly, B. M.; Lee, G. M.;
Farrell, D.; O'Meara, F.; Sndergaard, C. R.;
Teilum, K.; Hewage, C.; McIntosh, L. P. and
Nielsen, J. E., Proteins 2010, 79, 685-702.

The example system was prepared from PDB 4LZT
with our in-house program GCEM based on the MEAD
library (Ullmann, R. T. and Ullmann G. M., to
be published, see also the program manual of GMCT
for explanation of the energy terms).

The exceptional pH and stability of the structure
of HEWL makes is well suited for calculations with
a continuum electrostatics model. Very good agreement
with the experimental overall pH-titration curve
and site-specific pK(1/2) values from nuclear
magnetic resonance (NMR) experiments can be achieved.
(see the corresponding HEWL example provided with
the GMCT distribution). A single conformer was used
for each site as found in the crystal structure.
A detailed model was chosen for all titratable
sites with explicit positions for each of the
dissociable hydrogen atoms. All aspartate,
glutamate, and the C-terminus are modeled with
four tautomers for the protonated form
(cis and trans configuration of the titratable
proton at each carboxyl oxygen) and a single
deprotonated form. All lysine residues and the
N-terminus are modeled with a single instance
for the protonated form and three rotamers for
the deprotonated form (loss of each of the three
hydrogen atoms of the amino group). The argi-
nine residues are modeled with a single proto-
nated instance and five tautomers of the de-
protonated form (loss of each hydrogen atom of
the guanidino group). The tyrosine residues
are modeled with a single instance for the de-
protonated form and three rotamers for the
dissociable hydrogen atom of the hydroxyl
group. The histidine residue is modeled with
four instances corresponding to the fully de-
protonated imidazolate form, the singly proto-
nated delta and epsilon tautomers and the
fully protonated imidazolium form. This detailed
model yields pK1/2 values that compare very
well to the experimental ones (see for example
Nielsen2010). The increased detail relative
to traditional two-state models for proto-
nation state calculation is of course paid
with a somewhat higher computational cost.

The preprocessing run and the solver jobs are started
with the script run.csh (see comments in the script
for explanation of input files and program options).

./run.csh

The postprocessing run is again started with run.csh
after all continuum electrostatics jobs have finished.

./run.csh

The program output is written to d.gcem.out
for both runs.

With the standard settings, the example takes
less than 30 min on an Intel Core i7 860 using
8 threads.

bacteriorhodopsin from Halobacterium salinarum (bR)
===================================================

bacteriorhodopsin (bR) from Halobacterium salinarum
is a prototypical system of bioenergetics.

This example shows how to account for a membrane
environment. The influence of an electrochemical
transmembrane potential is considered (a concentration
difference and an electrostatic transmembrane potential).

See:
Calimet, N.; Ullmann, G. M., J. Mol. Biol.,
2004, 339, 571-589
and
Bombarda, E.; Becker, T.; Ullmann, G. M.,
J. Am. Chem. Soc., 2006, 128, 12129-12139.

Here, we consider a single monomer of bR in its membrane
environment. A simple model with two protonation forms
for each protonatable site (Asp, Glu, Lys, Arg, Tyr,
retinal Schiff-base / LYR) was used to comply with the
system setup used in the publications mentioned above.
The membrane is modeled with a three-slab continuum
dielectric model to mimic the hydrophilic and ion-accessible
head-group regions on each membrane side and the hydrophobic
and ion-inaccessible core region formed by the lipid tails.
Protein cavities are excluded from the membrane slabs and
filled with solvent dielectric. See also the manual of GMCT
for explanation of the energy terms and the corresponding
bR example for GMCT included included in the GMCT
distribution for the Monte Carlo simulation part.

The preprocessing run and the solver jobs are started
with the script run.csh (see comments in the script
for explanation of input files and program options).

./run.csh

The postprocessing run is again started with run.csh
after all continuum electrostatics jobs have finished.

./run.csh

The program output is written to d.gcem.out
for both runs.

With the standard settings, the example takes
about 45 min on an Intel Core i7 860 using
8 threads.
