<commit_before># coding: utf-8
# Copyright(c) Henniggroup.
# Distributed under the terms of the MIT License.

from __future__ import division, unicode_literals, print_function

"""
General module:

This module contains several classes central to the algorithm.

1. IDGenerator: generates consecutive integer ID numbers for organisms

2. Organism: represents an organism, which consists primarily of a
        structure and an energy

3. Cell: represents a structure, and extends pymatgen.core.structure.Structure

4. InitialPopulation: represents the initial population of organisms,
        produced from one or more of the organism creators

5. Pool: represents the population of organisms, after the initial population

6. OffspringGenerator: high-level class for making offspring organisms

7. SelectionProbDist: specifies the distribution from which selection
        probabilities are drawn

8. CompositionSpace: specifies the composition or composition range

9. CompositionFitnessWeight: specifies how the weight given to the composition
        fitness is to be computed

10. StoppingCriteria: specifies when a search should stop

11. DataWriter: writes information about the search to a file

"""

from pymatgen.core.structure import Structure, Molecule
from pymatgen.core.lattice import Lattice
from pymatgen.core.composition import Composition
from pymatgen.core.periodic_table import Element, DummySpecie
from pymatgen.phasediagram.maker import CompoundPhaseDiagram
from pymatgen.phasediagram.entries import PDEntry
from pymatgen.phasediagram.analyzer import PDAnalyzer
from pymatgen.transformations.standard_transformations import \
    RotationTransformation

import os
import copy
import math
import numpy as np
from collections import deque
from scipy.spatial.qhull import ConvexHull


class IDGenerator(object):
    """
    Generates successive integer ID numbers, starting from 1.
    """

    def __init__(self):
        """
        Makes an IDGenerator.
        """

        self.id = 0

    def make_id(self):
        """
        Returns the next id number.
        """

        self.id += 1
        return self.id


class Organism(object):
    """
    An organism, consisting primarily of a cell and an energy, as well as
    several derived quantities.
    """

    def __init__(self, cell, id_generator, maker, composition_space):
        """
        Makes an Organism.

        Args:
            cell: the cell of this organism, as a Cell object

            id_generator: the IDGenerator used to assign id numbers to all
                organisms

            maker: the name of algorithm that made the organism, as a string.
                Either a creator or a variation

            composition_space: the CompositionSpace of the search
        """

        self.cell = cell
        self.composition = self.cell.composition
        # position in composition space. Only used for phase diagram searches
        self.composition_vector = self.compute_composition_vector(
            composition_space)
        self.total_energy = None
        # the objective function value
        self.value = None
        # energy per atom
        self.epa = None
        # the fitness and relative fitness
        self.fitness = None
        self.relative_fitness = None
        # selection probability of this organism
        self.selection_prob = None
        # starting point of interval on number line of size selection_prob
        self.selection_loc = None
        # relative selection probability
        self.relative_selection_prob = None
        # starting point of interval of size relative_selection_prob
        self.relative_selection_loc = None
        # whether the organism is part of initial population or pool
        self.is_active = False
        # unique id number
        self._id = id_generator.make_id()
        # the name of the algorithm that created this organism
        self.made_by = maker

    # This keeps the id (sort of) immutable by causing an exception to be
    # raised if the id is attempted to be set with org.id = some_id.
    @property
    def id(self):
        return self._id

    def compute_composition_vector(self, composition_space):
        """
        Returns the composition vector of the organism, as a numpy array.

        Args:
            composition_space: the CompositionSpace of the search.
        """

        if composition_space.objective_function == 'epa':
            return None
        elif composition_space.objective_function == 'pd':
            # make CompoundPhaseDiagram and PDAnalyzer objects
            pdentries = []
            for endpoint in composition_space.endpoints:
                pdentries.append(PDEntry(endpoint, -10))
            compound_pd = CompoundPhaseDiagram(pdentries,
                                               composition_space.endpoints)
            pdanalyzer = PDAnalyzer(compound_pd)

            # transform the organism's composition
            transformed_entry = compound_pd.transform_entries(
                [PDEntry(self.composition, 10)], composition_space.endpoints)

            # get the transformed species and amounts
            if len(transformed_entry[0]) == 0:
                return None
            transformed_list = str(transformed_entry[0][0]).split()
            del transformed_list[0]
            popped = ''
            while popped != 'with':
                popped = transformed_list.pop()

            # separate the dummy species symbols from the amounts
            symbols = []
            amounts = []
            for entry in transformed_list:
                split_entry = entry.split('0+')
                symbols.append(split_entry[0])
                amounts.append(float(split_entry[1]))

            # make a dictionary mapping dummy species to amounts
            dummy_species_amounts = {}
            for i in range(len(symbols)):
                dummy_species_amounts[DummySpecie(symbol=symbols[i])] = \
                    amounts[i]

            # make Composition object with dummy species, get decomposition
            dummy_comp = Composition(dummy_species_amounts)
            decomp = pdanalyzer.get_decomposition(dummy_comp)

            # get amounts of the decomposition in terms of the (untransformed)
            # composition space endpoints
            formatted_decomp = {}
            for key in decomp:
                key_dict = key.as_dict()
                comp = Composition(key_dict['entry']['composition'])
                formatted_decomp[comp] = decomp[key]

            # make the composition vector
            composition_vector = []
            # because the random organism creator shuffles the endpoints
            composition_space.endpoints.sort()
            for endpoint in composition_space.endpoints:
                if endpoint in formatted_decomp:
                    composition_vector.append(formatted_decomp[endpoint])
                else:
                    composition_vector.append(0.0)
            return np.array(composition_vector)

    def is_at_endpoint(self, composition_space):
        """
        Returns a boolean indicating whether the organism is located at an
        endpoint of the composition space.

        Args:
            composition_space: the CompositionSpace of the search
        """

        for endpoint in composition_space.endpoints:
            if self.composition.reduced_composition.almost_equals(endpoint):
                return True
        return False


