Imaginary frequency

[mwchang]

Dear All

I try to calculate the frequencies of gas molecules.

e.g. CO, H2O, O2, H2, OH

I put a molecule into a cubic box with a 11.034*11.134*11.234 size , which has been broken the symmetry .

The size should enough avoid the interaction between images.





I used the INCAR and POSCAR as optimization input file to calculate the CO molecule .


SYSTEM = CO freq

ISTART = 1 !job : 0-new 1-cont 2-samecut
NWRITE = 3
ICHARG = 1 ! Read_CHGCAR

ENCUT = 500 eV
EDIFF = 1E-5
EDIFFG = -1E-2

NSW = 1
IBRION = 5
POTIM = 0.015
NFREE = 2

ISYM=0 ! 0 for broken Symm;
ISMEAR = 0 ! Gaussian
SIGMA = 0.1 eV ! DOS related value

GGA = 91 !GGA-PW91
VOSKOWN=1
ISPIN = 2 !open spin-polarized

PREC = High


POSCAR:

O C (CO (isolated))
1.000000000000000
11.0340000000000007 0.0000000000000000 0.0000000000000000
0.0000000000000000 11.1340000000000003 0.0000000000000000
0.0000000000000000 0.0000000000000000 11.2340000000000000
1 1
Selective dynamics
Direct
0.9999999897390340 0.0090479590898271 0.0199152517253264 T T T
0.0000000102609660 0.0009520409101709 0.1212583334671496 T T T


This is the calculation result.


1 f = 63.937991 THz 401.734248 2PiTHz 2132.741758 cm-1 264.426372 meV
X Y Z dx dy dz
11.034000 0.100740 0.223728 0.000033 0.013118 -0.163076
0.000000 0.010600 1.362216 0.000003 -0.017381 0.217485

2 f = 1.359656 THz 8.542970 2PiTHz 45.353238 cm-1 5.623087 meV
X Y Z dx dy dz
11.034000 0.100740 0.223728 0.190536 0.064953 0.028435
0.000000 0.010600 1.362216 -0.162346 -0.037917 0.020117

3 f = 1.180866 THz 7.419598 2PiTHz 39.389439 cm-1 4.883670 meV
X Y Z dx dy dz
11.034000 0.100740 0.223728 0.067358 -0.162797 -0.105973
0.000000 0.010600 1.362216 -0.047637 0.133948 -0.082079

4 f/i= 0.002957 THz 0.018580 2PiTHz 0.098639 cm-1 0.012230 meV
X Y Z dx dy dz
11.034000 0.100740 0.223728 -0.044688 0.169183 -0.071248
0.000000 0.010600 1.362216 -0.044671 0.169225 -0.071227

5 f/i= 0.321428 THz 2.019593 2PiTHz 10.721689 cm-1 1.329321 meV
X Y Z dx dy dz
11.034000 0.100740 0.223728 0.140220 0.043791 -0.010122
0.000000 0.010600 1.362216 0.229230 0.041485 -0.010345

6 f/i= 0.820376 THz 5.154576 2PiTHz 27.364805 cm-1 3.392805 meV
X Y Z dx dy dz
11.034000 0.100740 0.223728 0.000286 0.032644 -0.136706
0.000000 0.010600 1.362216 0.009465 -0.182256 -0.153738



As you can see, I got a few imaginary frequency modes.

These f/i modes are very small and they look like translations.

I try to remove them by adding dx ,dy ,dz to the original coor. but it still didn't work.

The f/i modes still exist and sometime become larger.


Is there any effect method to remove the imaginary frequency modes?

or I can ignore them because these modes are very small.

[alex]

Hi,

to answer your last question first: yes.
Check the book Wilson, Decius, Cross: 'Molecular Vibrations' for details.

Imaginary modes in the region a little smaller than 0 cm-1 often appear, because you are doing numerical 2nd derivatives. You might get rid of them doing the Hessian analytically.

Cheers,

alex

