// ////////////////////////////////////////////////////////
// # Title
// Cuboid layers
//
// # URL
// https://projecteuler.net/problem=126
// http://euler.stephan-brumme.com/126/
//
// # Problem
// The minimum number of cubes to cover every visible face on a cuboid measuring 3 x 2 x 1 is twenty-two.
//
// ![cuboids](p126.gif)
//
// If we then add a second layer to this solid it would require forty-six cubes to cover every visible face,
// the third layer would require seventy-eight cubes, and the fourth layer would require one-hundred and eighteen cubes to cover every visible face.
//
// However, the first layer on a cuboid measuring 5 x 1 x 1 also requires twenty-two cubes; similarly the first layer on cuboids measuring 5 x 3 x 1, 7 x 2 x 1, and 11 x 1 x 1 all contain forty-six cubes.
//
// We shall define `C(n)` to represent the number of cuboids that contain n cubes in one of its layers. So `C(22) = 2`, `C(46) = 4`, `C(78) = 5`, and `C(118) = 8`.
//
// It turns out that 154 is the least value of n for which `C(n) = 10`.
//
// Find the least value of n for which `C(n) = 1000`.
//
// # Solved by
// Stephan Brumme
// August 2017
//
// # Algorithm
// The initial cuboid (without any layers) has six sides:
// - front and back side, both with an area of `width * height`
// - left and right side, both with an area of `height * depth`
// - top and bottom side, both with an area of `depth * width`
// The total surface area is `2 * (width * height + height * depth + depth * width)`.
// Each part of the surface is covered by the first layer. That means:
// `C(1) = 2 * (width * height + height * depth + depth * width)`.
//
// For the 3 x 2 x 1 cuboid: `C(1) = 2 * (3 * 2 + 2 * 1 + 1 * 3) = 2 * (6 + 2 + 3) = 2 * 11 = 22`
//
// The computation of the second, third, fourth, ... layer wasn't obvious to me. My drawing skills are sub-par at best.
// So I decided to write some code: start with a shape and move that shape one unit left, right, up, down, forward and backward
// (start each move at the initial position). The difference between the old shape and the shape created by these six movements is the volume of the new layer.
// My function ''naive()'' (last lines of my code) performs exactly these steps. Starting with a specified cuboid it finds the volume of each additional layer and prints it.
//
// Playing around with several cuboid sizes I saw a pattern:
// `C(2) = C(1) + 4 * (width + height + depth)`
// `C(3) = C(2) + 4 * (width + height + depth) + 8`
// `C(4) = C(3) + 4 * (width + height + depth) + 16`
// `C(5) = C(4) + 4 * (width + height + depth) + 24`
//
// A bit more general:
// `C(1) = 2 * (width * height + height * depth + depth * width)`.
// `C(n) = C(n-1) + 4 * (width + height + depth) + 8 * (n - 2)`
//
// You can convert the second formula into a closed form but then the equation looks a bit messy.
// The C++ compiler won't generate much faster code, too. The function ''Cuboid::calculate'' evaluates these two equations.
//
// In order to solve `C(n) = x` for arbitrary `x` (it's `x = 1000` in the original problem) I wrote code based on ''std::map'':
// - each cuboids is represented by ''Cuboid'', it stores the size of the initial cuboid and the number of layers
// - ''Cuboid'' can be sorted by the size of the outermost layer (and then by it's sides)
// - since ''std::map'' is an ascendingly sorted container, it's first element is the cuboid with the smallest outermost layer
//
// Starting with the most basic cuboid (1,1,1) with only 1 layer I repeat the following:
// - pick and remove the first (= smallest) cuboid from the container
// - add a cuboid with ''width + 1'' to the container
// - add another cuboid with ''height + 1''
// - and a cuboid with ''depth + 1''
// - finally add a cuboid with ''layer + 1'' to that container
//
// Adding a cuboid is encapsulated in the ''add()'' function. It rejects cuboids where ''height > width'' and ''depth > height'' to avoid
// having the same (rotated) cuboid multiple times in the container. Adding exactly the same container twice is prevented, too.
// On top of that, ''add()'' keeps track of the volume of the most recent layer (in ''count[]'').
//
// That algorithm stops when the smallest cuboid has exactly 1000 "siblings" (==> the ''count'' of the current volume is 1000).
// It shows the correct result after about 4 seconds.
//
// However, my profiler showed that a few million cuboids where added (and removed) from the container.
// The overhead of memory allocation was responsible for the majority of the execution time.
//
// A much faster approach is to guess an upper limit of the result and then create all cuboids (and layers) that fit into that limit.
// Keeping track of just the various values of `C()` instead of the cuboids can be done with a small ''std::vector''.
// The program becomes about 50x faster and has almost no memory consumption. See my function ''fast()''.
//
// The main problem is to find a good guess for the upper limit.
// ''fast()'' returns 0 if no result was found, that's when the initial guess has to be incremented and another call to ''fast()'' is required.
// Starting with a very high guess will find the result but might waste lots of time computing cuboids with too large volumes.
// Currently, I use a ''StepSize = 10000'' ... yes, I know it's not optimal (it doesn't find the result in its first iteration)
// but I deliberately keep it this way.
// If ''StepSize = 50'' then ''fast()'' isn't faster than my original code.
// If ''StepSize = 150000'' then ''fast()'' isn't faster than my original code.
// That means for the original problem, ''fast()'' is only truly faster if ''50 < StepSize < 150000''.
//
// # Note
// Even though the code of ''fast()'' is shorter (and can be faster), I have a strange feeling about it:
// somehow it needs a good guess but I haven't found a way to come up with a straightforward way to generate that guess.
//
// I was stuck for a while because initially I stopped when I found a layer volume with __at least__ 1000 matching cuboids (instead of __exactly__ 1000).
// There are a few volumes below the correct solution with more than 1000 cuboids.
//
// Kristian Edlund visualized a nice way to enumerate the cubes added in each layer. His blog can found in my [links](#links).

