Take the 2-minute tour ×
Code Review Stack Exchange is a question and answer site for peer programmer code reviews. It's 100% free, no registration required.

This is my first go at genetic algorithms, using PyGame for the visual representation. Is there a way to improve the speed of the code? Click to place a point, press to start looking for the shortest route, or press c to clear all.

import pygame, time, random

#Size of next generation, change number at the end if needed
resize_rate = 1.005
#Controls crossover appearence
cross_rate = 0.7
#Controls number of chromosomes created at the start, exponential
exponent = 4

black = (0, 0, 0)
white = (255, 255, 255)
red = (255, 0, 0)

genes = []

def draw_background(screen):
    screen.fill(white)
def draw_connection(screen, pointlist):
    pygame.draw.lines(screen, black, True, pointlist, 1)
def draw_best_connection(screen, pointlist):
    pygame.draw.lines(screen, red, True, pointlist, 2)

def create_random_chromosomes(genes):
    temp_chromosomes = []
    chromosomes = []
    for i in range(len(genes)**exponent):
        genes_copy = genes[:]
        random.shuffle(genes_copy)
        temp_chromosomes.append(genes_copy)
    for element in temp_chromosomes:
        if element not in chromosomes:
            chromosomes.append(element)
    print(str(len(chromosomes)) + " random chromosomes created")
    return chromosomes

def calculate_distance(chromosome):
    distances = []
    total_distance = 0
    for i in range(len(chromosome)):
        if i < (len(chromosome)-1):
            distance_A_B = (((chromosome[i+1][0]-chromosome[i][0])**2)+(chromosome[i+1][1]-chromosome[i][1])**2)**0.5
            distances.append(distance_A_B)
        else:
            distance_A_B = (((chromosome[i][0]-chromosome[0][0])**2)+(chromosome[i][1]-chromosome[0][1])**2)**0.5
            distances.append(distance_A_B)
    for i in range(len(distances)):
        total_distance += distances[i]
    return total_distance

def calculate_fitness(chromosomes):
    fitness_scores = []
    for i in range(len(chromosomes)):
        fitness = 0
        fitness = 1/calculate_distance(chromosomes[i])
        fitness_scores.append(fitness)
    return fitness_scores

def build_chromosome_fitness_dict(chromosomes, fitness_scores):
    chrom_fit_dict = {}
    for i in range(len(chromosomes)):
        chrom_fit_dict[fitness_scores[i]] = chromosomes[i]
    return chrom_fit_dict

def roulette_selection(chromosomes, fitnesses, new_size):
    sum_fitness = sum(fitnesses)
    rel_fitness = [fitness/sum_fitness for fitness in fitnesses]
    #Generate probability intervals for each chromosome
    probabilities = [sum(rel_fitness[:i+1]) for i in range(len(rel_fitness))]
    #Select chromosomes
    selected_chromosomes = []
    for n in range(new_size):
        r = random.random()
        for (i, chromosome) in enumerate(chromosomes):
            if r <= probabilities[i]:
                selected_chromosomes.append(chromosome)
                break
    return selected_chromosomes

def find_closest(cities):
    dist_a = (((cities[0][0]-cities[1][0])**2)+(cities[0][1]-cities[1][1])**2)**0.5
    dist_b = (((cities[0][0]-cities[2][0])**2)+(cities[0][1]-cities[2][1])**2)**0.5
    if dist_a < dist_b:
        return cities[1]
    else:
        return cities[2]

