phonon_supercell:
          This supercell is used for harmonic phonons (frequencies,
          eigenvectors, group velocities, ...)
        phonon_supercell_matrix:
          Different supercell size can be specified.
196816
def _build_phonon_supercell(self):
        
        if self._phonon_supercell_matrix is None:
            self._phonon_supercell = self._supercell
        else:
            self._phonon_supercell = get_supercell(
                self._unitcell, self._phonon_supercell_matrix, self._symprec)
144171
def col_frequencies(col, weights=None, gap_chars='-.'):
    """"""
    counts = col_counts(col, weights, gap_chars)
    # Reduce to frequencies
    scale = 1.0 / sum(counts.values())
    return dict((aa, cnt * scale) for aa, cnt in counts.iteritems())
725003
def as_dict(self):
        
        d = {"vasp_version": self.vasp_version,
             "has_vasp_completed": self.converged,
             "nsites": len(self.final_structure)}
        comp = self.final_structure.composition
        d["unit_cell_formula"] = comp.as_dict()
        d["reduced_cell_formula"] = Composition(comp.reduced_formula).as_dict()
        d["pretty_formula"] = comp.reduced_formula
        symbols = [s.split()[1] for s in self.potcar_symbols]
        symbols = [re.split(r"_", s)[0] for s in symbols]
        d["is_hubbard"] = self.is_hubbard
        d["hubbards"] = self.hubbards

        unique_symbols = sorted(list(set(self.atomic_symbols)))
        d["elements"] = unique_symbols
        d["nelements"] = len(unique_symbols)

        d["run_type"] = self.run_type

        vin = {"incar": {k: v for k, v in self.incar.items()},
               "crystal": self.initial_structure.as_dict(),
               "kpoints": self.kpoints.as_dict()}
        actual_kpts = [{"abc": list(self.actual_kpoints[i]),
                        "weight": self.actual_kpoints_weights[i]}
                       for i in range(len(self.actual_kpoints))]
        vin["kpoints"]["actual_points"] = actual_kpts
        vin["nkpoints"] = len(actual_kpts)
        vin["potcar"] = [s.split(" ")[1] for s in self.potcar_symbols]
        vin["potcar_spec"] = self.potcar_spec
        vin["potcar_type"] = [s.split(" ")[0] for s in self.potcar_symbols]
        vin["parameters"] = {k: v for k, v in self.parameters.items()}
        vin["lattice_rec"] = self.lattice_rec.as_dict()
        d["input"] = vin

        nsites = len(self.final_structure)

        try:
            vout = {"ionic_steps": self.ionic_steps,
                    "final_energy": self.final_energy,
                    "final_energy_per_atom": self.final_energy / nsites,
                    "crystal": self.final_structure.as_dict(),
                    "efermi": self.efermi}
        except (ArithmeticError, TypeError):
            vout = {"ionic_steps": self.ionic_steps,
                    "final_energy": self.final_energy,
                    "final_energy_per_atom": None,
                    "crystal": self.final_structure.as_dict(),
                    "efermi": self.efermi}

        if self.eigenvalues:
            eigen = {str(spin): v.tolist()
                     for spin, v in self.eigenvalues.items()}
            vout["eigenvalues"] = eigen
            (gap, cbm, vbm, is_direct) = self.eigenvalue_band_properties
            vout.update(dict(bandgap=gap, cbm=cbm, vbm=vbm,
                             is_gap_direct=is_direct))

            if self.projected_eigenvalues:
                vout['projected_eigenvalues'] = {
                    str(spin): v.tolist()
                    for spin, v in self.projected_eigenvalues.items()}

        vout['epsilon_static'] = self.epsilon_static
        vout['epsilon_static_wolfe'] = self.epsilon_static_wolfe
        vout['epsilon_ionic'] = self.epsilon_ionic
        d['output'] = vout
        return jsanitize(d, strict=True)