#include <iostream>
#include <vector>
#include <set>
#include <map>

// a cuboid with 1+ added layers
struct Cuboid
{
  // cuboid size
  const unsigned short width;
  const unsigned short height;
  const unsigned short depth;
  // wrapped several times
  const unsigned short layers;
  // compute volume with calculated()
  const unsigned int   hullVolume;

  // create new cuboid
  Cuboid(unsigned short width_, unsigned short height_, unsigned short depth_, unsigned short layers_)
  : width(width_), height(height_), depth(depth_), layers(layers_),
    hullVolume(calculate())
  {}

  // calculate volume of the most recent / outermost layer
  unsigned int calculate() const
  {
    // innermost hull
    auto result = 2 * ((unsigned int)width  * height +
                       (unsigned int)height * depth  +
                       (unsigned int)depth  * width);

    // need just one layer ?
    if (layers == 1)
      return result;

    // add more layers
    auto increment = 4 * (width + height + depth);
    for (unsigned short current = 2; current <= layers; current++)
    {
      result    += increment;
      increment += 8;
    }

    return result;
  }

  // for std::set
  bool operator<(const Cuboid& other) const
  {
    // volume
    if (hullVolume < other.hullVolume) return true;
    if (hullVolume > other.hullVolume) return false;
    // the following lines are arbitrary => they only need to be consistent
    // sides
    if (width  < other.width)  return true;
    if (width  > other.width)  return false;
    if (height < other.height) return true;
    if (height > other.height) return false;
    if (depth  < other.depth)  return true;
    if (depth  > other.depth)  return false;
    // layers
    return layers < other.layers;
  }
};

// list of cuboids, sorted by volume of their outermost layer
std::set<Cuboid> todo;
// count how many cuboid layers have a certain volume (looking for 1000)
std::map<unsigned int, unsigned int> count; // could be replaced by std::vector if you can estimate the upper limit
// note: this information is implicitely stored in "todo" as well: count the consecutive cuboids with the same hullVolume

// add a cuboid with a certain amount of layers, return number of cuboids with the same hull volume
bool add(unsigned short width, unsigned short height, unsigned short depth, unsigned short layers)
{
  // accept only width >= height >= depth
  if (width  < height)
    return false;
  if (height < depth)
    return false;

  // create new cuboid
  Cuboid cuboid(width, height, depth, layers);

  // already queued ?
  if (todo.count(cuboid) != 0)
    return false;

  // add cuboid
  todo.insert(cuboid);
  count[cuboid.hullVolume]++;

  // successfully added
  return true;
}

const unsigned int Failed = 0;
// find the correct result if you know an upper bound maxVolume, if not found return Failed
unsigned int fast(unsigned int stopIf, unsigned int maxVolume)
{
  // allocate enough memory
  std::vector<unsigned int> count(maxVolume + 1, 0);

  // just four nested loops to generate all cuboids with hulls <= maxVolume
  // ensure that x >= y >= z to avoid enumerating some volume multiple times
  for (unsigned short x = 1; Cuboid(x,1,1,1).hullVolume <= maxVolume; x++)
    for (unsigned short y = 1; y <= x && Cuboid(x,y,1,1).hullVolume <= maxVolume; y++)
      for (unsigned short z = 1; z <= y && Cuboid(x,y,z,1).hullVolume <= maxVolume; z++)
        for (unsigned short layer = 1; ; layer++)
        {
          auto current = Cuboid(x,y,z,layer).hullVolume;
          if (current > maxVolume)
            break;

          count[current]++;
        }

  // find smallest hull volume that appears exactly 1000 times
  for (size_t i = 0; i < count.size(); i++)
    if (count[i] == stopIf)
      return i;

  // failed
  return Failed;
}