#Greedy crossover selects the first city of one parent,
#compares the cities leaving that city in both parents,
#and chooses the closer one to extend the tour.
#If one city has already appeared in the tour, we choose the other city.
#If both cities have already appeared, we randomly select a non-selected city.    
def crossover(chromosomes, rate):
    new_chromosomes = []
    #len -1 because index out of range error might occur otherwise
    for i in range(len(chromosomes)-1):
        #Crossover T/F
        if random.random() < rate:
            #Take 2 chromosomes
            a = chromosomes[i]
            b = chromosomes[i+1]
            #Randomly select crossover location
            cross_location = random.randint(1, len(chromosomes[i]))
            #Create parts
            c = a[:cross_location]
            for i in range((len(a)-len(c))):
                #City to start from
                start_pos = c[-1]
                #2 most plausible targets
                target_0 = a[cross_location+i]
                target_1 = b[cross_location+i]
                #Determine if targets already appear in c and add correct city
                if target_0 in c:
                    if target_1 in c:
                        not_yet_selected = []
                        for city in genes:
                            if city not in c:
                                not_yet_selected.append(city)
                        c.append(not_yet_selected[random.randint(0,(len(not_yet_selected)-1))])
                    else:
                        c.append(target_1)
                elif target_1 in c:
                    c.append(target_0)
                else:
                    cities = [start_pos, target_0, target_1]
                    closest = find_closest(cities)
                    c.append(closest)
            new_chromosomes.append(c)
        else:
            new_chromosomes.append(chromosomes[i])
    return new_chromosomes

pygame.init()
screen = pygame.display.set_mode((640,340))
clock = pygame.time.Clock()
done = False
first = True
while done == False:
    for event in pygame.event.get():
        if event.type == pygame.QUIT:
            done = True
        if event.type == pygame.KEYDOWN:
            if event.key == pygame.K_UP:
                if genes != []:
                    #Start Genetic Algorithm
                    print('Creating initial solutions')
                    chromosomes = create_random_chromosomes(genes)
                    #Draw most possible solutions
                    for i in range(len(chromosomes)):
                        draw_connection(screen, chromosomes[i])
                    #Main Logic Loop
                    print('Starting calculations')
                    while len(chromosomes) > 10:
                        print(str(len(chromosomes)) + ' chromosomes left')
                        new_size = int(len(chromosomes)/resize_rate)
                        #Calculate fitness scores
                        fitness_scores = calculate_fitness(chromosomes)
                        #Build dictionary to preserve order
                        chrom_fit_dict = build_chromosome_fitness_dict(chromosomes, fitness_scores)
                        #Roulette selection - select chromosomes for next generation
                        selected_chromosomes = roulette_selection(chromosomes, fitness_scores, new_size)
                        #Adopt newly selected chromosomes as population
                        chromosomes = selected_chromosomes
                        #Apply crossover mutation
                        mutated_chromosomes = crossover(chromosomes, cross_rate)
                        #Adopt mutated chromosomes as population
                        chromosomes = mutated_chromosomes
                    best_chrom = fitness_scores.index(max(chrom_fit_dict))
                    draw_best_connection(screen, chromosomes[best_chrom])
            if event.key == pygame.K_c:
                draw_background(screen)
                genes = []
    if first == True:
        draw_background(screen)
        first = False
    #Get points
    if pygame.mouse.get_pressed() == (1, 0, 0):
        mouse_pos = pygame.mouse.get_pos()
        genes.append(mouse_pos)
        pygame.draw.rect(screen, black, [mouse_pos[0], mouse_pos[1], 2, 2], 1)
        genes = list(set(genes))
    pygame.display.flip()
    clock.tick(30)
pygame.quit()
share|improve this question
    
How many vertexes (cities) do you want to process? Do you need a genetic algorithm only, not an exact one? –  cat_baxter Jul 3 '12 at 20:54
    
It isn't meant for real-world application so it doesn't really matter how many cities/points it can process, largest test I ran was with 30 different points, it executed in about 7-8 minutes and was most of the times I ran it 90% - 100% correct. I'm only interested in the genetic algorithm, this was a test project to learn how these work, so an exact result is not requiered. –  Daquicker Jul 4 '12 at 12:57
    
Ok, I see. If you need the exact TSP solution (may be in order to check you genetic algorithm solution) let me know. It solves 30 points TSP in about 0.5 second. –  cat_baxter Jul 4 '12 at 15:14
    
oh, to check solutions that would be interesting indeed. –  Daquicker Jul 4 '12 at 17:13

1 Answer 1

up vote 1 down vote accepted