566784
def types_of_specie(self):
        
        if not self.is_ordered:
            raise TypeError("""\
types_of_species cannot be used with disordered structures and partial occupancies.
Use OrderDisorderedStructureTransformation or EnumerateStructureTransformation
to build an appropriate supercell from partial occupancies.""")

        # Cannot use set since we want a deterministic algorithm.
        types = []
        for site in self:
            if site.specie not in types:
                types.append(site.specie)
        return types
57043
private void orthes() {
    // FIXME: does this fail on NaN/inf values?

    // This is derived from the Algol procedures orthes and ortran,
    // by Martin and Wilkinson, Handbook for Auto. Comp.,
    // Vol.ii-Linear Algebra, and the corresponding
    // Fortran subroutines in EISPACK.

    final int low = 0, high = n - 1;
    for(int m = low + 1; m <= high - 1; m++) {
      // Scale column.
      double scale = 0.0;
      for(int i = m; i <= high; i++) {
        scale += abs(H[i][m - 1]);
      }
      if(scale > 0.0 || scale < 0.0) {
        // Compute Householder transformation.
        double h = 0.0;
        for(int i = high; i >= m; i--) {
          double oi = ort[i] = H[i][m - 1] / scale;
          h += oi * oi;
        }
        double g = FastMath.sqrt(h), om = ort[m];
        g = (om > 0) ? -g : g;
        h -= om * g;
        ort[m] = om - g;

        // Apply Householder similarity transformation
        // H = (I-u*u'/h)*H*(I-u*u')/h)

        for(int j = m; j < n; j++) {
          double f = 0.0;
          for(int i = high; i >= m; i--) {
            f += ort[i] * H[i][j];
          }
          f /= h;
          for(int i = m; i <= high; i++) {
            H[i][j] -= f * ort[i];
          }
        }

        for(int i = 0; i <= high; i++) {
          final double[] Hi = H[i];
          double f = 0.0;
          for(int j = high; j >= m; j--) {
            f += ort[j] * Hi[j];
          }
          f /= h;
          for(int j = m; j <= high; j++) {
            Hi[j] -= f * ort[j];
          }
        }
        ort[m] *= scale;
        H[m][m - 1] = scale * g;
      }
    }

    // Accumulate transformations (Algol's ortran).
    for(int i = 0; i < n; i++) {
      Arrays.fill(V[i], 0.);
      V[i][i] = 1.;
    }

    for(int m = high - 1; m >= low + 1; m--) {
      final double[] H_m = H[m];
      if(H_m[m - 1] != 0.0) {
        for(int i = m + 1; i <= high; i++) {
          ort[i] = H[i][m - 1];
        }
        final double ort_m = ort[m];
        for(int j = m; j <= high; j++) {
          double g = 0.0;
          for(int i = m; i <= high; i++) {
            g += ort[i] * V[i][j];
          }
          // Double division avoids possible underflow
          g = (g / ort_m) / H_m[m - 1];
          for(int i = m; i <= high; i++) {
            V[i][j] += g * ort[i];
          }
        }
      }
    }
  }
694568
def lm_ffinal(freqs, damping_times, modes):
    

    f_max = {}
    for lmn in modes:
        l, m, nmodes = int(lmn[0]), int(lmn[1]), int(lmn[2])
        for n in range(nmodes):
            f_max['%d%d%d' %(l,m,n)] = qnm_freq_decay(freqs['%d%d%d' %(l,m,n)],
                                      damping_times['%d%d%d' %(l,m,n)], 1./1000)

    f_final = max(f_max.values())
    if f_final > max_freq:
        f_final = max_freq

    return f_final
969452
def multiplicity(self):
        
        sga = SpacegroupAnalyzer(self.bulk_structure)
        periodic_struc = sga.get_symmetrized_structure()
        poss_deflist = sorted(
            periodic_struc.get_sites_in_sphere(self.site.coords, 2, include_index=True), key=lambda x: x[1])
        defindex = poss_deflist[0][2]

        equivalent_sites = periodic_struc.find_equivalent_sites(self.bulk_structure[defindex])
        return len(equivalent_sites)
