This is an old revision of the document!
Calculate electrostatic potential of a bilayer with respect to z (position along bilayer normal)
Location of files on DT2: /a/fs-3/export/home/deepthought2/aleonar2/epp_bilayer
Location of files on MARCC:
1. Submit “chrg_dens.inp” using the submission file “chrg_dens.csh”. The input file <chrg_dens.inp> computes the total charge density along z (the bilayer normal) in slabs of thickness 0.1 Angstroms, from -32.95 to 32.95 Angstroms. It reads in the re-centered trajectories from the current directory with file names of the form “dyn@NFIL-recenter.dcd”. The script is set up to read inputs and trajectories from the current directory.
First update the inputs to suit your system:
- Topology file: toppar.str
- Protein structure file: step5_assembly.xplor.psf
- Coordinate file: step5_assembly.crd
- Thickness of box at start of simulation (“A”): xa
- First trajectory file number: nfil
- Last trajectory file number: In line 80, you can change “100”: “if nfil .le. 100 goto lp1”
And, in the file <chrg_dens.csh>, update the path of the CHARMM executable to your preferred CHARMM executable.
Data file: “chrg_dens.dat”
The data file produced by “chrg_dens.inp” will be called “chrg_dens.dat”. This file has two columns. The first column is the beginning of the z slab (of width 0.1 Angstrom) and the second column is the electron density computed in a frame. This file will be very large.
2. Use Matlab script “epp_bilayer.m” to process data file “chrg_dens.dat”.
The script “epp_bilayer.m” symmetrizes the electron density, integrates the electron density with respect to z, converts to volts, establishes electrostatic potential = 0 at the bilayer center, and plots the electrostatic potential as a function of z. It was written to process a pure DPPC bilayer, but will work with any system that is symmetric along z.
The vector containing the electrostatic potential relative to the bilayer center is “epp_dppc1”.