[mwchang]

Thanks for your prompt reply.

But I still have some questions.

How small imaginary modes can be ignored? If I use the EDIFF=1E-5 and EDIFFG=-0.01 as the criterion value.

In my opinion, the modes should be smaller than -80.65cm^-1.

Because 1eV equals 8065.5cm^-1 and the EDIFFG=-0.01eV/A,

8065.5*(-0.01)=-80.655cm^-1

Another question is what is the meaning of the Hessian analytically and how to do it?

I know the Hessian matrix is the second derivative of potential surface, but I don?t kwon how to do the Hessian analytically and use the method to remove the f/i mode.

Could you give me some suggestions?

[alex]

Hi mwchang,

surely your suggestion (v_i <= 80cm-1) is all right, but ideally it should became smaller, if a) the structure is well converged and b) the numerical 2nd derivative is good enough, e.g. small step width (like in your case) or well converged gradients (you might improve your criteria here).

If I need very accurate frequencies I'll converge with EDIFFG = -1.e-6, just to give you an idea. I'll use that also for the (numerical) frequency calculation.

Analytical 2nd derivate means you are solving the analytical form of the 2nd derivative instead of approaching it numerically by calculating gradients of displacements. You can access this type in VASP5 by using IBRION 7 or 8. Check the manual, please.

Cheers,

alex

[mwchang]

Thanks for your reply again, alex.

I do frequency calculations is in order to check optimized structures, which are a local minimum. So, I just care whether f/i modes are at present. I will increase the convergence criterion when calculate frequencies .

Actually, I didn't use a higher criterion for electronic and force because I found current criterion is enough and has a similar energy comparing with more highers(ΔE<1meV).

Take into consideration of future calculations , I will use a 3by3 with 9-atomic layers as models, so I choose an more appropriate criterion.



Unfortunately, we only have the licence of vasp4.6, but it don't support the Hessian calculation.

Maybe, I should suggest my boss to upgrade the vasp version. v.v


Best Regards.

[abmo]

Hi,

I am using VASP to find the TST for a number of reactions. I faced a problem when I performed frequency calculations for initial and final states. For some structure I got negative frequencies. Most of the time I re-relax those structures using more accurate criteria and it works then. But I could not solve the problem for some of my structures such as adsorbed CO2 and H over nickel. I was wondering if any body has a solution.
Here is my INCAR file for relaxation:

ENCUT = 600.000000
POTIM = 0.1
EDIFF = 1.00e-06
EDIFFG = -1.00e-03
ALGO = Fast
ISMEAR = 1
ISPIN = 2
IBRION = 2
NELM=100
NSW = 1000
LREAL = Auto
LWAVE = .FALSE.
LCHARG = .FALSE.
PREC = Accurate
ADDGRID = .TRUE.


Here is my INCAR file for frequency calcluations:

POTIM = 0.01
EDIFF = 1.00e-05
EDIFFG = -0.1E-03
ALGO = Fast
ISMEAR = 0
SIGMA = 0.01
ISPIN = 2
IBRION = 5
NELM = 200
NSW = 1
LREAL= Auto
NFREE = 2
LWAVE = FALSE
PREC = High
LWAVE = .FALSE.
LCHARG = .FALSE.


Best Regards
Abas

[alex]

Hi Abas, could you post the frequencies of the mentioned structure, please?

Cheers,

alex

[abmo]

Dear Alex,

Thanks for replying. I do not know how to upload files here. Would you please send an e-mail to : abas.mohsenzadeh@hb.se then I will forward the whole input and output files. If this is inconvenient to you I can just copy the frequency numbers here.

Regards
Abas

[alex]

just give it a

grep THz OUTCAR

...

Would be enough for starters.

Cheers,

alex

[abmo]

Dear Alex,

For the following INCAR of relaxation:

