From: jj mm (
Date: Tue Dec 14 2010 - 12:09:38 CST

Hi again:

I am sorry, but it seems that I am missing something here...

1) I have used the script below (first email) to calculate the phonon DOS in my system. My testing simulation lasts for 10 ps, and I am getting velocities every 0.1, 0.04, 0.02 and 0.01 ps. It means that I have 100, 250, 500 and 1000 frames, respectively. In every case I get a different DOS, the difference being up to a factor 3. I have tried to change the "delta t" value (written in the webpage for specden), but I didn't get any change. How can I get so large differences?

2) In any of the cases I have considered I got a normalized DOS, even thought I used the "fourier" switch, as stated below. Any light here?

Thanks in advance. Any help is very welcome...



P.S. This is the script I've used:

package require specden

    set mol [mol new {} waitfor all]

    set sel [atomselect $mol {all}]
    set nf [molinfo $mol get numframes]
    set na [$sel num]

    set reslist {}
    for {set a 0} {$a < $na} {incr a} {
        set dlist {}
        for {set f 0} {$f < $nf} {incr f} {
            $sel frame $f
            lappend dlist [lindex [$sel get {x y z}] $a]
        lassign [specden $dlist 413.414 2000.0 fourier] flist slist
        lappend reslist $slist

    # write out the result as: fequency, h1, h2, ...,
    set fp [open "dos-1500.dat" "w"]
    set ns [llength $flist]
    for {set i 0} {$i < $ns} {incr i} {
        puts -nonewline $fp "[lindex $flist $i] "
        set avg 0.0
        for {set a 0} {$a < $na} {incr a} {
            set val [lindex $reslist $a $i]
#           puts -nonewline $fp "$val "
            set avg [expr {$avg + $val}]
        puts $fp "[expr $avg / $na]"
    close $fp

--- El jue, 9/12/10, Axel Kohlmeyer <> escribió:

De: Axel Kohlmeyer <>
Asunto: Re: vmd-l: Phonon spectrum and VACF
Para: "jj mm" <>
Fecha: jueves, 9 de diciembre, 2010 12:26

> I have a few questions for you guys:
> 1) Is this procedure right? I mean: it is actually the phonon DOS what I get after the script is run?


> 2) I made the calculation setting k_max = 3000 cm^-1, and I got 90 frequency points. Is there any way to increase the number of wavenumber values without changing k_max? I am asking this because I calculated the integral below the DOS curve and I found it not to be normalized to 1.

the integral is normalized to 1.0, that is part of the process.
the only way to have more frequency points is to write out
velocities more frequently. the number of output frequencies
are selected based on what is possible with the spacing
provided from the input.

> 3) Is there any way to do the same thing directly from the LAMMPS trajectory file (i.e., a file which contains positions AND velocities of all the particles of the system)?

yes. if you use a recent version of VMD (you may need a alpha/beta version,
i don't remember exactly when i added this), the velocities will be
read in, too.


> Thanks a lot!
> Juanjo Melendez
> Associate Professor
> Department of Physics
> University of Extremadura
> Avda. de Elvas, s/n 06006 Badajoz (Spain)
> email:
> phone: +34 924 28 96 55
> fax: +34 924 28 96 51

Dr. Axel Kohlmeyer
Institute for Computational Molecular Science
Temple University, Philadelphia PA, USA.