265785
def frequency(self, context):
        
        channels = self._manager.spectral_window_table.getcol(MS.CHAN_FREQ)
        return channels.reshape(context.shape).astype(context.dtype)
1063673
def is_group_or_super_group(cls, obj) -> bool:
        
        return cls._check(obj, [cls.GROUP, cls.SUPER_GROUP])
1098289
def _transform(self, q):
        

        s2p = self._primitive.get_supercell_to_primitive_map()
        p2s = self._primitive.get_primitive_to_supercell_map()

        num_s = self._supercell.get_number_of_atoms()
        num_p = self._primitive.get_number_of_atoms()
        v = self._velocities

        q_array = np.reshape(q, (-1, 3))
        dtype = "c%d" % (np.dtype('double').itemsize * 2)
        v_q = np.zeros((v.shape[0], num_p, len(q_array), 3), dtype=dtype)

        for p_i, s_i in enumerate(p2s):
            for s_j, s2p_j in enumerate(s2p):
                if s2p_j == s_i:
                    for q_i, pf in enumerate(
                            self._get_phase_factor(p_i, s_j, q_array)):
                        v_q[:, p_i, q_i, :] += pf * v[:, s_j, :]
        return v_q
548195
def _init_tmatrix(self):
        

        if self.radius_type == Scatterer.RADIUS_MAXIMUM:
            # Maximum radius is not directly supported in the original
            # so we convert it to equal volume radius
            radius_type = Scatterer.RADIUS_EQUAL_VOLUME
            radius = self.equal_volume_from_maximum()
        else:
            radius_type = self.radius_type
            radius = self.radius

        self.nmax = pytmatrix.calctmat(radius, radius_type,
            self.wavelength, self.m.real, self.m.imag, self.axis_ratio,
            self.shape, self.ddelt, self.ndgs)
        self._tm_signature = (self.radius, self.radius_type, self.wavelength,
            self.m, self.axis_ratio, self.shape, self.ddelt, self.ndgs)
1218145
def as_phononwebsite(self):
        
        d = {}

        #define the lattice
        d["lattice"] = self.structure.lattice._matrix.tolist()
      
        #define atoms
        atom_pos_car = []
        atom_pos_red = []
        atom_types = []
        for site in self.structure.sites:
            atom_pos_car.append(site.coords.tolist())
            atom_pos_red.append(site.frac_coords.tolist())
            atom_types.append(site.species_string)

        #default for now
        d["repetitions"] = get_reasonable_repetitions(len(atom_pos_car))
    
        d["natoms"] = len(atom_pos_car)
        d["atom_pos_car"] = atom_pos_car
        d["atom_pos_red"] = atom_pos_red
        d["atom_types"] = atom_types
        d["atom_numbers"] = self.structure.atomic_numbers
        d["formula"] = self.structure.formula
        d["name"] = self.structure.formula

        #get qpoints
        qpoints = []
        for q in self.qpoints:
            qpoints.append(list(q.frac_coords))
        d["qpoints"] = qpoints

        # get labels
        hsq_dict = collections.OrderedDict()
        for nq,q in enumerate(self.qpoints):
            if q.label is not None:
                hsq_dict[nq] = q.label

        #get distances
        dist = 0
        nqstart = 0
        distances = [dist]
        line_breaks = []
        for nq in range(1,len(qpoints)):
            q1 = np.array(qpoints[nq])
            q2 = np.array(qpoints[nq-1])
            #detect jumps
            if ((nq in hsq_dict) and (nq-1 in hsq_dict)):
                if (hsq_dict[nq] != hsq_dict[nq-1]):
                    hsq_dict[nq-1] += "|"+hsq_dict[nq]
                del hsq_dict[nq]
                line_breaks.append((nqstart,nq))
                nqstart = nq
            else:
                dist += np.linalg.norm(q1-q2)
            distances.append(dist)
        line_breaks.append((nqstart,len(qpoints)))
        d["distances"] = distances
        d["line_breaks"] = line_breaks
        d["highsym_qpts"] = list(hsq_dict.items())

        #eigenvalues
        thz2cm1 = 33.35641
        bands = self.bands.copy()*thz2cm1
        d["eigenvalues"] = bands.T.tolist()

        #eigenvectors
        eigenvectors = self.eigendisplacements.copy()
        eigenvectors /= np.linalg.norm(eigenvectors[0,0])
        eigenvectors = eigenvectors.swapaxes(0,1)
        eigenvectors = np.array([eigenvectors.real, eigenvectors.imag])
        eigenvectors = np.rollaxis(eigenvectors,0,5)
        d["vectors"] = eigenvectors.tolist()
 
        return d
