Treatment of undefined atomic masses

phonopy/api_phonopy.py

@@ -419,9 +419,14 @@ class Phonopy:
     # Single q-point
     def get_dynamical_matrix_at_q(self, q):
         self._set_dynamical_matrix()
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            return None
+
         self._dynamical_matrix.set_dynamical_matrix(q)
         return self._dynamical_matrix.get_dynamical_matrix()
 
+
     def get_frequencies(self, q):
         """
         Calculate phonon frequencies at q
@@ -429,6 +434,10 @@ class Phonopy:
         q: q-vector in reduced coordinates of primitive cell
         """
         self._set_dynamical_matrix()
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            return None
+
         self._dynamical_matrix.set_dynamical_matrix(q)
         dm = self._dynamical_matrix.get_dynamical_matrix()
         frequencies = []
@@ -447,6 +456,10 @@ class Phonopy:
         q: q-vector in reduced coordinates of primitive cell
         """
         self._set_dynamical_matrix()
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            return None
+
         self._dynamical_matrix.set_dynamical_matrix(q)
         dm = self._dynamical_matrix.get_dynamical_matrix()
         frequencies = []
@@ -465,6 +478,11 @@ class Phonopy:
                            bands,
                            is_eigenvectors=False,
                            is_band_connection=False):
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            self._band_structure = None
+            return False
+
         self._band_structure = BandStructure(
             bands,
             self._dynamical_matrix,
@@ -472,6 +490,7 @@ class Phonopy:
             is_band_connection=is_band_connection,
             group_velocity=self._group_velocity,
             factor=self._factor)
+        return True
 
     def get_band_structure(self):
         band = self._band_structure
@@ -500,6 +519,11 @@ class Phonopy:
                  is_mesh_symmetry=True,
                  is_eigenvectors=False,
                  is_gamma_center=False):
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            self._mesh = None
+            return False
+
         self._mesh = Mesh(
             self._dynamical_matrix,
             mesh,
@@ -512,6 +536,7 @@ class Phonopy:
             rotations=self._primitive_symmetry.get_pointgroup_operations(),
             factor=self._factor,
             use_lapack_solver=self._use_lapack_solver)
+        return True
 
     def get_mesh(self):
         return (self._mesh.get_qpoints(),
@@ -827,6 +852,11 @@ class Phonopy:
                            is_eigenvectors=False,
                            write_dynamical_matrices=False,
                            factor=VaspToTHz):
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            self._qpoints_phonon = None
+            return False
+
         self._qpoints_phonon = QpointsPhonon(
             q_points,
             self._dynamical_matrix,
@@ -835,6 +865,7 @@ class Phonopy:
             group_velocity=self._group_velocity,
             write_dynamical_matrices=write_dynamical_matrices,
             factor=self._factor)
+        return True
         
     def get_qpoints_phonon(self):
         return (self._qpoints_phonon.get_frequencies(),
@@ -852,6 +883,10 @@ class Phonopy:
                         num_div=None,
                         shift=None,
                         filename=None):
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            return False
+
         if q_point is None:
             animation = Animation([0, 0, 0],
                                   self._dynamical_matrix,
@@ -926,6 +961,9 @@ class Phonopy:
                     animation.write_POSCAR(band_index,
                                            amplitude,
                                            num_div)
+
+        return True
+
     # Atomic modulation of normal mode
     def set_modulations(self,
                         dimension,
@@ -933,6 +971,11 @@ class Phonopy:
                         delta_q=None,
                         derivative_order=None,
                         nac_q_direction=None):
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            self._modulation = None
+            return False
+        
         self._modulation = Modulation(self._dynamical_matrix,
                                       dimension,
                                       phonon_modes,
@@ -941,6 +984,7 @@ class Phonopy:
                                       nac_q_direction=nac_q_direction,
                                       factor=self._factor)
         self._modulation.run()
+        return True
                     
     def get_modulated_supercells(self):
         """Returns cells with modulations as Atoms objects"""
@@ -970,6 +1014,11 @@ class Phonopy:
                    is_little_cogroup=False,
                    nac_q_direction=None,
                    degeneracy_tolerance=1e-4):
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            self._irreps = None
+            return None
+
         self._irreps = IrReps(
             self._dynamical_matrix,
             q,
@@ -993,11 +1042,17 @@ class Phonopy:
 
     # Group velocity
     def set_group_velocity(self, q_length=None):
+        if self._dynamical_matrix is None:
+            print "Warning: Dynamical matrix has not yet built."
+            self._group_velocity = None
+            return False
+
         self._group_velocity = GroupVelocity(
             self._dynamical_matrix,
             q_length=q_length,
             symmetry=self._primitive_symmetry,
             frequency_factor_to_THz=self._factor)
+        return True
 
     def get_group_velocity(self):
         return self._group_velocity.get_group_velocity()
@@ -1026,10 +1081,15 @@ class Phonopy:
                 computation_algorithm=computation_algorithm)
 
     def _set_dynamical_matrix(self):
-        if (self._supercell is None or
-            self._primitive is None or
-            self._force_constants is None):
-            return False
+        if (self._supercell is None or self._primitive is None):
+            print "Bug: Supercell or primitive is not created."
+            sys.exit(1)
+        elif self._force_constants is None:
+            print "Warning: Force constants are prepared."
+            self._dynamical_matrix = None
+        elif self._primitive.get_masses() is None:
+            print "Warning: Atomic masses are not correctly set."
+            self._dynamical_matrix = None
         else:
             if self._nac_params is None:
                 self._dynamical_matrix = DynamicalMatrix(
@@ -1046,7 +1106,6 @@ class Phonopy:
                     nac_params=self._nac_params,
                     decimals=self._dynamical_matrix_decimals,
                     symprec=self._symprec)
-            return True
 
     def _search_symmetry(self):
         self._symmetry = Symmetry(self._supercell,

phonopy/interface/phonopy_yaml.py

@@ -53,12 +53,32 @@ class phonopyYaml:
     def __init__(self, filename):
         self._data = yaml.load(open(filename), Loader=Loader)
         self._lattice = self._data['lattice']
-        self._points = [x['coordinates'] for x in self._data['points']]
-        self._symbols = [x['symbol'] for x in self._data['points']]
+
+        self._points = None
+        self._symbols = None
+        self._masses = None
+
+        points = []
+        symbols = []
+        masses = []
+        for x in self._data['points']:
+            if 'coordinates' in x:
+                points.append(x['coordinates'])
+            if 'symbol' in x:
+                symbols.append(x['symbol'])
+            if 'mass' in x:
+                masses.append(x['mass'])
+        if points:
+            self._points = points
+        if symbols:
+            self._symbols = symbols
+        if masses:
+            self._masses = masses
 
     def get_atoms(self):
         return Atoms(symbols=self._symbols,
                      cell=self._lattice,
+                     masses=self._masses,
                      scaled_positions=self._points)
 
     def __str__(self):

phonopy/structure/atoms.py

@@ -118,7 +118,10 @@ class Atoms:
             self.masses = np.array(masses, dtype='double')
 
     def get_masses(self):
-        return self.masses.copy()
+        if self.masses is None:
+            return None
+        else:
+            return self.masses.copy()
 
     def set_magnetic_moments(self, magmoms):
         if magmoms is None:
@@ -163,8 +166,11 @@ class Atoms:
                                  for s in self.symbols], dtype='intc')
         
     def _symbols_to_masses(self):
-        self.masses = np.array([atom_data[symbol_map[s]][3]
-                                for s in self.symbols], dtype='double')
+        masses = [atom_data[symbol_map[s]][3] for s in self.symbols]
+        if None in masses:
+            self.masses = None
+        else:
+            self.masses = np.array(masses, dtype='double')
 
 class PhonopyAtoms(Atoms):
     def __init__(self,
@@ -207,7 +213,7 @@ class PhonopyAtoms(Atoms):
         return "\n".join(lines)
 
 atom_data = [ 
-    [  0, "X", "X", 0], # 0
+    [  0, "X", "X", None], # 0
     [  1, "H", "Hydrogen", 1.00794], # 1
     [  2, "He", "Helium", 4.002602], # 2
     [  3, "Li", "Lithium", 6.941], # 3
@@ -250,7 +256,7 @@ atom_data = [
     [ 40, "Zr", "Zirconium", 91.224], # 40
     [ 41, "Nb", "Niobium", 92.90638], # 41
     [ 42, "Mo", "Molybdenum", 95.96], # 42
-    [ 43, "Tc", "Technetium", 0], # 43
+    [ 43, "Tc", "Technetium", None], # 43
     [ 44, "Ru", "Ruthenium", 101.07], # 44
     [ 45, "Rh", "Rhodium", 102.90550], # 45
     [ 46, "Pd", "Palladium", 106.42], # 46
@@ -268,7 +274,7 @@ atom_data = [
     [ 58, "Ce", "Cerium", 140.116], # 58
     [ 59, "Pr", "Praseodymium", 140.90765], # 59
     [ 60, "Nd", "Neodymium", 144.242], # 60
-    [ 61, "Pm", "Promethium", 0], # 61
+    [ 61, "Pm", "Promethium", None], # 61
     [ 62, "Sm", "Samarium", 150.36], # 62
     [ 63, "Eu", "Europium", 151.964], # 63
     [ 64, "Gd", "Gadolinium", 157.25], # 64
@@ -291,41 +297,41 @@ atom_data = [
     [ 81, "Tl", "Thallium", 204.3833], # 81
     [ 82, "Pb", "Lead", 207.2], # 82
     [ 83, "Bi", "Bismuth", 208.98040], # 83
-    [ 84, "Po", "Polonium", 0], # 84
-    [ 85, "At", "Astatine", 0], # 85
-    [ 86, "Rn", "Radon", 0], # 86
-    [ 87, "Fr", "Francium", 0], # 87
-    [ 88, "Ra", "Radium", 0], # 88
-    [ 89, "Ac", "Actinium", 0], # 89
+    [ 84, "Po", "Polonium", None], # 84
+    [ 85, "At", "Astatine", None], # 85
+    [ 86, "Rn", "Radon", None], # 86
+    [ 87, "Fr", "Francium", None], # 87
+    [ 88, "Ra", "Radium", None], # 88
+    [ 89, "Ac", "Actinium", None], # 89
     [ 90, "Th", "Thorium", 232.03806], # 90
     [ 91, "Pa", "Protactinium", 231.03588], # 91
     [ 92, "U", "Uranium", 238.02891], # 92
-    [ 93, "Np", "Neptunium", 0], # 93
-    [ 94, "Pu", "Plutonium", 0], # 94
-    [ 95, "Am", "Americium", 0], # 95
-    [ 96, "Cm", "Curium", 0], # 96
-    [ 97, "Bk", "Berkelium", 0], # 97
-    [ 98, "Cf", "Californium", 0], # 98
-    [ 99, "Es", "Einsteinium", 0], # 99
-    [100, "Fm", "Fermium", 0], # 100
-    [101, "Md", "Mendelevium", 0], # 101
-    [102, "No", "Nobelium", 0], # 102
-    [103, "Lr", "Lawrencium", 0], # 103
-    [104, "Rf", "Rutherfordium", 0], # 104
-    [105, "Db", "Dubnium", 0], # 105
-    [106, "Sg", "Seaborgium", 0], # 106
-    [107, "Bh", "Bohrium", 0], # 107
-    [108, "Hs", "Hassium", 0], # 108
-    [109, "Mt", "Meitnerium", 0], # 109
-    [110, "Ds", "Darmstadtium", 0], # 110
-    [111, "Rg", "Roentgenium", 0], # 111
-    [112, "Cn", "Copernicium", 0], # 112
-    [113, "Uut", "Ununtrium", 0], # 113
-    [114, "Uuq", "Ununquadium", 0], # 114
-    [115, "Uup", "Ununpentium", 0], # 115
-    [116, "Uuh", "Ununhexium", 0], # 116
-    [117, "Uus", "Ununseptium", 0], # 117
-    [118, "Uuo", "Ununoctium", 0], # 118
+    [ 93, "Np", "Neptunium", None], # 93
+    [ 94, "Pu", "Plutonium", None], # 94
+    [ 95, "Am", "Americium", None], # 95
+    [ 96, "Cm", "Curium", None], # 96
+    [ 97, "Bk", "Berkelium", None], # 97
+    [ 98, "Cf", "Californium", None], # 98
+    [ 99, "Es", "Einsteinium", None], # 99
+    [100, "Fm", "Fermium", None], # 100
+    [101, "Md", "Mendelevium", None], # 101
+    [102, "No", "Nobelium", None], # 102
+    [103, "Lr", "Lawrencium", None], # 103
+    [104, "Rf", "Rutherfordium", None], # 104
+    [105, "Db", "Dubnium", None], # 105
+    [106, "Sg", "Seaborgium", None], # 106
+    [107, "Bh", "Bohrium", None], # 107
+    [108, "Hs", "Hassium", None], # 108
+    [109, "Mt", "Meitnerium", None], # 109
+    [110, "Ds", "Darmstadtium", None], # 110
+    [111, "Rg", "Roentgenium", None], # 111
+    [112, "Cn", "Copernicium", None], # 112
+    [113, "Uut", "Ununtrium", None], # 113
+    [114, "Uuq", "Ununquadium", None], # 114
+    [115, "Uup", "Ununpentium", None], # 115
+    [116, "Uuh", "Ununhexium", None], # 116
+    [117, "Uus", "Ununseptium", None], # 117
+    [118, "Uuo", "Ununoctium", None], # 118
     ]
 
 symbol_map = {

phonopy/structure/cells.py

@@ -56,7 +56,10 @@ def trim_cell(relative_axes, cell, symprec):
 
     trimed_positions = []
     trimed_numbers = []
-    trimed_masses = []
+    if masses is None:
+        trimed_masses = None
+    else:
+        trimed_masses = []
     if magmoms is None:
         trimed_magmoms = None
     else:
@@ -84,7 +87,8 @@ def trim_cell(relative_axes, cell, symprec):
             trimed_positions[num_atom] = pos
             num_atom += 1
             trimed_numbers.append(numbers[i])
-            trimed_masses.append(masses[i])
+            if masses is not None:
+                trimed_masses.append(masses[i])
             if magmoms is not None:
                 trimed_magmoms.append(magmoms[i])
             atom_map.append(i)
@@ -114,17 +118,16 @@ def print_cell(cell, mapping=None, stars=None):
             if i in stars:
                 num = "*"
         num += "%d" % (i + 1)
-        if magmoms is None:
-            print "%5s %-2s%18.14f%18.14f%18.14f %7.3f" % \
-                (num, symbols[i], v[0], v[1], v[2], masses[i]),
-        else:
-            print "%5s %-2s%18.14f%18.14f%18.14f %7.3f  %5.3f" % \
-                (num, symbols[i], v[0], v[1], v[2], masses[i], magmoms[i]),
-
+        line = ("%5s %-2s%18.14f%18.14f%18.14f" % 
+                (num, symbols[i], v[0], v[1], v[2]))
+        if masses is not None:
+            line += " %7.3f" % masses[i]
+        if magmoms is not None:
+            line += "  %5.3f" % magmoms[i]
         if mapping is None:
-            print
+            print line
         else:
-            print ">", mapping[i]+1
+            print line + " >", mapping[i]+1
 
 class Supercell(Atoms):
     """Build supercell from supercell matrix
