Re: Pressure Discrepancy

From: Morad Alawneh (alawneh_at_chem.byu.edu)
Date: Wed Mar 28 2007 - 13:55:43 CDT

*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