src/make-ruletable.cpp

   1                         /*** /
   2
   3 [This file is part of Bedbugs. All I've done is made it implement
   4 the Conway transition function. Original copyright notice follows:]
   5
   6 This file is part of Golly, a Game of Life Simulator.
   7 Copyright (C) 2008 Andrew Trevorrow and Tomas Rokicki.
   8
   9 This program is free software; you can redistribute it and/or
  10 modify it under the terms of the GNU General Public License
  11 as published by the Free Software Foundation; either version 2
  12 of the License, or (at your option) any later version.
  13
  14 This program is distributed in the hope that it will be useful,
  15 but WITHOUT ANY WARRANTY; without even the implied warranty of
  16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  17 GNU General Public License for more details.
  18
  19 You should have received a copy of the GNU General Public License
  20 along with this program; if not, write to the Free Software
  21 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
  22
  23  Web site:  http://sourceforge.net/projects/golly
  24  Authors:   rokicki@gmail.com  andrew@trevorrow.com
  25
  26                         / ***/
  27 #include <math.h>
  28 #include <time.h>
  29
  30 #include <iostream>
  31 #include <set>
  32 #include <fstream>
  33 #include <vector>
  34 #include <map>
  35 #include <sstream>
  36 #include <algorithm>
  37 using namespace std;
  38
  39 /*
  40
  41 Makes a rule-table for your transition function.
  42
  43 To compile:
  44 g++ -O5 -o make-ruletable make-ruletable.cpp
  45 or in Microsoft Visual Studio, add to an empty CLR project.
  46
  47 To use:
  48 1) fill slowcalc with your own transition function.
  49 2) set the parameters in main() at the bottom.
  50 3) execute the program from the command line (no options taken)
  51
  52 For a 32-state 5-neighbourhood rule it took 16mins on a 2.2GHz machine.
  53 For a 4-state 9-neighbourhood rule it took 4s.
  54
  55 The merging is very simple - if a transition is compatible with the first rule,
  56 then merge the transition into the rule. If not, try the next rule. If there are
  57 no more rules, add a new rule containing just the transition.
  58
  59 === Worked example: ===
  60 Transitions:
  61 1,0,0->3
  62 1,2,0->3
  63 2,0,0->3
  64 2,2,0->1
  65 Start. First transition taken as first rule:
  66 rule 1: 1,0,0->3
  67 Next transition - is it compatible with rule 1?
  68 i.e. is 1,[0,2],0->3 valid for all permutations? Yes.
  69 rule 1 now: 1,[0,2],0->3
  70 Next transition - is it compatible with rule 1?
  71 i.e. is [1,2],[0,2],0->3 valid for all permutations?
  72 no - because 2,2,0 ->1, not ->3. so add the transition as a new rule:
  73 rule 1 still: 1,[0,2],0 -> 3
  74 rule 2 now : 2,0,0->3
  75 Next transition - is it compatible with rule 1? no - output is different.
  76 Is it compatible with rule 2? no - output is different.
  77 Final output:
  78 1,[0,2],0 -> 3
  79 2,0,0 -> 3
  80 2,2,0 -> 1
  81 Written with variables:
  82 var a={0,2}
  83 1,a,0,3
  84 2,0,0,3
  85 2,2,0,1
  86 ===============
  87
  88 In the function produce_rule_table, the state space is exhaustively traversed.
  89 If your transition function consists of transition rules already then you can
  90 optimise by running through the transitions instead. You might also want to
  91 turn off the optimisation in rule::can_merge, to see if it gives you better
  92 compression.
  93
  94 Also note: I feel sure there are better ways to compress rule tables than this...
  95
  96 Contact: Tim Hutton <tim.hutton@gmail.com>
  97
  98 */
  99
 100 // some typedefs and compile-time constants
 101 typedef unsigned short state;
 102 enum TSymm { none, rotate4, rotate8, reflect, rotate4reflect, rotate8reflect };
 103 static const string symmetry_strings[] = {"none","rotate4","rotate8","reflect","rotate4reflect","rotate8reflect"};
 104
 105 // fill in this function with your desired transition rules
 106 // (for von Neumann neighbourhoods, just ignore the nw,se,sw,ne inputs)
 107 state slowcalc(state nw,state n,state ne,state w,state c,state e,state sw,state s,state se)
 108 {
 109    int neighbours = nw+n+ne+e+se+s+sw+w;
 110    // Conway:
 111    switch (c)
 112    {
 113      case 0:
 114         if (neighbours == 3) {
 115             return 1;
 116         } else {
 117             return 0;
 118         }
 119      case 1:
 120         if (neighbours == 2 || neighbours == 3) {
 121             return 1;
 122         } else {
 123             return 0;
 124         }
 125      default:
 126         return 0 ; // should throw an error here
 127    }
 128 }
 129
 130 vector<state> rotate_inputs(const vector<state>& inputs,int rot)
 131 {
 132    vector<state> rotinp(inputs);
 133    rotate_copy(inputs.begin()+1,inputs.begin()+1+rot,inputs.end(),rotinp.begin()+1);
 134    return rotinp;
 135 }
 136
 137 vector<state> reflect_inputs(const vector<state>& inputs,int neighbourhood_size)
 138 {
 139    vector<state> refinp(inputs);
 140    if(neighbourhood_size==5) // CNESW
 141    {
 142       refinp[2]=inputs[4]; // swap E and W
 143       refinp[4]=inputs[2];
 144    }
 145    else // neighbourhood_size==9 (C,N,NE,E,SE,S,SW,W,NW)
 146    {
 147       refinp[2]=inputs[8];
 148       refinp[8]=inputs[2];
 149       refinp[3]=inputs[7];
 150       refinp[7]=inputs[3];
 151       refinp[4]=inputs[6];
 152       refinp[6]=inputs[4]; // swap all E and W
 153    }
 154    return refinp;
 155 }
 156
 157 // simple rule structure, e.g. 1,2,[4,5],8,2 -> 0
 158 class rule {
 159
 160 public:
 161    set<state> inputs[9]; // c,n,ne,e,se,s,sw,w,nw  or  c,n,e,s,w
 162    state ns; // new state
 163
 164    int n_inputs; // 5: von Neumann; 9: Moore
 165    TSymm symm;
 166
 167 public:
 168    // constructors
 169    rule(const rule& r) : ns(r.ns),n_inputs(r.n_inputs),symm(r.symm)
 170    {
 171       for(int i=0;i<n_inputs;i++)
 172          inputs[i]=r.inputs[i];
 173    }
 174    rule& operator=(const rule& r)
 175    {
 176       n_inputs=r.n_inputs;
 177       symm = r.symm;
 178       ns = r.ns;
 179       for(int i=0;i<n_inputs;i++)
 180          inputs[i]=r.inputs[i];
 181       return *this;
 182    }
 183    rule(const vector<state>& inputs,int n_inputs,state ns1,TSymm symm1)
 184       : ns(ns1),n_inputs(n_inputs),symm(symm1)
 185    {
 186       merge(inputs);
 187    }
 188
 189    // if we merge the rule and the supplied transition, will the rule remain true for all cases?
 190    bool can_merge(const vector<state>& test_inputs,state ns1) const
 191    {
 192       if(ns1!=ns) return false; // can't merge if the outputs are different
 193
 194       // If you are running through your own transitions, or for small state spaces,
 195       // you might want to turn off this optimisation, to get better compression.
 196       // On JvN29 it doesn't make any difference but on Nobili32 it does.
 197       const bool forbid_multiple_input_differences = true;
 198
 199       if(forbid_multiple_input_differences)
 200       {
 201          // optimisation: we skip this rule if more than 2 entries are different, we
 202          // assume we will have considered the relevant one-change rules before this.
 203          int n_different=0;
 204          for(int i=0;i<n_inputs;i++)
 205             if(inputs[i].find(test_inputs[i])==inputs[i].end())
 206                if(++n_different>1)
 207                   return false;
 208          // just check the new permutations
 209          for(int i=0;i<n_inputs;i++)
 210          {
 211             if(inputs[i].find(test_inputs[i])==inputs[i].end())
 212             {
 213                rule r1(*this);
 214                r1.inputs[i].clear(); // (since we don't need to re-test all the other permutations)
 215                r1.inputs[i].insert(test_inputs[i]);
 216                if(!r1.all_true()) return false;
 217             }
 218          }
 219       }
 220       else
 221       {
 222          // need to check all combinations - this can be very slow for large state spaces
 223          for(int i=0;i<n_inputs;i++)
 224          {
 225             if(inputs[i].find(test_inputs[i])==inputs[i].end())
 226             {
 227                rule r1(*this);
 228                r1.merge(test_inputs); // this is what makes it slow, we may introduce many new permutations
 229                r1.inputs[i].clear(); // (since we don't need to re-test all the old permutations)
 230                r1.inputs[i].insert(test_inputs[i]);
 231                if(!r1.all_true()) return false;
 232             }
 233          }
 234       }
 235       return true;
 236    }
 237    // merge the inputs with this rule
 238    void merge(const vector<state>& new_inputs)
 239    {
 240       for(int i=0;i<n_inputs;i++)
 241          inputs[i].insert(new_inputs[i]); // may already exist, set ignores if so
 242    }
 243
 244    // is this set of inputs a match for the rule, for the given symmetry?
 245    bool matches(const vector<state>& test_inputs) const
 246    {
 247       int n_rotations,rotation_skip;
 248       bool do_reflect;
 249       switch(symm)
 250       {
 251          default:
 252          case none: n_rotations=1; rotation_skip=1; do_reflect=false; break;
 253          case rotate4:
 254             if(n_inputs==5)
 255             {
 256                n_rotations=4; rotation_skip=1; do_reflect=false;
 257             }
 258             else
 259             {
 260                n_rotations=4; rotation_skip=2; do_reflect=false;
 261             }
 262             break;
 263          case rotate8: n_rotations=8; rotation_skip=1; do_reflect=false; break;
 264          case reflect: n_rotations=1; rotation_skip=1; do_reflect=true; break;
 265          case rotate4reflect:
 266             if(n_inputs==5)
 267             {
 268                n_rotations=4; rotation_skip=1; do_reflect=true;
 269             }
 270             else
 271             {
 272                n_rotations=4; rotation_skip=2; do_reflect=true;
 273             }
 274             break;
 275          case rotate8reflect: n_rotations=8; rotation_skip=1; do_reflect=true; break;
 276       }
 277       for(int iRot=0;iRot<n_rotations;iRot++)
 278       {
 279          if(nosymm_matches(rotate_inputs(test_inputs,iRot*rotation_skip)))
 280             return true;
 281          if(do_reflect && nosymm_matches(reflect_inputs(rotate_inputs(test_inputs,iRot*rotation_skip),n_inputs)))
 282             return true;
 283       }
 284       return false; // no match found
 285    }
 286
 287 protected:
 288
 289    // ignoring symmetry, does this set of inputs match the rule?
 290    bool nosymm_matches(const vector<state>& test_inputs) const
 291    {
 292       for(int i=0;i<n_inputs;i++)
 293          if(inputs[i].find(test_inputs[i])==inputs[i].end())
 294             return false;
 295       return true;
 296    }
 297
 298    // is the rule true in all permutations?
 299    bool all_true() const
 300    {
 301       set<state>::const_iterator c_it,n_it,ne_it,e_it,se_it,s_it,sw_it,w_it,nw_it;
 302       if(n_inputs==9)
 303       {
 304          for(c_it = inputs[0].begin();c_it!=inputs[0].end();c_it++)
 305             for(n_it = inputs[1].begin();n_it!=inputs[1].end();n_it++)
 306                for(ne_it = inputs[2].begin();ne_it!=inputs[2].end();ne_it++)
 307                   for(e_it = inputs[3].begin();e_it!=inputs[3].end();e_it++)
 308                      for(se_it = inputs[4].begin();se_it!=inputs[4].end();se_it++)
 309                         for(s_it = inputs[5].begin();s_it!=inputs[5].end();s_it++)
 310                            for(sw_it = inputs[6].begin();sw_it!=inputs[6].end();sw_it++)
 311                               for(w_it = inputs[7].begin();w_it!=inputs[7].end();w_it++)
 312                                  for(nw_it = inputs[8].begin();nw_it!=inputs[8].end();nw_it++)
 313                                     if(slowcalc(*nw_it,*n_it,*ne_it,*w_it,*c_it,*e_it,*sw_it,*s_it,*se_it)!=ns)
 314                                        return false;
 315       }
 316       else
 317       {
 318          for(c_it = inputs[0].begin();c_it!=inputs[0].end();c_it++)
 319             for(n_it = inputs[1].begin();n_it!=inputs[1].end();n_it++)
 320                for(e_it = inputs[2].begin();e_it!=inputs[2].end();e_it++)
 321                   for(s_it = inputs[3].begin();s_it!=inputs[3].end();s_it++)
 322                      for(w_it = inputs[4].begin();w_it!=inputs[4].end();w_it++)
 323                         if(slowcalc(0,*n_it,0,*w_it,*c_it,*e_it,0,*s_it,0)!=ns)
 324                            return false;
 325       }
 326       return true;
 327    }
 328 };
 329
 330 // makes a unique variable name for a given value
 331 string get_variable_name(unsigned int iVar)
 332 {
 333    const char VARS[52]={'a','b','c','d','e','f','g','h','i','j',
 334       'k','l','m','n','o','p','q','r','s','t','u','v','w','x',
 335       'y','z','A','B','C','D','E','F','G','H','I','J','K','L',
 336       'M','N','O','P','Q','R','S','T','U','V','W','X','Y','Z'};
 337    ostringstream oss;
 338    if(iVar<52)
 339       oss << VARS[iVar];
 340    else if(iVar<52*52)
 341       oss << VARS[(iVar-(iVar%52))/52 - 1] << VARS[iVar%52];
 342    else
 343       oss << "!"; // we have a 52*52 limit ("should be enough for anyone")
 344    return oss.str();
 345 }
 346
 347 void print_rules(const vector<rule>& rules,ostream& out)
 348 {
 349    // first collect all variables (makes reading easier)
 350    map< string , set<state> > vars;
 351    ostringstream rules_out;
 352    for(vector<rule>::const_iterator r_it=rules.begin();r_it!=rules.end();r_it++)
 353    {
 354       vector<string> variables_used;
 355       for(int i=0;i<r_it->n_inputs;i++)
 356       {
 357          // if more than one state for this input, we need a variable
 358          if(r_it->inputs[i].size()>1)
 359          {
 360             string var;
 361             // is there a variable that matches these inputs, and that we haven't used?
 362             bool found_unused_var=false;
 363             for(map<string, set<state> >::const_iterator v_it=vars.begin();v_it!=vars.end();v_it++)
 364             {
 365                if(v_it->second==r_it->inputs[i] && find(variables_used.begin(),variables_used.end(),v_it->first)==variables_used.end())
 366                {
 367                   found_unused_var = true;
 368                   var = v_it->first;
 369                   break;
 370                }
 371             }
 372             if(!found_unused_var)
 373             {
 374                // we need to make a new one for this set of inputs
 375                var = get_variable_name(vars.size());
 376                // add it to the list of made variables
 377                vars[var] = r_it->inputs[i];
 378                // print it
 379                out << "var " << var << "={";
 380                set<state>::const_iterator it=r_it->inputs[i].begin();
 381                while(true)
 382                {
 383                   out << (int)*it;
 384                   it++;
 385                   if(it!=r_it->inputs[i].end()) out << ",";
 386                   else break;
 387                }
 388                out << "}\n";
 389             }
 390             // add the variable to the list of used ones
 391             variables_used.push_back(var);
 392             rules_out << var << ",";
 393          }
 394          else
 395          {
 396             // just a state, output it
 397             rules_out << (int)*r_it->inputs[i].begin() << ",";
 398          }
 399       }
 400       rules_out << (int)r_it->ns << endl;
 401    }
 402    out << rules_out.str();
 403 }
 404
 405 void produce_rule_table(vector<rule>& rules,int N,int nhood_size,TSymm symm,bool remove_stasis)
 406 {
 407    int n_rotations,rotation_skip;
 408    bool do_reflect;
 409    switch(symm)
 410    {
 411       default:
 412       case none: n_rotations=1; rotation_skip=1; do_reflect=false; break;
 413       case rotate4:
 414          if(nhood_size==5)
 415          {
 416             n_rotations=4; rotation_skip=1; do_reflect=false;
 417          }
 418          else
 419          {
 420             n_rotations=4; rotation_skip=2; do_reflect=false;
 421          }
 422          break;
 423       case rotate8: n_rotations=8; rotation_skip=1; do_reflect=false; break;
 424       case reflect: n_rotations=1; rotation_skip=1; do_reflect=true; break;
 425       case rotate4reflect:
 426          if(nhood_size==5)
 427          {
 428             n_rotations=4; rotation_skip=1; do_reflect=true;
 429          }
 430          else
 431          {
 432             n_rotations=4; rotation_skip=2; do_reflect=true;
 433          }
 434          break;
 435       case rotate8reflect: n_rotations=8; rotation_skip=1; do_reflect=true; break;
 436    }
 437
 438    state c,n,ne,nw,sw,s,se,e,w,ns;
 439    vector<rule>::iterator it;
 440    bool merged;
 441    for(c=0;c<N;c++)
 442    {
 443       cout << "\nProcessing for c=" << (int)c << ", " << rules.size() << " rules so far." << endl;
 444
 445       if(nhood_size==9)
 446       {
 447          vector<state> inputs(9);
 448          inputs[0]=c;
 449          for(n=0;n<N;n++)
 450          {
 451             cout << ".";
 452             cout.flush();
 453             inputs[1]=n;
 454             for(ne=0;ne<N;ne++)
 455             {
 456                inputs[2]=ne;
 457                for(e=0;e<N;e++)
 458                {
 459                   inputs[3]=e;
 460                   for(se=0;se<N;se++)
 461                   {
 462                      inputs[4]=se;
 463                      for(s=0;s<N;s++)
 464                      {
 465                         inputs[5]=s;
 466                         for(sw=0;sw<N;sw++)
 467                         {
 468                            inputs[6]=sw;
 469                            for(w=0;w<N;w++)
 470                            {
 471                               inputs[7]=w;
 472                               for(nw=0;nw<N;nw++)
 473                               {
 474                                  ns = slowcalc(nw,n,ne,w,c,e,sw,s,se);
 475                                  if(remove_stasis && ns == c)
 476                                     continue; // we can ignore stasis transitions
 477                                  // can we merge this transition with any existing rule?
 478                                  inputs[8]=nw;
 479                                  merged = false;
 480                                  for(it=rules.begin();!merged && it!=rules.end();it++)
 481                                  {
 482                                     rule &r = *it;
 483                                     for(int iRot=0;!merged && iRot<n_rotations;iRot++)
 484                                     {
 485                                        if(r.can_merge(rotate_inputs(inputs,iRot*rotation_skip),ns))
 486                                        {
 487                                           r.merge(rotate_inputs(inputs,iRot*rotation_skip));
 488                                           merged = true;
 489                                        }
 490                                        else if(do_reflect && r.can_merge(reflect_inputs(rotate_inputs(inputs,iRot*rotation_skip),nhood_size),ns))
 491                                        {
 492                                           r.merge(reflect_inputs(rotate_inputs(inputs,iRot*rotation_skip),nhood_size));
 493                                           merged = true;
 494                                        }
 495                                     }
 496                                  }
 497                                  if(!merged)
 498                                  {
 499                                     // need to make a new rule starting with this transition
 500                                     rule r(inputs,nhood_size,ns,symm);
 501                                     rules.push_back(r);
 502                                  }
 503                               }
 504                            }
 505                         }
 506                      }
 507                   }
 508                }
 509             }
 510          }
 511       }
 512       else // nhood_size==5
 513       {
 514          vector<state> inputs(5);
 515          inputs[0]=c;
 516          for(n=0;n<N;n++)
 517          {
 518             cout << ".";
 519             cout.flush();
 520             inputs[1]=n;
 521             for(e=0;e<N;e++)
 522             {
 523                inputs[2]=e;
 524                for(s=0;s<N;s++)
 525                {
 526                   inputs[3]=s;
 527                   for(w=0;w<N;w++)
 528                   {
 529                      ns = slowcalc(0,n,0,w,c,e,0,s,0);
 530                      if(remove_stasis && ns == c)
 531                         continue; // we can ignore stasis transitions
 532
 533                      // can we merge this transition with any existing rule?
 534                      inputs[4]=w;
 535                      merged = false;
 536                      for(it=rules.begin();!merged && it!=rules.end();it++)
 537                      {
 538                         rule &r = *it;
 539                         for(int iRot=0;!merged && iRot<n_rotations;iRot++)
 540                         {
 541                            if(r.can_merge(rotate_inputs(inputs,iRot*rotation_skip),ns))
 542                            {
 543                               r.merge(rotate_inputs(inputs,iRot*rotation_skip));
 544                               merged = true;
 545                            }
 546                            else if(do_reflect && r.can_merge(reflect_inputs(rotate_inputs(inputs,iRot*rotation_skip),nhood_size),ns))
 547                            {
 548                               r.merge(reflect_inputs(rotate_inputs(inputs,iRot*rotation_skip),nhood_size));
 549                               merged = true;
 550                            }
 551                         }
 552                      }
 553                      if(!merged)
 554                      {
 555                         // need to make a new rule starting with this transition
 556                         rule r(inputs,nhood_size,ns,symm);
 557                         rules.push_back(r);
 558                      }
 559                   }
 560                }
 561             }
 562          }
 563       }
 564    }
 565 }
 566
 567 // here we use the computed rule table as a replacement slowcalc, for checking correctness
 568 state new_slowcalc(const vector<rule>& rules,const vector<state>& inputs)
 569 {
 570    for(vector<rule>::const_iterator it=rules.begin();it!=rules.end();it++)
 571       if(it->matches(inputs))
 572          return it->ns;
 573     return inputs[0]; // default: no change
 574 }
 575
 576 bool is_correct(const vector<rule>&rules,int N,int neighbourhood_size)
 577 {
 578    // exhaustive check
 579    state c,n,ne,nw,sw,s,se,e,w;
 580    if(neighbourhood_size==9)
 581    {
 582       vector<state> inputs(9);
 583       for(c=0;c<N;c++)
 584       {
 585          inputs[0]=c;
 586          for(n=0;n<N;n++)
 587          {
 588             inputs[1]=n;
 589             for(ne=0;ne<N;ne++)
 590             {
 591                inputs[2]=ne;
 592                for(e=0;e<N;e++)
 593                {
 594                   inputs[3]=e;
 595                   for(se=0;se<N;se++)
 596                   {
 597                      inputs[4]=se;
 598                      for(s=0;s<N;s++)
 599                      {
 600                         inputs[5]=s;
 601                         for(sw=0;sw<N;sw++)
 602                         {
 603                            inputs[6]=sw;
 604                            for(w=0;w<N;w++)
 605                            {
 606                               inputs[7]=w;
 607                               for(nw=0;nw<N;nw++)
 608                               {
 609                                  inputs[8]=nw;
 610                                  if(new_slowcalc(rules,inputs)
 611                                     != slowcalc(nw,n,ne,w,c,e,sw,s,se))
 612                                        return false;
 613                               }
 614                            }
 615                         }
 616                      }
 617                   }
 618                }
 619             }
 620          }
 621       }
 622    }
 623    else
 624    {
 625       vector<state> inputs(5);
 626       for(c=0;c<N;c++)
 627       {
 628          inputs[0]=c;
 629          for(n=0;n<N;n++)
 630          {
 631             inputs[1]=n;
 632             for(e=0;e<N;e++)
 633             {
 634                inputs[2]=e;
 635                for(s=0;s<N;s++)
 636                {
 637                   inputs[3]=s;
 638                   for(w=0;w<N;w++)
 639                   {
 640                      inputs[4]=w;
 641                      if(new_slowcalc(rules,inputs)
 642                         != slowcalc(0,n,0,w,c,e,0,s,0))
 643                         return false;
 644                   }
 645                }
 646             }
 647          }
 648       }
 649    }
 650    return true;
 651 }
 652
 653 int main()
 654 {
 655    // parameters for use:
 656    const int N_STATES = 2;
 657    // Bedbugs: FIXME: output would probably be more efficient if we used
 658    // this, but also complicated
 659    const TSymm symmetry = none;
 660    const int nhood_size = 9;
 661    const string output_filename = "life.table";
 662    // Bedbugs: we need these because they will be munged into real
 663    // transitions:
 664    const bool remove_stasis_transitions = false;
 665
 666    vector<rule> rules;
 667    time_t t1,t2;
 668    time(&t1);
 669    produce_rule_table(rules,N_STATES,nhood_size,symmetry,remove_stasis_transitions);
 670    time(&t2);
 671    int n_secs = (int)difftime(t2,t1);
 672    cout << "\nProcessing took: " << n_secs << " seconds." << endl;
 673
 674    {
 675       ofstream out(output_filename.c_str());
 676       out << "# rules: " << rules.size() << "\n#";
 677       out << "\n# Golly rule-table format.\n# Each rule: C,";
 678       if(nhood_size==5)
 679          out << "N,E,S,W";
 680       else
 681          out << "N,NE,E,SE,S,SW,W,NW";
 682       out << ",C'";
 683       out << "\n# N.B. Where the same variable appears multiple times in a transition,\n# it takes the same value each time.\n#";
 684       if(remove_stasis_transitions)
 685          out << "\n# Default for transitions not listed: no change\n#";
 686       else
 687          out << "\n# All transitions should be included below, including no-change ones.\n#";
 688       out << "\nn_states:" << N_STATES;
 689       out << "\nneighborhood:" << ((nhood_size==5)?("vonNeumann"):("Moore"));
 690       out << "\nsymmetries:" << symmetry_strings[symmetry] << "\n";
 691       print_rules(rules,out);
 692    }
 693    cout << rules.size() << " rules written." << endl;
 694
 695    // optional: run through the entire state space, checking that new_slowcalc matches slowcalc
 696    cout << "Verifying is correct (can abort if you're confident)...";
 697    cout.flush();
 698    if(is_correct(rules,N_STATES,nhood_size))
 699       cout << "yes." << endl;
 700    else
 701       cout << "no! Either there's a bug in the code, or the transition function does not have the symmetry you selected." << endl;
 702 }