Dear GULP community,
I am working on obtaining the phonon dispersion for a carbon nanotube with chirality (10,10). I relaxed the structure under constant pressure using the optimized Tersoff potential (Lindsay 2010). Afterward, I calculated the phonon dispersion but encountered a negative frequency at the Gamma point.
Has my system truly reached its energy minimum? If not, how can I achieve that?
if the problem isn’t related to the optimization procedure, how can I address this issue?
optimise conp comp phonon dispersion eigenvector
name cnt
vectors
100 0 0
0 100 0
0 0 2.456
# ************** system ************* #
cartesian region 1
C core 21.78000000 15.00000000 0.00000000 0.0 1.0
C core 21.63200000 16.41000000 0.00000000 0.0 1.0
C core 21.44800000 17.09500000 1.228 0.0 1.0
C core 20.87100000 18.39000000 1.228 0.0 1.0
C core 20.48500000 18.98500000 0.00000000 0.0 1.0
C core 19.53700000 20.03800000 0.00000000 0.0 1.0
C core 18.98500000 20.48500000 1.228 0.0 1.0
C core 17.75800000 21.19400000 1.228 0.0 1.0
C core 17.09500000 21.44800000 0.00000000 0.0 1.0
C core 15.70900000 21.74300000 0.00000000 0.0 1.0
C core 15.00000000 21.78000000 1.228 0.0 1.0
C core 13.59000000 21.63200000 1.228 0.0 1.0
C core 12.90500000 21.44800000 0.00000000 0.0 1.0
C core 11.61000000 20.87100000 0.00000000 0.0 1.0
C core 11.01500000 20.48500000 1.228 0.0 1.0
C core 9.96200000 19.53700000 1.228 0.0 1.0
C core 9.51500000 18.98500000 0.00000000 0.0 1.0
C core 8.80600000 17.75800000 0.00000000 0.0 1.0
C core 8.55200000 17.09500000 1.228 0.0 1.0
C core 8.25700000 15.70900000 1.228 0.0 1.0
C core 8.22000000 15.00000000 0.00000000 0.0 1.0
C core 8.36800000 13.59000000 0.00000000 0.0 1.0
C core 8.55200000 12.90500000 1.228 0.0 1.0
C core 9.12900000 11.61000000 1.228 0.0 1.0
C core 9.51500000 11.01500000 0.00000000 0.0 1.0
C core 10.46300000 9.96200000 0.00000000 0.0 1.0
C core 11.01500000 9.51500000 1.228 0.0 1.0
C core 12.24200000 8.80600000 1.228 0.0 1.0
C core 12.90500000 8.55200000 0.00000000 0.0 1.0
C core 14.29100000 8.25700000 0.00000000 0.0 1.0
C core 15.00000000 8.22000000 1.228 0.0 1.0
C core 16.41000000 8.36800000 1.228 0.0 1.0
C core 17.09500000 8.55200000 0.00000000 0.0 1.0
C core 18.39000000 9.12900000 0.00000000 0.0 1.0
C core 18.98500000 9.51500000 1.228 0.0 1.0
C core 20.03800000 10.46300000 1.228 0.0 1.0
C core 20.48500000 11.01500000 0.00000000 0.0 1.0
C core 21.19400000 12.24200000 0.00000000 0.0 1.0
C core 21.44800000 12.90500000 1.228 0.0 1.0
C core 21.74300000 14.29100000 1.228 0.0 1.0
element
mass C 12.01
end
# ************** Optimized tersoff potential ************* #
species
C core 0.0
boonebody
C core 0.0 0.00000015724 0.72751 0.72751
botwobody
C core C core 1393.600 430.0 3.48790 2.2119 1.8 2.1
borepulsive
C core 3 1.0 2.2119
boattractive theta
C core 3 1.0 2.2119 38049.0 4.3484 -0.930
shrink 0 0 1
dispersion 1 10
0.0 0.0 0.0 to 0.0 0.0 0.5
# ************** output ************* #
output xyz opt_sys
output xr opt_sys
output phonon dos
output eig eigen_vec
# ********************************************** #
Hi Vivekkumar
An imaginary frequency (which GULP writes as negative for simplicity) means your structure is unstable and wants to lower its symmetry along the eigenvector direction. You can use the “lower” keyword to displace the atoms in a way to break the symmetry & then restart the optimisation using “rfo” which checks the Hessian structure for stability.
Regards,
Julian
Dear Julian,
Thank you for your response. I used the lower_symmetry command along with RFO, and as expected, the negative frequency disappeared. Instead of breaking the symmetry, I modified the lattice vectors from (60, 60, 2.456) to (30, 30, 2.456). This change in the lattice vectors influence the frequency. My CNT’s z-axis is approximately aligned at (14 Å, 14 Å). I systematically calculated the final energy of my system for different lattice vectors, and noticed that the energy increases and then stabilizes at a constant value.
Here, it can be observed that the final enegry for the RFO method and lattice vectors above 60Å (without breaking symmetry )are identical, except for the frequency. Based on the above data, it can be concluded that the minimized structure was achieved using the RFO optimization.
However, my question is: despite lattice vectors larger than 60 Å having the minimum energy, why does the structure remain unstable?
My second question is related to the supercell method. Suppose I have a 1 x 1 x N supercell (where N represents the number of unit cells). During the dispersion calculation, which option should I enable to obtain 3n phonon modes (where n is the number of atoms per unit cell) instead of 3Nn modes?
Hi Vivek,
The problem with what you’re doing is varying the cell parameters for a model that is short-range only & so there is no interaction across the vacuum part between images. Therefore it’s not sensible to optimise these directions. It would be smarter to use 1-D boundary boundary conditions for a 1-D tube model.
Using a N times supercell you will always get 3Nn modes. If you want 3n modes then use a single cell & k points.
Regards,
Julian
Dear Julian,
Thank you for your suggestion. Following your advice, I conducted the 1D relaxation for my quasi-1D system. After adjusting the ftol to 1e-6, the negative frequency was eliminated. The Total lattice energy is now -316.98988968 eV, which aligns with my earlier analysis using the RFO method.
Dear Julian,
From my above 1D dispersion calculation (with the CNT being periodic along the x-direction), I noticed that if the eigenvectors along the x-direction are non-zero, then the eigenvectors along the y and z directions are zero, or vice versa. This pattern holds for all k-points except at k=0 and k=0.5. However, when I treat the system as 3D, this behavior doesn’t occur.
Could you please suggest which approach is more suitable for determining the dispersion and eigenvectors at various k-points for a quasi-1D system?
Dear Vivek,
It makes perfect sense that the motions along the tube axis and normal to it are relatively decoupled because one direction is bond stretching/compression & the other angle bending/out of plane motion. It’s also important to remember that when phonon modes are degenerate (typically happens at high symmetry points such as gamma) then the mixing of the phonon eigenvectors becomes arbitrary. Hopefully this will help you interpret the modes.
Regards,
Julian