def new_entry(self, structure, active=True):
     properties = {'forces': None, 'stress': None, 'energy': None}
     status = {self.tag: active}
     entry = {'structure': structure.to_dict, 'properties': properties, 'status': status}
     entry_id = self.insert_entry(entry)
     pcm_log.debug('Added new entry: %s with tag=%s: %s' % (str(entry_id), self.tag, str(active)))
     return entry_id

 def get2_duplicates(self, ids, fast=True):
     dupes_dict = {}
     dupes_list = []
     selection = self.ids_sorted(ids)
     pcm_log.debug('Searching duplicates in %d structures' % len(selection))
     for i in range(len(selection) - 1):
         ncomps = 0
         entry_id = selection[i]
         if fast and entry_id in dupes_list:
             continue
         sys.stdout.write(" %5d of %5d: " % (i, len(selection)))
         value_i = self.value(entry_id)
         for j in range(i + 1, len(selection)):
             entry_jd = selection[j]
             if fast and entry_jd in dupes_list:
                 continue
             value_j = self.value(entry_jd)
             if abs(value_i - value_j) < self.value_tol:
                 ncomps += 1
                 distance = self.distance(entry_id, entry_jd)
                 if distance < self.distance_tolerance:
                     if entry_id in dupes_dict:
                         dupes_dict[entry_id].append(entry_jd)
                     else:
                         if fast:
                             dupes_dict[entry_id] = [entry_jd]
                         else:
                             dupes_dict[entry_id] = entry_jd
                     dupes_list.append(entry_jd)
         sys.stdout.write(' comparisons: %d\n' % ncomps)
     return dupes_dict, [x for x in selection if x in dupes_list]

 def refine(self, entry_id, gtol=None):
     if self.refine:
         pcm_log.debug('Evaluating %s target density= %7.3F' % (entry_id, self.minimal_density))
         structure = self.get_structure(entry_id)
         structure, properties = self.evaluate(structure, gtol=gtol)
         self.update_properties(entry_id=entry_id, new_properties=properties)
         return self.set_structure(entry_id, structure)