Ok, so here is exact ATSP solver, just in order to check your genetic algorithm solution. It's based on pascal code implementing the branch-and- bound algorithm by Little, Murty, Sweeney, and Karel. It's very very fast (30 points in 0.5 sec!).

'''

Branch and bound Algorithm for Travelling Salesman Problem

based on

http://www.cs.sunysb.edu/~algorith/implement/syslo/distrib/processed/babtsp.p

'''

import time

INF = 100000000

best_cost = 0

def reduce(size, w, row, col, rowred, colred):
    rvalue = 0
    for i in xrange(size):
        temp = INF
        for j in xrange(size):
            temp = min(temp, w[row[i]][col[j]])
        if temp > 0:
            for j in xrange(size):
                if w[row[i]][col[j]] < INF:
                    w[row[i]][col[j]] -= temp
            rvalue += temp
        rowred[i] = temp
    for j in xrange(size):
        temp = INF
        for i in xrange(size):
            temp = min(temp, w[row[i]][col[j]])
        if temp > 0:
            for i in xrange(size):
                if w[row[i]][col[j]] < INF:
                    w[row[i]][col[j]] -= temp
            rvalue += temp
        colred[j] = temp
    return rvalue


def bestEdge(size, w, row, col):
    mosti = -INF
    ri = 0
    ci = 0
    for i in xrange(size):
        for j in xrange(size):
            if not w[row[i]][col[j]]:
                minrowwelt = INF
                zeroes = 0
                for k in xrange(size):
                    if not w[row[i]][col[k]]:
                        zeroes += 1
                    else:
                        minrowwelt = min(minrowwelt, w[row[i]][col[k]])
                if zeroes > 1: minrowwelt = 0
                mincolwelt = INF
                zeroes = 0
                for k in xrange(size):
                    if not w[row[k]][col[j]]:
                        zeroes += 1
                    else:
                        mincolwelt = min(mincolwelt, w[row[k]][col[j]])
                if zeroes > 1: mincolwelt = 0
                if minrowwelt + mincolwelt > mosti:
                    mosti = minrowwelt + mincolwelt
                    ri = i
                    ci = j
    return mosti, ri, ci


def explore(n, w, edges, cost, row, col, best, fwdptr, backptr):
    global best_cost

    colred = [0 for _ in xrange(n)]
    rowred = [0 for _ in xrange(n)]
    size = n - edges
    cost += reduce(size, w, row, col, rowred, colred)
    if cost < best_cost:
        if edges == n - 2:
            for i in xrange(n): best[i] = fwdptr[i]
            if w[row[0]][col[0]] >= INF:
                avoid = 0
            else:
                avoid = 1
            best[row[0]] = col[1 - avoid]
            best[row[1]] = col[avoid]
            best_cost = cost
        else:
            mostv, rv, cv = bestEdge(size, w, row, col)
            lowerbound = cost + mostv
            fwdptr[row[rv]] = col[cv]
            backptr[col[cv]] = row[rv]
            last = col[cv]
            while fwdptr[last] != INF: last = fwdptr[last]
            first = row[rv]
            while backptr[first] != INF: first = backptr[first]
            colrowval = w[last][first]
            w[last][first] = INF
            newcol = [INF for _ in xrange(size)]
            newrow = [INF for _ in xrange(size)]
            for i in xrange(rv): newrow[i] = row[i]
            for i in xrange(rv, size - 1): newrow[i] = row[i + 1]
            for i in xrange(cv): newcol[i] = col[i]
            for i in xrange(cv, size - 1): newcol[i] = col[i + 1]
            explore(n, w, edges + 1, cost, newrow, newcol, best, fwdptr, backptr)
            w[last][first] = colrowval
            backptr[col[cv]] = INF
            fwdptr[row[rv]] = INF
            if lowerbound < best_cost:
                w[row[rv]][col[cv]] = INF
                explore(n, w, edges, cost, row, col, best, fwdptr, backptr)
                w[row[rv]][col[cv]] = 0

    for i in xrange(size):
        for j in xrange(size):
            w[row[i]][col[j]] = w[row[i]][col[j]] + rowred[i] + colred[j]