class Cell(Structure):
    """
    Represents a cell. Provides additional functionality to the
    pymatgen.core.structure.Structure class.
    """

    def rotate_to_principal_directions(self):
        """
        Rotates the cell into the principal directions. That is, lattice vector
        a is parallel to the Cartesian x-axis, lattice vector b lies in the
        Cartesian x-y plane and the z-component of lattice vector c is
        positive.

        Note: this method doesn't change the fractional coordinates of the
        sites. However, the Cartesian coordinates may be changed.
        """

        # rotate about the z-axis to align a vertically with the x-axis
        rotation = RotationTransformation(
            [0, 0, 1], 180 - (180/np.pi)*np.arctan2(
                    self.lattice.matrix[0][1],
                    self.lattice.matrix[0][0]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # rotate about the y-axis to make a parallel to the x-axis
        rotation = RotationTransformation(
                [0, 1, 0], (180/np.pi)*np.arctan2(
                    self.lattice.matrix[0][2],
                    self.lattice.matrix[0][0]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # rotate about the x-axis to make b lie in the x-y plane
        rotation = RotationTransformation(
                [1, 0, 0], 180 - (180/np.pi)*np.arctan2(
                    self.lattice.matrix[1][2],
                    self.lattice.matrix[1][1]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # make sure they are all pointing in positive directions
        if self.lattice.matrix[0][0] < 0:
            # rotate about y-axis to make a positive
            rotation = RotationTransformation([0, 1, 0], 180)
            new_structure = rotation.apply_transformation(self)
            self.modify_lattice(new_structure.lattice)
        if self.lattice.matrix[1][1] < 0:
            # rotate about x-axis to make b positive
            rotation = RotationTransformation([1, 0, 0], 180)
            new_structure = rotation.apply_transformation(self)
            self.modify_lattice(new_structure.lattice)
        if self.lattice.matrix[2][2] < 0:
            # mirror c across the x-y plane to make it positive
            # a and b
            a = self.lattice.matrix[0]
            b = self.lattice.matrix[1]
            # the components of c
            cx = self.lattice.matrix[2][0]
            cy = self.lattice.matrix[2][1]
            cz = -1*self.lattice.matrix[2][2]
            self.modify_lattice(Lattice([a, b, [cx, cy, cz]]))

    def rotate_c_parallel_to_z(self):
        """
        Rotates the cell such that lattice vector c is parallel to the
        Cartesian z-axis.

        Note: this method doesn't change the fractional coordinates of the
        sites. However, the Cartesian coordinates may be changed.
        """

        # rotate about the z-axis until c lies in the x-z plane
        rotation = RotationTransformation(
                [0, 0, 1], 180 - (180/np.pi)*np.arctan2(
                    self.lattice.matrix[2][1],
                    self.lattice.matrix[2][0]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # rotate about the y-axis to make c parallel to the z-axis
        rotation = RotationTransformation(
            [0, 1, 0], 180 - (180/np.pi)*np.arctan2(
                self.lattice.matrix[2][0],
                self.lattice.matrix[2][2]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # make sure c is pointing along the positive z-axis
        if self.lattice.matrix[2][2] < 0:
            # rotate 180 degrees about the x-axis
            rotation = RotationTransformation([1, 0, 0], 180)
            new_structure = rotation.apply_transformation(self)
            self.modify_lattice(new_structure.lattice)

    def translate_atoms_into_cell(self):
        """
        Translates all the atoms into the cell, so that their fractional
        coordinates are between 0 and 1.

        Precondition: all the atoms can fit inside the structure's lattice,
            and there is at least 0.001 of (fractional) extra space along
            each lattice vector.
        """

        # get the bounding box of the atoms, in fractional coordinates
        bounding_box = self.get_bounding_box(cart_coords=False)
        translation_vector = []

        # determine the needed shift along each lattice vector
        for i in range(3):
            if bounding_box[i][0] < 0.0:
                translation_vector.append(-1*(bounding_box[i][0]) + 0.001)
            elif bounding_box[i][1] >= 1.0:
                translation_vector.append(-1*(bounding_box[i][1] + 0.001) +
                                          1.0)
            else:
                translation_vector.append(0.0)

        # do the shift
        site_indices = [i for i in range(len(self.sites))]
        self.translate_sites(site_indices, translation_vector,
                             frac_coords=True, to_unit_cell=False)

    def reduce_cell(self):
        """
        Applies Niggli cell reduction to cell.

        Returns a boolean indicated whether cell reduction was successful.

        Note: pymatgen's Niggli cell reduction algorithm sometimes moves the
            atoms' relative positions a little (I've seen up to 0.5 A...)
        """
        # get the reduced structure
        try:  # sometimes pymatgen's reduction routine fails
            reduced_structure = self.get_reduced_structure()
        except:
            return False

        # modify the cell to correspond to the reduced structure
        rcartesian_coords = reduced_structure.cart_coords
        rspecies = reduced_structure.species
        self.modify_lattice(reduced_structure.lattice)
        site_indices = []
        for i in range(len(self.sites)):
            site_indices.append(i)
        self.remove_sites(site_indices)
        for i in range(len(rcartesian_coords)):
            self.append(rspecies[i], rcartesian_coords[i],
                        coords_are_cartesian=True)

        # rotate the cell into the principal directions
        self.rotate_to_principal_directions()
        return True

    def reduce_sheet_cell(self, geometry, constraints):
        """
        Applies Niggli cell reduction to a sheet cell. The idea is to make
        lattice vector c vertical and add lots of vertical vacuum so that the
        standard reduction algorithm only changes the a and b lattice vectors.

        Returns a boolean indicating whether cell reduction was successful.

        Args:
            geometry: the Geometry of the search

            constraints: the Constraints of the search
        """
        # pad the cell with 100 Angstroms of vertical vacuum
        geometry.pad(self, padding=100)

        # reduce the padded cell
        successful_reduction = self.reduce_cell()

        # unpad the reduced cell
        geometry.unpad(self, constraints)

        return successful_reduction

    def get_bounding_box(self, cart_coords=True):
        """
        Returns the smallest and largest coordinates in each dimension of all
        the sites in a cell, as a list of three lists:

            [[min, max], [min, max], [min, max]]

        where the first inner list contains data for the Cartesian x-coordinate
        or lattice vector a, the second inner list contains data for the
        Cartesian y-coordinate or lattice vector b, and the third inner list
        contains data for the Cartesian z-coordinate or lattice vector c.

        Args:
            cart_coords: whether to give the result in Cartesian or fractional
                coordinates
        """

        if cart_coords:
            coords = self.cart_coords
        else:
            coords = self.frac_coords

        # find the largest and smallest coordinates in each dimension
        minx = np.inf
        maxx = -np.inf
        miny = np.inf
        maxy = -np.inf
        minz = np.inf
        maxz = -np.inf
        for coord in coords:
            if coord[0] < minx:
                minx = coord[0]
            if coord[0] > maxx:
                maxx = coord[0]
            if coord[1] < miny:
                miny = coord[1]
            if coord[1] > maxy:
                maxy = coord[1]
            if coord[2] < minz:
                minz = coord[2]
            if coord[2] > maxz:
                maxz = coord[2]
        return [[minx, maxx], [miny, maxy], [minz, maxz]]


class InitialPopulation():
    """
    The initial population of organisms.
    """

    def __init__(self, run_dir_name):
        """
        Makes an InitialPopulation.

        Args:
            run_dir_name: the name (not path) of the garun directory where the
                search will be done
        """

        self.run_dir_name = run_dir_name
        self.initial_population = []

    def add_organism(self, organism_to_add, composition_space):
        """
        Adds a relaxed organism to the initial population.

        Args:
            organism_to_add: the Organism to add to the initial population
        """

        organism_to_add.cell.sort()
        organism_to_add.cell.to('poscar', os.getcwd() + '/POSCAR.' +
                                str(organism_to_add.id))
        print('Adding organism {} to the initial population'.format(
            organism_to_add.id))
        self.initial_population.append(organism_to_add)
        organism_to_add.is_active = True

    def replace_organism(self, old_org, new_org, composition_space):
        """
        Replaces an organism in the initial population with a new organism.

        Precondition: old_org is a current member of the initial population

        Args:
            old_org: the Organism in the initial population to replace

            new_org: the new Organism to replace the old one
        """

        new_org.cell.sort()
        new_org.cell.to('poscar', os.getcwd() + '/POSCAR.' +
                        str(new_org.id))
        print('Replacing organism {} with organism {} in the initial '
              'population'.format(old_org.id, new_org.id))

        # remove the redundant organism and add the new one
        for org in self.initial_population:
            if org.id == old_org.id:
                self.initial_population.remove(org)
                org.is_active = False
        self.initial_population.append(new_org)
        new_org.is_active = True

    def get_progress(self, composition_space):
        """
        Returns either the best energy per atom (for fixed-composition
        search) or the area/volume of the convex hull (for phase diagram
        search).

        Args:
            composition_space: the CompositionSpace of the search
        """

        if composition_space.objective_function == 'epa':
            return self.get_best_epa()
        elif composition_space.objective_function == 'pd':
            return self.get_convex_hull_area(composition_space)

    def get_best_epa(self):
        """
        Returns the epa of the best organism in the initial population.
        """

        sorted_list = copy.deepcopy(self.initial_population)
        sorted_list.sort(key=lambda x: x.epa, reverse=False)
        return sorted_list[0].epa

    def get_convex_hull_area(self, composition_space):
        """
        Returns the area/volume of the current convex hull.

        Args:
            composition_space: the CompositionSpace of the search
        """

        # check if the initial population contains organisms at all the
        # endpoints of the composition space
        if self.has_endpoints(composition_space) and \
                self.has_non_endpoint(composition_space):

            # compute and print the area or volume of the convex hull
            pdentries = []
            for organism in self.initial_population:
                pdentries.append(PDEntry(organism.composition,
                                         organism.total_energy))
            compound_pd = CompoundPhaseDiagram(pdentries,
                                               composition_space.endpoints)

            # get the data for the convex hull
            qhull_data = compound_pd.qhull_data
            # for some reason, the last point is positive, so remove it
            hull_data = np.delete(qhull_data, -1, 0)

            # make a ConvexHull object from the hull data
            # Sometime this fails, saying that only two points were given to
            # construct the convex hull, even though the if statement above
            # checks that there are enough points.
            try:
                convex_hull = ConvexHull(hull_data)
            except:
                return None
            if len(composition_space.endpoints) == 2:
                return convex_hull.area
            else:
                return convex_hull.volume

    def has_endpoints(self, composition_space):
        """
        Checks if there are organisms in the initial population at all of
        endpoints of the composition space.

        Returns a boolean indicating whether or not there are organisms at all
        the endpoints.

        Args:
            composition_space: the CompositionSpace of the search
        """

        for endpoint in composition_space.endpoints:
            has_endpoint = False
            for organism in self.initial_population:
                if endpoint.almost_equals(
                        organism.composition.reduced_composition):
                    has_endpoint = True
            if not has_endpoint:
                return False
        return True

    def has_non_endpoint(self, composition_space):
        """
        Checks that the initial population contains at least one organism not
        at one of the endpoint compositions.

        Returns a boolean indicating whether or not there is an organism not
        at one of the endpoints.

        Args:
            composition_space: the CompositionSpace of the search
        """

        for organism in self.initial_population:
            not_endpoint = True
            for endpoint in composition_space.endpoints:
                if endpoint.almost_equals(
                        organism.composition.reduced_composition):
                    not_endpoint = False
            if not_endpoint:
                return True
        return False


class Pool(object):
    """
    Pool to hold all the organisms that are currently candidates for selection
    to become parents.

    Composed of two parts: a promotion set containing the best few organisms,
    and a queue containing the rest of the organisms in the pool.
    """
    def __init__(self, pool_params, composition_space, run_dir_name):
        """
        Makes a Pool, and sets default parameter values if necessary.

        Args:
            pool_params: the parameters of the pool, as a dictionary

            composition_space: the CompositionSpace of the search

            run_dir_name: the name (not path) of the garun directory
        """

        # the default values
        if len(composition_space.endpoints) == 1:
            # then we're doing a fixed composition search
            self.default_size = 20
        elif len(composition_space.endpoints) == 2:
            # then we're doing a binary phase diagram search search
            self.default_size = 25
        elif len(composition_space.endpoints) == 3:
            # then we're doing a ternary phase diagram search
            self.default_size = 75
        elif len(composition_space.endpoints) == 4:
            # then we're doing a quaternary phase diagram search
            self.default_size = 150

        # this only gets used for fixed composition searches
        self.default_num_promoted = 3

        # if entire Pool block was set to default or left blank
        if pool_params in (None, 'default'):
            self.size = self.default_size
            self.num_promoted = self.default_num_promoted
        else:
            # check if size parameter was used
            if 'size' not in pool_params:
                self.size = self.default_size
            elif pool_params['size'] in (None, 'default'):
                self.size = self.default_size
            else:
                self.size = pool_params['size']

            # check if num_promoted parameter was used
            if 'num_promoted' not in pool_params:
                self.num_promoted = self.default_num_promoted
            elif pool_params['num_promoted'] in (None, 'default'):
                self.num_promoted = self.default_num_promoted
            else:
                self.num_promoted = pool_params['num_promoted']

        # the best few organisms in the pool
        self.promotion_set = []
        # the rest of the organisms in the pool
        self.queue = deque()
        # the parameters for the selection distribution and composition fitness
        # weight
        # note: these are just placeholders - they get set in
        # objects_maker.make_objects()
        self.selection = []
        self.comp_fitness_weight = []
        # the number of organisms added to the pool (excluding the initial
        # population)
        self.num_adds = 0
        # the name (not  path) of the garun directory
        self.run_dir_name = run_dir_name

    def add_initial_population(self, initial_population, composition_space):
        """
        Adds the organisms of the initial population to the pool.

        Args:
            initial_population: the InitialPopulation containing the relaxed
                organisms of the initial population

            composition_space: the CompositionSpace of the search

        Precondition: the pool is empty (no organisms have been previously
            added to it)
        """

        print('Populating the pool with the initial population...')
        organisms_list = initial_population.initial_population

        # populate the promotion set and queue
        if composition_space.objective_function == 'epa':
            for organism in organisms_list:
                organism.value = organism.epa
            organisms_list.sort(key=lambda x: x.value, reverse=False)

            # assign them to either the promotion set or the queue
            for i in range(len(organisms_list)):<commit_msg>Added check to ensure the initial population contains at least three organisms.<commit_after># coding: utf-8
# Copyright(c) Henniggroup.
# Distributed under the terms of the MIT License.

from __future__ import division, unicode_literals, print_function

"""
General module:

This module contains several classes central to the algorithm.

1. IDGenerator: generates consecutive integer ID numbers for organisms

2. Organism: represents an organism, which consists primarily of a
        structure and an energy

3. Cell: represents a structure, and extends pymatgen.core.structure.Structure

4. InitialPopulation: represents the initial population of organisms,
        produced from one or more of the organism creators

5. Pool: represents the population of organisms, after the initial population

6. OffspringGenerator: high-level class for making offspring organisms

7. SelectionProbDist: specifies the distribution from which selection
        probabilities are drawn

8. CompositionSpace: specifies the composition or composition range

9. CompositionFitnessWeight: specifies how the weight given to the composition
        fitness is to be computed

10. StoppingCriteria: specifies when a search should stop

11. DataWriter: writes information about the search to a file

"""

from pymatgen.core.structure import Structure, Molecule
from pymatgen.core.lattice import Lattice
from pymatgen.core.composition import Composition
from pymatgen.core.periodic_table import Element, DummySpecie
from pymatgen.phasediagram.maker import CompoundPhaseDiagram
from pymatgen.phasediagram.entries import PDEntry
from pymatgen.phasediagram.analyzer import PDAnalyzer
from pymatgen.transformations.standard_transformations import \
    RotationTransformation

import os
import copy
import math
import numpy as np
from collections import deque
from scipy.spatial.qhull import ConvexHull


class IDGenerator(object):
    """
    Generates successive integer ID numbers, starting from 1.
    """

    def __init__(self):
        """
        Makes an IDGenerator.
        """

        self.id = 0

    def make_id(self):
        """
        Returns the next id number.
        """

        self.id += 1
        return self.id


class Organism(object):
    """
    An organism, consisting primarily of a cell and an energy, as well as
    several derived quantities.
    """

    def __init__(self, cell, id_generator, maker, composition_space):
        """
        Makes an Organism.

        Args:
            cell: the cell of this organism, as a Cell object

            id_generator: the IDGenerator used to assign id numbers to all
                organisms

            maker: the name of algorithm that made the organism, as a string.
                Either a creator or a variation

            composition_space: the CompositionSpace of the search
        """

        self.cell = cell
        self.composition = self.cell.composition
        # position in composition space. Only used for phase diagram searches
        self.composition_vector = self.compute_composition_vector(
            composition_space)
        self.total_energy = None
        # the objective function value
        self.value = None
        # energy per atom
        self.epa = None
        # the fitness and relative fitness
        self.fitness = None
        self.relative_fitness = None
        # selection probability of this organism
        self.selection_prob = None
        # starting point of interval on number line of size selection_prob
        self.selection_loc = None
        # relative selection probability
        self.relative_selection_prob = None
        # starting point of interval of size relative_selection_prob
        self.relative_selection_loc = None
        # whether the organism is part of initial population or pool
        self.is_active = False
        # unique id number
        self._id = id_generator.make_id()
        # the name of the algorithm that created this organism
        self.made_by = maker

    # This keeps the id (sort of) immutable by causing an exception to be
    # raised if the id is attempted to be set with org.id = some_id.
    @property
    def id(self):
        return self._id

    def compute_composition_vector(self, composition_space):
        """
        Returns the composition vector of the organism, as a numpy array.

        Args:
            composition_space: the CompositionSpace of the search.
        """

        if composition_space.objective_function == 'epa':
            return None
        elif composition_space.objective_function == 'pd':
            # make CompoundPhaseDiagram and PDAnalyzer objects
            pdentries = []
            for endpoint in composition_space.endpoints:
                pdentries.append(PDEntry(endpoint, -10))
            compound_pd = CompoundPhaseDiagram(pdentries,
                                               composition_space.endpoints)
            pdanalyzer = PDAnalyzer(compound_pd)

            # transform the organism's composition
            transformed_entry = compound_pd.transform_entries(
                [PDEntry(self.composition, 10)], composition_space.endpoints)

            # get the transformed species and amounts
            if len(transformed_entry[0]) == 0:
                return None
            transformed_list = str(transformed_entry[0][0]).split()
            del transformed_list[0]
            popped = ''
            while popped != 'with':
                popped = transformed_list.pop()

            # separate the dummy species symbols from the amounts
            symbols = []
            amounts = []
            for entry in transformed_list:
                split_entry = entry.split('0+')
                symbols.append(split_entry[0])
                amounts.append(float(split_entry[1]))

            # make a dictionary mapping dummy species to amounts
            dummy_species_amounts = {}
            for i in range(len(symbols)):
                dummy_species_amounts[DummySpecie(symbol=symbols[i])] = \
                    amounts[i]

            # make Composition object with dummy species, get decomposition
            dummy_comp = Composition(dummy_species_amounts)
            decomp = pdanalyzer.get_decomposition(dummy_comp)

            # get amounts of the decomposition in terms of the (untransformed)
            # composition space endpoints
            formatted_decomp = {}
            for key in decomp:
                key_dict = key.as_dict()
                comp = Composition(key_dict['entry']['composition'])
                formatted_decomp[comp] = decomp[key]

            # make the composition vector
            composition_vector = []
            # because the random organism creator shuffles the endpoints
            composition_space.endpoints.sort()
            for endpoint in composition_space.endpoints:
                if endpoint in formatted_decomp:
                    composition_vector.append(formatted_decomp[endpoint])
                else:
                    composition_vector.append(0.0)
            return np.array(composition_vector)

    def is_at_endpoint(self, composition_space):
        """
        Returns a boolean indicating whether the organism is located at an
        endpoint of the composition space.

        Args:
            composition_space: the CompositionSpace of the search
        """

        for endpoint in composition_space.endpoints:
            if self.composition.reduced_composition.almost_equals(endpoint):
                return True
        return False


class Cell(Structure):
    """
    Represents a cell. Provides additional functionality to the
    pymatgen.core.structure.Structure class.
    """

    def rotate_to_principal_directions(self):
        """
        Rotates the cell into the principal directions. That is, lattice vector
        a is parallel to the Cartesian x-axis, lattice vector b lies in the
        Cartesian x-y plane and the z-component of lattice vector c is
        positive.

        Note: this method doesn't change the fractional coordinates of the
        sites. However, the Cartesian coordinates may be changed.
        """

        # rotate about the z-axis to align a vertically with the x-axis
        rotation = RotationTransformation(
            [0, 0, 1], 180 - (180/np.pi)*np.arctan2(
                    self.lattice.matrix[0][1],
                    self.lattice.matrix[0][0]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # rotate about the y-axis to make a parallel to the x-axis
        rotation = RotationTransformation(
                [0, 1, 0], (180/np.pi)*np.arctan2(
                    self.lattice.matrix[0][2],
                    self.lattice.matrix[0][0]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # rotate about the x-axis to make b lie in the x-y plane
        rotation = RotationTransformation(
                [1, 0, 0], 180 - (180/np.pi)*np.arctan2(
                    self.lattice.matrix[1][2],
                    self.lattice.matrix[1][1]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # make sure they are all pointing in positive directions
        if self.lattice.matrix[0][0] < 0:
            # rotate about y-axis to make a positive
            rotation = RotationTransformation([0, 1, 0], 180)
            new_structure = rotation.apply_transformation(self)
            self.modify_lattice(new_structure.lattice)
        if self.lattice.matrix[1][1] < 0:
            # rotate about x-axis to make b positive
            rotation = RotationTransformation([1, 0, 0], 180)
            new_structure = rotation.apply_transformation(self)
            self.modify_lattice(new_structure.lattice)
        if self.lattice.matrix[2][2] < 0:
            # mirror c across the x-y plane to make it positive
            # a and b
            a = self.lattice.matrix[0]
            b = self.lattice.matrix[1]
            # the components of c
            cx = self.lattice.matrix[2][0]
            cy = self.lattice.matrix[2][1]
            cz = -1*self.lattice.matrix[2][2]
            self.modify_lattice(Lattice([a, b, [cx, cy, cz]]))

    def rotate_c_parallel_to_z(self):
        """
        Rotates the cell such that lattice vector c is parallel to the
        Cartesian z-axis.

        Note: this method doesn't change the fractional coordinates of the
        sites. However, the Cartesian coordinates may be changed.
        """

        # rotate about the z-axis until c lies in the x-z plane
        rotation = RotationTransformation(
                [0, 0, 1], 180 - (180/np.pi)*np.arctan2(
                    self.lattice.matrix[2][1],
                    self.lattice.matrix[2][0]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # rotate about the y-axis to make c parallel to the z-axis
        rotation = RotationTransformation(
            [0, 1, 0], 180 - (180/np.pi)*np.arctan2(
                self.lattice.matrix[2][0],
                self.lattice.matrix[2][2]))
        new_structure = rotation.apply_transformation(self)
        self.modify_lattice(new_structure.lattice)

        # make sure c is pointing along the positive z-axis
        if self.lattice.matrix[2][2] < 0:
            # rotate 180 degrees about the x-axis
            rotation = RotationTransformation([1, 0, 0], 180)
            new_structure = rotation.apply_transformation(self)
            self.modify_lattice(new_structure.lattice)

    def translate_atoms_into_cell(self):
        """
        Translates all the atoms into the cell, so that their fractional
        coordinates are between 0 and 1.

        Precondition: all the atoms can fit inside the structure's lattice,
            and there is at least 0.001 of (fractional) extra space along
            each lattice vector.
        """

        # get the bounding box of the atoms, in fractional coordinates
        bounding_box = self.get_bounding_box(cart_coords=False)
        translation_vector = []

        # determine the needed shift along each lattice vector
        for i in range(3):
            if bounding_box[i][0] < 0.0:
                translation_vector.append(-1*(bounding_box[i][0]) + 0.001)
            elif bounding_box[i][1] >= 1.0:
                translation_vector.append(-1*(bounding_box[i][1] + 0.001) +
                                          1.0)
            else:
                translation_vector.append(0.0)

        # do the shift
        site_indices = [i for i in range(len(self.sites))]
        self.translate_sites(site_indices, translation_vector,
                             frac_coords=True, to_unit_cell=False)

    def reduce_cell(self):
        """
        Applies Niggli cell reduction to cell.

        Returns a boolean indicated whether cell reduction was successful.

        Note: pymatgen's Niggli cell reduction algorithm sometimes moves the
            atoms' relative positions a little (I've seen up to 0.5 A...)
        """
        # get the reduced structure
        try:  # sometimes pymatgen's reduction routine fails
            reduced_structure = self.get_reduced_structure()
        except:
            return False

        # modify the cell to correspond to the reduced structure
        rcartesian_coords = reduced_structure.cart_coords
        rspecies = reduced_structure.species
        self.modify_lattice(reduced_structure.lattice)
        site_indices = []
        for i in range(len(self.sites)):
            site_indices.append(i)
        self.remove_sites(site_indices)
        for i in range(len(rcartesian_coords)):
            self.append(rspecies[i], rcartesian_coords[i],
                        coords_are_cartesian=True)

        # rotate the cell into the principal directions
        self.rotate_to_principal_directions()
        return True

    def reduce_sheet_cell(self, geometry, constraints):
        """
        Applies Niggli cell reduction to a sheet cell. The idea is to make
        lattice vector c vertical and add lots of vertical vacuum so that the
        standard reduction algorithm only changes the a and b lattice vectors.

        Returns a boolean indicating whether cell reduction was successful.

        Args:
            geometry: the Geometry of the search

            constraints: the Constraints of the search
        """
        # pad the cell with 100 Angstroms of vertical vacuum
        geometry.pad(self, padding=100)

        # reduce the padded cell
        successful_reduction = self.reduce_cell()

        # unpad the reduced cell
        geometry.unpad(self, constraints)

        return successful_reduction

    def get_bounding_box(self, cart_coords=True):
        """
        Returns the smallest and largest coordinates in each dimension of all
        the sites in a cell, as a list of three lists:

            [[min, max], [min, max], [min, max]]

        where the first inner list contains data for the Cartesian x-coordinate
        or lattice vector a, the second inner list contains data for the
        Cartesian y-coordinate or lattice vector b, and the third inner list
        contains data for the Cartesian z-coordinate or lattice vector c.

        Args:
            cart_coords: whether to give the result in Cartesian or fractional
                coordinates
        """

        if cart_coords:
            coords = self.cart_coords
        else:
            coords = self.frac_coords

        # find the largest and smallest coordinates in each dimension
        minx = np.inf
        maxx = -np.inf
        miny = np.inf
        maxy = -np.inf
        minz = np.inf
        maxz = -np.inf
        for coord in coords:
            if coord[0] < minx:
                minx = coord[0]
            if coord[0] > maxx:
                maxx = coord[0]
            if coord[1] < miny:
                miny = coord[1]
            if coord[1] > maxy:
                maxy = coord[1]
            if coord[2] < minz:
                minz = coord[2]
            if coord[2] > maxz:
                maxz = coord[2]
        return [[minx, maxx], [miny, maxy], [minz, maxz]]


class InitialPopulation():
    """
    The initial population of organisms.
    """

    def __init__(self, run_dir_name):
        """
        Makes an InitialPopulation.

        Args:
            run_dir_name: the name (not path) of the garun directory where the
                search will be done
        """

        self.run_dir_name = run_dir_name
        self.initial_population = []

    def add_organism(self, organism_to_add, composition_space):
        """
        Adds a relaxed organism to the initial population.

        Args:
            organism_to_add: the Organism to add to the initial population
        """

        organism_to_add.cell.sort()
        organism_to_add.cell.to('poscar', os.getcwd() + '/POSCAR.' +
                                str(organism_to_add.id))
        print('Adding organism {} to the initial population'.format(
            organism_to_add.id))
        self.initial_population.append(organism_to_add)
        organism_to_add.is_active = True

    def replace_organism(self, old_org, new_org, composition_space):
        """
        Replaces an organism in the initial population with a new organism.

        Precondition: old_org is a current member of the initial population

        Args:
            old_org: the Organism in the initial population to replace

            new_org: the new Organism to replace the old one
        """

        new_org.cell.sort()
        new_org.cell.to('poscar', os.getcwd() + '/POSCAR.' +
                        str(new_org.id))
        print('Replacing organism {} with organism {} in the initial '
              'population'.format(old_org.id, new_org.id))

        # remove the redundant organism and add the new one
        for org in self.initial_population:
            if org.id == old_org.id:
                self.initial_population.remove(org)
                org.is_active = False
        self.initial_population.append(new_org)
        new_org.is_active = True

    def get_progress(self, composition_space):
        """
        Returns either the best energy per atom (for fixed-composition
        search) or the area/volume of the convex hull (for phase diagram
        search).

        Args:
            composition_space: the CompositionSpace of the search
        """

        if composition_space.objective_function == 'epa':
            return self.get_best_epa()
        elif composition_space.objective_function == 'pd':
            return self.get_convex_hull_area(composition_space)

    def get_best_epa(self):
        """
        Returns the epa of the best organism in the initial population.
        """

        sorted_list = copy.deepcopy(self.initial_population)
        sorted_list.sort(key=lambda x: x.epa, reverse=False)
        return sorted_list[0].epa

    def get_convex_hull_area(self, composition_space):
        """
        Returns the area/volume of the current convex hull.

        Args:
            composition_space: the CompositionSpace of the search
        """

        # check if the initial population contains organisms at all the
        # endpoints of the composition space
        if self.has_endpoints(composition_space) and \
                self.has_non_endpoint(composition_space):

            # compute and print the area or volume of the convex hull
            pdentries = []
            for organism in self.initial_population:
                pdentries.append(PDEntry(organism.composition,
                                         organism.total_energy))
            compound_pd = CompoundPhaseDiagram(pdentries,
                                               composition_space.endpoints)

            # get the data for the convex hull
            qhull_data = compound_pd.qhull_data
            # for some reason, the last point is positive, so remove it
            hull_data = np.delete(qhull_data, -1, 0)

            # make a ConvexHull object from the hull data
            # Sometime this fails, saying that only two points were given to
            # construct the convex hull, even though the if statement above
            # checks that there are enough points.
            try:
                convex_hull = ConvexHull(hull_data)
            except:
                return None
            if len(composition_space.endpoints) == 2:
                return convex_hull.area
            else:
                return convex_hull.volume

    def has_endpoints(self, composition_space):
        """
        Checks if there are organisms in the initial population at all of
        endpoints of the composition space.

        Returns a boolean indicating whether or not there are organisms at all
        the endpoints.

        Args:
            composition_space: the CompositionSpace of the search
        """

        for endpoint in composition_space.endpoints:
            has_endpoint = False
            for organism in self.initial_population:
                if endpoint.almost_equals(
                        organism.composition.reduced_composition):
                    has_endpoint = True
            if not has_endpoint:
                return False
        return True

    def has_non_endpoint(self, composition_space):
        """
        Checks that the initial population contains at least one organism not
        at one of the endpoint compositions.

        Returns a boolean indicating whether or not there is an organism not
        at one of the endpoints.

        Args:
            composition_space: the CompositionSpace of the search
        """

        for organism in self.initial_population:
            not_endpoint = True
            for endpoint in composition_space.endpoints:
                if endpoint.almost_equals(
                        organism.composition.reduced_composition):
                    not_endpoint = False
            if not_endpoint:
                return True
        return False


class Pool(object):
    """
    Pool to hold all the organisms that are currently candidates for selection
    to become parents.

    Composed of two parts: a promotion set containing the best few organisms,
    and a queue containing the rest of the organisms in the pool.
    """
    def __init__(self, pool_params, composition_space, run_dir_name):
        """
        Makes a Pool, and sets default parameter values if necessary.

        Args:
            pool_params: the parameters of the pool, as a dictionary

            composition_space: the CompositionSpace of the search

            run_dir_name: the name (not path) of the garun directory
        """

        # the default values
        if len(composition_space.endpoints) == 1:
            # then we're doing a fixed composition search
            self.default_size = 20
        elif len(composition_space.endpoints) == 2:
            # then we're doing a binary phase diagram search search
            self.default_size = 25
        elif len(composition_space.endpoints) == 3:
            # then we're doing a ternary phase diagram search
            self.default_size = 75
        elif len(composition_space.endpoints) == 4:
            # then we're doing a quaternary phase diagram search
            self.default_size = 150

        # this only gets used for fixed composition searches
        self.default_num_promoted = 3

        # if entire Pool block was set to default or left blank
        if pool_params in (None, 'default'):
            self.size = self.default_size
            self.num_promoted = self.default_num_promoted
        else:
            # check if size parameter was used
            if 'size' not in pool_params:
                self.size = self.default_size
            elif pool_params['size'] in (None, 'default'):
                self.size = self.default_size
            else:
                self.size = pool_params['size']

            # check if num_promoted parameter was used
            if 'num_promoted' not in pool_params:
                self.num_promoted = self.default_num_promoted
            elif pool_params['num_promoted'] in (None, 'default'):
                self.num_promoted = self.default_num_promoted
            else:
                self.num_promoted = pool_params['num_promoted']

        # the best few organisms in the pool
        self.promotion_set = []
        # the rest of the organisms in the pool
        self.queue = deque()
        # the parameters for the selection distribution and composition fitness
        # weight
        # note: these are just placeholders - they get set in
        # objects_maker.make_objects()
        self.selection = []
        self.comp_fitness_weight = []
        # the number of organisms added to the pool (excluding the initial
        # population)
        self.num_adds = 0
        # the name (not  path) of the garun directory
        self.run_dir_name = run_dir_name

    def add_initial_population(self, initial_population, composition_space):
        """
        Adds the organisms of the initial population to the pool.

        Args:
            initial_population: the InitialPopulation containing the relaxed
                organisms of the initial population

            composition_space: the CompositionSpace of the search

        Precondition: the pool is empty (no organisms have been previously
            added to it)
        """

        print('Populating the pool with the initial population...')
        organisms_list = initial_population.initial_population

        # check that the initial population contains at least three organisms
        if len(organisms_list) < 3:
            print('The initial population must contain at least three '
                  'organisms.')
            print('Quitting...')
            quit()

        # populate the promotion set and queue
        if composition_space.objective_function == 'epa':
            for organism in organisms_list:
                organism.value = organism.epa
            organisms_list.sort(key=lambda x: x.value, reverse=False)

            # assign them to either the promotion set or the queue
            for i in range(len(organisms_list)):
{'fraction_non_alphanumeric': 0.06610816663817176, 'fraction_numerical': 0.00792158203681541, 'mean_word_length': 2.0367464666858957, 'pattern_counts': {'":': 0, '<': 20, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 16, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '3986583', 'n_tokens_mistral': 14043, 'n_tokens_neox': 13071, 'n_words': 4486}