ENCUT = 600.000000
POTIM = 0.5
EDIFF = 1.00e-05
EDIFFG = -1.00e-02
ALGO = VeryFast
ISMEAR = 1
ISPIN = 2
IBRION = 2
NELM=100
NSW = 1000
LREAL = Auto
LWAVE = .FALSE.
LCHARG = .FALSE.


I got the following frequencies:


1 f = 52.066804 THz 327.145377 2PiTHz 1736.761579 cm-1 215.331070 meV
2 f = 34.351786 THz 215.838638 2PiTHz 1145.852212 cm-1 142.067619 meV
3 f = 33.072701 THz 207.801910 2PiTHz 1103.186533 cm-1 136.777747 meV
4 f = 27.602455 THz 173.431337 2PiTHz 920.718753 cm-1 114.154617 meV
5 f = 25.680770 THz 161.357038 2PiTHz 856.618261 cm-1 106.207166 meV
6 f = 19.108837 THz 120.064367 2PiTHz 637.402189 cm-1 79.027828 meV
7 f = 14.613836 THz 91.821440 2PiTHz 487.465084 cm-1 60.437989 meV
8 f = 8.583086 THz 53.929122 2PiTHz 286.300934 cm-1 35.496805 meV
9 f = 5.500522 THz 34.560797 2PiTHz 183.477648 cm-1 22.748337 meV
10 f = 4.271285 THz 26.837275 2PiTHz 142.474729 cm-1 17.664621 meV
11 f/i= 1.715960 THz 10.781693 2PiTHz 57.238254 cm-1 7.096642 meV
12 f/i= 2.559729 THz 16.083249 2PiTHz 85.383350 cm-1 10.586190 meV

Then I increased the accuracy to:

ENCUT = 600.000000
POTIM = 0.1
EDIFF = 1.00e-06
EDIFFG = -1.00e-03
ALGO = Fast
ISMEAR = 1
ISPIN = 2
IBRION = 2
NELM=100
NSW = 1000
LREAL = Auto
LWAVE = .FALSE.
LCHARG = .FALSE.
PREC = Accurate
ADDGRID = .TRUE.

And I got the frequencies as:

1 f = 52.068336 THz 327.155003 2PiTHz 1736.812683 cm-1 215.337406 meV
2 f = 34.529316 THz 216.954088 2PiTHz 1151.773954 cm-1 142.801822 meV
3 f = 33.094628 THz 207.939681 2PiTHz 1103.917934 cm-1 136.868429 meV
4 f = 27.543985 THz 173.063961 2PiTHz 918.768412 cm-1 113.912806 meV
5 f = 25.404552 THz 159.621511 2PiTHz 847.404628 cm-1 105.064821 meV
6 f = 19.235130 THz 120.857884 2PiTHz 641.614842 cm-1 79.550130 meV
7 f = 14.868838 THz 93.423665 2PiTHz 495.971036 cm-1 61.492593 meV
8 f = 8.529719 THz 53.593806 2PiTHz 284.520797 cm-1 35.276096 meV
9 f = 5.052488 THz 31.745721 2PiTHz 168.532868 cm-1 20.895420 meV
10 f = 3.916881 THz 24.610489 2PiTHz 130.653081 cm-1 16.198923 meV
11 f/i= 1.357475 THz 8.529265 2PiTHz 45.280482 cm-1 5.614066 meV
12 f/i= 2.635993 THz 16.562434 2PiTHz 87.927265 cm-1 10.901595 meV


I also tried other settings for frequency calculations. I set Ediff =10^ -8 and POTIM= 0.005 and I got following results respectively:

1 f = 52.093810 THz 327.315064 2PiTHz 1737.662422 cm-1 215.442760 meV
2 f = 34.345045 THz 215.796282 2PiTHz 1145.627352 cm-1 142.039740 meV
3 f = 33.124068 THz 208.124658 2PiTHz 1104.899948 cm-1 136.990184 meV
4 f = 27.279611 THz 171.402851 2PiTHz 909.949849 cm-1 112.819443 meV
5 f = 25.217507 THz 158.446270 2PiTHz 841.165468 cm-1 104.291264 meV
6 f = 19.348855 THz 121.572441 2PiTHz 645.408310 cm-1 80.020461 meV
7 f = 14.911182 THz 93.689717 2PiTHz 497.383466 cm-1 61.667712 meV
8 f = 8.681225 THz 54.545748 2PiTHz 289.574502 cm-1 35.902675 meV
9 f = 5.195181 THz 32.642283 2PiTHz 173.292568 cm-1 21.485548 meV
10 f = 4.376014 THz 27.495304 2PiTHz 145.968097 cm-1 18.097744 meV
11 f/i= 1.135514 THz 7.134644 2PiTHz 37.876663 cm-1 4.696109 meV
12 f/i= 2.223951 THz 13.973497 2PiTHz 74.183024 cm-1 9.197526 meV

and

1 f = 51.917831 THz 326.209352 2PiTHz 1731.792373 cm-1 214.714966 meV
2 f = 34.080322 THz 214.132978 2PiTHz 1136.797139 cm-1 140.944933 meV
3 f = 32.646803 THz 205.125915 2PiTHz 1088.980106 cm-1 135.016374 meV
4 f = 26.325181 THz 165.405993 2PiTHz 878.113503 cm-1 108.872238 meV
5 f = 24.957581 THz 156.813109 2PiTHz 832.495282 cm-1 103.216297 meV
6 f = 18.995056 THz 119.349455 2PiTHz 633.606837 cm-1 78.557264 meV
7 f = 14.820772 THz 93.121659 2PiTHz 494.367736 cm-1 61.293809 meV
8 f = 8.387895 THz 52.702701 2PiTHz 279.790064 cm-1 34.689559 meV
9 f = 5.548401 THz 34.861632 2PiTHz 185.074733 cm-1 22.946351 meV
10 f/i= 0.757611 THz 4.760208 2PiTHz 25.271169 cm-1 3.133227 meV
11 f/i= 1.960391 THz 12.317497 2PiTHz 65.391589 cm-1 8.107527 meV
12 f/i= 3.724753 THz 23.403311 2PiTHz 124.244368 cm-1 15.404344 meV


I hope this is what you wanted.

Regards
Abas

[alex]

Hi Abas,

yes, it is.

I do not understand the very last set. Is it a different system?

The set(s) for the other structure looks like you are trapped in some local minimum. Try to figure out which atom(s) are involved in vib 11 and 12 (by looking at the modes printed above or below the frequencies 11 and 12). Then move the atom a little and reoptimise again. Freq calc. afterwards.

Hth,

alex

[abmo]

Dear Alex,

Thanks for the tips. I will try. The last numbers are related to frequency calculations where I tried to check if the problem is related to freq calc seeting or not i.e. I used Ediff =10^ -8 and POTIM= 0.005 and in two different runs.
POTIM = 0.01
EDIFF = 1.00e-05
EDIFFG = -0.1E-03
ALGO = Fast
ISMEAR = 0
SIGMA = 0.01
ISPIN = 2
IBRION = 5
NELM = 200
NSW = 1
LREAL= Auto
NFREE = 2
LWAVE = FALSE
PREC = High
LWAVE = .FALSE.
LCHARG = .FALSE.

Best Regards
Abas

[abmo]

Dear Alex,

Thanks for the tips. I will try what you suggested. In the last two set I used Ediff =10^ -8 and POTIM= 0.005 to check whether I can improve the frequency calculations.

Regards
Abas

[abmo]

Dear Alex,

Thanks for the tips. I will try what you suggested. In the last two set I used Ediff =10^ -8 and POTIM= 0.005 to check whether I can improve the frequency calculations.

Regards
Abas

[abmo]

Dear Alex,

Thanks for the tips. I will try what you suggested. In the last two set I used Ediff =10^ -8 and POTIM= 0.005 to check whether I can improve the frequency calculations.

Regards
Abas