generate pics of generalized chaos sierpinski-like figures on polygons

chaos-polygons.cxx

#include <iostream>
#include <tuple>
#include <cmath>
#include <cstdlib>
#include <iterator>
#include <vector>
#include <algorithm>
#include <string>
// STEPS=100;for i in $(seq $((3 * $STEPS)) $((8 * $STEPS))); do ./chaos-polygons $(bc -l <<< "scale=2;$i / $STEPS") > generated-cmax-100-$i.pgm; done
// for i in generated-cmax-100-*.pgm; do convert $i ${i%.pgm}.png; done

template<typename T, typename Op>
std::tuple<T,T> op(std::tuple<T,T> const& p0, std::tuple<T,T> const& p1, Op o){
  return std::make_tuple(o(std::get<0>(p0), std::get<0>(p1)), o(std::get<1>(p0),std::get<1>(p1)));
}

template<typename T, typename Op>
std::tuple<T,T> op(std::tuple<T,T> const& p, Op o){
  return std::tuple<T,T>(static_cast<T>(o(std::get<0>(p))), static_cast<T>(o(std::get<1>(p))));
}

template<typename T>
std::tuple<T,T> average(std::tuple<T,T> const& p0, std::tuple<T,T> const& p1){
  return op(op(p0, p1,[](auto a, auto b){return a+b;}), [](auto a){return a/2;});
}

template<typename T>
std::tuple<T,T> rotate(float a, std::tuple<T,T> const& p){
  return std::tuple<T,T>(std::cos(a)*std::get<0>(p) - std::sin(a)*std::get<1>(p)
                         , std::sin(a)*std::get<0>(p) + std::cos(a)*std::get<1>(p));
}

constexpr double pi() { return std::atan(1)*4; }

template<typename C, typename T, typename F>
void plot(std::vector<std::vector<C> >& screen, std::tuple<T,T> p, F f){
  f(screen[(std::get<1>(p)+1.)*screen.size()/2][(std::get<0>(p)+1.)*screen[0].size()/2]);
}

template<typename T, typename N, typename Out>
Out polygon(std::tuple<T,T> p, N n, Out o){
  for(int i= 0; i < std::ceil(n); ++i, ++o){
    *o= p;
    p= rotate(2*pi()/n, p);
  }
  return o;
}

template<typename C> std::vector<std::vector<C> > create_screen(C init, std::size_t w, std::size_t h){
  std::vector<std::vector<C> > screen(h, std::vector<C>(w, init));
  return screen;
}

template<typename C> std::ostream& operator<<(std::ostream& os, std::vector<std::vector<C> > const& screen){
  std::ostream_iterator<C> out(os, " ");
  std::for_each(screen.rbegin(), screen.rend()
                , [&out](std::vector<C>const& row){ std::copy(row.begin(), row.end(), out);});
  return os;
}

template<typename C> typename C::value_type random_elt(C const& c){
  return *(c.cbegin() + std::rand()% c.size());
}

int main(int argc, char* argv[]){
  float n_edges= argc > 1 ? std::atof(argv[1]) : 3.f;
  int constexpr w= 1024;
  int constexpr h= 1024;
  int constexpr white= 255;
  std::size_t n= (argc > 2 ? std::atol(argv[2]) : 1L<<25)*n_edges;
  
  typedef std::tuple<float, float> p_t;
  std::vector<p_t> poly;
  polygon(p_t(0.f, 0.75f), n_edges, std::back_inserter(poly));
  auto sc= create_screen(0, w, h);
  auto to_white= [white](int& c){ c= white;};
  std::for_each(poly.begin(), poly.end(), [&sc, &to_white](p_t const&p){plot(sc, p, to_white);});
  std::vector<p_t> points(poly.begin(), poly.end());
  auto brighten=[white](int& c){ c= std::min((c+1), white);};
  while(points.size()< n){
    auto p= average(random_elt(poly), random_elt(points));
    plot(sc, p, brighten);
    points.push_back(p);
  }
  std::cout << "P2" << std::endl << w << ' ' << h << std::endl << white << std::endl;
  std::cout<< sc << std::endl;
  return 0;
}