338540
function DFT(bufferSize, sampleRate) {
  FourierTransform.call(this, bufferSize, sampleRate);

  var N = bufferSize/2 * bufferSize;
  var TWO_PI = 2 * Math.PI;

  this.sinTable = new Float32Array(N);
  this.cosTable = new Float32Array(N);

  for (var i = 0; i < N; i++) {
    this.sinTable[i] = Math.sin(i * TWO_PI / bufferSize);
    this.cosTable[i] = Math.cos(i * TWO_PI / bufferSize);
  }
}
662843
def create_supercell(system, supercell=[[-1, 1], [-1, 1], [-1, 1]]):
    """"""
    if 'lattice' not in system.keys():
        matrix = unit_cell_to_lattice_array(system['unit_cell'])
    else:
        matrix = system['lattice']
    coordinates = deepcopy(system['coordinates'])
    multiplication_matrices = []
    for a_ in range(supercell[0][0], supercell[0][1] + 1):
        for b_ in range(supercell[1][0], supercell[1][1] + 1):
            for c_ in range(supercell[2][0], supercell[2][1] + 1):
                mult_matrix = np.array([[a_, b_, c_]])
                mult_matrix = np.repeat(
                    mult_matrix, coordinates.shape[0], axis=0)
                multiplication_matrices.append(mult_matrix)
    frac_coordinates = cart2frac_all(coordinates, matrix)
    updated_coordinates = []
    for mat in multiplication_matrices:
        updated_coor = frac_coordinates + mat
        updated_coordinates.append(updated_coor)
    supercell_frac_coordinates = np.concatenate(updated_coordinates, axis=0)
    supercell_coordinates = frac2cart_all(supercell_frac_coordinates, matrix)
    # Now for each new cell in the supercell we need to repeat the
    # elements array so that it maches
    new_elements = deepcopy(system['elements'])
    new_ids = deepcopy(system['atom_ids'])
    for i in range(len(updated_coordinates) - 1):
        new_elements = np.concatenate((new_elements, system['elements']))
        new_ids = np.concatenate((new_ids, system['atom_ids']))
    cryst = lattice_array_to_unit_cell(matrix)
    supercell_system = {
        'elements': new_elements,
        'atom_ids': new_ids,
        'coordinates': supercell_coordinates,
        'unit_cell': cryst,
        'lattice': matrix,
    }
    return supercell_system
