# Simulation of a short coarse-grained peptide.
# (details omitted. someone else's research in progress)
# ----------------- Init Section -----------------

units           real
atom_style      full

bond_style      harmonic
angle_style     harmonic
dihedral_style  charmm

pair_style hybrid hbond/dreiding/morse 2 9 11 90 lj/charmm/coul/charmm 9 11

special_bonds   lj 0.0 0.0 1.0

# ----------------- Atom Definition Section -----------------

read_data "hybrid_pair_example.data"

# ----------------- Settings Section -----------------

pair_coeff 1 1   lj/charmm/coul/charmm  0.6 2.0
pair_coeff 2 2   lj/charmm/coul/charmm  0.6 2.0
pair_coeff 3 3   lj/charmm/coul/charmm  0.6 2.0
pair_coeff 4 4   lj/charmm/coul/charmm  3.0 3.6
pair_coeff 5 5   lj/charmm/coul/charmm  0.6 3.6
pair_coeff 6 6   lj/charmm/coul/charmm  0.6 2.0
pair_coeff 7 7   lj/charmm/coul/charmm  0.6 2.0
pair_coeff 8 8   lj/charmm/coul/charmm  0.6 2.0
pair_coeff 9 9   lj/charmm/coul/charmm  0.6 2.0
pair_coeff 10 10 lj/charmm/coul/charmm  0.6 2.0
pair_coeff 11 11 lj/charmm/coul/charmm  0.6 2.0
pair_coeff 12 12 lj/charmm/coul/charmm  0.6 3.6
pair_coeff 13 13 lj/charmm/coul/charmm  0.6 3.6

# Mixing rules are used for interactions between (nearly) all pairs of unlike 
# atoms.  However we want to override these mixing rules for interactions 
# between atoms of type 1 and 3:

pair_coeff 1 3   hbond/dreiding/morse 2 i 3.880 1.7241 2.90 2 9 11.0 90

bond_coeff  1     60.0      3.944616585
bond_coeff  2     60.0      4.422668877
bond_coeff  3     200.0     1.0
bond_coeff  4     60.0      2.828427125
bond_coeff  5     60.0      6.0
angle_coeff 1     60.0     64.32004681
angle_coeff 2     60.0     90.0
angle_coeff 3     60.0    111.00917051
angle_coeff 4     200.0   135.0
angle_coeff 5     200.0    90.0
angle_coeff 6     0.0     140.243257095
angle_coeff 7     0.0     149.534455081
dihedral_coeff 1   -12.0 1  113  0
dihedral_coeff 2   -12.0 1 -162  0
dihedral_coeff 3   -12.0 1 -161  0
dihedral_coeff 4   -12.0 1  60   0
dihedral_coeff 5   -12.0 1 -50   0

group ghydrophobic type 4
group gdonors      type 1
group gacceptors   type 3


# ----------------- Run Section -----------------

timestep        2.0
dump            1 all custom 10000 traj_nvt.lammpstrj id mol type x y z ix iy iz

# To use Langevin dynamics in LAMMPS you need both "fix langevin" and "fix nve".
# (See http://lammps.sandia.gov/doc/fix_langevin.html for details.)

fix fxlan all langevin 300.0 300.0 5000.0 48279
fix fxnve all nve

# Calculate the Hbond energy (HBOND) and the Hydrophobic sidechain energy (HYD)
compute HBOND gdonors group/group gacceptors
compute HYD ghydrophobic group/group ghydrophobic
variable vHBOND equal c_HBOND
variable vHYD equal c_HYD
variable vtime equal step

fix fPrintE all print 10000 "${vtime} ${vHBOND} ${vHYD}" file E_hbond_hyd.dat


thermo_style    custom step temp pe etotal press vol epair ebond eangle edihed
thermo          10000  # time interval for printing out "thermo" data

run		20000000

write_restart  system_after_nvt.rst