[QE-users] Phonon: negative frequencies for adsorbate molecule using "nat-todo" option

Tamas Karpati tkarpati at gmail.com
Thu Nov 19 18:44:07 CET 2020


Dear Omer,

Sorry, I'm not familiar with the EMT calculator, please consult the
ASE or the EMT forum.
Just a hint: the Python World is undergoing greater earthquakes, watch out for
Python, ASE etc... versions (you need to match all).

As for your theoretical problem to solve, I believe that you have a reaction

   SO2  +  slab  -->   TS   -->   slab~SO2

In such a reaction you need the partition ftn. (ie. phonon or
frequency calc.) for
four species: SO2, slab, TS and slab~SO2 (the latter being the adsorbed form).
Why would you calculate such a weird thing?
I mean it is easy to do with QE/PH.x but for what aim?
   t

On Wed, Nov 18, 2020 at 4:51 PM Omer Mutasim <omermutasim at ymail.com> wrote:
>
> Dear Dr. Tamas
>
> I tried to use ASE base on your earlier suggestion. I want to calculate the vibration frequencies for slab with SO2 adsorbate , i want to vibrate only the adsorbate, while fixing the surface using finite difference method. I used the following script, similar to the examples but it doesn't work, i got the following error :
>
>
>   File "D:\Anaconda\lib\site-packages\ase\formula.py", line 402, in parse2
>     raise ValueError
>
> ValueError
>
> Please tell me what is wrong with this script shown below. or please can you write the correct script for this calculation
>
>
> import ase.io
> from ase import Atoms
> from ase.io import read, write
> from ase.calculators.emt import EMT
> from ase.optimize import BFGS
> from ase.vibrations import Vibrations
> slab = read("CuO_slab.cif")
> slab = Atoms('slab', calculator=EMT())
> BFGS(slab).run(fmax=0.01)
> vib = Vibrations(slab)
> vib.run()
> vib.summary()
>
>
> what is the right procedure for reading cif file of the optimized structure
>
>
> On Thursday, November 5, 2020, 03:57:18 PM GMT+4, Tamas Karpati <tkarpati at gmail.com> wrote:
>
>
> Dear Omer,
>
> Yes, i meant SO2 gas phase sim. This is an alternative to using the physisorbed
> slab+SO2 complex as "reactant", R. Question of methodology and the nature of
> materials. I cannot recall whether S+2O were together (as SO2) or decomposed in
> your starting structure but in the second case you may use this approach.
> Also, kinetic parameters are probably easier to derive for the gas
> phase reactant (?)
> If your SO2 is in a single piece on the slab in R, then you can drop this idea.
>
> With XCrySDen, to visualize vibrations, you need to use dynmat.x and specify
> filxsf, then open *.xsf with XCrySDen. Select to visualize forces and use the
> animation controls to choose which normal mode is to be shown.
> Expect arrows, not animated vibrations. Such arrows indicate the atomic
> replacements to do in order to get 1 negative freq. better: move them
> where the arrows
> point to, or the exact opposite direction. One by one you can find a
> real minimum str.
> Note that unless you do an all-atom phonon, you can never trust those arrows
> (neither direction nor frequency, even sign).
>
> As for H+O (if distant enough, yet in the same simbox), you formally expect
> the same as for the HO radical -a single bond with a single freq but since
> this is not a real bond (strength is negligible), freq will be very
> small (and yes, S=0).
>   t
>
> On Thu, Nov 5, 2020 at 5:15 AM Omer Mutasim <omermutasim at ymail.com> wrote:
> >
> > Dear Dr. Tamas
> >
> > Sorry , what do you mean by “ you most probably need an SO2 simulation (optimization+phonons)
> > rather than the same for a surface attached SO2 or SO+O. Big difference! “  do you mean i should do phonon for SO2 in gas phase  ?
> >
> > I do agree with you , 3-atoms phonon is non-physical, i will include the top layer also.
> >
> > Regarding Xcrysden , I don’t find any axsf file in my output files , how do you visualise it ?
> > Another question:
> > For the reaction : H+O = OH , if i did phonon for the initial state( reactant ) , I should expect to get no vibrational modes at all , right ? I will get only 6 translational modes . So the vibrational entropy will be zero,  Please correct me if i am wrong.
> >
> > By the way , i tried ASE for rate constants, as you recommended, it is really helpful.
> >
> > Thanks a lot for you unwavering help.
> >
> > Sent from Yahoo Mail for iPhone
> >
> > On Wednesday, November 4, 2020, 10:46 PM, Tamas Karpati <tkarpati at gmail.com> wrote:
> >
> > In addition to my earlier comments, i'd like to mention that for kinetics
> > you most probably need an SO2 simulation (optimization+phonons)
> > rather than the same for a surface attached SO2 or SO+O. Big difference!
> >
> > Back to the earliers, 3-atoms phonon is so unphysical that it is
> > recommended to do an all-atom one, or add at least the directly
> > bonded surface atoms (and extend towards all-atoms if you can).
> >
> > In addition, only all-atom phonon will show you whether your big negative
> > freqs. indicate a non-minimum structure (ie. "freqs" ~ 2nd derivatives
> > of the PEHS).
> > With somewhat less atoms you can be lucky, and by visualizing normal
> > modes (phonons) by eg. XCrysDen will show you where to move atoms
> > to get into the local minimum. Kinetics theory builds on minima (and TS-es).
> > Anyways, 3 atoms are too few (also see first section above).
> >
> > On Mon, Nov 2, 2020 at 11:29 AM Tamas Karpati <tkarpati at gmail.com> wrote:
> > >
> > > Dear Omer,
> > >
> > > I guess that your input is fine, your structure is not.
> > > (By the way, tr2_ph could be lower.)
> > >
> > > You woud expect 6 near zero and 3N-6 positive freqs. for a completely
> > > relaxed minimum structure (and 5 + 3N-5 for a TS).
> > > This is ruined if you run PH.x on a different potential energy
> > > hypersurface, PEHS.
> > >
> > > It is really very easy to spoil: use a different no. of k point or functional
> > > for PW/vc-relax and PH and you're there. Another temptation is to
> > > use experimental crystal structure and fix some/most atoms as such.
> > > These all mean different PEHS'.
> > > In addition, unconverged relaxation (or too loose convergence),
> > > while moves on the same PEHS, provides you with inappropriate
> > > freqs, as it does not bring your structure close enough to the local minimum.
> > >
> > > I would recommend to reconsider the "life" of your structure
> > > (origin, optimization method, other parameters) and adjust if necessary.
> > >    t
> > >
> > > On Mon, Nov 2, 2020 at 7:43 AM Omer Mutasim <omermutasim at ymail.com> wrote:
> > > >
> > > >  Dear all
> > > >
> > > >  I'm doing phonon calculation at Gamma point (q)  in order to estimate the reaction rate constants for a micro-kinetic model.  I have perturbed only the adsorbate molecule with the 3 surface atoms, connected to adsorbate, using "nat-todo" option. However, i got 15 negative frequencies (should be 6 as i know ) ,with high absolute value.
> > > >
> > > > Can you please help me to know what is wrong with my input files ?
> > > >
> > > > Below are the output & input files:
> > > >
> > > >
> > > >
> > > >      Mode symmetry, C_1 (1)    point group:
> > > >
> > > >      freq (  1 -  1) =      -2762.6  [cm-1]  --> A              I+R
> > > >      freq (  2 -  2) =      -2570.3  [cm-1]  --> A              I+R
> > > >      freq (  3 -  3) =      -2460.4  [cm-1]  --> A              I+R
> > > >      freq (  4 -  4) =      -2423.6  [cm-1]  --> A              I+R
> > > >      freq (  5 -  5) =      -2356.3  [cm-1]  --> A              I+R
> > > >      freq (  6 -  6) =      -2158.0  [cm-1]  --> A              I+R
> > > >      freq (  7 -  7) =      -2151.1  [cm-1]  --> A              I+R
> > > >      freq (  8 -  8) =      -2067.5  [cm-1]  --> A              I+R
> > > >      freq (  9 -  9) =      -2034.8  [cm-1]  --> A              I+R
> > > >      freq ( 10 - 10) =      -2025.2  [cm-1]  --> A              I+R
> > > >      freq ( 11 - 11) =      -1864.3  [cm-1]  --> A              I+R
> > > >      freq ( 12 - 12) =      -1804.5  [cm-1]  --> A              I+R
> > > >      freq ( 13 - 13) =      -1099.4  [cm-1]  --> A              I+R
> > > >      freq ( 14 - 14) =      -947.6  [cm-1]  --> A              I+R
> > > >      freq ( 15 - 15) =      -912.5  [cm-1]  --> A              I+R
> > > >      freq (316 -316) =        179.3  [cm-1]  --> A              I+R
> > > >      freq (317 -317) =        193.0  [cm-1]  --> A              I+R
> > > >      freq (318 -318) =        215.8  [cm-1]  --> A              I+R
> > > >      freq (319 -319) =        240.2  [cm-1]  --> A              I+R
> > > >      freq (320 -320) =        270.4  [cm-1]  --> A              I+R
> > > >      freq (321 -321) =        317.0  [cm-1]  --> A              I+R
> > > >      freq (322 -322) =        370.8  [cm-1]  --> A              I+R
> > > >      freq (323 -323) =        377.3  [cm-1]  --> A              I+R
> > > >      freq (324 -324) =        398.3  [cm-1]  --> A              I+R
> > > >      freq (325 -325) =        417.8  [cm-1]  --> A              I+R
> > > >      freq (326 -326) =        468.2  [cm-1]  --> A              I+R
> > > >      freq (327 -327) =        659.0  [cm-1]  --> A              I+R
> > > >      freq (328 -328) =      1096.6  [cm-1]  --> A              I+R
> > > >      freq (329 -329) =      1795.5  [cm-1]  --> A              I+R
> > > >      freq (330 -330) =      2199.3  [cm-1]  --> A              I+R
> > > >
> > > >
> > > > ph.x input file:
> > > >
> > > > phonon calculation at Gamma point.
> > > > &inputph
> > > >  outdir = './outdir'
> > > >  prefix = 'HS'
> > > >  tr2_ph = 1.0d-09
> > > >  epsil = .false.
> > > >  amass(1) = 58.69340
> > > >  amass(2) = 30.97376
> > > >  amass(3) = 1.00784
> > > >  amass(4) = 32.065
> > > >  fildyn = 'HS.dyn'
> > > >
> > > > alpha_mix(1)=0.1
> > > >
> > > >  recover=.true
> > > >  nogg = .true
> > > >  nat_todo = 5
> > > >
> > > > /
> > > > 0.0 0.0 0.0
> > > >
> > > > 1 2 37 46 54
> > > >
> > > >
> > > >
> > > > scf input file:
> > > >
> > > > &CONTROL
> > > >    calculation  = "scf"
> > > > prefix = 'HS'
> > > >    outdir = './outdir'
> > > >    pseudo_dir = '/home/'
> > > > restart_mode = 'from_scratch'
> > > >    forc_conv_thr =  1.0e-03
> > > > etot_conv_thr = 1e-04
> > > >    nstep        = 999
> > > > /
> > > > &SYSTEM
> > > > ibrav  =  0
> > > >    ecutrho                  =  200
> > > >    ecutwfc                  =  25
> > > >    nat                      = 110
> > > >    ntyp                      = 4
> > > > occupations='smearing',smearing='gaussian',degauss=0.005
> > > > vdw_corr = 'DFT-D2'
> > > >      nspin = 2
> > > >  starting_magnetization(1)=  0.01
> > > > /
> > > > &ELECTRONS
> > > >    conv_thr        = 1e-8
> > > >    electron_maxstep = 200
> > > > mixing_mode ='local-TF'
> > > >    mixing_beta      =  0.3
> > > > /
> > > > &IONS
> > > > /
> > > > K_POINTS {automatic}
> > > > 3 3 1    0 0 0
> > > > ATOMIC_SPECIES
> > > > Ni 58.69340 Ni.pbe-n-rrkjus_psl.0.1.UPF
> > > > P 30.97376 P.pbe-n-rrkjus_psl.1.0.0.UPF
> > > > H 1.00784 H.pbe-rrkjus_psl.0.1.UPF
> > > > S  32.065      S.pbe-n-rrkjus_psl.1.0.0.UPF
> > > > CELL_PARAMETERS {angstrom}
> > > >        11.76538354          0.0000000000        0.0000000000
> > > >        -5.8826917709        10.189121032        0.0000000000
> > > >        0.0000000000        0.0000000000        30.993869056
> > > > ATOMIC_POSITIONS (angstrom)
> > > > H        0.879694621  3.392266427  10.708999692
> > > > S        2.266698845  3.396363162  10.560733430
> > > > Ni      -2.744571590  4.755054131  0.244939179
> > > > Ni      3.134031329  1.363792691  0.248008546
> > > > .
> > > > .
> > > > .
> > > > P      -1.060403962  1.841094610  1.604930623
> > > > P      -3.921453199  6.792156181  0.000000000    0  0  0
> > > > P        1.960697149  3.396027080  0.000000000    0  0  0
> > > > P        7.842906399  0.000000000  0.000000000    0  0  0
> > > >
> > > >
> > > > regards
> > > >
> > > >
> > > > _______________________________________________
> > > > Quantum ESPRESSO is supported by MaX (www.max-centre.eu)
> > > > users mailing list users at lists.quantum-espresso.org
> > > > https://lists.quantum-espresso.org/mailman/listinfo/users
> > _______________________________________________
> > Quantum ESPRESSO is supported by MaX (www.max-centre.eu)
> > users mailing list users at lists.quantum-espresso.org
> > https://lists.quantum-espresso.org/mailman/listinfo/users


More information about the users mailing list