1320422
function computeHOAMatrices(matrix) {
  // We start by computing the 2nd-order matrix from the 1st-order matrix.
  for (let i = 2; i <= matrix.length; i++) {
    computeBandRotation(matrix, i);
  }
}
1325181
private void calcInterfacesCrystal(StructureInterfaceList set, double cutoff) {


		// initialising debugging vars
		long start = -1;
		long end = -1;
		int trialCount = 0;
		int skippedRedundant = 0;
		int skippedAUsNoOverlap = 0;
		int skippedChainsNoOverlap = 0;
		int skippedSelfEquivalent = 0;

		// The bounding boxes of all AUs of the unit cell
		UnitCellBoundingBox bbGrid = new UnitCellBoundingBox(numOperatorsSg, numPolyChainsAu);;
		// we calculate all the bounds of each of the asym units, those will then be reused and translated
		bbGrid.setBbs(structure, ops, INCLUDE_HETATOMS);


		// if not crystallographic there's no search to do in other cells, only chains within "AU" will be checked
		if (!searchBeyondAU) numCells = 0;

		boolean verbose = logger.isDebugEnabled();

		if (verbose) {
			trialCount = 0;
			start= System.currentTimeMillis();
			int neighbors = (2*numCells+1)*(2*numCells+1)*(2*numCells+1)-1;
			int auTrials = (numPolyChainsAu*(numPolyChainsAu-1))/2;
			int trials = numPolyChainsAu*numOperatorsSg*numPolyChainsAu*neighbors;
			logger.debug("Chain clash trials within original AU: "+auTrials);
			logger.debug(
					"Chain clash trials between the original AU and the neighbouring "+neighbors+
					" whole unit cells ("+numCells+" neighbours)" +
					"(2x"+numPolyChainsAu+"chains x "+numOperatorsSg+"AUs x "+neighbors+"cells) : "+trials);
			logger.debug("Total trials: "+(auTrials+trials));
		}

		List<Chain> polyChains = structure.getPolyChains();

		for (int a=-numCells;a<=numCells;a++) {
			for (int b=-numCells;b<=numCells;b++) {
				for (int c=-numCells;c<=numCells;c++) {

					Point3i trans = new Point3i(a,b,c);
					Vector3d transOrth = new Vector3d(a,b,c);
					if (a!=0 || b!=0 || c!=0) {
						// we avoid doing the transformation for 0,0,0 (in case it's not crystallographic)
						this.crystallographicInfo.getCrystalCell().transfToOrthonormal(transOrth);
					}

					UnitCellBoundingBox bbGridTrans = bbGrid.getTranslatedBbs(transOrth);

					for (int n=0;n<numOperatorsSg;n++) {

						// short-cut strategies
						// 1) we skip first of all if the bounding boxes of the AUs don't overlap
						if (!bbGrid.getAuBoundingBox(0).overlaps(bbGridTrans.getAuBoundingBox(n), cutoff)) {
							skippedAUsNoOverlap++;
							continue;
						}

						// 2) we check if we didn't already see its equivalent symmetry operator partner
						CrystalTransform tt = new CrystalTransform(this.crystallographicInfo.getSpaceGroup(), n);
						tt.translate(trans);
						if (isRedundantTransform(tt)) {
							skippedRedundant++;
							continue;
						}
						addVisitedTransform(tt);


						boolean selfEquivalent = false;

						// 3) an operator can be "self redundant" if it is the inverse of itself (involutory, e.g. all pure 2-folds with no translation)
						if (tt.isEquivalent(tt)) {
							logger.debug("Transform {} is equivalent to itself, will skip half of i-chains to j-chains comparisons", tt.toString());
							// in this case we can't skip the operator, but we can skip half of the matrix comparisons e.g. j>i
							// we set a flag and do that within the loop below
							selfEquivalent = true;
						}

						StringBuilder builder = null;
						if (verbose) builder = new StringBuilder(String.valueOf(tt)).append(" ");

						// Now that we know that boxes overlap and operator is not redundant, we have to go to the details
						int contactsFound = 0;

						for (int j=0;j<numPolyChainsAu;j++) {

							for (int i=0;i<numPolyChainsAu;i++) { // we only have to compare the original asymmetric unit to every full cell around

								if(selfEquivalent && (j>i)) {
									// in case of self equivalency of the operator we can safely skip half of the matrix
									skippedSelfEquivalent++;
									continue;
								}
								// special case of original AU, we don't compare a chain to itself
								if (n==0 && a==0 && b==0 && c==0 && i==j) continue;

								// before calculating the AtomContactSet we check for overlap, then we save putting atoms into the grid
								if (!bbGrid.getChainBoundingBox(0,i).overlaps(bbGridTrans.getChainBoundingBox(n,j),cutoff)) {
									skippedChainsNoOverlap++;
									if (verbose) {
										builder.append(".");
									}
									continue;
								}

								trialCount++;

								// finally we've gone through all short-cuts and the 2 chains seem to be close enough:
								// we do the calculation of contacts
								Chain chaini = polyChains.get(i);
								Chain chainj = polyChains.get(j);

								if (n!=0 || a!=0 || b!=0 || c!=0) {
									Matrix4d mJCryst = new Matrix4d(ops[n]);
									translate(mJCryst, transOrth);
									chainj = (Chain)chainj.clone();
									Calc.transform(chainj,mJCryst);
								}

								StructureInterface interf = calcContacts(chaini, chainj, cutoff, tt, builder);
								if (interf == null) {
									continue;
								}

								contactsFound++;
								if(this.hasNcsOps()) {
									StructureInterface interfNcsRef = findNcsRef(interf);
									set.addNcsEquivalent(interf,interfNcsRef);
								} else {
									set.add(interf);
								}
							}
						}

						if( verbose ) {
							if (a==0 && b==0 && c==0 && n==0)
								builder.append(" "+contactsFound+"("+(numPolyChainsAu*(numPolyChainsAu-1))/2+")");
							else if (selfEquivalent)
								builder.append(" "+contactsFound+"("+(numPolyChainsAu*(numPolyChainsAu+1))/2+")");
							else
								builder.append(" "+contactsFound+"("+numPolyChainsAu*numPolyChainsAu+")");

							logger.debug(builder.toString());
						}
					}
				}
			}
		}

		end = System.currentTimeMillis();
		logger.debug("\n"+trialCount+" chain-chain clash trials done. Time "+(end-start)/1000+"s");
		logger.debug("  skipped (not overlapping AUs)       : "+skippedAUsNoOverlap);
		logger.debug("  skipped (not overlapping chains)    : "+skippedChainsNoOverlap);
		logger.debug("  skipped (sym redundant op pairs)    : "+skippedRedundant);
		logger.debug("  skipped (sym redundant self op)     : "+skippedSelfEquivalent);
		logger.debug("Found "+set.size()+" interfaces.");
	}
