non_dominated_sort.jl

{"nbformat_minor": 0, "cells": [{"execution_count": 1, "cell_type": "code", "source": "# ------------------------------------------------------------------------------\n# A generic type of a evaluated member of a population\nimmutable Member\n    objectives::Vector{Float64}\n    variables::Vector{Float64}\nend\n\ntypealias Population Vector{Member}\n# ------------------------------------------------------------------------------", "outputs": [{"execution_count": 1, "output_type": "execute_result", "data": {"text/plain": "Array{Member,1}"}, "metadata": {}}], "metadata": {"collapsed": false, "trusted": true}}, {"execution_count": 2, "cell_type": "code", "source": "# ==============================================================================\n# Sorting Algorithms\n#\n# Overloading the sorting algorithms to work with the population type.\n# ------------------------------------------------------------------------------\nfunction Base.sortperm(p::Population, obj_types::Vector{Float64})\n    pop_size = length(p)\n    dim_size = length(p[1].objectives)\n    pop_arrays = SubArray[sub(p[i].objectives .* obj_types, 1:dim_size)\n                          for i in 1:pop_size]\n    return sortperm(pop_arrays; order=Base.Sort.Lexicographic)\nend\n\n# ==============================================================================\n# Non dominated sorting\n#\n# This is an implementation of the algorithm purposed by Zhang et al. (2014)\n# \"An Efficient Approach to Nondominated Sorting for Evolutionary Multiobjective\n# Optimization\". Evolutionary Computation, IEEE Transactions on , vol.19, no.2\n# ------------------------------------------------------------------------------\nfunction is_dominated(member::Vector{Float64}, frontier::Vector{Float64},\n                      start::Int, n_vars::Int)\n    # This function checks if the current member is dominated by the current\n    # member of the current pareto frontier. Note the frontier is held in\n    # column (member) major order in the format of a vector.\n\n    for i in 1:n_vars\n        if isless(member[i], frontier[(i + start) - 1])\n            return false\n        end\n    end\n\n    return true\nend\n\n# ------------------------------------------------------------------------------\nfunction on_frontier(member::Vector{Float64}, frontier::Vector{Float64},\n                     n_vars::Int)\n    # This funciton checks if the current_member is dominated by any members of\n    # current pareto frontier. If no member dominates the current member then\n    # it is considered to be on the current frontier.\n\n    for i in 1:n_vars:length(frontier)\n        if is_dominated(member, frontier, i, n_vars)\n            return false\n        end\n    end\n\n    return true\nend\n\n# ------------------------------------------------------------------------------\nfunction get_frontier_index(member::Vector{Float64}, frontiers::Vector{Vector{Float64}},\n                            n_vars::Int)\n    # This function gets the index value of the pareto frontier of the current\n    # member. To do this, it checks if the current member should be placed on\n    # any of the initialized frontiers, otherwise it returns\n    # 1 + number of initialized frontiers.\n    #\n    # Args:\n    #    member: These are the values for the current member being indexed.\n    #    frontiers: A vector of vectors that holds the values of each member\n    #         on the existing frontiers.\n    #    n_vars: The number of variables for each member (this equals\n    #        length(member.field)).\n    #\n    # Returns\n    #    The index value of the pareto frontier the current member belongs to.\n\n    for i in 1:length(frontiers)\n        if on_frontier(member, frontiers[i], n_vars)\n            return i\n        end\n    end\n\n    return length(frontiers) + 1\nend\n\n# ------------------------------------------------------------------------------\nfunction non_dominated_sort(p::Population, obj_types::Vector{Float64}, n_vars::Int)\n    # This algorithm determines the pareto frontier for each member in a population.\n    #\n    # Args:\n    #    p: The population.\n    #    field: The field to use to determine the sorting (i.e. objectives or\n    #        constraints).\n    #    obj_types: An identifier for each value to determine if the comparison\n    #        should be a maximum or a minimum (1.0 = minimum, -1.0 = maximum).\n    #    n_vars: The number of variables (length) of the Member.field vector.\n    #\n    # Returns:\n    #    A vector of indices for each Member of a population that corresponds\n    #    to frontier which it belongs to.\n\n    # Get the sort permutation the population and retain it for later reordering\n    sort_index = sortperm(p, obj_types)\n    sorted_p = p[sort_index]\n\n    # Initialize the frontier index vector\n    frontier_index = zeros(Int, length(p))\n    frontier_index[1]  = 1\n\n    # Initialize the frontiers archive, with the first frontier and the first member\n    # of the sorted population (this is always on the first frontier).\n    frontiers = Array(Vector{Float64}, 1)\n    member_values = sorted_p[1].objectives\n    frontiers[1] =  Float64[member_values[i] * obj_types[i] for i in 1:n_vars]\n\n    for i in 2:length(p)\n        member_values = sorted_p[i].objectives\n        member_values =  Float64[member_values[i] * obj_types[i] for i in 1:n_vars]\n        frontier_index[i] = get_frontier_index(member_values, frontiers, n_vars)\n\n        if frontier_index[i] > length(frontiers)\n            # If the current member is on a new frontier, create a new frontier.\n            push!(frontiers, [j for j in member_values])\n        else\n            # If the current member is on an existing frontier prepend its values\n            # to the frontier. This is done for efficiency as comparisons between\n            # a member and pareto frontier are done from newest added to oldest.\n            # Also not that holding the frontier values in vector is done for\n            # performance reason (as it is faster than holding a vector of\n            # vectors).\n            prepend!(frontiers[frontier_index[i]], member_values)\n        end\n    end\n\n    return frontier_index[sortperm(sort_index)]\nend\n\n# ------------------------------------------------------------------------------\nfunction non_dominated_sort(p::Population, obj_types::Vector{Float64})\n    return non_dominated_sort(p, obj_types, length(p[1].objectives))\nend", "outputs": [{"execution_count": 2, "output_type": "execute_result", "data": {"text/plain": "non_dominated_sort (generic function with 2 methods)"}, "metadata": {}}], "metadata": {"collapsed": false, "trusted": true}}, {"execution_count": 3, "cell_type": "code", "source": "# ------------------------------------------------------------------------------\nfunction Base.searchsortedfirst(v::Population, m::Member)\n    pop = Vector{Float64}[i.objectives for i in v]\n    return searchsortedfirst(pop, m.objectives; order= Base.Sort.Lexicographic) - 1\nend\n\n# ------------------------------------------------------------------------------\nfunction update_frontier!(v::Population, m::Member, n_vars::Int)\n    split = searchsortedfirst(v, m)\n    \n    for i in split:-1:1\n        if is_dominated(m.objectives, v[i].objectives, 1, n_vars)\n            return nothing\n        end\n    end\n    \n    \n    for i in (split+1):length(v)\n        if !is_dominated(v[i].objectives, m.objectives, 1, n_vars)\n            splice!(v, (split+1):(i - 1), [m])\n            return nothing\n        end\n    end\n    \n    splice!(v, (split+1):length(v), [m])\n    return nothing\nend\n\n# ------------------------------------------------------------------------------\nfunction update_frontier!(v::Population, m::Member)\n    update_frontier!(v, m, length(m.objectives))\nend", "outputs": [{"execution_count": 3, "output_type": "execute_result", "data": {"text/plain": "update_frontier! (generic function with 2 methods)"}, "metadata": {}}], "metadata": {"collapsed": false, "trusted": true}}, {"execution_count": 4, "cell_type": "code", "source": "function add_points(pop::Population,m::Int)\n     for i in 1:m\n    a = rand() * 20  - 10\n    x = a^2\n    y = (a-2)^2 * -1.0\n    \n    mem = Member([x,y], rand(3))\n    update_frontier!(pop, mem)\n    end\nend", "outputs": [{"execution_count": 4, "output_type": "execute_result", "data": {"text/plain": "add_points (generic function with 1 method)"}, "metadata": {}}], "metadata": {"collapsed": false, "trusted": true}}, {"execution_count": 5, "cell_type": "code", "source": "# Warm up\nn = 10\nto_add = 10\n\npop = Array(Member, n)\n\nfor i in 1:n\n    a = rand() * 10  - 10\n    x = a^2\n    y = (a-2)^2 * -1.0\n    \n    pop[i] = Member([x,y], rand(3))\nend\n\npop = pop[sortperm(pop, [1.0, 1.0])]\n\n@time add_points(pop, to_add)", "outputs": [{"output_type": "stream", "name": "stdout", "text": "elapsed time: 0.121690709 seconds (3210788 bytes allocated)\n"}], "metadata": {"collapsed": false, "trusted": true}}, {"execution_count": 7, "cell_type": "code", "source": "n = 6000\npoints_to_add = 16000\n\npop = Array(Member, n)\n\nfor i in 1:n\n    a = rand() * 10  - 10\n    x = a^2\n    y = (a-2)^2 * -1.0\n    \n    pop[i] = Member([x,y], rand(3))\nend\n\npop = pop[sortperm(pop, [1.0, 1.0])]\n\n@time add_points(pop, points_to_add)", "outputs": [{"output_type": "stream", "name": "stdout", "text": "elapsed time: 2.73976081 seconds (1288508000 bytes allocated, 25.96% gc time)\n"}], "metadata": {"collapsed": false, "trusted": true}}], "nbformat": 4, "metadata": {"kernelspec": {"display_name": "Python 2", "name": "python2", "language": "python"}, "language_info": {"version": "0.3.6", "name": "julia"}}}