def get_current_status(pcdb, entry_id, relaxator_params, verbose=False):
    max_force = 1
    max_diag_stress = 1
    max_nondiag_stress = 1
    entry = pcdb.get_entry(entry_id)

    if entry is not None and 'properties' in entry and entry['properties'] is not None:
        if 'forces' in entry['properties'] and entry['properties']['forces'] is not None:
            forces = np.array(entry['properties']['forces']).reshape((-1, 3))
            max_force = np.max(np.apply_along_axis(np.linalg.norm, 1, forces))
        if 'stress' in entry['properties'] and entry['properties']['stress'] is not None:
            stress = np.array(entry['properties']['stress']).reshape((-1, 3))
            assert (stress.shape == (2, 3))
            max_diag_stress = np.abs(np.max(np.abs(stress[0])) - (relaxator_params['pressure'] / 1602.1766208))
            max_nondiag_stress = np.max(np.abs(stress[1]))
    else:
        pcm_log.debug('Bad entry')
        print(entry)
    if max(max_force, max_diag_stress, max_nondiag_stress) > relaxator_params['target_forces'] and verbose:
        if (max_force * max_diag_stress * max_nondiag_stress) == 1.0:
            print('No forces/stress information for entry: %s' % entry['_id'])
        else:
            print('Convergence status for entry: %s' % entry['_id'])
            print('Max Interatomic Force: %9.2E [eV/Ang]' % max_force)
            print('Max Diag stress      : %9.2E with pressure: %9.2E [eV/Ang^3]' % (max_diag_stress,
                                                                                    relaxator_params['pressure'] /
                                                                                    1602.1766208))
            print('Max NonDiag stress   : %9.2E [eV/Ang^3]' % max_nondiag_stress)
    return max(max_force, max_diag_stress, max_nondiag_stress)

 def replace_by_changed(self, entry_id_old, reason=None):
     entry_id_new = self.population.move_random(entry_id_old, factor=self.delta_change, in_place=False,
                                                kind='change')
     change = {'change': 'modified', 'to': entry_id_new, 'reason': reason}
     pcm_log.debug('Modified  %s -> %s' % (str(entry_id_old), str(entry_id_new)))
     self.population.disable(entry_id_old)
     self.advance(entry_id_old, entry_id_new, change)

 def set_slater_koster(self, search_paths):
     ret = True
     if isinstance(search_paths, str):
         search_paths = [search_paths]
     elif not isinstance(search_paths, list):
         raise ValueError('search_path is not an string or list')
     self.slater_koster = []
     for ispecie in self.structure.species:
         for jspecie in self.structure.species:
             pair = ispecie + '-' + jspecie
             pair_found = False
             path = self.workdir + os.sep + pair + '.skf'
             if os.path.exists(path):
                 self.slater_koster.append(path)
                 pair_found = True
             else:
                 for ipath in search_paths:
                     path = ipath + os.sep + pair + '.skf'
                     if os.path.exists(ipath + os.sep + ispecie + '-' + jspecie + '.skf'):
                         self.slater_koster.append(path)
                         pair_found = True
             if pair_found:
                 # print 'Slater-Koster %7s : %s' % (pair, path)
                 pcm_log.debug('Slater-Koster %7s : %s' % (pair, path))
             else:
                 pcm_log.debug('ERROR: Slater_Koster for ' + pair + ' not found')
                 ret = False
     return ret

    def print_status(self):
        pcm_log.info(' %s (tag: %s)' % (self.population.name, self.population.tag))
        pcm_log.info(' Current Generation             : %4d' % self.current_generation)
        pcm_log.info(' Population (evaluated/total)   : %4d /%4d' % (len(self.population.evaluated),
                                                                     len(self.population.members)))
        pcm_log.info(' Actives (evaluated/total)      : %4d /%4d' % (len(self.population.actives_evaluated),
                                                                     len(self.population.actives)))
        pcm_log.info(' Size of Generation (this/next) : %4d /%4d\n' % (len(self.get_generation()),
                                                                       len(self.get_generation(
                                                                         self.current_generation + 1))))

        if len(self.get_generation(self.current_generation + 1)) + len(self.population.actives) != self.generation_size:
            pcm_log.debug('Change in generations')
            for i in self.old_actives:
                if i not in self.population.actives:
                    pcm_log.debug('This active disappeared: %s' % str(i))
            for i in self.population.actives:
                if i not in self.old_actives:
                    pcm_log.debug('This active appeared: %s' % str(i))
            for i in self.old_nextgen:
                if i not in self.get_generation(self.current_generation + 1):
                    pcm_log.debug('This nextgen disappeared: %s' % str(i))
            for i in self.get_generation(self.current_generation + 1):
                if i not in self.old_nextgen:
                    pcm_log.debug('This nextgen appeared: %s' % str(i))
        self.old_actives = self.population.actives
        self.old_nextgen = self.get_generation(self.current_generation + 1)

 def new_entry(self, data, active=True):
     properties = {'eigvec': list(data['eigvec'].flatten()),
                   'D': list(data['D'].flatten()),
                   'I': list(data['I'].flatten())}
     status = {self.tag: active}
     entry = {'structure': self.structure.to_dict, 'properties': properties, 'status': status}
     entry_id = self.insert_entry(entry)
     pcm_log.debug('Added new entry: %s with tag=%s: %s' % (str(entry_id), self.tag, str(active)))
     return entry_id

    def run(self):
        """
        Continuously search for suitable candidates to evaluation among a list of databases.

        :return:
        """

        procs = []
        ids_running = []

        # Creating a list to store the 'nconcurrent' running jobs
        for i in range(self.nconcurrent):
            procs.append(None)
            ids_running.append(None)

        # Main loop looking permanently for candidates for evaluation
        while True:

            to_evaluate = self.get_list_candidates()

            index = 0
            currently_evaluating = 0
            for j in range(self.nconcurrent):
                if procs[j] is not None and procs[j].is_alive():
                    currently_evaluating += 1
            print('Candidates to evaluate: %d  Candidates in evaluation: %d' % (len(to_evaluate), currently_evaluating))

            while index < len(to_evaluate):

                db_settings = dict(self.db_settings)
                # The first component of each pair in to_evaluate is the name of the database
                db_settings['name'] = self.dbname
                pcdb = get_database(db_settings)

                slot = None
                while True:
                    for j in range(self.nconcurrent):
                        if procs[j] is None or not procs[j].is_alive():
                            slot = j
                    if slot is None:
                        time.sleep(self.sleeping_time)
                    else:
                        break

                # The function is_evaluated needs two arguments, the database object and entry identifier and
                # must return a boolean to decide if the candidate should be evaluated.
                if self.is_evaluated(pcdb, to_evaluate[index]):
                    ids_running[slot] = to_evaluate[index]

                    # This is the actual call to the worker, it must be a function with 4 arguments:
                    # The database settings, the entry identifier, the working directory and arguments for the worker
                    procs[j] = Process(target=self.worker, args=(db_settings, to_evaluate[index]))
                    procs[j].start()
                else:
                    pcm_log.debug('Not evaluable: %s' % to_evaluate[index])
                index += 1
            time.sleep(self.sleeping_time)

 def new_entry(self, data, active=True):
     data = np.array(data)
     # Magnetic moments are stored in spherical coordinates
     properties = {'magmom': list(data.flatten())}
     status = {self.tag: active}
     entry={'structure': self.structure.to_dict, 'properties': properties, 'status': status}
     entry_id = self.insert_entry(entry)
     pcm_log.debug('Added new entry: %s with tag=%s: %s' % (str(entry_id), self.tag, str(active)))
     return entry_id

    def run(self):

        n = 10
        dftb = DFTBplus()
        kpoints = KPoints.optimized_grid(self.structure.lattice, kp_density=10000, force_odd=True)
        dftb.initialize(workdir=self.workdir, structure=self.structure, kpoints=kpoints)
        ans = dftb.set_slater_koster(search_paths=self.slater_path)
        if not ans:
            print('Slater-Koster files not complete')
            return

        grid = None
        energies = []

        while True:
            density = n ** 3
            kpoints = KPoints.optimized_grid(self.structure.lattice, kp_density=density, force_odd=True)
            if np.sum(grid) != np.sum(kpoints.grid):
                pcm_log.debug('Trial density: %d  Grid: %s' % (density, kpoints.grid))
                grid = list(kpoints.grid)
                dftb.kpoints = kpoints

                dftb.basic_input()
                dftb.hamiltonian['MaxSCCIterations'] = 50
                if os.path.isfile('charges.bin'):
                    dftb.hamiltonian['ReadInitialCharges'] = True
                dftb.hamiltonian['Mixer'] = {'name': 'DIIS'}
                dftb.set_static()
                dftb.set_inputs()
                dftb.run()
                if self.waiting:
                    dftb.runner.wait()
                while True:
                    if dftb.runner is not None and dftb.runner.poll() is not None:
                        pcm_log.info('Execution completed. Return code %d' % dftb.runner.returncode)
                        filename = dftb.workdir + os.sep + 'detailed.out'
                        if os.path.exists(filename):
                            ret = read_detailed_out(filename)
                            line = 'KPoint_grid= %15s  iSCC= %4d  Total_energy= %10.4f  SCC_error= %9.3E'
                            print(line % (grid, ret['SCC']['iSCC'], ret['total_energy'], ret['SCC']['SCC_error']))
                        else:
                            print('detailed.out could not be found, exiting...')
                            return
                        n += 2
                        energies.append(ret['total_energy'])
                        break
                    time.sleep(10)

                self.results.append({'kp_grid': grid,
                                     'iSCC': ret['SCC']['iSCC'],
                                     'Total_energy': ret['total_energy'],
                                     'SCC_error': ret['SCC']['SCC_error']})
            else:
                n += 2

            if len(energies) > 2 and abs(max(energies[-3:]) - min(energies[-3:])) < self.energy_tolerance:
                break

 def basic_input(self):
     self.basic_driver()
     if not self.slater_koster:
         pcm_log.debug('The Slater-Koster files were not selected')
         pcm_log.debug('The Hamiltonian could not be setup')
     else:
         self.basic_hamiltonian()
     self.basic_options()
     self.basic_analysis()
     self.basic_parser_options()

    def close_distances(self):
        """
        Computes the closest distances for all the atoms

        :return: (tuple) Return a bond's dictionary and distance's list
        """
        if self._pairs is None or self._distances is None:

            if self.structure.is_periodic:
                pcm_log.debug('Computing distances from scratch...')
                pairs_dict = {}
                distances_list = []
                index = 0
                for i, j in itertools.combinations(range(self.structure.natom), 2):
                    if index % 100 == 0:
                        pcm_log.debug('Computing distance between atoms %d and %d' % (i, j))
                    ret = self.structure.lattice.distance2(self.structure.reduced[i], self.structure.reduced[j],
                                                           radius=self.radius)
                    for k in ret:
                        if str(i) not in pairs_dict:
                            pairs_dict[str(i)] = [index]
                        else:
                            pairs_dict[str(i)].append(index)
                        if str(j) not in pairs_dict:
                            pairs_dict[str(j)] = [index]
                        else:
                            pairs_dict[str(j)].append(index)
                        ret[k]['pair'] = (i, j)
                        distances_list.append(ret[k])
                        index += 1
                for i in range(self.structure.natom):
                    ret = self.structure.lattice.distance2(self.structure.reduced[i], self.structure.reduced[i])
                    for k in ret:
                        if str(i) not in pairs_dict:
                            pairs_dict[str(i)] = [index]
                        else:
                            pairs_dict[str(i)].append(index)
                        ret[k]['pair'] = (i, i)
                        distances_list.append(ret[k])
                        index += 1
                self._pairs = pairs_dict
                self._distances = distances_list
            else:
                dm = self.structure.distance_matrix()
                dm += np.eye(len(dm)) * max(dm.flatten())
                pairs_dict = {}
                distances_list = []
                for i in range(self.structure.natom):
                    index = dm[:, i].argmin()
                    pairs_dict[str(i)] = [index]
                    distances_list.append(dm[index, i])
                self._pairs = pairs_dict
                self._distances = distances_list

        return self._pairs, self._distances

    def test_swarm(self):
        """
        Tests (pychemia.searcher.swarm) with LJ Clusters             :
        """

        pcm_log.debug('ParticleSwarm')
        popu = LJCluster('test', composition='Xe13', refine=False, direct_evaluation=True)
        popu.pcdb.clean()
        searcher = ParticleSwarm(popu, generation_size=8, stabilization_limit=3)
        searcher.run()
        popu.pcdb.clean()

    def new_entry(self, structure, active=True):

        if active and self.direct_evaluation:
            structure, properties = self.evaluate(structure, gtol=self.target_forces)
        else:
            properties = {}
        status = {self.tag: active}
        entry = {'structure': structure.to_dict, 'properties': properties, 'status': status}
        entry_id = self.set_entry(entry)
        pcm_log.debug('Added new entry: %s with tag=%s: %s' % (str(entry_id), self.tag, str(active)))
        return entry_id

    def test_harmony(self):
        """
        Tests (pychemia.searcher.harmony) with LJ Clusters           :
        """

        pcm_log.debug('HarmonySearch')
        popu = LJCluster('test', composition='Xe13', refine=False, direct_evaluation=True)
        popu.pcdb.clean()
        searcher = HarmonySearch(popu, generation_size=8, stabilization_limit=3)
        searcher.run()
        popu.pcdb.clean()

    def test_genetic(self):
        """

        Tests (pychemia.searcher.genetic) with LJ Clusters           :
        """

        pcm_log.debug('GeneticAlgorithm')
        popu = LJCluster('test', composition='Xe13', refine=False, direct_evaluation=True)
        popu.pcdb.clean()
        searcher = GeneticAlgorithm(popu, generation_size=8, stabilization_limit=5)
        searcher.run()
        popu.pcdb.clean()

    def steepest_descent(self, dt=1, tolerance=1E-3):
        forces = self.get_forces()
        while np.max(self.get_magnitude_forces()) > tolerance:
            for i in range(self.structure.natom):
                imass = mass(self.structure.symbols[i])
                self.structure.positions[i] += 0.5 * forces[i] / imass * dt ** 2

            forces = self.get_forces()
            maxforce = np.max(self.get_magnitude_forces())
            dt = max(100.0 / maxforce, 1.0)
            pcm_log.debug('Steepest Descent maxforce= %9.3E   density=%7.4f dt=%7.3f' % (maxforce,
                                                                                         self.structure.density, dt))

    def move(self, entry_id, entry_jd, factor=0.2, in_place=False):

        st_orig = self.get_structure(entry_id)
        st_dest = self.get_structure(entry_jd)

        cm = ClusterMatch(st_orig, st_dest)
        cm.match()

        # pos_orig = np.array(entry_orig['structure']['positions']).reshape((-1, 3))
        # pos_dest = np.array(entry_dest['structure']['positions']).reshape((-1, 3))
        pos_orig = cm.structure1.positions
        pos_dest = cm.structure2.positions

        # Move to a position with negative energy
        reduc = 1
        new_positions = np.array(pos_orig)
        while True:
            new_positions = rotation_move(pos_orig, pos_dest, fraction=reduc * factor)
            new_structure = Structure(positions=new_positions, symbols=st_orig.symbols, periodicity=False)
            lj = LennardJones(new_structure)
            if lj.get_energy() < 0.0:
                break
            reduc -= 0.05
            pcm_log.debug('Effective factor will be reduced to %7.3f, original factor %7.3f' % (reduc * factor, factor))
            if reduc <= 0.0:
                # print 'No movement effective'
                break

        # Avoid condition with atoms too close
        distance_matrix = scipy.spatial.distance_matrix(new_positions, new_positions)
        tmp = np.max(distance_matrix.flatten())
        # print 'Scaling by', tmp
        minimal_distance = np.min((distance_matrix + tmp * np.eye(len(new_positions))).flatten())

        if minimal_distance < 1E-8:
            pcm_log.debug("Null distance between different atoms, no moving")
            new_positions = pos_orig

        if tmp > 5:
            # print 'Big scaling, better not to move'
            new_positions = pos_orig
        else:
            max_cov = np.max(covalent_radius(st_orig.symbols))
            new_positions *= max_cov / minimal_distance

        new_structure = Structure(positions=new_positions, symbols=st_orig.symbols, periodicity=False)
        # print 'Density of cluster', new_structure.density

        if in_place:
            self.unset_properties(entry_id)
            return self.set_structure(entry_id, new_structure)
        else:
            return self.new_entry(new_structure, active=False)

 def value(self, entry_id):
     entry = self.get_entry(entry_id)
     structure = self.get_structure(entry_id)
     if 'properties' not in entry:
         pcm_log.debug('This entry has no properties %s' % str(entry['_id']))
         return None
     elif entry['properties'] is None:
         return None
     elif 'energy' not in entry['properties']:
         pcm_log.debug('This entry has no energy in properties %s' % str(entry['_id']))
         return None
     else:
         return entry['properties']['energy'] / structure.get_composition().gcd