// forward declaration of my old helper code
//unsigned int naive(unsigned int x, unsigned int y, unsigned int z, unsigned int layers);

int main()
{
  // start old code to analyze a single cuboid
  //naive(3,2,1, 5); return 0;

  // look for that number of cuboids with the same volume
  unsigned int stopIf = 1000;
  std::cin >> stopIf;

  // ---------- fast enumeration ----------

  // find solution if an upper limit of the result is known
#define FAST_ENUMERATION
#ifdef  FAST_ENUMERATION
  const unsigned int StepSize = 10000; // I deliberately choose not the optimal value !
  // assume that the result is below 10000, then 20000, ... until fast() returns > 0
  for (unsigned int maxVolume = StepSize; ; maxVolume += StepSize)
  {
    auto current = fast(stopIf, maxVolume); // returns 0 if not found
    // found ?
    if (current != Failed)
    {
      std::cout << current << std::endl;
      return 0;
    }
  }
#endif

  // ---------- slower enumeration ----------

  // note: I found the correct result with this approach but later created the faster enumeration you find above

  // initial cuboid: one unit wide, high and deep with one layer around it
  add(1,1,1,1);
  unsigned int volume = 0;
  while (count[volume] != stopIf)
  {
    // pick cuboid with smallest volume and remove it
    auto current = *(todo.begin());
    todo.erase(todo.begin());

    volume = current.hullVolume;

    // enlarge cuboid
    add(current.width + 1, current.height, current.depth, current.layers);
    add(current.width, current.height + 1, current.depth, current.layers);
    add(current.width, current.height, current.depth + 1, current.layers);
    // add one layer
    add(current.width, current.height, current.depth, current.layers + 1);
  }

  // display result
  std::cout << volume << std::endl;
  return 0;
}

// ---------- compute hull volume, slow approach ----------

// a single cube that is part of an cuboid (more or less a coordinate)
struct Cube
{
  // coordinates
  const int x, y, z;

  // create new cube
  Cube(int xx, int yy, int zz)
  : x(xx), y(yy), z(zz)
  {}

  // for std::set
  bool operator<(const Cube& other) const
  {
    if (x < other.x) return true;
    if (x > other.x) return false;
    if (y < other.y) return true;
    if (y > other.y) return false;
    return z < other.z;
  }
};

unsigned int naive(unsigned int x, unsigned int y, unsigned int z, unsigned int layers)
{
  // keep track of all cubes belonging to the shape
  std::set<Cube> cuboid;
  // cubes added in the most recent layer
  std::vector<Cube> hull;

  // initial shape
  for (int xx = 0; xx < (int)x; xx++)
    for (int yy = 0; yy < (int)y; yy++)
      for (int zz = 0; zz < (int)z; zz++)
      {
        Cube current(xx, yy, zz);
        cuboid.insert(current);
        hull.push_back(current);
      }

  // add desired numbers of layers
  for (unsigned int layer = 1; layer <= layers; layer++)
  {
    std::vector<Cube> next;
    for (auto current : hull)
    {
      // grow in each direction x,y,z
      Cube left (current.x - 1, current.y, current.z);
      Cube right(current.x + 1, current.y, current.z);
      Cube down (current.x, current.y - 1, current.z);
      Cube up   (current.x, current.y + 1, current.z);
      Cube back (current.x, current.y, current.z - 1);
      Cube front(current.x, current.y, current.z + 1);

      // add cube (if not already existing)
      if (cuboid.count(left)  == 0)
      {
        cuboid.insert(left);
        next.push_back(left);
      }
      if (cuboid.count(right) == 0)
      {
        cuboid.insert(right);
        next.push_back(right);
      }
      if (cuboid.count(down)  == 0)
      {
        cuboid.insert(down);
        next.push_back(down);
      }
      if (cuboid.count(up)    == 0)
      {
        cuboid.insert(up);
        next.push_back(up);
      }
      if (cuboid.count(back)  == 0)
      {
        cuboid.insert(back);
        next.push_back(back);
      }
      if (cuboid.count(front) == 0)
      {
        cuboid.insert(front);
        next.push_back(front);
      }
    }

    std::cout << "layer " << layer << ": "
              << "C(" << x << "," << y << "," << z << ")=" << next.size();
    if (layer > 1)
      std::cout << " delta=" << next.size() - hull.size();
    std::cout << std::endl;

    // next iteration
    hull = std::move(next);
  }

  return hull.size();
}