← Back to branch summary

~ubuntu-branches/debian/sid/lammps/sid

« back to all changes in this revision

Viewing changes to tools/moltemplate/examples/coarse_grained_examples/vesicle_Brannigan2005+Bellesia2010/moltemplate_files/1beadProtSci2010.lt

  • Committer: Package Import Robot
  • Author(s): Anton Gladky
  • Date: 2015-04-29 23:44:49 UTC
  • mfrom: (5.1.3 experimental)
  • Revision ID: package-import@ubuntu.com-20150429234449-mbhy9utku6hp6oq8
Tags: 0~20150313.gitfa668e1-1
Upload into unstable.

Show diffs side-by-side

added added

removed removed

Lines of Context:
tools/moltemplate/examples/coarse_grained_examples/vesicle_Brannigan2005+Bellesia2010/moltemplate_files/1beadProtSci2010.lt
 
1
# This file defines a 4-helix bundle coarse-grained protein model (AUF2) used in
 
2
#      G. Bellesia, AI Jewett, and J-E Shea, 
 
3
#      Protein Science, Vol19 141-154 (2010)
 
4
#
 
5
# Strategy:
 
6
#
 
7
#1) First I'll define some building blocks (A16, B16, T3)
 
8
#   which are helices, sheets and turns of a predetermined length)
 
9
#
 
10
#2)  Then I'll copy and paste them together to build
 
11
#    a 4-helix bundle (or a 4-strand beta-barrel).
 
12
#    This approach is optional.  If your protein has helices which are not
 
13
#    identical, you should probably just include all 4 helices in a single
 
14
#    "Data Atoms" section and don't try to subdivide the protein into pieces.)
 
