Differences

This shows you the differences between two versions of the page.

— random:datacenter_cooling_or_counting_hamiltonian_paths_in_a_grid [2012/03/23 15:15] (current)
grant created
Line 1: Line 1:
 +====== Datacenter Cooling (or Counting Hamiltonian Paths in a Grid) ======
 +
 +<code cpp>
 +// Datacenter Cooling (or Counting Hamiltonian Paths in a Grid)
 +//
 +// Written by Grant Jenks
 +// http://www.grantjenks.com/
 +//
 +// DISCLAIMER
 +// THERE IS NO WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE
 +// LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR
 +// OTHER PARTIES PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY OF ANY KIND,
 +// EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
 +// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE
 +// ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU.
 +// SHOULD THE PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY
 +// SERVICING, REPAIR OR CORRECTION.
 +//
 +// Copyright 2012
 +// This work is licensed under the
 +// Creative Commons Attribution-NonCommercial-NoDerivs 3.0 Unported License
 +// To view a copy of this license, visit
 +//    http://creativecommons.org/licenses/by-nc-nd/3.0/
 +// or send a letter to
 +//    Creative Commons
 +//    171 Second Street, Suite 300
 +//    San Francisco, California, 94105, USA
 +//
 +// I originally saw this problem at http://www.quora.com/challenges, but it
 +// looks like an ACM competition / AI textbook problem.
 +//
 +// Build with:
 +//    cl.exe /EHsc /O2 [/GL] [/DMEASURETIME] datacenter.cpp
 +//
 +//                             Datacenter Cooling
 +// We have some rooms in our datacenter, and we need to connect them all with
 +// a single cooling duct.
 +// Here are the rules:
 +//     * The datacenter is represented by a 2D grid.
 +//     * Rooms we own are represented by a 0.
 +//     * Rooms we do not own are represented by a 1.
 +//     * The duct has to start at the air intake valve, which is represented
 +//       by a 2.
 +//     * The duct has to end at the air conditioner, which is represented
 +//       by a 3.
 +//     * The duct cannot go in multiple directions out of the intake or the
 +//       AC - they must be the two endpoints of the duct.
 +//     * The duct must pass through each room exactly once.
 +//     * The duct cannot pass through rooms we do not own.
 +//     * The duct can connect between rooms horizontally or vertically but
 +//       not diagonally.
 +// Here is an example datacenter:
 +//     2 0 0 0
 +//     0 0 0 0
 +//     0 0 3 1
 +// There are two possible ways to run the duct here:
 +//     2--0--0--0
 +//              |
 +//     0--0--0--0
 +//     |
 +//     0--0--3  1
 +// or
 +//     2  0--0--0
 +//     |  |     |
 +//     0  0  0--0
 +//     |  |  |
 +//     0--0  3  1
 +// Write a program to compute the number of possible ways to run the duct. For
 +// the above example, the correct answer is 2.
 +//
 +// Input format:
 +// Your program should read from stdin. The input will be a series of
 +// integers separated by whitespace.
 +// The first two integers will be W and H, with width and height of the
 +// datacenter.
 +// These will be followed by W*H more integers, specifying the 2D grid
 +// representing the datacenter.
 +//
 +// Output format:
 +// Your program should write a single integer out to stdout: the number of
 +// possible ways to run the duct.
 +//
 +//                                 Test Inputs
 +// Small (result: 2)
 +//    4 3
 +//    2 0 0 0
 +//    0 0 0 0
 +//    0 0 3 1
 +// Medium (result: 20)
 +//    5 4
 +//    2 0 0 0 0
 +//    0 0 0 0 0
 +//    0 0 0 0 0
 +//    0 0 0 0 3
 +// Benchmarking (result: 301716)
 +//    7 8
 +//    2 0 0 0 0 0 0
 +//    0 0 0 0 0 0 0
 +//    0 0 0 0 0 0 0
 +//    0 0 0 0 0 0 0
 +//    0 0 0 0 0 0 0
 +//    0 0 0 0 0 0 0
 +//    0 0 0 0 0 0 0
 +//    3 0 0 0 0 1 1
 +//
 +//                                  Changelog
 +// 2012-03-15 Grant Jenks
 +// * Did a little research online and learned a new trick. If, while computing
 +//   reachability (the inner-dfs loop), a dead-end is discovered and the
 +//   dead-end is not the end position, then there are no possible paths. This
 +//   is a cheaper and narrower test for the observation in 2010-07-24.
 +//   Processing time: ~2.2 seconds
 +//
 +// 2012-03-14 Grant Jenks
 +// * Realized that I could memoize the tape. What if you occupy the same rooms
 +//   in two different configurations? Wouldn't the result be the same after
 +//   that point? For example:
 +//   If the starting board is:
 +//   S 0 0 0 0
 +//   0 0 0 0 0
 +//   0 0 0 0 0
 +//   0 0 0 0 0
 +//   0 0 0 0 E
 +//   And we reach this point (numbers indicate traversal pattern):
 +//   S 1 8 0 0
 +//   3 2 7 0 0
 +//   4 5 6 0 0
 +//   0 0 0 0 0
 +//   0 0 0 0 E
 +//   Then we will compute how many paths are reachable.
 +//   And if later we reach this point:
 +//   S 7 8 0 0
 +//   1 6 5 0 0
 +//   2 3 4 0 0
 +//   0 0 0 0 0
 +//   0 0 0 0 E
 +//   Then we know already (from memoization of the results) that this will
 +//   yield the same number of paths.
 +//   I prototyped these changes in a Python version of this program. While
 +//   the memoized results were frequently hit, cache maintenance and size
 +//   were impractical.
 +//
 +// 2010-07-24 Grant Jenks
 +// * Came up with an idea to approximate this problem as a max-flow problem.
 +//   This led me to consider what would happen if there existed an open
 +//   room, r, not adjacent to the end room, e, such that closing the room
 +//   would make the DFS test fail. Before this is done, I think it would
 +//   be good to add more comments and refactor the code some.
 +//   I got distracted from completing the above and didn't come back to the
 +//   code until a different idea come to me two years later.
 +//
 +// 2010-04-03 Grant Jenks
 +// * Converted conditionals in loop to straight line code.
 +//   Processing time: ~40 seconds
 +//
 +// 2010-03-27 Grant Jenks
 +// * Added BOUNDARY to datacenter array to eliminate bounds checking.
 +//   Processing time: ~60 seconds
 +//
 +// 2010-03-20 Grant Jenks
 +// * Added post-inner-DFS loop check that all open rooms were reachable.
 +//   Processing time: ~2 minutes
 +//
 +// 2010-03-13 Grant Jenks
 +// * Added the inner DFS loop but only checked that end was reachable.
 +//   Processing time: ~2 hours
 +//
 +// 2010-03-06 Grant Jenks
 +// * Created the basic outer loop with tape.
 +//   Processing time: Several hours.
 +
 +#define OPEN 1
 +#define CLOSED 0
 +#define VISITED 0
 +#define BOUNDARY 0
 +#define DFS_VISIT 8
 +#define VISIT_MASK 0x80000000
 +
 +// For data input/output
 +
 +#include <iostream>
 +#include <iterator>
 +
 +// For measuring time accurately
 +
 +#ifdef MEASURETIME
 +#include <windows.h>
 +#endif MEASURETIME
 +
 +// To coerce the compiler to use certain forms
 +
 +#include <intrin.h>
 +
 +// Unsigned int sizes
 +
 +typedef unsigned __int8 __uint8;
 +typedef unsigned __int16 __uint16;
 +typedef unsigned __int32 __uint32;
 +
 +using namespace std;
 +
 +int main()
 +{
 +#ifdef MEASURETIME
 +  LARGE_INTEGER start_time;
 +  QueryPerformanceCounter(&start_time);
 +#endif
 +
 +  __uint32 o_width, o_height, a_width, a_height;
 +  __uint32 datacenter_size, start, end, not_owned;
 +
 +  // Use __uint16 though the backing array will be __uint8 so that input is not
 +  // interpreted as characters.
 +
 +  typedef istream_iterator<__uint16> iis;
 +
 +  // Read in the initial width and height and initialize the datacenter array.
 +
 +  iis input = iis(std::cin);
 +
 +  o_width = *input++;
 +  o_height = *input++;
 +  a_width = o_width + 2;
 +  a_height = o_height + 2;
 +  datacenter_size = a_width * a_height;
 +  __uint8 * datacenter = new __uint8[datacenter_size];
 +
 +  // The perimeter of the datacenter is denoted by the BOUNDARY. This allows us
 +  // to skip bounds checking when looking at adjacent rooms to visit.
 +
 +  for (__uint32 pos = 0; pos < a_width; pos++)
 +    {
 +      datacenter[pos] = BOUNDARY;
 +      datacenter[pos + (a_width * (a_height - 1))] = BOUNDARY;
 +    }
 +  for (__uint32 pos = 0; pos < a_height; pos++)
 +    {
 +      datacenter[pos * a_width] = BOUNDARY;
 +      datacenter[(pos * a_width) + (a_width - 1)] = BOUNDARY;
 +    }
 +
 +  // Read in the rest of the datacenter.
 +
 +  start = 0;
 +  end = 0;
 +  not_owned = 0;
 +
 +  for(__uint32 o_pos = 0; input != iis(); input++, o_pos++)
 +    {
 +      __uint8 value = *input;
 +      __uint32 o_pos_w = o_pos % o_width;
 +      __uint32 o_pos_h = o_pos / o_width;
 +      __uint32 pos = (o_pos_w + 1) + ((o_pos_h + 1) * a_width);
 +
 +      switch (value)
 + {
 + case 0:
 +  datacenter[pos] = OPEN;
 +  break;
 + case 1:
 +  not_owned++;
 +  datacenter[pos] = CLOSED;
 +  break;
 + case 2:
 +  if (start != 0)
 +    {
 +      cerr<<"Start defined twice!"<<endl;
 +              return 1;
 +    }
 +  start = pos;
 +  datacenter[pos] = OPEN;
 +  break;
 + case 3:
 +  if (end != 0)
 +    {
 +      cerr<<"End defined twice!"<<endl;
 +      return 2;
 +    }
 +  end = pos;
 +  datacenter[pos] = OPEN;
 +  break;
 + default:
 +  cerr<<"Invalid input: "<<value<<endl;
 +  return 3;
 + }
 +    }
 +
 +  // Define the tape.
 +
 +  // Note that we're playing a game with the visited bit by putting
 +  // it in the high-bit of the tape position. The positions stored
 +  // on the tape are most likely less than 2^31-1 so this won't be
 +  // an issue.
 +
 +  __uint32 * tape = new __uint32[4 * datacenter_size];
 +  __uint32 * tape_pos = tape;
 +  __uint32 * dfs = new __uint32[datacenter_size + 1];
 +
 +  // Because we've recorded the 'start' and 'end' positions, we can
 +  // open the rooms.
 +
 +  datacenter[start] = OPEN;
 +  datacenter[end] = OPEN;
 +
 +  *tape_pos = 0x7FFFFFFF;
 +
 +  // Rather than recursing, the tape is used to solve the problem
 +  // iteratively. This gives finer-grain control over the assembly.
 +
 +  tape_pos++;
 +  *tape_pos = start;
 +
 +  __uint32 steps = 0;
 +  __uint32 routes = 0;
 +  __uint32 max_steps = o_width * o_height - not_owned;
 +
 +  while (tape_pos != tape)
 +    {
 +      __uint32 pos = *tape_pos;
 +
 +      if (pos == end)
 + {
 +          steps++;
 +  *tape_pos |= VISIT_MASK;
 +  if (steps == max_steps)
 +    {
 +      routes++;
 +    }
 + }
 +      else
 + {
 +  __uint32 * stack_pos = tape_pos;
 +  __uint32 * dfs_pos = dfs;
 +
 +  int visited_count = 0;
 +
 +          // Do a DFS through the open rooms to make sure all are
 +          // reachable and none (but the end room) are a dead-end.
 +
 +  while (stack_pos >= tape_pos)
 +    {
 +      __uint32 cur = *stack_pos;
 +      stack_pos--;
 +
 +      if (datacenter[cur] != OPEN)
 +                  continue;
 +
 +      visited_count++;
 +
 +      datacenter[cur] = DFS_VISIT;
 +      dfs_pos++;
 +      *dfs_pos = cur;
 +
 +      __uint32 above = cur - a_width;
 +      __uint32 left = cur - 1;
 +      __uint32 right = cur + 1;
 +      __uint32 below = cur + a_width;
 +
 +      __uint8 value, sum = 0;
 +
 +      stack_pos++;
 +
 +              // A trick is being played here with 'value &= 1' so that
 +              // the following is straight-line code. (Conditional branches
 +              // will cause hiccups in the processor.) Since only open
 +              // rooms have value 1, they alone will be stored on the stack.
 +              // This avoids visiting closed or visited rooms.
 +
 +      value = datacenter[above];
 +              sum += value;
 +              value &= 1;
 +      *stack_pos = above;
 +      stack_pos += value;
 +
 +      value = datacenter[right];
 +              sum += value;
 +              value &= 1;
 +      *stack_pos = right;
 +      stack_pos += value;
 +
 +      value = datacenter[below];
 +              sum += value;
 +              value &= 1;
 +      *stack_pos = below;
 +      stack_pos += value;
 +
 +      value = datacenter[left];
 +              sum += value;
 +              value &= 1;
 +      *stack_pos = left;
 +      stack_pos += value;
 +
 +              if (sum == DFS_VISIT && cur != end)
 +                 {
 +                   // If the sum of the values of the surrounding rooms
 +                   // is DFS_VISIT then it must be that all but one are
 +                   // closed - e.g. with p denoting the current position:
 +                   // 0 0 0
 +                   // 0 p 8
 +                   // 0 0 0
 +                   // This is a dead-end. Unless this position is also
 +                   // end, there will be no path through here that connects
 +                   // start and end.
 +                   // Set visited_count to max_steps only so the
 +                   // reachability check will fail.
 +
 +                   visited_count = max_steps;
 +                   break;
 +                 }
 +
 +      stack_pos--;
 +    }
 +
 +          // Re-open those rooms which were marked DFS_VISIT within the
 +          // above loop.
 +
 +  for(; dfs_pos > dfs; dfs_pos--)
 +    {
 +      datacenter[*dfs_pos] = OPEN;
 +    }
 +
 +  *tape_pos = (pos | VISIT_MASK);
 +  datacenter[pos] = VISITED;
 +
 +          // This is the reachability check. We will only append new
 +          // positions to visit on the tape if
 +
 +  if ((visited_count + steps) == max_steps)
 +    {
 +      __uint32 above = pos - a_width;
 +      __uint32 left = pos - 1;
 +      __uint32 right = pos + 1;
 +      __uint32 below = pos + a_width;
 +
 +      __uint8 value;
 +
 +      tape_pos++;
 +
 +              // This is using the same trick as the DFS loop for
 +              // reachability. There's no need to 'value &= 1'
 +              // because every room has either a 1 or 0 in it. If
 +              // we add a new position to the tape it is initially
 +              // added un-visited. We will mark it visited later
 +              // with VISIT_MASK.
 +
 +      value = datacenter[above];
 +      *tape_pos = above;
 +      tape_pos += value;
 +
 +      value = datacenter[right];
 +      *tape_pos = right;
 +      tape_pos += value;
 +
 +      value = datacenter[below];
 +      *tape_pos = below;
 +      tape_pos += value;
 +
 +      value = datacenter[left];
 +      *tape_pos = left;
 +      tape_pos += value;
 +
 +      tape_pos--;
 +    }
 +
 +          steps++;
 + }
 +
 +      // This is where the tape is used to rewind rooms occupied in the
 +      // datacenter. Above, on either sides of the if/else, was:
 +      //    *tape_pos = (pos | VISIT_MASK);
 +      // which set the high-bit in the tape position. Now, _bittestandreset,
 +      // returns that bit and sets it to 0. So the pseudocode reads:
 +      //    while the high bit at tape_pos is set
 +      //       clear the high bit
 +      //       mark the room at tape_pos as OPEN
 +      //       decrement tape_pos and steps
 +
 +      while(_bittestandreset((long *)tape_pos, 31))
 + {
 +  datacenter[*tape_pos] = OPEN;
 +  tape_pos--;
 +  steps--;
 + }
 +    }
 +
 +  cout<<routes<<endl;
 +
 +#ifdef MEASURETIME
 +  LARGE_INTEGER stop_time;
 +  QueryPerformanceCounter(&stop_time);
 +  LARGE_INTEGER frequency_time;
 +  QueryPerformanceFrequency(&frequency_time);
 +  double seconds = (double)(stop_time.QuadPart - start_time.QuadPart) /
 +    (double)(frequency_time.QuadPart);
 +  cout<<"Computed in "<<seconds<<" seconds."<<endl;
 +#endif MEASURETIME
 +}
 +</code>
Trace:
random/datacenter_cooling_or_counting_hamiltonian_paths_in_a_grid.txt · Last modified: 2012/03/23 15:15 by grant