Re: Pressure Discrepancy

From: Peter Freddolino (petefred_at_ks.uiuc.edu)
Date: Wed Mar 28 2007 - 14:00:56 CDT

Since you're measuring something that's strongly reliant on the pressure
coupling, yes, you should...
Peter

Morad Alawneh wrote:
> *Hi Peter,
>
> Thanks for your suggestion, I build my system and followed the
> following tutorial instruction:
> http://www.ks.uiuc.edu/Research/namd/tutorial/NCSA2002/hands-on/
> In that tutorial they suggest for equilibration to use
> *
> *langevinPistonPeriod 200
> langevinPistonDecay 100*
> *
> **whereas for NPT production to use the following with a commented line
> **# use lighter damping now that system is equilibrated*
> *langevinPistonPeriod 200
> ***langevinPistonDecay 500***
> *
>
> *Should I follow your suggestion even I have a very long equilibration
> of about ~ 25 ns?
>
> Thanks
>
> Morad Alawneh
> *
>
>
>
> Peter Freddolino wrote:
>> Hi Morad,
>> I saw that you had a piston period of 200 fs, but decay time of 500 fs.
>> As noted at http://www.ks.uiuc.edu/Research/namd/2.6/ug/node32.html,
>> this will result in underdamped dynamics, and thus pressure control of
>> your system. You should try a decay time of 50 or 100 fs and see if
>> that stabilizes your pressure.
>> Peter
>>
>> Morad Alawneh wrote:
>>
>>> *No, I do not.
>>> I have attached the NAMD script file to those who want to check my the
>>> input parameters.
>>>
>>> Thanks
>>>
>>> Morad Alawneh
>>> *
>>> Marcos Sotomayor wrote:
>>>
>>>> Are you using multiple time stepping?
>>>>
>>>> Marcos
>>>>
>>>> On Wed, 28 Mar 2007, Morad Alawneh wrote:
>>>>
>>>>
>>>>> *There is one suggestion, which is to use NVT ensemble, by having
>>>>> different box sizes and choose the box size that gives you the right
>>>>> pressure. This suggestion is fine if you can use NVT ensemble, which
>>>>> varies the pressure and keeps the box size constant.
>>>>>
>>>>> Unfortunately, this suggestion is not suitable for my case, since I
>>>>> have
>>>>> to use NPAT ensemble.
>>>>>
>>>>> I think there are some bugs in the NAMD algorithms, since I found other
>>>>> two problems. I was not able to constrain the com of my protein around
>>>>> the origin for whatever the force constant I use in both SMD or Free
>>>>> Energy of Conformational Change Calculations.
>>>>>
>>>>> I sent a copy of my notes to James Phillips, NAMD developer, hopping to
>>>>> here back from him soon.
>>>>>
>>>>> Thanks
>>>>>
>>>>> Morad Alawneh
>>>>> **
>>>>> *
>>>>> tamal_at_iitk.ac.in wrote:
>>>>>
>>>>>> I am also finding the same problem.It never relates with our input
>>>>>> pressure.Have you got any reply?.
>>>>>> Thanks for your referred paper in your mail.It is really a
>>>>>> fundamental paper.
>>>>>>
>>>>>> thanks
>>>>>> Tamal Banerjee
>>>>>>
>>>>>>
>>>>>>
>>>>>>> *Dear NAMD Developers and users,
>>>>>>>
>>>>>>> I am trying to evaluate the surface tension of my system (a protein
>>>>>>> surrounded by a lipid bilayer in water and salt) using this equation:
>>>>>>> 0.5*Lz*(Pzz - Pt), where
>>>>>>>
>>>>>>> Lz: the box length in z direction, which is the reaction coordinate.
>>>>>>> Pzz: the pressure sensor in z direction
>>>>>>> pt = 0.5*(Pxx + Pyy)
>>>>>>>
>>>>>>> I am using the NPAT ensemble with the following settings in NAMD2.6:
>>>>>>>
>>>>>>> # Langevin Dynamics Parameters
>>>>>>> langevin on
>>>>>>> langevinTemp 303.15
>>>>>>> langevinDamping 5.0
>>>>>>> langevinHydrogen on
>>>>>>> # Pressure Control
>>>>>>> useGroupPressure yes
>>>>>>> useFlexibleCell yes
>>>>>>> useConstantRatio no
>>>>>>> useConstantArea yes
>>>>>>> # Nose-Hoover Langevin Piston Pressure Control
>>>>>>> LangevinPiston on
>>>>>>> LangevinPistonTarget 1.01325 ;# 1 atm 1.01325 bar
>>>>>>> LangevinPistonPeriod 200.0
>>>>>>> LangevinPistonDecay 500.0
>>>>>>> LangevinPistonTemp 303.15
>>>>>>> # Periodic Boundary Conditions
>>>>>>> cellBasisVector1 40.0 00.0 00.0
>>>>>>> cellBasisVector2 00.0 40.0 00.0
>>>>>>> cellBasisVector3 00.0 00.0 90.0
>>>>>>> cellOrigin 0.0 0.0 0.0
>>>>>>> wrapWater on
>>>>>>> wrapAll on
>>>>>>> wrapNearest off
>>>>>>>
>>>>>>> According to the article
>>>>>>> ** ZHANG YH, FELLER SE, BROOKS BR, et al.
>>>>>>> COMPUTER-SIMULATION OF LIQUID/LIQUID INTERFACES .1. THEORY AND
>>>>>>> APPLICATION TO OCTANE/WATER
>>>>>>> <http://apps.isiknowledge.com/WoS/CIW.cgi?SID=2EPNc3mJfLCGDA6hnmm&Func=Abstract&doc=9/5>
>>>>>>>
>>>>>>> JOURNAL OF CHEMICAL PHYSICS 103 (23): 10252-10266 DEC 15 1995 *
>>>>>>> *Pzz should take the value of Pn, which is ~ 1 bar, but when I tested
>>>>>>> its values it had different values.
>>>>>>> If I take the average of those pressure values, they are far away
>>>>>>> from 1
>>>>>>> bar, here is an example:
>>>>>>>
>>>>>>> Pxx Pyy Pzz Pt
>>>>>>> P <P>
>>>>>>> 720.2230 415.1900 -427.7030 567.7070 235.9037
>>>>>>> 235.9037
>>>>>>> 196.3230 577.8960 266.0390 387.1090 346.7525
>>>>>>> 107.5640
>>>>>>> ............... ............... ...............
>>>>>>> ................ ............... .............
>>>>>>> **simulation average: -16.4359 -112.3742 -80.3948
>>>>>>> -54.0773
>>>>>>>
>>>>>>> I know there is a huge fluctuations in pressure, but the average
>>>>>>> values
>>>>>>> should take care of that and be close to 1 bar for Pzz, P, or <P>.
>>>>>>>
>>>>>>> Would you explain to me why do I have this kind of pressure
>>>>>>> discrepancy
>>>>>>> for Pzz, P, or <P>?
>>>>>>>
>>>>>>> Thanks
>>>>>>>
>>>>>>> Morad Alawneh
>>>>>>>
>>>>>>>
>>>>>>> **
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> *
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>> *****************************************
>>>>>> Tamal Banerjee
>>>>>> Research Scholar
>>>>>> Department of Chemical Engineering
>>>>>> I.I.T Kanpur
>>>>>> Phone: +91-512-2597533 (Lab)
>>>>>> +91-512-2595433 (Residence)
>>>>>> URL: home.iitk.ac.in/~tamal
>>>>>> *****************************************
>>>>>>
>>>>>>
>>>>>>
>>> ------------------------------------------------------------------------
>>>
>>> # Equilibrate the system under NPAT
>>>
>>> set InputName sys
>>> set Process0 equil10
>>> set Process prod
>>> set LastRun 0
>>> set OutputName [concat "$InputName\_$Process"]
>>> #set FixedAtoms fixed.pdb
>>> #set Restraints [concat "restraints_$Process0\.pdb"]
>>> set ESF [concat "../prep/$InputName\_$Process0\.xsc"]
>>> set Numer_of_Steps 5000000
>>> set FirstTimeStep [expr {$Numer_of_Steps*$LastRun}]period
>>> set Temp 303.15
>>> set Dist_Temp 303.15
>>> #set Gamma ??.0 ;# dyn/cm
>>> ####################################################################################
>>> # Input Files
>>>
>>> structure "../psfgen/$InputName.psf"
>>> coordinates "../psfgen/$InputName.pdb"
>>> bincoordinates [concat "../prep/$InputName\_$Process0\.coor"]
>>> binvelocities [concat "../prep/$InputName\_$Process0\.vel"];# use this or
>>> ;# temperature
>>> parameters "../psfgen/par_all27_prot_lipid.prm"
>>> paraTypeCharmm on
>>> ####################################################################################
>>> # Output Files
>>>
>>> outputname $OutputName
>>> binaryoutput yes ;# default
>>> #restartname res_$OutputName ;# produces pdb restart file
>>> #restartfreq 100 ;# generate restart file every N timesteps
>>> DCDfile $OutputName.dcd
>>> DCDfreq 500
>>> ####################################################################################
>>> # Standard Output
>>>
>>> outputEnergies 100 ;# timesteps between energy output
>>> outputMomenta 100 ;# timesteps between momenta output
>>> outputPressure 100 ;# timesteps between pressure output
>>> outputTiming 100 ;# timesteps between timing output
>>> ####################################################################################
>>> # Timestep Parameters
>>>
>>> numsteps $Numer_of_Steps ;# run stops when this step is reached
>>> timestep 2.0 ;# fs
>>> firsttimestep $FirstTimeStep ;# use this instead of restart from 0
>>> stepspercycle 20 ;# default
>>> ####################################################################################
>>> # Simulation Space Partitioning
>>>
>>> switching on
>>> switchDist 8.0
>>> cutoff 12.0
>>> pairlistdist 14.0
>>> margin 2.0
>>> pairlistsPerCycle 4
>>> ####################################################################################
>>> # Basic Dynamics
>>>
>>> #temperature $Temp ;# use this or binvelocities
>>> exclude scaled1-4 ;# specified by CHARMM
>>> 1-4scaling 1.0 ;# specified by CHARMM
>>> seed 2006
>>> rigidBonds all ;# activate ShakeH
>>> rigidTolerance 0.00000001 ;# default
>>> rigidIterations 100 ;# default
>>> useSettle on ;# default
>>> dielectric 1.0 ;# default
>>> COMmotion no ;# default
>>> zeroMomentum yes
>>> ####################################################################################
>>> # PME Parameters
>>>
>>> PME yes
>>> PMEGridSizeX 64 ;# Should be > 1.5 the box width
>>> PMEGridSizeY 64 ;# Should be > 1.5 the box height
>>> PMEGridSizeZ 256 ;# Should be > 1.5 the box length
>>> PMEInterpOrder 6
>>> ####################################################################################
>>> # Multiple Timestep Parameters
>>>
>>> fullElectFrequency 2
>>> nonbondedFreq 1
>>> MTSAlgorithm impulse ;# default
>>> longSplitting c1 ;# default
>>> ####################################################################################
>>> # Harmonic Constraint Parameters
>>>
>>> #constraints on
>>> #consexp 2 ;# default, use only even integers
>>> #consref $Restraints
>>> #conskfile $Restraints
>>> #conskcol B
>>> #selectConstraints on ;# Planer restraint
>>> #selectConstrZ on
>>> ####################################################################################
>>> # Fixed Atoms Parameters
>>>
>>> #fixedAtoms on
>>> #fixedAtomsForces on ;# keep it off at CP
>>> #fixedAtomsFile $FixedAtoms
>>> #fixedAtomsCol B
>>> ####################################################################################
>>> # Temperature Control
>>>
>>> # Langevin Dynamics Parameters
>>> langevin on
>>> langevinTemp $Temp
>>> langevinDamping 5.0 ;# /ps
>>> langevinHydrogen on ;# default
>>>
>>> # Temperature Coupling Parameters
>>> #tCouple on
>>> #tCoupleTemp $Temp
>>> #tCoupleFile
>>> #tCoupleCol B
>>>
>>> # Temperature Rescalling Parameters
>>> #rescaleFreq 50
>>> #rescaleTemp $Temp
>>>
>>> # Temperature Reassignment Parameters, usefull for annealing
>>> #reassignFreq 50
>>> #reassignTemp $Temp
>>> #reassignIncr 5.0
>>> #reassignHold $Dist_Temp
>>> ####################################################################################
>>> # Boundary Conditions
>>>
>>> # Spherical Harmonic Conditions
>>>
>>> # Cylindrical Harmonic Conditions
>>>
>>> # Periodic Boundary Conditions
>>> cellBasisVector1 40.0 00.0 00.0
>>> cellBasisVector2 00.0 40.0 00.0
>>> cellBasisVector3 00.0 00.0 90.0
>>> cellOrigin 0.0 0.0 0.0 ;# default
>>> wrapWater on
>>> wrapAll on
>>> wrapNearest off ;# use for non-rectangular cells
>>> extendedSystem $ESF
>>> ####################################################################################
>>> # Pressure Control
>>>
>>> useGroupPressure yes
>>> useFlexibleCell yes ;# anisotropic cell fluctuations
>>> useConstantRatio no ;# constant shape in x-y plane
>>> useConstantArea yes ;# constant area; constant dimentions in x-y
>>> ;#plane
>>>
>>> # Nose-Hoover Langevin Piston Pressure Control
>>> LangevinPiston on
>>> LangevinPistonTarget 1.01325 ;# 1 atm 1.01325 bar
>>> LangevinPistonPeriod 200.0 ;# fs
>>> LangevinPistonDecay 500.0 ;# fs
>>> LangevinPistonTemp $Temp
>>> #SurfaceTensionTarget $Gamma
>>> ####################################################################################
>>> # Interactive Molecular Dynamics
>>>
>>> IMDon on
>>> IMDport 2006
>>> IMDfreq 1
>>> IMDignore yes
>>> ####################################################################################
>>> # Energy Minimization: Conjugate Graidient Parameters
>>>
>>> #minimization on
>>> #minTinyStep 0.000001 ;# defalut
>>> #minBabyStep 0.0100 ;# default
>>> ####################################################################################
>>> # # Tcl Forces and Analysis
>>> #
>>> # set PDB [concat "../psfgen/$InputName\.pdb"]
>>> # ###----------------
>>> # set SelSegID1 {GA1}
>>> # set SelSegID2 {GA2}
>>> # set SelSegID3 {}
>>> # set SelSegID4 {}
>>> # ###----------------
>>> # set Protein_Center [list 0.0 0.0 0.0]
>>> # ###----------------
>>> # tclForces on
>>> # tclForcesScript "../psfgen/restraints.tcl"
>>> ####################################################################################
>>>
>>> #minimize 5000
>>>
>>> run $Numer_of_Steps
>>>
>>
>>

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:44:31 CST