Re: ERROR: CALCULATED VDWA FORCE MAY NOT MATCH ENERGY! POSSIBLE BUG!

From: Axel Kohlmeyer (akohlmey_at_gmail.com)
Date: Mon Feb 17 2014 - 21:39:00 CST

On Mon, Feb 17, 2014 at 10:00 PM, Shalton Evans <swe31_at_msstate.edu> wrote:
> Axel Kohlmeyer,
>
>
>
> I really appreciate the response. I apologize for the lack of information. I
> am trying to simulate a half-cell lithium battery.This consists of Graphite
> sheets/Carbon(Electrode) Lithium ions, propylene carbonate molecules and
> ethylene carbonate molecules(electrolyte).The only force field data I have
> is within my 'prmtop' file. I used the GAFF force field to build this file
> (attached). Below is my configuration file of how I produced the results of
> which I am inquiring. Any further help would be appreciated.

since this is an amber input, i would suggest to first test this with
amber's MD codes and see if that works at all. but before doing that,
i would made some tests with individual components of your system
independently. you have a rather complex and complicated setup, this
practically never works immediately, particularly if you have limited
experience with the code. some of your parameter choices look odd,
e.g. your cutoffs seem awfully short. is there a reason for that?

apart from that, i would worry that your choice of force field may
lead to significant errors. have you made a survey of the published
literature to find out what kind of force field people have used for
similar studies?

axel.

>
>
> Thank you,
>
> Shalton Evans
>
>
> #############################################################
>
> ## JOB DESCRIPTION ##
>
> #############################################################
>
>
> # Minimization and Equilibration of
>
> # Ubiquitin in a Water Box
>
>
>
> #############################################################
>
> ## ADJUSTABLE PARAMETERS ##
>
> #############################################################
>
>
>
> set temperature 300
>
> set outputname PD.10000skip.25GHz
>
> set inputname all.big.fixed.atoms.consP.sym3
>
>
> firsttimestep 0
>
>
>
> #############################################################
>
> ## SIMULATION PARAMETERS ##
>
> #############################################################
>
>
> # Input
>
> amber on
>
> #ambercoor all.big.fixed.atoms.consP.Z.rst
>
> coordinates min1.pdb
>
> bincoordinates $inputname.restart.coor
>
> #binvelocities $inputname.restart.vel
>
> parmfile all.big.10Li2PF6.prmtop
>
> temperature $temperature
>
> ExtendedSystem $inputname.restart.xsc
>
>
> # Force-Field Parameters
>
> exclude scaled1-4
>
> 1-4scaling 1.0
>
> cutoff 6.1
>
> switching on
>
> switchdist 6.
>
> pairlistdist 8.5
>
>
>
> # Integrator Parameters
>
> timestep 2.0 ;# 2fs/step
>
> rigidBonds all ;# needed for 2fs steps
>
> nonbondedFreq 1
>
> fullElectFrequency 2
>
> stepspercycle 10
>
>
>
> # Constant Temperature Control
>
> langevin on ;# do langevin dynamics
>
> langevinDamping 5 ;# damping coefficient (gamma) of 5/ps
>
> langevinTemp $temperature
>
> langevinHydrogen off ;# don't couple langevin bath to hydrogens
>
>
>
> # Periodic Boundary Conditions
>
> #cellBasisVector1 15.6893 0. 0.
>
> #cellBasisVector2 0. 29.01327 0.
>
> #cellBasisVector3 0. 0. 81.4132
>
> #cellOrigin 6.94532728195 12.9714078903 55.0701789856
>
>
> wrapAll on
>
>
> # PME (for full-system periodic electrostatics)
>
> PME yes
>
> PMEGridSizeX 16
>
> PMEGridSizeY 48
>
> PMEGridSizeZ 48
>
>
>
> # Constant Pressure Control (variable volume)
>
> #useGroupPressure yes ;# needed for rigidBonds
>
> #useFlexibleCell yes
>
> #useConstantArea yes
>
>
> #langevinPiston on
>
> #langevinPistonTarget 1.01325 ;# in bar -> 1 atm
>
> #langevinPistonPeriod 100.
>
> #langevinPistonDecay 50.
>
> #langevinPistonTemp $temperature
>
>
> # Electric Field
>
> eFieldon yes
>
> eField 0.0 0.0 0.0
>
>
> #Output Name
>
> outputName $outputname
>
>
> restartfreq 500 ;# 500steps = every 1ps
>
> dcdfreq 250
>
> xstFreq 250
>
> outputEnergies 100
>
> outputPressure 100
>
>
>
> #############################################################
>
> ## EXTRA PARAMETERS ##
>
> #############################################################
>
> #Constrained atoms
>
> constraints on
>
> consref min1.fixed.pdb
>
> conskfile min1.fixed.pdb
>
> conskcol O
>
> # Fixed Atoms
>
> #fixedAtoms on
>
> #fixedAtomsFile min1.fixed.pdb
>
> #fixedAtomsCol O
>
> #fixedAtomsForces on
>
>
> #############################################################
>
> ## EXECUTION SCRIPT ##
>
> #############################################################
>
>
>
> #minimize 100
>
> reinitvels $temperature
>
>
> set period [expr 2*$halfperiod]
>
> for { set t 1 } { $t <= 1000000000 } { incr t 1 } {
>
> if { ([getstep] - $timer1) % $period == 0 } { eField $A 0.0 0.0 }
>
> if { ([getstep] - $timer1) % $period == $halfperiod } { eField -$A 0.0
> 0.0 }
>
> run 10
>
> }
>
>

-- 
Dr. Axel Kohlmeyer  akohlmey_at_gmail.com  http://goo.gl/1wk0
College of Science & Technology, Temple University, Philadelphia PA, USA
International Centre for Theoretical Physics, Trieste. Italy.

This archive was generated by hypermail 2.1.6 : Wed Dec 31 2014 - 23:22:08 CST