import sys
import numpy as np
from phonopy import Phonopy
from phonopy.structure.symmetry import Symmetry
from phonopy.interface.vasp import read_vasp
from phonopy.file_IO import parse_FORCE_SETS, parse_BORN
from phonopy.structure.spglib import get_stabilized_reciprocal_mesh
from phonopy.structure.tetrahedron_method import TetrahedronMethod
from phonopy.phonon.tetrahedron_mesh import TetrahedronMesh

cell = read_vasp(sys.argv[1])
phonon = Phonopy(cell, [[2, 0, 0], [0, 2, 0], [0, 0, 2]],
                 is_auto_displacements=False)
force_sets = parse_FORCE_SETS()
phonon.set_force_sets(force_sets)
phonon.set_post_process([[0, 0.5, 0.5], [0.5, 0, 0.5], [0.5, 0.5, 0]])
primitive = phonon.get_primitive()
born = parse_BORN(primitive)
phonon.set_nac_params(born)
symmetry = phonon.get_primitive_symmetry()
mesh = [4, 4, 4]
# phonon.set_mesh(mesh)
# phonon.set_total_DOS(sigma=0.1)
# phonon.plot_total_DOS().show()

rotations = symmetry.get_pointgroup_operations()
thm = TetrahedronMesh(phonon.get_dynamical_matrix(),
                      mesh,
                      rotations,
                      is_gamma_center=True)
                      
thm.run()

@@ -484,33 +484,16 @@
  0, 2%5D%5D,
-%0A                
  is_auto
@@ -514,16 +514,16 @@
 =False)%0A
+
 force_se
@@ -798,15 +798,18 @@
  = %5B
-4, 4, 4
+20, 20, 20
 %5D%0A# 
@@ -1135,11 +1135,108 @@
 %0Athm.run
-()
+_dos()%0Afor f, iw in zip(thm.get_frequency_points(), thm.get_integration_weights()):%0A    print f, iw
 %0A

Update tetrahedron method test