pwtools.pydos.pdos

pwtools.pydos.pdos(vel, dt=1.0, m=None, full_out=False, area=1.0, window=True, npad=None, tonext=False, mirr=False, method='direct')[source]

Phonon DOS by FFT of the VACF or direct FFT of atomic velocities.

Integral area is normalized to area. It is possible (and recommended) to zero-padd the velocities (see npad).

Parameters:
  • vel (3d array (nstep, natoms, 3)) – atomic velocities

  • dt (time step)

  • m (1d array (natoms,),) – atomic mass array, if None then mass=1.0 for all atoms is used

  • full_out (bool)

  • area (float) – normalize area under frequency-PDOS curve to this value

  • window (bool) – use Welch windowing on data before FFT (reduces leaking effect, recommended)

  • npad ({None, int}) – method=’direct’ only: Length of zero padding along axis. npad=None = no padding, npad > 0 = pad by a length of (nstep-1)*npad. npad > 5 usually results in sufficient interpolation.

  • tonext (bool) – method=’direct’ only: Pad vel with zeros along axis up to the next power of two after the array length determined by npad. This gives you speed, but variable (better) frequency resolution.

  • mirr (bool) – method=’vacf’ only: mirror one-sided VACF at t=0 before fft

Returns:

  • if full_out = False – | (faxis, pdos) | faxis : 1d array [1/unit(dt)] | pdos : 1d array, the phonon DOS, normalized to area

  • if full_out = True – | if method == ‘direct’: | (faxis, pdos, (full_faxis, full_pdos, split_idx)) | if method == ‘vavcf’: | (faxis, pdos, (full_faxis, full_pdos, split_idx, vacf, fft_vacf)) | fft_vacf : 1d complex array, result of fft(vacf) or fft(mirror(vacf)) | vacf : 1d array, the VACF

Examples

>>> from pwtools.constants import fs,rcm_to_Hz
>>> tr = Trajectory(...)
>>> # freq in [Hz] if timestep in [s]
>>> freq,dos = pdos(tr.velocity, m=tr.mass, dt=tr.timestep*fs,
>>>                 method='direct', npad=1)
>>> # frequency in [1/cm]
>>> plot(freq/rcm_to_Hz, dos)

Notes

padding (only method=’direct’): With npad we pad the velocities vel with npad*(nstep-1) zeros along axis (the time axis) before FFT b/c the signal is not periodic. For npad=1, this gives us the exact same spectrum and frequency resolution as with pdos(..., method='vacf',mirr=True) b/c the array to be fft’ed has length 2*nstep-1 along the time axis in both cases (remember that the array length = length of the time axis influences the freq. resolution). FFT is only fast for arrays with length = a power of two. Therefore, you may get very different fft speeds depending on whether 2*nstep-1 is a power of two or not (in most cases it won’t). Try using tonext but remember that you get another (better) frequency resolution.

References

[1] Phys Rev B 47(9) 4863, 1993