A (one dimensional) cellular automaton is a function1 F : Σ → Σ with the property that there is a K > 0 such that F (x)i depends only on the 2K + 1 coordinates xi−K , xi−K+1, . . . , xi−1, xi, xi+1, . . . , xi+K . A periodic point of σ is any x such that σ^p (x) = x for some p ∈ N, and a periodic point of F is any x such that F^q (x) = x for some q ∈ N. Given a cellular automaton F, a point x ∈ Σ is jointly periodic if there are p, q ∈ N such that σ^p (x) = F^q (x) = x, that is, it is a periodic point under both functions.
This project aims to explore the nature of one-dimensional Cellular Automata, in the hope of finding the structure of cellular automata through its periodic points.
License: MIT
ubuntu2004
/*
* Copyright (C) 2004 Bryant Lee
*
* This file is part of FPeriod.
*
* FPeriod is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* FPeriod is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with FPeriod; if not, write to the Free Software
* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
/* Main function
* FPeriod program finds points that are periodic under the block code
* induced by F
*
* Written by: Bryant Lee
* Date: 10/29/04
*/
#define USEARRAY 0
#define DEBUG 0
//Standard C++
#include <string.h>
#include <iostream>
#include <iomanip>
#include <fstream> //for file ops
//C
#include <stdlib.h> //for atoi
#include <math.h> //for pow
//STL
#include <vector>
#include <list>
#include <map>
#include <hash_map>
//local
#include "StringOps.h"
#include "Comp.h"
#include "StorageKey.h"
#include "StorageVal.h"
#include "RecNode.h"
//function definitions
void convertPeriodList(const string &strInput, vector<int> &periodList);
void exhaustiveExp(unsigned int shiftSize, const vector<int> &periodList,
string &prefix, Comp *comp, bool btrunc,
list<RecNode*> &recList,
list<RecNode*> &leastRecList,
vector<string> &commandbuffer);
void exhaustiveExpTwoShift(const vector<int> &periodList,
string &prefix, Comp *comp, bool btrunc,
list<RecNode*> &recList,
list<RecNode*> &leastRecList,
vector<string> &commandbuffer);
void clearStorageNodes(vector<StorageVal*> &table);
void clearRecList(list<RecNode*> &recList);
void output(ofstream & fout, ofstream &fout_preperiod, Comp * comp,
list<RecNode*> &recList, list<RecNode*> &leastRecList);
void truncoutput(ofstream & fout, Comp * comp, list<RecNode*> &leastRecList);
void findDivisors(unsigned int x, vector<int> &divisorsList);
bool leastPeriod(byte *word, unsigned int currPeriod,
const vector<int> &divisorsList);
bool leastPeriod(unsigned int word, unsigned int currPeriod,
const vector<int> &divisorsList);
void outcontrol(string &prefix,
Comp * comp,
bool btrunc,
list<RecNode*> &recList,
list<RecNode*> &leastRecList,
vector<string> &commandbuffer);
/*
* GLOBAL variables
*/
int UINTSIZE = 32;
//boolean that determines whether output is given after each k
//or only at end
bool OUTPUT_AFTER_EACH = false;
//main
int main() {
unsigned int shiftSize; //# symbols in shift
string prefix; //prefix for output files
string trunc; //"trunc" or "no trunc": determine whether to truncate output
bool btrunc; //true if truncate output is on
string strInput, strInput2, strInput3; //temp variables
vector<int> periodList; //list of periods
unsigned int i = 0, temp; //counter
list<RecNode*> recList; //list of records
list<RecNode*> leastRecList; //list of records for points of least period k
Comp *comp = NULL; //composition to apply
string outoption;
//hold the command strings
vector<string> commandbuffer;
int commandbufferit;
//grab all the input from cin and put in an array of strings
while(!cin.eof()) {
std::getline(cin, strInput);
commandbuffer.push_back(strInput);
}
commandbufferit = 0;
//get "output after each" or "output at end"
outoption = commandbuffer[commandbufferit++];
trim(outoption);
transform(outoption.begin(), outoption.end(), outoption.begin(), tolower);
if(outoption == "output after each")
OUTPUT_AFTER_EACH = true;
else if(outoption == "output at end")
OUTPUT_AFTER_EACH = false;
else {
cout << "Error on line 1: invalid value for output mode\n";
goto cleanup;
}
//get truncate value
trunc = commandbuffer[commandbufferit++];
trim(trunc);
transform(trunc.begin(), trunc.end(), trunc.begin(), tolower);
if(trunc == "trunc")
btrunc = true;
else if(trunc == "no trunc")
btrunc = false;
else {
cout << "Error on line 2: invalid value for truncate\n";
goto cleanup;
}
//get filename prefix
prefix = commandbuffer[commandbufferit++];
trim(prefix);
//get N, the number of symbols in the full shift
strInput = commandbuffer[commandbufferit++];
trim(strInput);
shiftSize = atoi(strInput.c_str());
//get set of K, the periodicities to be tested
strInput = commandbuffer[commandbufferit++];
convertPeriodList(strInput, periodList);
//Allocate composition object
comp = new Comp(shiftSize);
//get "Begin Definitions"
strInput = commandbuffer[commandbufferit++];
trim(strInput);
//make it lowercase
transform(strInput.begin(), strInput.end(), strInput.begin(), tolower);
if(strInput.compare("begin definitions") != 0) {
cout << "ERROR: Expecting \"BEGIN DEFINITIONS\"\n";
goto cleanup;
}
//read definitions
while(1) {
strInput = commandbuffer[commandbufferit++];
//strInput2 is lowercased and trimmed
strInput2.resize(strInput.length());
transform(strInput.begin(), strInput.end(), strInput2.begin(), tolower);
trim(strInput2);
//strInput3 is lowercased and NOT TRIMMED
strInput3.resize(strInput.length());
transform(strInput.begin(), strInput.end(), strInput3.begin(), tolower);
if(strInput2 == "" || (strInput2[0] == '/' && strInput2[1] == '/')) {
//ignore blank lines and allow comment out
}
else if(strInput2 == "begin commands") {
break; //BREAK
}
else {
i = 0;
while(i < strInput.length()) {
if(strInput[i] == '=') {
break;
}
i++;
}
if(strInput2[0] == 'o' && strInput2[1] == 'p' && strInput2[2] == 't') {
temp = strInput3.find("opt");
temp += 3;
comp->addFunc(strInput.substr(i+1, strInput.length() - (i+1)),
strInput.substr(temp, i-temp), true);
}
else {
comp->addFunc(strInput.substr(i+1, strInput.length() - (i+1)),
strInput.substr(0,i), false);
}
}
}
//read commands
while(commandbufferit < (int) commandbuffer.size()) { //read until end
strInput = commandbuffer[commandbufferit++];
trim(strInput);
//ignore blank lines and allow comment out with "//"
if(strInput != "" && !(strInput[0] == '/' && strInput[1] == '/')) {
//set composition
if(!comp->setComposition(strInput)) {
cout << "ERROR: Command \"" << strInput
<< "\" references undefined function\n";
goto cleanup;
}
//check maxshiftperiod
int maxshiftperiod = 0;
for(i = 0; i < periodList.size(); i++)
maxshiftperiod = (maxshiftperiod < periodList[i] ?
periodList[i] : maxshiftperiod);
//do experiment
if(shiftSize == 2
&& sizeof(unsigned int) == 4
&& !USEARRAY
&& maxshiftperiod <= 32) //if shift period goes > 32 use array
exhaustiveExpTwoShift(periodList, prefix, comp, btrunc,
recList, leastRecList,
commandbuffer);
else
exhaustiveExp(shiftSize, periodList,
prefix, comp, btrunc, recList, leastRecList,
commandbuffer);
if(!OUTPUT_AFTER_EACH)
outcontrol(prefix, comp, btrunc,
recList, leastRecList, commandbuffer);
//clear records (delete the nodes with ptrs in lists)
clearRecList(recList);
clearRecList(leastRecList);
recList.clear();
leastRecList.clear();
}
}
//cleanup and exit
cleanup:
//delete the nodes with ptrs in lists
clearRecList(recList);
clearRecList(leastRecList);
//delete composition (and the functions in the Comp object's storage)
delete comp;
return 0;
}
/*
* This function does the output. It chooses between full output
* and truncated output accordingly.
*/
void outcontrol(string &prefix,
Comp * comp,
bool btrunc,
list<RecNode*> &recList,
list<RecNode*> &leastRecList,
vector<string> &commandbuffer) {
unsigned int commandbufferit=0;
string test;
ofstream fout;
string fname = prefix + comp->returnName() + ".txt";
fout.open(fname.c_str(), ios::trunc | ios::out);
if(fout.is_open()) {
//HEADER
fout << "FPeriod 3.67\n\n";
//copy of the input
fout << "INPUT\n";
fout << "--------------------\n";
for(commandbufferit = 0; commandbufferit < commandbuffer.size(); commandbufferit++) {
test = commandbuffer[commandbufferit];
trim(test);
if(test.compare("") != 0) {
fout << commandbuffer[commandbufferit] << "\n";
}
}
fout << "--------------------\n";
fout << "\n";
if(btrunc)
truncoutput(fout, comp, leastRecList);
else {
ofstream fout_preperiod; //file for preperiods
string fname_preperiod = prefix + comp->returnName() +
"_preperiods.txt";
fout_preperiod.open(fname_preperiod.c_str(), ios::trunc | ios::out);
if(fout_preperiod.is_open()) {
//HEADER
fout_preperiod << "FPeriod 3.67\n\n";
//copy of the input
fout_preperiod << "INPUT\n";
fout_preperiod << "--------------------\n";
for(commandbufferit = 0; commandbufferit < commandbuffer.size(); commandbufferit++) {
test = commandbuffer[commandbufferit];
trim(test);
if(test.compare("") != 0) {
fout_preperiod << commandbuffer[commandbufferit] << "\n";
}
}
fout_preperiod << "--------------------\n";
fout_preperiod << "\n";
output(fout, fout_preperiod, comp, recList, leastRecList);
fout_preperiod.close();
}
else {
cout << "ERROR: Could not create file " << fname_preperiod << "\n";
}
}
fout.close();
}
else {
cout << "ERROR: Could not create file " << fname << "\n";
}
}
/*
* Truncated output
*
* Outputs only a table for least period k with the column P^(1/k) i.e.
*
* K P^(1/k)
* 1 1.000
* 2 1.4112
* 3 1.4432
* etc.
*/
void truncoutput(ofstream & fout, Comp * comp, list<RecNode*> &leastRecList) {
list<RecNode*>::const_iterator it;
fout << "Truncated output\n\n";
//Output
fout << "COMMAND:\n";
fout << comp->returnName() << "\n";
fout << "WITH DEFINITIONS:\n";
comp->printDefinitions(fout);
fout << "\n";
//Summary info
fout << "LEAST PERIOD K\n\n";
fout << "K P^(1/k)\n";
//fill on right side, show in fixed notation
fout << setiosflags(ios::left | ios::fixed);
for(it = leastRecList.begin(); it != leastRecList.end(); it++) {
fout << setw(3) << (*it)->currPeriod << " ";
if((*it)->pending) {
fout << "PENDING" << "\n";
}
else {
fout << setw(8) << setprecision(4) << pow((*it)->numPeriodic, ((double)1)/((double)(*it)->currPeriod)) << "\n";
}
}
}
/*
* Full output
*/
void output(ofstream & fout, ofstream & fout_preperiod,
Comp * comp, list<RecNode*> &recList,
list<RecNode*> &leastRecList) {
list<RecNode*>::const_iterator it;
/* Output file */
//fill on right side, show in fixed notation
fout << setiosflags(ios::left | ios::fixed);
fout << "Full output\n\n";
fout << "COMMAND:\n";
fout << comp->returnName() << "\n";
fout << "WITH DEFINITIONS:\n";
comp->printDefinitions(fout);
fout << "\n";
/* Preperiod file */
//fill on right side, show in fixed notation
fout_preperiod << setiosflags(ios::left | ios::fixed);
fout_preperiod << "\nPreperiod output\n\n";
fout_preperiod << "COMMAND:\n";
fout_preperiod << comp->returnName() << "\n";
fout_preperiod << "WITH DEFINITIONS:\n";
comp->printDefinitions(fout_preperiod);
fout_preperiod << "\n";
/* output file */
fout << "PERIOD K (NOT LEAST)\n\n";
/* preperiod file */
fout_preperiod << "PERIOD K (NOT LEAST)\n\n";
/* output file */
//Summary info
fout << "K FracPdc P^(1/k) L^(1/k) P L Not-P AvgPd AvgPrepd MaxPrepd\n";
for(it = recList.begin(); it != recList.end(); it++) {
fout << setw(3) << (*it)->currPeriod << " ";
if((*it)->pending) {
fout << "PENDING" << "\n";
}
else {
fout << setw(8) << setprecision(6) << (*it)->fractionPeriodic() << " ";
fout << setw(8) << setprecision(4) << pow((*it)->numPeriodic, ((double)1)/((double)(*it)->currPeriod)) << " ";
fout << setw(8) << setprecision(4) << pow((*it)->largestorbit, ((double)1)/((double)(*it)->currPeriod)) << " ";
fout << setw(9) << (*it)->numPeriodic << " ";
fout << setw(9) << (*it)->largestorbit << " ";
fout << setw(9) << (*it)->numNonPeriodic() << " ";
fout << setw(12) << setprecision(2) << (*it)->avgPeriod() << " ";
fout << setw(12) << setprecision(2) << (*it)->avgPreperiod() << " ";
fout << setw(9) << (*it)->maxPreperiod() <<"\n";
}
}
//orbit and period multiplicity
fout << "\n";
fout << "K Period Mult(#orbits) Mult(#perpoints) Mult(#points)\n";
for(it = recList.begin(); it != recList.end(); it++) {
fout << setw(3) << (*it)->currPeriod << " ";
if((*it)->pending) {
fout << "PENDING" << "\n";
}
else {
map<unsigned long long, unsigned long long>::const_iterator pitorbit, pitpoint, pitperpoint;
bool pfirst = true;
//note that the period and orbit maps have the same periods in the same orders
// so we can iterate through both at once
for(pitorbit = (*it)->orbitMap.begin(),pitpoint = (*it)->periodMap.begin();
pitorbit != (*it)->orbitMap.end();
pitorbit++, pitpoint++) {
if(!pfirst) {
fout << " ";
}
else {
pfirst = false;
}
fout << setw(11) << pitorbit->first << " ";
fout << setw(13) << pitorbit->second << " "; //orbits
if((pitperpoint = (*it)->periodicPeriodMap.find(pitorbit->first))
!= (map<unsigned long long, unsigned long long>::const_iterator) (*it)->periodicPeriodMap.end()) {
fout << setw(16) << pitperpoint->second << " "; //periodic points
}
else {
fout << setw(16) << 0 << " "; //periodic points
}
fout << setw(13) << pitpoint->second << "\n"; //points
}
}
}
/* Preperiod file */
//Preperiod multiplicity
fout_preperiod << "K Preperiod Mult\n";
for(it = recList.begin(); it != recList.end(); it++) {
fout_preperiod << setw(3) << (*it)->currPeriod << " ";
if((*it)->pending) {
fout_preperiod << "PENDING" << "\n";
}
else {
map<unsigned long long, unsigned long long>::const_iterator pit;
bool pfirst = true;
for(pit = (*it)->preperiodMap.begin(); pit != (*it)->preperiodMap.end();
pit++) {
if(!pfirst) {
fout_preperiod << " ";
}
else {
pfirst = false;
}
fout_preperiod << setw(11) << pit->first << " ";
fout_preperiod << setw(13) << pit->second << "\n";
}
}
}
/* output file */
fout << "\n";
fout << "LEAST PERIOD K\n\n";
/* preperiod file */
fout_preperiod << "\n";
fout_preperiod << "LEAST PERIOD K\n\n";
/* output file */
//Summary info
fout << "K FracPdc P^(1/k) L^(1/k) P L Not-P AvgPd AvgPrepd MaxPrepd\n";
for(it = leastRecList.begin(); it != leastRecList.end(); it++) {
fout << setw(3) << (*it)->currPeriod << " ";
if((*it)->pending) {
fout << "PENDING" << "\n";
}
else {
fout << setw(8) << setprecision(6) << (*it)->fractionPeriodic() << " ";
fout << setw(8) << setprecision(4) << pow((*it)->numPeriodic, ((double)1)/((double)(*it)->currPeriod)) << " ";
fout << setw(8) << setprecision(4) << pow((*it)->largestorbit, ((double)1)/((double)(*it)->currPeriod)) << " ";
fout << setw(9) << (*it)->numPeriodic << " ";
fout << setw(9) << (*it)->largestorbit << " ";
fout << setw(9) << (*it)->numNonPeriodic() << " ";
fout << setw(12) << setprecision(2) << (*it)->avgPeriod() << " ";
fout << setw(12) << setprecision(2) << (*it)->avgPreperiod() << " ";
fout << setw(9) << (*it)->maxPreperiod() <<"\n";
}
}
//orbit and period multiplicity
fout << "\n";
fout << "K Period Mult(#perpoints) Mult(#points)\n";
for(it = leastRecList.begin(); it != leastRecList.end(); it++) {
fout << setw(3) << (*it)->currPeriod << " ";
if((*it)->pending) {
fout << "PENDING" << "\n";
}
else {
map<unsigned long long, unsigned long long>::const_iterator pit, pitperpoint;
bool pfirst = true;
for(pit = (*it)->periodMap.begin();
pit != (*it)->periodMap.end();
pit++) {
if(!pfirst) {
fout << " ";
}
else {
pfirst = false;
}
fout << setw(11) << pit->first << " ";
if((pitperpoint = (*it)->periodicPeriodMap.find(pit->first))
!= (map<unsigned long long, unsigned long long>::const_iterator) (*it)->periodicPeriodMap.end()) {
fout << setw(16) << pitperpoint->second << " "; //periodic points
}
else {
fout << setw(16) << 0 << " "; //periodic points
}
fout << setw(13) << pit->second << "\n"; //points
}
}
}
/* Preperiod file */
//Preperiod multiplicity
fout_preperiod << "K Preperiod Mult\n";
for(it = leastRecList.begin(); it != leastRecList.end(); it++) {
fout_preperiod << setw(3) << (*it)->currPeriod << " ";
if((*it)->pending) {
fout_preperiod << "PENDING" << "\n";
}
else {
map<unsigned long long, unsigned long long>::const_iterator pit;
bool pfirst = true;
for(pit = (*it)->preperiodMap.begin(); pit != (*it)->preperiodMap.end();
pit++) {
if(!pfirst) {
fout_preperiod << " ";
}
else {
pfirst = false;
}
fout_preperiod << setw(11) << pit->first << " ";
fout_preperiod << setw(13) << pit->second << "\n";
}
}
}
//DONE
}
void exhaustiveExp(unsigned int shiftSize, const vector<int> &periodList,
string &prefix, Comp *comp, bool btrunc,
list<RecNode*> &recList,
list<RecNode*> &leastRecList,
vector<string> &commandbuffer) {
unsigned int i = 0, j = 0; //k = 0; //counters
unsigned int currPeriod;
list<StorageVal*> orbitList;
list<bool> leastPeriodList;
list<RecNode*>::const_iterator recIterator, leastRecIterator;
vector<StorageVal*> table;
vector<int> divisorsList;
if(DEBUG)
cout << "Array" << endl;
//make a list of recNodes and leastRecNodes
//each list has one node per K value
for(i = 0; i < periodList.size(); i++) {
recList.push_back(new RecNode(periodList[i], pow(shiftSize,periodList[i])));
leastRecList.push_back(new RecNode(periodList[i], 0));
}
//initialize recIterator and leastRecIterator
recIterator = recList.begin();
leastRecIterator = leastRecList.begin();
for(i = 0; i < periodList.size(); i++) {
currPeriod = periodList[i];
byte itWord[currPeriod]; //iterating word
//set table size
table.resize(pow(shiftSize,currPeriod));
//record keeping
RecNode *record = *recIterator;
RecNode *leastRecord = *leastRecIterator;
//zero the itWord
for(j = 0; j < currPeriod; j++) {
itWord[j] = 0;
}
//Find divisors
findDivisors(currPeriod, divisorsList);
while(true) {
byte cWord[currPeriod]; //curr word
byte output[currPeriod]; //helper word
unsigned long long period = 0, preperiod = 0;
StorageVal *cStore;
StorageVal *rt;
bool halt = false;
while(halt != true &&
table[StorageKey(itWord,currPeriod).hash(shiftSize)] != NULL)
{
halt = iterateWord(itWord, currPeriod, shiftSize);
}
if(halt == true)
break; //exit loop
//Set current word
copyWord(cWord, itWord, currPeriod);
cStore = new StorageVal(0);
StorageKey sk1 = StorageKey(cWord, currPeriod);
table[sk1.hash(shiftSize)] = cStore;
orbitList.push_back(cStore);
if(leastPeriod(cWord, currPeriod, divisorsList))
leastPeriodList.push_back(true);
else
leastPeriodList.push_back(false);
//DEBUG
//cout << "-- Start orbit --\n";
//printArray(cWord, currPeriod);
//Perform experiment
j = 1;
while(true) {
comp->image(cWord, currPeriod, output);
copyWord(cWord, output, currPeriod);
//DEBUG
//printArray(cWord, currPeriod);
StorageKey sk2 = StorageKey(cWord, currPeriod);
rt = table[sk2.hash(shiftSize)];
if(rt != NULL) {
break;
}
else {
cStore = new StorageVal(j);
table[sk2.hash(shiftSize)] = cStore;
orbitList.push_back(cStore);
if(leastPeriod(cWord, currPeriod, divisorsList))
leastPeriodList.push_back(true);
else
leastPeriodList.push_back(false);
}
j++;
}
//DEBUG
//cout << "-- Returned to seen --\n";
if(rt->imageNum == -1) {
//we hit a previously found point
period = rt->period;
preperiod = rt->preperiod + j;
list<StorageVal*>::iterator it;
list<bool>::iterator lpit = leastPeriodList.begin();
for(it = orbitList.begin(); it != orbitList.end(); it++) {
(*it)->imageNum = -1;
(*it)->period = period;
(*it)->preperiod = preperiod;
if(*lpit == true) {
leastRecord->recordPeriod(period, preperiod);
leastRecord->numTrials++;
}
lpit++;
record->recordPeriod(period,preperiod);
preperiod--;
//NOTE: none of these points can be periodic
}
}
else {
//we circled back to an active point
period = j - rt->imageNum;
preperiod = rt->imageNum;
record->recordOrbit(period);
list<StorageVal*>::iterator it;
list<bool>::iterator lpit = leastPeriodList.begin();
for(it = orbitList.begin(); it != orbitList.end(); it++) {
(*it)->imageNum = -1;
(*it)->period = period;
(*it)->preperiod = preperiod;
if(*lpit == true) {
leastRecord->recordPeriod(period, preperiod);
leastRecord->numTrials++;
if(preperiod == 0) {
leastRecord->recordPeriodicPeriod(period);
leastRecord->numPeriodic++;
}
}
lpit++;
record->recordPeriod(period, preperiod);
if(preperiod > 0) {
preperiod--;
}
else {
record->recordPeriodicPeriod(period);
record->numPeriodic++;
}
}
}
//cout << "Clearing list\n";
orbitList.clear();
leastPeriodList.clear();
}
record->pending = false;
leastRecord->pending = false;
//advance recIterator and leastRecIterator
recIterator++;
leastRecIterator++;
//cout << "Clearing table\n";
clearStorageNodes(table);
table.clear();
if(OUTPUT_AFTER_EACH) {
outcontrol(prefix, comp, btrunc,
recList, leastRecList, commandbuffer);
}
}
}
void exhaustiveExpTwoShift(const vector<int> &periodList,
string &prefix, Comp *comp, bool btrunc,
list<RecNode*> &recList,
list<RecNode*> &leastRecList,
vector<string> &commandbuffer) {
unsigned int i = 0, j = 0; //k = 0; //counters
unsigned int currPeriod;
unsigned long long numTrials;
list<StorageVal*> orbitList;
list<bool> leastPeriodList;
vector<StorageVal*> table;
vector<int> divisorsList;
list<RecNode*>::const_iterator recIterator, leastRecIterator;
//DEBUG
if(DEBUG)
cout << "Special case 2-shift" << endl;
//make a list of recNodes and leastRecNodes
for(i = 0; i < periodList.size(); i++) {
recList.push_back(new RecNode(periodList[i], pow(2, periodList[i])));
leastRecList.push_back(new RecNode(periodList[i], 0));
}
//initialize recIterator
recIterator = recList.begin();
leastRecIterator = leastRecList.begin();
for(i = 0; i < periodList.size(); i++) {
currPeriod = periodList[i];
unsigned int itWord = 0; //iterating word
numTrials = (unsigned long long) pow(2, periodList[i]);
table.resize(numTrials);
RecNode *record = *recIterator;
RecNode *leastRecord = *leastRecIterator;
//Find divisors
findDivisors(currPeriod, divisorsList);
while(true) {
unsigned int cWord = 0; //curr word
unsigned int output = 0; //helper word
unsigned long long period = 0, preperiod = 0;
StorageVal *cStore;
StorageVal* rt;
while(itWord != numTrials &&
table[itWord] != NULL) {
itWord++;
}
if(itWord == numTrials)
break; //exit loop
//Set current word
cWord = itWord;
cStore = new StorageVal(0);
table[cWord] = cStore;
orbitList.push_back(cStore);
if(leastPeriod(cWord, currPeriod, divisorsList))
leastPeriodList.push_back(true);
else
leastPeriodList.push_back(false);
//DEBUG
//cout << "-- start orbit --\n";
//printUInt(cWord, 4);
//Perform experiment
j = 1;
while(true) {
comp->image(&cWord, currPeriod, &output, 1);
cWord = output;
//DEBUG
//printUInt(cWord, 4);
rt = table[cWord];
if(rt != NULL) {
break;
}
else {
cStore = new StorageVal(j);
table[cWord] = cStore;
orbitList.push_back(cStore);
if(leastPeriod(cWord, currPeriod, divisorsList))
leastPeriodList.push_back(true);
else
leastPeriodList.push_back(false);
}
j++;
}
//DEBUG
//cout << "-- end orbit --\n";
if(rt->imageNum == -1) {
//we hit a previously found point
period = rt->period;
preperiod = rt->preperiod + j;
list<StorageVal*>::iterator it;
list<bool>::iterator lpit = leastPeriodList.begin();
for(it = orbitList.begin(); it != orbitList.end(); it++) {
(*it)->imageNum = -1;
(*it)->period = period;
(*it)->preperiod = preperiod;
if(*lpit == true) {
leastRecord->recordPeriod(period, preperiod);
leastRecord->numTrials++;
}
lpit++;
record->recordPeriod(period,preperiod);
preperiod--;
//NOTE: none of these points can be periodic
}
}
else {
//we circled back to an active point
period = j - rt->imageNum;
preperiod = rt->imageNum;
record->recordOrbit(period);
list<StorageVal*>::iterator it;
list<bool>::iterator lpit = leastPeriodList.begin();
for(it = orbitList.begin(); it != orbitList.end(); it++) {
(*it)->imageNum = -1;
(*it)->period = period;
(*it)->preperiod = preperiod;
if(*lpit == true) {
leastRecord->recordPeriod(period, preperiod);
leastRecord->numTrials++;
if(preperiod == 0) {
leastRecord->recordPeriodicPeriod(period);
leastRecord->numPeriodic++;
}
}
lpit++;
record->recordPeriod(period, preperiod);
if(preperiod > 0) {
preperiod--;
}
else {
record->recordPeriodicPeriod(period);
record->numPeriodic++;
}
}
}
//cout << "Clearing list\n";
orbitList.clear();
leastPeriodList.clear();
}
record->pending = false;
leastRecord->pending = false;
//advance recIterator and leastRecIterator
recIterator++;
leastRecIterator++;
//cout << "Clearing table\n";
clearStorageNodes(table);
table.clear();
if(OUTPUT_AFTER_EACH) {
outcontrol(prefix, comp, btrunc,
recList, leastRecList, commandbuffer);
}
}
}
void convertPeriodList(const string &strInput,
vector<int> &periodList) {
vector<string> strList; //all terms
vector<string> subList; //vector representation of an "[a,b]" term
unsigned int strListLength = 0, i = 0;
split(strInput,';', strList);
strListLength = strList.size();
for(i = 0; i < strListLength; i++) {
string *temp = &(strList[i]);
trim(*temp);
if((*temp)[0] == '[') {
int j, low, high;
temp->erase(0,1); //remove '['
temp->erase(temp->length() - 1, 1); //remove ']'
split(*temp, ',', subList);
trim(subList[0]);
trim(subList[1]);
low = atoi(subList[0].c_str());
high = atoi(subList[1].c_str());
for(j = low; j <= high; j++) {
periodList.push_back(j);
}
subList.clear();
}
else {
periodList.push_back(atoi(temp->c_str()));
}
}
}
/*
* Delete the nodes pointed to by ptrs in recList
*/
void clearRecList(list<RecNode*> &recList) {
list<RecNode*>::iterator it;
for(it = recList.begin(); it != recList.end(); it++) {
delete (*it);
}
}
/*
* Delete the nodes pointed to by the values in the table
*/
void clearStorageNodes(vector<StorageVal*> &table) {
vector<StorageVal*>::iterator it;
for(it = table.begin(); it != table.end(); it++) {
delete (*it);
}
}
/*
* Find divisors of X, excluding X itself
*/
void findDivisors(unsigned int x, vector<int> &divisorsList) {
unsigned int i = 2;
unsigned int halfX = x/2;
divisorsList.clear();
if(x > 1) {
divisorsList.push_back(1);
}
for(i = 2; i <= halfX; i++) {
if(x % i == 0) {
divisorsList.push_back(i);
}
}
}
/*
* returns true if this word has least period p
*/
bool leastPeriod(byte *word, unsigned int currPeriod,
const vector<int> &divisorsList) {
bool result = false;
unsigned int i, j, k;
if(currPeriod == 1)
return true;
//for every divisor d
for(i = 0; i < divisorsList.size(); i++) {
unsigned int d = divisorsList[i];
result = false;
//check in chunks of length d
for(j = 0; j < currPeriod; j+= d) {
//check that this chunk is same as other chunks
for(k = 0; k < d; k++) {
if(word[k + j] != word[k]) {
result = true;
break;
}
}
if(result == true)
break;
}
if(result == false)
break;
}
return result;
}
/*
* returns true if this word has least period p
*/
bool leastPeriod(unsigned int word, unsigned int currPeriod,
const vector<int> &divisorsList) {
bool result = false;
unsigned int i, j, k;
if(currPeriod == 1)
return true;
//for every divisor d
for(i = 0; i < divisorsList.size(); i++) {
unsigned int d = divisorsList[i];
result = false;
//check in chunks of length d
for(j = 0; j < currPeriod; j+= d) {
//check that this chunk is same as other chunks
for(k = 0; k < d; k++) {
if(((word >> (k + j)) & 1) != ((word >> k) & 1)) {
result = true;
break;
}
}
if(result == true)
break;
}
if(result == false)
break;
}
return result;
}