Colvars and COM restraints for protein+membrane; plus a minor bug

From: Marc Baaden (baaden_at_smplinux.de)
Date: Thu Sep 23 2010 - 09:24:59 CDT

Hi,

Following some exchange with a colleague using geometrical plane and cylinder
restraints in Charmm to keep membrane and protein from drifting too much, I
am trying to implement something similar within NAMD. Some discussion on
the mailing list suggested that the colvars module could be used for
this. Unfortunately I couldn't find any specific example, and so I'd
like to cross-check whether I went for the right implementation, and
also post my implementation if anybody wants to use it.

In the NAMD configuration file I added the following two lines:

colvars on
colvarsConfig colvars.namd

And the colvars.namd file contains:

colvar {
  name memb_plane
  distanceZ {
    axis (0.0, 0.0, 1.0)
    main { # add all the P's within residues 1 to 253 of segment "POPC"
                        atomNameResidueRange P 1-253
                        psfSegID POPC }
    ref { dummyAtom ( 0.0, 0.0, -30.0 ) }
  }
}

colvar {
  name prot_cyl
  distanceXY {
    axis (0.0, 0.0, 1.0)
    main { atomNumbers { 5 21 34 .. }
psfSegID PEPA
 }
    ref { dummyAtom ( 0.0, 0.0, 0.0 ) }
  }
}

harmonic {
  colvars memb_plane prot_cyl
  centers 0.0 0.0
  forceConstant 5.0
}

The purpose of the "memb_plane" restraint is to keep the membrane
centered at -30 on the z-axis. The purpose of the "prot_cyl" restraint
is to keep the protein around (0,0) on the xy axis. Is this the way to
do it?

I also noticed a bug in the definition of the atom groups. I have a
multimeric protein, with each chain called "PEPA", "PEPB",.. ; I want
all Calpha atoms of all chains, and don't want to use the "psfSegID"
specification. But if I don't include it, NAMD complains and doesn't
start. On the other hand, the atom numbers seem to be properly taken
into account (eg also those not in segID PEPA are counted). THis seems
to be a but to me.

It's also a bit unfortunate that one cannot use several segment
identifiers as in all Calphas of PEPA, PEPB, etc.

I also wonder whether it would make a difference to include all lipid
(respective protein) atoms in the atom groups (this is how it was setup
in Charmm), or only the P and Calpha atoms as done here. In particular,
does the forceConstant have to be scaled according to the number of
particles (I don't think so)? Also, it seems that using too large an
atom group might hurt performance. Any experience/numbers available?

Thanks in advance for any suggestions,
Marc Baaden

-- 
 Dr. Marc Baaden  - Institut de Biologie Physico-Chimique, Paris
 mailto:baaden_at_smplinux.de      -      http://www.baaden.ibpc.fr
 FAX: +33 15841 5026  -  Tel: +33 15841 5176  ou  +33 609 843217

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:54:32 CST