[SIESTA-L] Structure does not converge
Dear SIESTA users, I am very new to the SIESTA code, so to start with I am familiarising myself with the protocols involved to run calculations with the code. I am trying to reproduce some numbers from a paper, where a crystal structure (triclinic system) of a Ru-allyl complex was relaxed (Organometallics, 2013, 32, 3784). I have set up a SIESTA input for relaxing this structure, using PBE with the default DZP basis sets for all atoms, along with pseudopotentials from the SIESTA database (I am aware that this might not be the best choice, but to start with it might be a fair option). The calculation seems to fail to converge to a minimum structure and stops after the default 500 opt cycles. This makes me suspicious, since my feeling is that the structure should have converged within 500 cycles. The SCF cycles converge fine. I have attached my input file below, the structure is read in from an external xyz file (attached). I have used as many default settings for the moment as possible. I would appreciate if you could share your experience and advice on this subject. Is such an input in principle a good starting point for the task at hand? What could I improve? Would the use of a Monkhorst-Pack grid be a better choice? There are some other Questions regarding dispersion corrections but first I need to play around with what I understand from the manual. Kind regards, thanks for everybody's help. Really trying to learn how to get reliable results. Tobias SystemName [Ru(AcCN)2(Cp*)(eta3-phenylallyl) unit cell SystemLabel ru_allyl NumberOfAtoms 160 NumberOfSpecies 7 WriteMullikenPop1 %block ChemicalSpeciesLabel 1 44 Ru # Species index, atomic number, species label 2 15 P 3 9 F 4 8 O 5 7 N 6 6 C 7 1 H %endblock ChemicalSpeciesLabel MeshCutoff 400.0 Ry NetCharge 0.00 XC.functional GGA XC.authors PBE SpinPolarized.false. DM.NumberPulay4 DM.MixingWeight 0.3 DM.UseSaveDM .true. NeglNonOverlapInt False AtomicCoordinatesFormat NotScaledCartesianAng AtomicCoordinatesAndAtomicSpecies ru_allylinput.xyz LatticeConstant 1.0 Ang %block LatticeParameters 9.5347 11.7495 15.8787 109.395 96.719 100.965 %endblock LatticeParameters WriteCoorXmol True WriteForces .true. WriteMullikenPop 1 SaveDeltaRho .true. MD.TypeOfRun cg # Type of dynamics: Conjugate gradients MD.NumCGsteps 500 # number of CG steps MD.MaxCGDispl 0.15 Ang MD.MaxForceTol 0.04 eV/Ang MD.UseSaveXV yes MD.VariableCell.true. Dr. Tobias Kraemer Research Associate Institute of Chemical Sciences School of Engineering Physical Sciences Heriot-Watt University Edinburgh EH14 4AS United Kingdom * t.krae...@hw.ac.uk * +44 (0)131 451 3259 - Sunday Times Scottish University of the Year 2011-2013 Top in the UK for student experience Fourth university in the UK and top in Scotland (National Student Survey 2012) We invite research leaders and ambitious early career researchers to join us in leading and driving research in key inter-disciplinary themes. Please see www.hw.ac.uk/researchleaders for further information and how to apply. Heriot-Watt University is a Scottish charity registered under charity number SC000278. ru_allylinput.xyz Description: ru_allylinput.xyz ru_allyl.fdf Description: ru_allyl.fdf
Re: [SIESTA-L] Structure does not converge
Tobias, I'm curious, what do the results of the following commands look like? 1) fgrep constrained 1.out 2) fgrep outcell: Cell vector modules (Ang) 1.out *where 1.out is the name of the output log file You don't have to print them and send them our way, as they'll be very large, but what can you say about them in general? Are the forces fluctuating about a small amount that isn't quite as small as the threshold you've chosen? Are the lattice parameters very different from what you started with? Is everything blowing up? In general, this is what I do. Instead of doing one relaxation with a very large number of CG steps, perform several in steps, making sure the UseSaveData flag is set to true so that you don't have to come up with a new density matrix each time and can use the previous CG trajectory information as well: 1) Do one with MD.NumCGSteps ~ 10 and MD.MaxCGDispl ~ 0.15 Ang, to perform a very coarse relaxation (it is unlikely you will meet the threshold, and that's okay. You're not done). Call the output file 1.out or something like that. 2) Do one with MD.NumCGSteps ~ 20 and MD.MaxCGDispl ~ 0.05 Ang, to perform a coarse relaxation. Call the output file 2.out 3) Do one with MD.NumCGSteps ~ 20 and MD.MaxCGDispl ~ 0.01 Ang, to perform a finer relaxation. Call the output file 3.out ... as many successive steps as you need with as much fineness as desired or until you meet your force convergence criteria, at which point you will have a relaxed structure. If 'everything is blowing up' then you may have to take another look at your simulation parameters. Sometimes choosing the wrong Pulay mixing parameters will prevent your simulation from falling into a convergence basin, so another look at that may be helpful. Transition metals sometimes call for using small (~0.01 or smaller) DM.MixingWeight with a kick every several steps (see the manual). Warm regards, -- *Abraham Hmiel* Katherine Belz Groves Fellow in Nanoscience Xue Group, SUNY College of Nanoscale Science and Engineering http://abehmiel.net/about On Fri, Jan 24, 2014 at 1:09 PM, Kraemer, Tobias t.krae...@hw.ac.uk wrote: Dear SIESTA users, I am very new to the SIESTA code, so to start with I am familiarising myself with the protocols involved to run calculations with the code. I am trying to reproduce some numbers from a paper, where a crystal structure (triclinic system) of a Ru-allyl complex was relaxed (Organometallics, 2013, 32, 3784). I have set up a SIESTA input for relaxing this structure, using PBE with the default DZP basis sets for all atoms, along with pseudopotentials from the SIESTA database (I am aware that this might not be the best choice, but to start with it might be a fair option). The calculation seems to fail to converge to a minimum structure and stops after the default 500 opt cycles. This makes me suspicious, since my feeling is that the structure should have converged within 500 cycles. The SCF cycles converge fine. I have attached my input file below, the structure is read in from an external xyz file (attached). I have used as many default settings for the moment as possible. I would appreciate if you could share your experience and advice on this subject. Is such an input in principle a good starting point for the task at hand? What could I improve? Would the use of a Monkhorst-Pack grid be a better choice? There are some other Questions regarding dispersion corrections but first I need to play around with what I understand from the manual. Kind regards, thanks for everybody’s help. Really trying to learn how to get reliable results. Tobias SystemName [Ru(AcCN)2(Cp*)(eta3-phenylallyl) unit cell SystemLabel ru_allyl NumberOfAtoms 160 NumberOfSpecies 7 WriteMullikenPop1 %block ChemicalSpeciesLabel 1 44 Ru # Species index, atomic number, species label 2 15 P 3 9 F 4 8 O 5 7 N 6 6 C 7 1 H %endblock ChemicalSpeciesLabel MeshCutoff 400.0 Ry NetCharge 0.00 XC.functional GGA XC.authors PBE SpinPolarized.false. DM.NumberPulay4 DM.MixingWeight 0.3 DM.UseSaveDM .true. NeglNonOverlapInt False AtomicCoordinatesFormat NotScaledCartesianAng AtomicCoordinatesAndAtomicSpecies ru_allylinput.xyz LatticeConstant 1.0 Ang %block LatticeParameters 9.5347 11.7495 15.8787 109.395 96.719 100.965 %endblock LatticeParameters WriteCoorXmol True WriteForces .true. WriteMullikenPop 1 SaveDeltaRho .true. MD.TypeOfRun cg # Type of dynamics: Conjugate gradients MD.NumCGsteps 500 # number of CG steps MD.MaxCGDispl 0.15 Ang MD.MaxForceTol 0.04 eV/Ang MD.UseSaveXV yes MD.VariableCell.true. Dr. Tobias Kraemer Research Associate Institute of
Re: [SIESTA-L] atoms
Suman,
I think the answer you're looking for is 42. Only 42-atom simulations are
considered by peer review ;)
But joking aside, you should look to previously-published ab-initio results
in your system of interest or something close to it. Even if such results
exist, you should perform methodological checks on system size, i.e.
size-convergence simulations with respect to: the number of layers in a
surface slab, or supercell dimensions if you're simulating a bulk compound
(perhaps with a defect or other non-homogeneity), different sizes of
nanowire, etc. Some of the most important current problems in nanoscience
involve quantum size-effects, and you don't really know how many atoms a
good simulation requires until you try it- the physics of a particular
solid-state system can vary dramatically based on small size increments!
So, to answer the question you posed, no real number of atoms guideline
exists. Use what is feasible and timely for your computational setup while
giving trustworthy, meaningful, repeatable results! By the way, the T in
siesta stands for Thousands (but this only really makes sense with the
order-N functional). The number of atoms in a SIESTA simulation could be
anywhere in the range {10^0, 10^4}.
--
*Abraham Hmiel*
Katherine Belz Groves Fellow in Nanoscience
Xue Group, SUNY College of Nanoscale Science and Engineering
http://abehmiel.net/about
On Fri, Jan 24, 2014 at 8:21 AM, Suman Chowdhury sumanchowdhur...@gmail.com
wrote:
Dear all,
what is the typical number of atoms one should use for publishable results
using SIESTA?
Re: [SIESTA-L] atoms
thanks a lot
On Sat, Jan 25, 2014 at 2:50 AM, Abraham Hmiel abehm...@gmail.com wrote:
Suman,
I think the answer you're looking for is 42. Only 42-atom simulations are
considered by peer review ;)
But joking aside, you should look to previously-published ab-initio
results in your system of interest or something close to it. Even if such
results exist, you should perform methodological checks on system size,
i.e. size-convergence simulations with respect to: the number of layers in
a surface slab, or supercell dimensions if you're simulating a bulk
compound (perhaps with a defect or other non-homogeneity), different sizes
of nanowire, etc. Some of the most important current problems in
nanoscience involve quantum size-effects, and you don't really know how
many atoms a good simulation requires until you try it- the physics of a
particular solid-state system can vary dramatically based on small size
increments!
So, to answer the question you posed, no real number of atoms guideline
exists. Use what is feasible and timely for your computational setup while
giving trustworthy, meaningful, repeatable results! By the way, the T in
siesta stands for Thousands (but this only really makes sense with the
order-N functional). The number of atoms in a SIESTA simulation could be
anywhere in the range {10^0, 10^4}.
--
*Abraham Hmiel*
Katherine Belz Groves Fellow in Nanoscience
Xue Group, SUNY College of Nanoscale Science and Engineering
http://abehmiel.net/about
On Fri, Jan 24, 2014 at 8:21 AM, Suman Chowdhury
sumanchowdhur...@gmail.com wrote:
Dear all,
what is the typical number of atoms one should use for publishable
results using SIESTA?
