Source Discussion


Let's now look at the source that implements the simulated annealing algorithm to solve the N-queens problem.

First, we'll look at the symbolic constants and types that are used by the algorithm. See Listing 2.1.

Listing 2.1: Types and Symbolic Constants.
start example
 #define MAX_LENGTH 30 typedef int solutionType[MAX_LENGTH]; typedef struct {   solutionType solution;   float energy; } memberType; /* Annealing Schedule */ #define INITIAL_TEMPERATURE      30.0 #define FINAL_TEMPERATURE        0.5 #define ALPHA                    0.99 #define STEPS_PER_CHANGE         100 
end example
 

The solutionType is our encoding for the N-Queens problem. The MAX_LENGTH symbolic defines the size of the board (in this case, we're solving the 30-Queens problem). The MAX_LENGTH can be increased (to 50 or more), though some changes to the cooling schedule may be necessary for larger solutions.

We encapsulate the solution within another type called memberType , which also includes the energy assessed for the solution.

The remainder of Listing 2.1 defines the cooling schedule. The INITIAL_TEMPERATURE and FINAL_TEMPERATURE define the bounds of the schedule and ALPHA is the constant used for geometric cooling. STEPS_PER_CHANGE defines the number of iterations that we'll perform at each temperature change (plateau).

Next , we'll look at the support functions for the algorithm. Listing 2.2 contains the encoding initialization and tweaking functions. These are used to create an initial solution and to randomly perturb it.

Listing 2.2: Initialization and Tweaking Functions.
start example
 void tweakSolution( memberType *member ) {   int temp, x, y;   x = getRand(MAX_LENGTH);   do {     y = getRand(MAX_LENGTH);   } while (x == y);   temp = member->solution[x];   member->solution[x] = member->solution[y];   member->solution[y] = temp; } void initializeSolution( memberType *member ) {   int i;   /* Initial setup of the solution */   for (i = 0 ; i < MAX_LENGTH ; i++) {     member->solution[i] = i;   }   /* Randomly perturb the solution */   for (i = 0 ; i < MAX_LENGTH ; i++) {     tweakSolution( member );   } } 
end example
 

In initializeSolution , we create a solution where each of the N queens is placed on the board. Each queen has identical row and column indices, thereby ensuring that no horizontal or vertical conflicts exist. The solution is then perturbed using the tweakSolution function. The tweakSolution function is used later in the algorithm to perturb the working solution as derived from the current solution.

Assessing a solution is the job of computeEnergy . This function identifies any conflicts that exist with the current solution (see Listing 2.3).

Listing 2.3: Assessing the Solution.
start example
 void computeEnergy( memberType *member ) {   int i, j, x, y, tempx, tempy;   char board[MAX_LENGTH][MAX_LENGTH];   int conflicts;   const int dx[4] = {-1,  1, -1,  1};   const int dy[4] = {-1,  1,  1, -1};   /* Standard library function to clear memory*/   bzero( (void *)board, MAX_LENGTH * MAX_LENGTH );   for (i = 0 ; i < MAX_LENGTH ; i++) {     board[i][member->solution[i]] = 'Q';   }   /* Walk through each of the Queens, and compute the number   * of conflicts   */   conflicts = 0;   for (i = 0 ; i < MAX_LENGTH ; i++) {     x = i; y = member->solution[i];     /* NOTE: Based upon the encoding, horizontal and vertical      * conflicts will never occur!      */     /* Check diagonals */     for (j = 0 ; j < 4 ; j++) {       tempx = x ; tempy = y;       while(1) {         tempx += dx[j]; tempy += dy[j];         if ((tempx < 0)  (tempx >= MAX_LENGTH)                 (tempy < 0)  (tempy >= MAX_LENGTH)) break;         if (board[tempx][tempy] == 'Q') conflicts++;       }     }   }   member->energy = (float)conflicts; } 
end example
 

For demonstration purposes, we'll actually build the chessboard. This isn't necessary, but it makes the resulting solution much simpler to visualize. Note the dx and dy arrays. These delta arrays are used to calculate the next position on the board for each queen that we've found. In the first case, dx = -1 and dy = -1 . This corresponds to movement in the northwest direction. The final case, dx = 1 , dy = -1 , corresponds to movement in the northeast direction.

We select each queen on the board (denoted by the x and y pair) and then walk each of the four diagonals looking for conflicts (other queens in the path ). Each time one is found, we increment the conflicts variable. At the end of this routine, we load the conflicts into the member structure as the energy component.

The next support function is used to copy one solution to another. Recall from Listing 2.1 that a solution is encoded and stored within memberType . The copySolution routine copies the contents of one memberType to another (see Listing 2.4).

Listing 2.4: Copying One Solution to Another.
start example
 void copySolution( memberType *dest, memberType *src ) {   int i;   for (i = 0 ; i < MAX_LENGTH ; i++) {     dest->solution[i] = src->solution[i];   }   dest->energy = src->energy; } 
end example
 

The final support function that we'll investigate is emitSolution . This function quite simply creates a chessboard representation from an encoded solution and emits it to standard out. The solution to be emitted is passed as an argument to the function. See Listing 2.5.

Listing 2.5: Emitting a Solution in Chessboard Format.
start example
 void emitSolution( memberType *member ) {   char board[MAX_LENGTH][MAX_LENGTH];   int x, y;   bzero( (void *)board, MAX_LENGTH * MAX_LENGTH );   for (x = 0 ; x < MAX_LENGTH ; x++) {     board[x][member->solution[x]] = 'Q';   }   printf("board:\n");   for (y = 0 ; y < MAX_LENGTH ; y++) {     for (x = 0 ; x < MAX_LENGTH ; x++) {       if (board[x][y] == 'Q') printf("Q ");       else printf(". ");     }     printf("\n");   }   printf("\n\n"); } 
end example
 

Within emitSolution , the board is constructed based upon the encoding (row = index, contents = column). The board is then emitted with "Q" as a location with a queen and a "." as an empty location.

Now that we've looked at all of the support functions, we'll look at the actual simulated annealing algorithm (coded within the main() function in Listing 2.6).

Listing 2.6: Simulated Annealing Algorithm.
start example
 int main() {   int timer=0, step, solution=0, useNew, accepted;   float temperature = INITIAL_TEMPERATURE;   memberType current, working, best;   FILE *fp;   fp = fopen("stats.txt", "w");   srand(time(NULL));   initializeSolution( &current );   computeEnergy( &current );   best.energy = 100.0;   copySolution( &working, &current);   while (temperature > FINAL_TEMPERATURE) {     printf("Temperature : %f\n", temperature);     accepted = 0;     /* Monte Carlo Step */     for (step = 0 ; step < STEPS_PER_CHANGE ; step++) {       useNew = 0;       tweakSolution( &working );       computeEnergy( &working );       if (working.energy <= current.energy) {         useNew = 1;       } else {         float test = getSRand();         float delta = working.energy - current.energy;         float calc = exp(-delta/temperature);         if (calc > test) {           accepted++;           useNew = 1;         }       }       if (useNew){         useNew = 0;         copySolution( &current, &working );         if (current.energy < best.energy) {           copySolution( &best, &current );           solution = 1;         }       } else {         copySolution( &working, &current);       }     }     fprintf(fp, "%d %f %f %d\n",              timer++, temperature, best.energy, accepted);     printf("Best energy = %f\n", best.energy);     temperature *= ALPHA;   }   fclose(fp);   if (solution) {     emitSolution( &best );   }   return 0; } 
end example
 

The simulated annealing algorithm shown in Listing 2.6 follows very closely the summarized algorithm in Figure 2.1. After two ancillary steps (opening a file for plot data and seeding the random number generator), we initialize our current solution named current and assess the energy of the solution with computeEnergy . We then copy the current solution to the working solution and start the algorithm.

The outer loop of the algorithm tests the current temperature against the final temperature. Note that we exit the loop once the current temperature is less than or equal to the FINAL_TEMPERATURE . This avoids us having to use a zero temperature in the acceptance probability equation.

The inner loop of the algorithm is called the Metropolis Monte Carlo simulation [Metropolis et al. 1953]. This loop is performed for a number of iterations at the current temperature to fully exploit the search capabilities of the temperature.

The first step in the Monte Carlo simulation is to tweak the working solution using the tweakSolution function. We then compute the energy of the working solution and compare it against our current solution. If the energy of the new working solution is less than (or equal to) the current solution, we accept it by default. Otherwise, we perform Equation 2.1 (acceptance probability equation) to determine if this worse solution should propagate. The delta is calculated as the difference of the working energy and the current energy. Note that to get here, the working energy is always greater than the current . We calculate a random number between 0 and 1 for our probability and then compare it with the result of Equation 2.1. If we've met the acceptance criteria (the result of Equation 2.1 is greater than our random number), we accept the working solution. This entails copying the working solution to current , since the working solution will again be tweaked at the next iteration of the Monte Carlo inner loop.

If the working solution is not accepted, we copy the current solution over the working solution. At the next iteration, the older working solution is lost and we tweak the existing current solution and try again.

After emitting some information to our statistics file, we reduce the temperature after performing the required number of iterations of the inner loop. Recall from Equation 2.2 that our temperature schedule is a simple geometric function. We multiply our current temperature by ALPHA , and repeat our outer loop.

At the end of the algorithm, we emit the best solution that was found (as long as one was actually found, as indicated by the solution variable). This variable is set within the inner loop after determining that a solution with less energy was found than the current best solution.




Visual Basic Developer
Visual Basic Developers Guide to ASP and IIS: Build Powerful Server-Side Web Applications with Visual Basic. (Visual Basic Developers Guides)
ISBN: 0782125573
EAN: 2147483647
Year: 1999
Pages: 175

flylib.com © 2008-2017.
If you may any questions please contact us: flylib@qtcs.net