[gmx-users] Re: angle constrain, constrained PF6 anion
Vitaly Chaban
vvchaban at gmail.com
Thu Feb 17 19:12:59 CET 2011
Hey, Gyorgy -
The problem is definitely with the incorrectly defined constraints as
I mentioned yesterday. You may find Chapter 5 of the gromacs manual
useful to repair this.
I hope to be able to look through your attached system on the weekend,
if you don't locate the problem before yourself...
Regards and no runtime errors in your simulations,
Vitaly
On Thu, Feb 17, 2011 at 11:48 AM, <gyorgy.hantal at fc.up.pt> wrote:
>
> Hi Vitaly,
>
> Thanks for your answer.
> Actually, we tried in a lot of ways to run the simulations. This was also
> one way, and in fact, if you use a much smaller time step, the simulation
> doesn't crash (I'm not saying it gives a physical result).
>
> If I remove the [bonds] part from the topology, the simulation runs and the
> geometry is conserved, but produces positive total energy. Maybe this is
> related to the unphysically high force constants of the angle bendings, but
> I think, this doesn't matter if the angles are contrained.
> (Maybe the problem is somewhere here, but I don't know how to constrain
> properly the angles if not like this (and nobody can tell me in my
> proximity))
>
> If I remove the [angles] section, de system's total energy is negative but
> the geometry of the anion is completely wrong.
>
> I also have problems when running on more than 1 processor, even if I use
> -pd for mdrun (to switch to particle decomposition instead of domain
> decomposition), the simulation crashes immediately:
>
> Range checking error:
> Variable n has value 7. It should have been within [ 0 .. 7 ]
>
> If you want to test the topology on one single ion pair, you can find a
> configuration in the attachement.
>
> Thanks in advance,
> Best,
>
> Gyorgy
>
>
> Quoting Vitaly Chaban <vvchaban at gmail.com>:
>
>> Hello,
>>
>> I got your files. However, the error is completely another than you
>> mention and it occurs directly at the first time-step of the
>> simulation.
>>
>> ============================================
>> Inner product between old and new vector <= 0.0!
>> constraint #1 atoms 1377 and 1381
>> Wrote pdb files with previous and current coordinates
>> step 0Inner product between old and new vector <= 0.0!
>> constraint #1 atoms 2177 and 2181
>> Wrote pdb files with previous and current coordinates
>> Inner product between old and new vector <= 0.0!
>> constraint #1 atoms 2 and 6
>> Wrote pdb files with previous and current coordinates
>> Segmentation fault
>> ============================================
>>
>> There are a number of strange things in your configuration with the
>> unhealthy defined constraints in topology as probably the most severe
>> one. If you have constraints, why do you define bonds force constants
>> at the same time?
>>
>> Good luck with your study,
>> Vitaly
>>
>>
>>
>>
>> On Wed, Feb 16, 2011 at 1:58 PM, <gyorgy.hantal at fc.up.pt> wrote:
>>>
>>> Hi Vitaly,
>>>
>>> Thanks for your answer. I was busy with other stuff these last two weeks,
>>> that's why I have time to reply only now.
>>>
>>> So the system is the BMIM PF6 ionic liquid.
>>> We chose two potentials to describe the system, and both treat the anion
>>> as
>>> rigid. We cannot use LINCS because it shouldn't be used if one constrains
>>> not only bonds but angles too, so we have to use SHAKE.
>>>
>>> But with SHAKE I encountered some strange things as I wrote in my first
>>> thread. In the attachement I am sending you all the files needed for the
>>> simulation, maybe you'll find the problem.
>>>
>>> I haven"t mentioned but we try to constrain the angles by adding fictious
>>> bonds. Maybe there is something wrong with this...
>>>
>>> Thank you in advance for your help.
>>>
>>> All the best,
>>> Gyorgy
>>>
>>>
>>> Quoting Vitaly Chaban <vvchaban at gmail.com>:
>>>
>>>> Hey Gyorgy,
>>>>
>>>> Your current topology file(s) is(are) also important to analyze the
>>>> situation and fix the problem.
>>>>
>>>> Pass my kind regards to Julia...
>>>>
>>>> ========================
>>>> Dr. Vitaly V. Chaban, Ph.D.
>>>> Department of Chemistry
>>>> University of Rochester
>>>> Rochester, New York 14627-0216
>>>> The United States of America
>>>> =========================
>>>>
>>>>
>>>>> The mdp file is attached.
>>>>> Best,
>>>>> Gyorgy
>>>>>
>>>>> Quoting "Justin A. Lemkul" <jalemkul at vt.edu>:
>>>>>
>>>>>>
>>>>>>
>>>>>> gyorgy.hantal at fc.up.pt wrote:
>>>>>>>
>>>>>>> Dear all,
>>>>>>>
>>>>>>> I am setting up a simulation of ionic liquids with the PF6 anion.
>>>>>>> According to the potential, the anion should be kept rigid, wich
>>>>>>> obviously means that bond lengths and angles have to be
>>>>>>> constrained. LINCS doesn't work with angle constraints (i.e.
>>>>>>> constraing a triangle), so we decided to use SHAKE. However, SHAKE
>>>>>>> seems to work a bit strangely: I know SHAKE mustn't be used with
>>>>>>> domain decomposition, but even if I set the corresponding variable
>>>>>>> to NO in the mdp file, the simulation crashes on 8 procs and gives
>>>>>>> the following error message:
>>>>>>>
>>>>>>> Fatal error:
>>>>>>> 1 particles communicated to PME node 7 are more than a cell length
>>>>>>> out of the domain decomposition cell of their charge group.
>>>>>>>
>>>>>>> If I try to run mdrun with -pd (to 'really' switch off domain
>>>>>>> decomposition), the simulation doesn't chrash but gives nonsense
>>>>>>> (the energy seems to increase constantly).
>>>>>>>
>>>>>>> I am not an expert user so maybe I do something wrong but, anyway,
>>>>>>> does anyone have an idea how to constrain this anion with Gromacs?
>>>>>>> I checked mailing list archive but couldn't find any answer
>>>>>>> corresponding to my question.
>>>>>>>
>>>>>>
>>>>>> Without seeing a complete .mdp file, it's not possible to fully
>>>>>> diagnose this problem. The combination of SHAKE + particle
>>>>>> decomposition should be stable, but there are a whole host of
>>>>>> different
>>>>>> things that can go wrong.
>>>>>>
>>>>>> -Justin
>>>>>>
>>>>>>> Thanks in advance.
>>>>>>>
>>>>>>> Gyorgy
>>>>>>>
>>>>>>>
>>>>>>
>>>>>> --
>>>>>> ========================================
>>>>>>
>>>>>> Justin A. Lemkul
>>>>>> Ph.D. Candidate
>>>>>> ICTAS Doctoral Scholar
>>>>>> MILES-IGERT Trainee
>>>>>> Department of Biochemistry
>>>>>> Virginia Tech
>>>>>> Blacksburg, VA
>>>>>> jalemkul[at]vt.edu | (540) 231-9080
>>>>>> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
>>>>>>
>>>>>> ========================================
>>>>>> --
>>>>>> gmx-users mailing list gmx-users at gromacs.org
>>>>>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>>>>>> Please search the archive at
>>>>>> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>>>>>> Please don't post (un)subscribe requests to the list. Use the www
>>>>>> interface or send it to gmx-users-request at gromacs.org.
>>>>>> Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>>>>
>>>>>
>>>>> -------------- next part --------------
>>>>> title = BMIM PF6 bulk simulation
>>>>> cpp = /usr/bin/cpp
>>>>>
>>>>> ; RUN CONTROL PARAMETERS ;l_bfgs
>>>>> integrator = md;steep
>>>>> ; Start time and timestep in ps
>>>>> tinit = 0
>>>>> dt = 0.00000001 ! just to see if it starts
>>>>> nsteps = 25000
>>>>> ; mode for center of mass motion removal
>>>>> comm-mode = linear
>>>>> ; number of steps for center of mass motion removal
>>>>> nstcomm = 10
>>>>> ; group(s) for center of mass motion removal
>>>>>
>>>>> ; ENERGY MINIMIZATION OPTIONS
>>>>> ; Force tolerance and initial step-size
>>>>> ;emtol = 50
>>>>> ;emstep = 0.5
>>>>> ;lincs_iter = 3
>>>>> ;lincs_warnangle = 50
>>>>> ; Frequency of steepest descents steps when doing CG
>>>>> ;nstcgsteep = 1000
>>>>> ;nbfgscorr = 10
>>>>>
>>>>> ; OUTPUT CONTROL OPTIONS
>>>>> ; Output frequency for coords (x), velocities (v) and forces (f)
>>>>> nstxout = 0
>>>>> nstvout = 0
>>>>> nstfout = 0
>>>>> ; Checkpointing helps you continue after crashes
>>>>> nstcheckpoint = 10000
>>>>> ; Output frequency for energies to log file and energy file
>>>>> nstlog = 10000
>>>>> nstenergy = 10000
>>>>> ; Output frequency and precision for xtc file
>>>>> nstxtcout = 10000
>>>>> xtc-precision = 1000
>>>>> ; This selects the subset of atoms for the xtc file. You can
>>>>> ; select multiple groups. By default all atoms will be written.
>>>>> xtc-grps =
>>>>> ; Selection of energy groups
>>>>> energygrps =
>>>>>
>>>>> ; NEIGHBORSEARCHING PARAMETERS
>>>>> ; nblist update frequency
>>>>> nstlist = 25
>>>>> ; ns algorithm (simple or grid)
>>>>> ns_type = grid
>>>>> ; Periodic boundary conditions: xyz (default), no (vacuum)
>>>>> ; or full (infinite systems only)
>>>>> pbc = xyz
>>>>> ; nblist cut-off
>>>>> rlist = 1.5
>>>>> domain-decomposition = no
>>>>>
>>>>> ; OPTIONS FOR ELECTROSTATICS AND VDW
>>>>> ; Method for doing electrostatics
>>>>> coulombtype = PME ; can RF also be used?
>>>>> rcoulomb-switch = 0
>>>>> rcoulomb = 1.5
>>>>> ; Dielectric constant (DC) for cut-off or DC of reaction field
>>>>> epsilon_rf =
>>>>> ; Method for doing Van der Waals
>>>>> vdw-type = Shift
>>>>> ; cut-off lengths
>>>>> rvdw-switch = 1.1
>>>>> rvdw = 1.2
>>>>> ; Apply long range dispersion corrections for Energy and Pressure
>>>>> DispCorr = EnerPres
>>>>> ; Extension of the potential lookup tables beyond the cut-off
>>>>> table-extension = 1
>>>>> ; Spacing for the PME/PPPM FFT grid
>>>>> fourierspacing = 0.12
>>>>> ; FFT grid size, when a value is 0 fourierspacing will be used
>>>>> fourier_nx = 0
>>>>> fourier_ny = 0
>>>>> fourier_nz = 0
>>>>> ; EWALD/PME/PPPM parameters
>>>>> pme_order = 4
>>>>> ewald_rtol = 1e-05
>>>>> ewald_geometry = 3d
>>>>> epsilon_surface = 0
>>>>> optimize_fft = no
>>>>>
>>>>> ; OPTIONS FOR WEAK COUPLING ALGORITHMS
>>>>> ; Temperature coupling
>>>>> Tcoupl = V-rescale
>>>>> ; Groups to couple separately
>>>>> tc-grps = System
>>>>> ; Time constant (ps) and reference temperature (K)
>>>>> tau_t = 1.0
>>>>> ref_t = 298.15
>>>>> ; Pressure coupling
>>>>> Pcoupl = Parrinello-Rahman
>>>>> Pcoupltype = isotropic
>>>>> ; Time constant (ps), compressibility (1/bar) and reference P (bar)
>>>>> tau_p = 4.0
>>>>> compressibility = 4.5e-5 ; proper value for IL???
>>>>> ref_p = 1.00
>>>>>
>>>>> ; SIMULATED ANNEALING
>>>>> ; Type of annealing for each temperature group (no/single/periodic)
>>>>> ; annealing = no
>>>>>
>>>>> ; GENERATE VELOCITIES FOR STARTUP RUN
>>>>> gen_vel = yes
>>>>> gen_temp = 298.15
>>>>> gen_seed = 1993
>>>>>
>>>>> ; OPTIONS FOR BONDS
>>>>> constraints = hbonds
>>>>> constraint_algorithm = SHAKE
>>>>
More information about the gromacs.org_gmx-users
mailing list