15
 
 
16
 
 
17
 
 
18
1beadProtSci2010 {  # <-- enclose definitions in a namespace for portability
 
19
 
 
20
  # A16 is a coarse-grained alpha-helix containing 16 residues (one "atom" each)
 
21
 
 
22
  A16 {
 
23
 
 
24
    #  AtomID   MoleculeID AtomType  Charge       X        Y         Z
 
25
 
 
26
    write('Data Atoms') {
 
27
      $atom:a1  $mol:... @atom:../sL  0.0      -2.4      -2.4       0.0
 
28
      $atom:a2  $mol:... @atom:../sL  0.0       2.4      -2.4       3.6
 
29
      $atom:a3  $mol:... @atom:../sH  0.0       2.4       2.4       7.2
 
30
      $atom:a4  $mol:... @atom:../sH  0.0      -2.4       2.4       10.8
 
31
      $atom:a5  $mol:... @atom:../sL  0.0      -2.4      -2.4       14.4
 
32
      $atom:a6  $mol:... @atom:../sL  0.0       2.4      -2.4       18.0
 
33
      $atom:a7  $mol:... @atom:../sH  0.0       2.4       2.4       21.6
 
34
      $atom:a8  $mol:... @atom:../sH  0.0      -2.4       2.4       25.2
 
35
      $atom:a9  $mol:... @atom:../sL  0.0      -2.4      -2.4       28.8
 
36
      $atom:a10 $mol:... @atom:../sL  0.0       2.4      -2.4       32.4
 
37
      $atom:a11 $mol:... @atom:../sH  0.0       2.4       2.4       36.0
 
38
      $atom:a12 $mol:... @atom:../sH  0.0      -2.4       2.4       39.6
 
39
      $atom:a13 $mol:... @atom:../sL  0.0      -2.4      -2.4       43.2
 
40
      $atom:a14 $mol:... @atom:../sL  0.0       2.4      -2.4       46.8
 
41
      $atom:a15 $mol:... @atom:../sH  0.0       2.4       2.4       50.4
 
42
      $atom:a16 $mol:... @atom:../sH  0.0      -2.4       2.4       54.0
 
43
    }
 
44
 
 
45
    write('Data Bonds') {
 
46
      $bond:b1 @bond:../backbone $atom:a1 $atom:a2
 
47
      $bond:b2 @bond:../backbone $atom:a2 $atom:a3
 
48
      $bond:b3 @bond:../backbone $atom:a3 $atom:a4
 
49
      $bond:b4 @bond:../backbone $atom:a4 $atom:a5
 
50
      $bond:b5 @bond:../backbone $atom:a5 $atom:a6
 
51
      $bond:b6 @bond:../backbone $atom:a6 $atom:a7
 
52
      $bond:b7 @bond:../backbone $atom:a7 $atom:a8
 
53
      $bond:b8 @bond:../backbone $atom:a8 $atom:a9
 
54
      $bond:b9 @bond:../backbone $atom:a9 $atom:a10
 
55
      $bond:b10 @bond:../backbone $atom:a10 $atom:a11
 
56
      $bond:b11 @bond:../backbone $atom:a11 $atom:a12
 
57
      $bond:b12 @bond:../backbone $atom:a12 $atom:a13
 
58
      $bond:b13 @bond:../backbone $atom:a13 $atom:a14
 
59
      $bond:b14 @bond:../backbone $atom:a14 $atom:a15
 
60
      $bond:b15 @bond:../backbone $atom:a15 $atom:a16
 
61
    }
 
62
 
 
63
  } # A16
 
64
 
 
65
 
 
66
  T3 { # T3 is a "turn" region consisting of 3 beads
 
67
 
 
68
    #  AtomID   MoleculeID AtomType  Charge       X        Y         Z
 
69
 
 
70
    write('Data Atoms') {
 
71
      $atom:a1  $mol:... @atom:../tN  0.0      -4.8       0.0       0.0
 
72
      $atom:a2  $mol:... @atom:../tN  0.0       0.0       3.3     -1.44
 
73
      $atom:a3  $mol:... @atom:../tN  0.0       4.8       0.0       0.0
 
74
    }
 
75
 
 
76
    write('Data Bonds') {
 
77
      $bond:b1 @bond:../backbone $atom:a1 $atom:a2
 
78
      $bond:b2 @bond:../backbone $atom:a2 $atom:a3
 
79
    }
 
80
 
 
81
  } # T3
 
82
 
 
83
 
 
84
  # ----- Now build a larger molecule using A16 and T3 -------
 
85
 
 
86
  # Create a 4-Helix bundle.
 
87
  # In this version, the hydrophobic beads are poing outward.
 
88
  # I oriented them this way because I want to place this protein in a membrane.
 
89
  # (There is another file in this directory containing alternate version
 
90
  #  of this same molecule with the hydrophobic beads pointing inward.)
 
91
 
 
92
  4HelixInsideOut {
 
93
    helix1 = new A16.rot(-225, 0,0,1).move(-5.70,-5.70,-32.4)
 
94
    helix2 = new A16.rot(-135, 0,0,1).move( 5.70,-5.70,-28.8)
 
95
    helix3 = new A16.rot( -45, 0,0,1).move( 5.70, 5.70,-25.2)
 
96
    helix4 = new A16.rot(  45, 0,0,1).move(-5.70, 5.70,-21.6)
 
97
 
 
98
    turn1  = new T3.rot(180,1,0,0).rot(-20,0,1,0).rot( 10,0,0,1).move(0.78,-4.2, 27.9)
 
99
    turn2  = new T3.rot(-10,1,0,0).rot( 20,0,1,0).rot(-70,0,0,1).move(4.55, 2.4,-33.0)
 
100
    turn3  = new T3.rot(180,1,0,0).rot(-20,0,1,0).rot(190,0,0,1).move(-0.78,4.2, 34.2)
 
101
 
 
102
    write('Data Bonds') {
 
103
      $bond:turn1a  @bond:../backbone   $atom:turn1/a1   $atom:helix1/a16
 
104
      $bond:turn1b  @bond:../backbone   $atom:turn1/a3   $atom:helix2/a16
 
105
      $bond:turn2a  @bond:../backbone   $atom:turn2/a1   $atom:helix3/a1
 
106
      $bond:turn2b  @bond:../backbone   $atom:turn2/a3   $atom:helix2/a1
 
107
      $bond:turn3a  @bond:../backbone   $atom:turn3/a1   $atom:helix3/a16
 
108
      $bond:turn3b  @bond:../backbone   $atom:turn3/a3   $atom:helix4/a16
 
109
    }
 
110
    create_var { $mol } # molecule ID number shared by all atoms in this protein
 
111
 
 
112
  } # 4HelixInsideOut
 
113
 
 
114
 
 
115
  # -------- Minor coordinates adjustment: -----------
 
116
 
 
117
  # Those coordinates in the commands above are a little too large.
 
118
  # To make it easier to type them in, I was using sigma=6.0 Angstroms.
 
119
  # Instead, here I'll try using sigma=5.5 Angstroms.   5.5/6 = 0.916667)
 
120
 
 
121
  4HelixInsideOut.scale(0.9166666666666666)
 
122
 
 
123
  # Note: "scale()" only effects the initial coordinates of
 
124
  #       the molecule, not the force field parameters.
 
125
  #       (If you plan to minimize the molecule, you don't need to 
 
126
  #        be so careful about the initial coordinates.  In that case,
 
127
  #        you don't have worry about "scale()".  Feel free to remove.)
 
128
 
 
129
 
 
130
 
 
131
  #  -------------- Force-Field Parameters ------------
 
132
 
 
133
  #  Units and force-field styles for this protein model
 
134
  #  (These can be overridden later.)
 
135
 
 
136
  write_once("In Init") {
 
137
    units           real
 
138
    atom_style      full
 
139
    bond_style      hybrid harmonic
 
140
    angle_style     hybrid harmonic
 
141
    dihedral_style  hybrid fourier
 
142
    pair_style hybrid lj/charmm/coul/charmm/inter es4k4l maxmax 21.0 24.0
 
143
    pair_modify     mix arithmetic
 
144
    special_bonds   lj 0.0 0.0 1.0   #(turn on "1-4" interactions)
 
145
  }
 
146
 
 
147
  #  --- Distance Units ---
 
148
  # In this version of the model, sigma (the bond-length 
 
149
  # and particle diameter) is rounded to 5.5 Angstroms.
 
150
  #
 
151
  #  --- Energy & Temperature Units ---
 
152
  #   In this protein model, "epsilon" represents the free energy 
 
153
  # bonus for bringing two hydrophobic amino acids together.  
 
154
  # Here I choose to set epsilon to 1.806551818181818 kCal/mole.
 
155
  # This value was chosen so that a temperature of 300 Kelvin lies at 
 
156
  # 0.33 epsilon, which is the unfolding temperature of the marginally stable 
 
157
  # "ASF1" protein model from the Bellesia et al 2010 paper.
 
158
  #   This choice insures that both the "ASF1" model from that paper,
 
159
  # as well as the much more stable "AUF2" protein we use here (which
 
160
  # unfolds at 0.42*eps) should definitely remain stable at 300 degrees Kelvin,
 
161
  # in the bulk at least.  (However it's not clear that these energy
 
162
  # parameters will work well for a protein in membrane.  Perhaps I'll
 
163
  # run some tests and fine tune these parameters for this scenario.)
 
164
 
 
165
 
 
166
  # 2-body (non-bonded) interactions:
 
167
  #
 
168
  #   Uij(r) = 4*eps_ij * (K*(sig_ij/r)^12 + L*(sig_ij/r)^6)
 
169
  #
 
170
  #              i       j            pairstylename            eps sig K L
 
171
  #
 
172
  write_once("In Settings") {
 
173
    pair_coeff @atom:sH @atom:sH lj/charmm/coul/charmm/inter  1.8065518 5.5 1 -1
 
174
    pair_coeff @atom:sL @atom:sL lj/charmm/coul/charmm/inter  1.8065518 5.5 1 0
 
175
    pair_coeff @atom:tN @atom:tN lj/charmm/coul/charmm/inter  1.8065518 5.5 1 0
 
176
  }
 
177
 
 
178
  # The exact value of the bond_coeff does not matter too much as long as 
 
179
  # it is "stiff enough".  Here I use a softer bond spring than the one 
 
180
  # used in the paper so that I can increase the time step.
 
181
  # I also use a relatively soft spring to constrain the bond angles.
 
182
 
 
183
  #   bond_coeff     bondType            bondstylename    k    r0
 
184
 
 
185
  write_once("In Settings") {
 
186
    bond_coeff @bond:1beadProtSci2010/backbone harmonic  10.0  5.5
 
187
  }
 
188
 
 
189
 
 
190
  # angleType       atomtypes1 2 3   bondtypes1 2
 
191
 
 
192
  write_once("Data Angles By Type") {
 
193
    @angle:backbone    @atom:*  @atom:*  @atom:*   @bond:*  @bond:*
 
194
  }
 
195
 
 
196
  #   angle_coeff  angleType  anglestylename   k   theta0
 
197
  write_once("In Settings") {
 
198
    angle_coeff @angle:backbone   harmonic   100.0  105.0
 
199
  }
 
200
 
 
201
 
 
202
  #   dihedralType         atomtypes1 2 3 4                  bondtypes1 2 3
 
203
 
 
204
  write_once("Data Dihedrals By Type") {
 
205
    # For a chain of sH and sL atoms, use the @dihedral:delta65_0 
 
206
    # parameters.  (This corresponds to the "AUF2" model from the
 
207
    # Bellesia et. al 2010 paper.)
 
208
 
 
209
    @dihedral:delta65_0  @atom:s* @atom:s* @atom:s* @atom:s*   *  *  *
 
210
 
 
211
    # If "tN" (turn) atoms are present, use the @dihedral:turn parameters
 
212
 
 
213
    @dihedral:turn       @atom:tN @atom:*  @atom:*  @atom:*    *  *  *
 
214
  }
 
215
 
 
216
  write_once("In Settings") {
 
217
    dihedral_coeff @dihedral:delta60_0 fourier 2 2.167862 3 0 2.167862 1 -60.0
 
218
    dihedral_coeff @dihedral:delta65_0 fourier 2 2.167862 3 0 2.167862 1 -65.0
 
219
    dihedral_coeff @dihedral:turn      fourier 1 0.361310 3 0
 
220
    # Note: 2.167862=1.2*epsilon and 0.361310=0.2*epsilon.
 
221
  }
 
222
 
 
223
 
 
224
  #  --- Mass Units ---
 
225
  # Typical amino acids weigh approximately 110.0 grams/mole.  (Rounding down):
 
226
  write_once("Data Masses") {
 
227
    @atom:1beadProtSci2010/sH 100.0
 
228
    @atom:1beadProtSci2010/sL 100.0
 
229
    @atom:1beadProtSci2010/tN 100.0
 
230
  }
 
231
 
 
232
} # 1beadProtSci2010 (namespace)
 
233