1124921
def repeat(nx=1, ny=1, nz=1)
      
      raise "Not a periodic system." if self.lattice_vectors.nil?
      
      u = self.copy
      v1 = self.lattice_vectors[0]
      v2 = self.lattice_vectors[1]
      v3 = self.lattice_vectors[2]

      nx_sign = (0 < nx) ? 1 : -1
      ny_sign = (0 < ny) ? 1 : -1
      nz_sign = (0 < nz) ? 1 : -1
      
      new_atoms = []
      nx.to_i.abs.times do |i|        
        ny.to_i.abs.times do |j|
          nz.to_i.abs.times do |k|
            new_atoms << self.displace(nx_sign*i*v1[0] + ny_sign*j*v2[0] + nz_sign*k*v3[0], 
                                       nx_sign*i*v1[1] + ny_sign*j*v2[1] + nz_sign*k*v3[1], 
                                       nx_sign*i*v1[2] + ny_sign*j*v2[2] + nz_sign*k*v3[2]).atoms
          end
        end
      end

      u.atoms = new_atoms.flatten
      u.lattice_vectors = [Vector[nx.abs*v1[0], nx.abs*v1[1], nx.abs*v1[2]], 
                           Vector[ny.abs*v2[0], ny.abs*v2[1], ny.abs*v2[2]], 
                           Vector[nz.abs*v3[0], nz.abs*v3[1], nz.abs*v3[2]]]
      # u.make_bonds 
      return u     
    end