@@ -221,7 +224,10 @@ class Supercell(Atoms):
         atom_map = []
         positions_multi = []
         numbers_multi = []
-        masses_multi = []
+        if masses is None:
+            masses_multi = None
+        else:
+            masses_multi = []
         if magmoms is None:
             magmoms_multi = None
         else:
@@ -234,9 +240,10 @@ class Supercell(Atoms):
                                                 (pos[1] + j) / multi[1],
                                                 (pos[2] + i) / multi[2]])
                         numbers_multi.append(numbers[l])
-                        masses_multi.append(masses[l])
+                        if masses is not None:
+                            masses_multi.append(masses[l])
                         atom_map.append(l)
-                        if not magmoms is None:
+                        if magmoms is not None:
                             magmoms_multi.append(magmoms[l])
 
         simple_supercell = Atoms(numbers=numbers_multi,

scripts/phonopy

@@ -45,6 +45,7 @@ from phonopy.cui.show_symmetry import check_symmetry
 from phonopy.units import *
 from phonopy.version import phonopy_version
 from phonopy.structure.cells import print_cell, determinant
+from phonopy.structure.atoms import atom_data, symbol_map
 from phonopy.interface import create_FORCE_SETS, read_crystal_structure
 
 
@@ -856,12 +857,13 @@ if settings.get_is_rotational_invariance():
 
 # Atomic species without mass case
 symbols_with_no_mass = []
-for m, s in zip(primitive.get_masses(), 
-                primitive.get_chemical_symbols()):
-    if m < 0.1:
-        if symbols_with_no_mass.count(s) == 0:
+if primitive.get_masses() is None:
+    for s in primitive.get_chemical_symbols():
+        if (atom_data[symbol_map[s]][3] is None and
+            s not in symbols_with_no_mass):
             symbols_with_no_mass.append(s)
-            print_error_message("Atomic mass of \'%s\' is not implemented." % s)
+            print_error_message(
+                "Atomic mass of \'%s\' is not implemented in phonopy." % s)
             print_error_message("MASS tag can be used to set atomic masses.")
 
 if len(symbols_with_no_mass) > 0: