MBG wiki
|
RecentChanges
|
Blog
|
2024-03-29
|
2024-03-28
Editing Running a molecular dynamics simulation of a tetrapeptide in a truncated octahedral cell
----------------- ----------------- Go read & do the previous tutorials on Rop and protein G, then come back ... (read also the page on truncated octahedra) ----------------- ----------------- If you are seriously into peptides, do read [[Spoel's thesis]] ----------------- ----------------- === Peptide PDB generation === The program you need is *ribosome*. Prepare a script file with your sequence : <code> # cat > ribosome.script title FPGF default helix res phe res pro res gly res phe <CTRL-D> </code> Run the program : <code> /usr/local/bin/ribosome ribosome.script FPGF.pdb /usr/local/bin/ribosome.dat </code> This should create a PDB file with the name FPGF.pdb containing your tetrapeptide. The backbone conformation is helical (as requested in the script file), the side chains' conformations are random. === Continuing ... === Use /moleman/ to align the axes of inertia : <code> #!/bin/tcsh -f # # This will read a PDB file and rotate/translate it so that # 1. the centre of gravity will be at the origin # 2. the axes of inertia will be aligned with the orthogonal frame # lx_moleman2 >& moleman.log << eof /usr/local/xutil/moleman2.lib REad FPGF.pdb XYz ALign_inertia_axes WRite aligned.pdb quit eof exit </code> Edit the PDB file created by moleman, remove the header, give all atoms a chain identifier (e.g. 'A'). Use /psfgen/ to generate PSF and PDB files : <code> #!/bin/tcsh -f /usr/local/NAMD_2.5/psfgen >& psfgen.log << END topology top_all27_prot_na.inp segment A { pdb aligned.pdb } alias atom ILE CD1 CD coordpdb aligned.pdb A guesscoord writepsf psfgen.psf writepdb psfgen.pdb END exit </code> Use /pdbset/ to determine limits (along x, y and z) of the aligned molecule. Decide the size of the cube of water that has to be added. Use /vmd/ to add waters and ions for a truncated octahedral cell. Your vmd script should look similar to <code> # # Make water box # package require vexpr package require toctsolvate toctsolvate psfgen.psf psfgen.pdb -o hydrated -minmax {{-17 -17 -17} {17 17 17}} -b 1.80 # # Add ions to neutralise charge # package require autoionize autoionize -psf hydrated.psf -pdb hydrated.pdb -is 0.20 # # Prepare restraints files # mol load psf ionized.psf pdb ionized.pdb set all [atomselect top all] set sel [atomselect top "protein and name CA"] $all set beta 0 $sel set beta 0.5 $all writepdb restrain_ca.pdb set all [atomselect top all] set to_fix [atomselect top "protein and backbone"] $all set beta 0 $to_fix set beta 1 $all writepdb fix_backbone.pdb </code> Examine the system (ionised.pdb). Heat it up (for a real run you should increase the time spend in the heating phase) : <code> # # Input files # structure ionized.psf coordinates ionized.pdb parameters par_all27_prot_na.inp paraTypeCharmm on # # Output files & writing frequency for DCD # and restart files # outputname output/heat_out binaryoutput off restartname output/restart restartfreq 1000 binaryrestart yes dcdFile output/heat_out.dcd dcdFreq 200 DCDunitcell on # # Frequencies for logs and the xst file # outputEnergies 20 outputTiming 200 xstFreq 200 # # Timestep & friends # timestep 2.0 stepsPerCycle 8 nonBondedFreq 2 fullElectFrequency 4 # # Simulation space partitioning # switching on switchDist 8 cutoff 9 pairlistdist 11 # # Basic dynamics # temperature 0 COMmotion no dielectric 1.0 exclude scaled1-4 1-4scaling 1.0 rigidbonds all # # Particle Mesh Ewald parameters. # Pme on PmeGridsizeX 32 # <===== CHANGE ME PmeGridsizeY 32 # <===== CHANGE ME PmeGridsizeZ 27 # <===== CHANGE ME # # Periodic boundary things # wrapWater on wrapNearest on cellBasisVector1 32.00 0.00 0.00 # <===== CHANGE ME cellBasisVector2 0.00 32.00 0.00 # <===== CHANGE ME cellBasisVector3 16.00 16.00 16.00 # <===== CHANGE ME cellOrigin 0.00 0.00 0.00 # <===== CHANGE ME # # Fixed atoms for initial heating-up steps # fixedAtoms on fixedAtomsForces on fixedAtomsFile fix_backbone.pdb fixedAtomsCol B # # Restrained atoms for initial heating-up steps # constraints on consRef restrain_ca.pdb consKFile restrain_ca.pdb consKCol B # # Langevin dynamics parameters # langevin on langevinDamping 10 langevinTemp 320 # <===== Check me langevinHydrogen on langevinPiston on langevinPistonTarget 1.01325 langevinPistonPeriod 200 langevinPistonDecay 100 langevinPistonTemp 320 # <===== Check me useGroupPressure yes ########################################## # The actual minimisation and heating-up # protocol follows. The number of steps # shown below are too small for a real run ########################################## # # run one step to get into scripting mode # minimize 0 # # turn off pressure control until later # langevinPiston off # # minimize nonbackbone atoms # minimize 400 ;# <===== CHANGE ME output output/min_fix # # min all atoms # fixedAtoms off minimize 400 ;# <===== CHANGE ME output output/min_all # # heat with CAs restrained # set temp 20; while { $temp < 321 } { ;# <===== Check me langevinTemp $temp run 400 ;# <===== CHANGE ME output output/heat_ca set temp [expr $temp + 20] } # # equilibrate volume with CAs restrained # langevinPiston on run 400 ;# <===== CHANGE ME output output/equil_ca # # equilibrate volume without restraints # constraintScaling 0 run 2000 ;# <===== CHANGE ME </code> Production run : <code> # # Input files # structure ionized.psf coordinates heat_out.coor velocities heat_out.vel extendedSystem heat_out.xsc parameters par_all27_prot_na.inp paraTypeCharmm on # # Output files & writing frequency for DCD # and restart files # outputname output/equi_out binaryoutput off restartname output/restart restartfreq 1000 binaryrestart yes dcdFile output/equi_out.dcd dcdFreq 200 DCDunitcell on # # Frequencies for logs and the xst file # outputEnergies 20 outputTiming 200 xstFreq 200 # # Timestep & friends # timestep 2.0 stepsPerCycle 8 nonBondedFreq 2 fullElectFrequency 4 # # Simulation space partitioning # switching on switchDist 8 cutoff 9 pairlistdist 11 # # Basic dynamics # COMmotion no dielectric 1.0 exclude scaled1-4 1-4scaling 1.0 rigidbonds all # # Particle Mesh Ewald parameters. # Pme on PmeGridsizeX 32 # <===== CHANGE ME PmeGridsizeY 32 # <===== CHANGE ME PmeGridsizeZ 27 # <===== CHANGE ME # # Periodic boundary things # wrapWater on wrapNearest on wrapAll on # # Langevin dynamics parameters # langevin on langevinDamping 1 langevinTemp 320 # <===== Check me langevinHydrogen on langevinPiston on langevinPistonTarget 1.01325 langevinPistonPeriod 200 langevinPistonDecay 100 langevinPistonTemp 320 # <===== Check me useGroupPressure yes firsttimestep 9600 # <===== CHANGE ME run 5000000 ;# <===== CHANGE ME </code>
Summary:
This change is a minor edit.
Username:
Replace this text with a file.