def atsp(w):
    global best_cost
    size = len(w)
    col = [i for i in xrange(size)]
    row = [i for i in xrange(size)]
    best = [0 for _ in xrange(size)]
    route = [0 for _ in xrange(size)]
    fwdptr = [INF for _ in xrange(size)]
    backptr = [INF for _ in xrange(size)]
    best_cost = INF

    explore(size, w, 0, 0, row, col, best, fwdptr, backptr)

    index = 0
    for i in xrange(size):
        route[i] = index
        index = best[index]
    index = []
    cost = 0

    for i in xrange(size):
        if i != size - 1:
            src = route[i]
            dst = route[i + 1]
        else:
            src = route[i]
            dst = 0
        cost += w[src][dst]
        index.append([src, dst])
    return cost, index


def main():
    # adjasted matrix

    m = [
        [INF, 514, 230, 92, 172, 201, 320, 205, 329, 285, 404, 128, 352, 357, 190, 269, 484, 237, 567, 209, 413, 404,
         333, 480, 371, 175, 470, 281, 423, 328],
        [452, INF, 359, 410, 450, 401, 315, 293, 244, 200, 275, 446, 469, 143, 278, 362, 399, 344, 544, 644, 335, 381,
         248, 457, 626, 375, 551, 258, 397, 305],
        [414, 360, INF, 90, 170, 42, 90, 46, 282, 238, 186, 126, 320, 198, 31, 267, 437, 190, 358, 623, 246, 345, 286,
         271, 267, 16, 326, 222, 421, 119],
        [324, 427, 269, INF, 80, 109, 359, 113, 263, 233, 312, 36, 260, 265, 98, 177, 418, 257, 475, 533, 321, 312, 281,
         388, 279, 83, 378, 189, 331, 236],
        [244, 382, 312, 226, INF, 335, 402, 246, 197, 153, 334, 262, 486, 466, 231, 315, 352, 297, 670, 453, 241, 538,
         201, 583, 505, 309, 298, 415, 350, 431],
        [448, 318, 355, 216, 268, INF, 445, 111, 240, 196, 288, 252, 385, 156, 96, 358, 395, 148, 691, 640, 348, 528,
         244, 455, 373, 299, 566, 405, 393, 452],
        [330, 437, 126, 216, 280, 119, INF, 123, 206, 315, 312, 252, 348, 275, 108, 324, 361, 267, 268, 539, 372, 255,
         242, 181, 393, 93, 236, 132, 478, 29],
        [495, 557, 521, 539, 313, 563, 469, INF, 406, 444, 543, 575, 274, 700, 440, 524, 442, 506, 729, 633, 554, 716,
         442, 642, 788, 537, 611, 593, 641, 490],
        [208, 256, 115, 168, 206, 157, 205, 49, INF, 148, 137, 204, 323, 313, 34, 118, 155, 100, 473, 417, 197, 460, 36,
         386, 382, 131, 441, 337, 272, 234],
        [252, 229, 159, 212, 250, 201, 249, 93, 44, INF, 181, 248, 367, 357, 78, 162, 199, 144, 517, 444, 241, 398, 48,
         430, 426, 175, 485, 275, 197, 278],
        [534, 568, 363, 453, 411, 250, 237, 254, 326, 344, INF, 225, 449, 406, 239, 366, 481, 398, 505, 547, 60, 375,
         362, 418, 623, 224, 473, 252, 520, 266],
        [445, 507, 334, 389, 469, 332, 393, 336, 413, 369, 276, INF, 224, 488, 321, 141, 392, 321, 659, 583, 336, 496,
         417, 572, 601, 306, 629, 373, 295, 420],
        [221, 283, 313, 313, 393, 314, 195, 263, 214, 170, 351, 349, INF, 426, 248, 325, 168, 314, 455, 359, 411, 442,
         218, 368, 580, 288, 431, 319, 367, 216],
        [309, 286, 216, 269, 307, 258, 306, 150, 101, 57, 132, 305, 326, INF, 135, 219, 256, 201, 574, 501, 192, 455,
         105, 487, 483, 232, 542, 332, 254, 335],
        [416, 554, 484, 398, 172, 507, 484, 15, 369, 325, 506, 434, 289, 638, INF, 487, 457, 469, 717, 625, 413, 554,
         373, 359, 677, 481, 470, 431, 522, 478],
        [514, 509, 193, 283, 352, 191, 252, 195, 306, 387, 379, 319, 469, 347, 180, INF, 461, 339, 518, 723, 249, 355,
         342, 431, 460, 165, 488, 232, 154, 279],
        [53, 115, 278, 145, 225, 254, 368, 212, 163, 119, 300, 181, 405, 258, 197, 281, INF, 263, 620, 262, 360, 457,
         167, 533, 424, 228, 523, 334, 316, 381],
        [300, 277, 207, 68, 148, 177, 297, 141, 92, 48, 229, 104, 328, 333, 126, 210, 247, INF, 543, 492, 289, 380, 96,
         456, 347, 151, 446, 257, 245, 304],
        [758, 471, 620, 710, 754, 593, 710, 597, 550, 671, 687, 746, 571, 614, 582, 668, 705, 650, INF, 930, 747, 599,
         586, 404, 887, 567, 915, 476, 822, 523],
        [39, 434, 269, 131, 211, 240, 359, 227, 178, 324, 315, 167, 391, 396, 212, 296, 333, 276, 452, INF, 375, 439,
         214, 365, 410, 214, 509, 316, 450, 213],
        [526, 513, 443, 496, 534, 485, 533, 377, 328, 284, 441, 165, 389, 641, 362, 306, 483, 428, 698, 487, INF, 535,
         332, 611, 710, 459, 769, 412, 460, 459],
        [784, 654, 691, 552, 604, 336, 781, 447, 576, 532, 624, 588, 721, 492, 432, 694, 731, 484, 1027, 976, 684, INF,
         580, 791, 709, 635, 902, 741, 729, 788],
        [435, 220, 576, 527, 607, 569, 450, 513, 464, 420, 495, 563, 689, 363, 498, 582, 619, 564, 718, 396, 555, 601,
         INF, 631, 806, 543, 686, 478, 617, 479],
        [354, 402, 216, 306, 350, 189, 306, 193, 146, 294, 283, 342, 167, 345, 178, 264, 301, 246, 358, 526, 343, 195,
         182, INF, 483, 163, 511, 72, 418, 119],
        [311, 288, 218, 79, 159, 188, 308, 152, 103, 59, 240, 115, 113, 344, 137, 221, 258, 11, 554, 472, 300, 391, 107,
         462, INF, 162, 457, 268, 256, 315],
        [431, 344, 186, 242, 187, 26, 276, 30, 266, 222, 314, 278, 304, 182, 15, 384, 421, 174, 544, 640, 374, 531, 270,
         374, 399, INF, 348, 408, 419, 305],
        [477, 454, 384, 437, 475, 426, 474, 318, 269, 225, 300, 473, 494, 168, 303, 387, 424, 369, 742, 669, 360, 623,
         273, 655, 651, 400, INF, 500, 422, 503],
        [282, 330, 144, 234, 280, 186, 234, 123, 74, 222, 211, 270, 366, 342, 108, 192, 229, 174, 286, 491, 271, 123,
         110, 199, 411, 160, 470, INF, 346, 47],
        [360, 399, 39, 129, 209, 81, 98, 85, 152, 277, 225, 165, 359, 237, 70, 270, 307, 229, 364, 569, 95, 201, 188,
         277, 306, 55, 334, 78, INF, 125],
        [301, 433, 97, 187, 267, 139, 187, 143, 177, 325, 283, 223, 319, 295, 128, 295, 332, 277, 239, 510, 343, 226,
         213, 152, 364, 113, 423, 103, 449, INF]
    ]

    start_time = time.time()
    cost, path = atsp(m)
    print "Cost = ", cost
    print "Path = ", path
    print "Time (s)", time.time() - start_time

if __name__ == "__main__":
    main()
share|improve this answer

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.