|
|
|
#include "Attractor.hpp"
|
|
|
|
|
|
|
|
#include <iostream>
|
|
|
|
#include <fstream>
|
|
|
|
#include <cstdlib>
|
|
|
|
#include <string>
|
|
|
|
#include "stfu/stf.hpp"
|
|
|
|
|
|
|
|
#include "Projector.hpp"
|
|
|
|
|
|
|
|
#include "AttractorKernel.hpp"
|
|
|
|
#include "kernels/Logistic.hpp"
|
|
|
|
#include "kernels/Lorenz3D.hpp"
|
|
|
|
#include "kernels/Polynomial.hpp"
|
|
|
|
#include "kernels/PolynomialA3D.hpp"
|
|
|
|
#include "kernels/Unravel3D.hpp"
|
|
|
|
|
|
|
|
Attractor::Attractor() {
|
|
|
|
myAttractor = new Lorenz3D();
|
|
|
|
}
|
|
|
|
|
|
|
|
Attractor::Attractor(const char* const fileName) {
|
|
|
|
// opening file
|
|
|
|
std::cout << "Reading file " << fileName << "..." << std::endl;
|
|
|
|
|
|
|
|
stfu::node system;
|
|
|
|
system.read(fileName);
|
|
|
|
stfu::node attractor = system.getChild("attractor");
|
|
|
|
|
|
|
|
|
|
|
|
// reading basic stuff
|
|
|
|
std::string attractorType = attractor.getValue("type");
|
|
|
|
const std::string attractorDimension = attractor.getValue("dimensions");
|
|
|
|
|
|
|
|
for ( unsigned int i = 0; attractorType[i] != '\0'; i++ ) {
|
|
|
|
attractorType[i] = tolower(attractorType[i]);
|
|
|
|
}
|
|
|
|
const unsigned int dimension = atoi(attractorDimension.c_str());
|
|
|
|
|
|
|
|
std::cout << " Formula: " << attractorType << std::endl;
|
|
|
|
std::cout << " Dimensions: " << dimension << std::endl;
|
|
|
|
|
|
|
|
|
|
|
|
// depending on type, make the formula object
|
|
|
|
if ( attractorType == "lorenz" ){
|
|
|
|
if ( dimension == 3 ) {
|
|
|
|
myAttractor = new Lorenz3D();
|
|
|
|
} else {
|
|
|
|
std::cerr << "something wrong";
|
|
|
|
exit(37);
|
|
|
|
}
|
|
|
|
} else if ( attractorType == "polynomial" ) {
|
|
|
|
const std::string attractorOrde = attractor.getValue("orde");
|
|
|
|
const unsigned int orde = atoi(attractorOrde.c_str());
|
|
|
|
std::cout << " Orde: " << orde << std::endl;
|
|
|
|
myAttractor = new Polynomial(dimension, orde);
|
|
|
|
|
|
|
|
} else if ( attractorType == "polynomial a" ) {
|
|
|
|
if ( dimension == 3 ) {
|
|
|
|
myAttractor = new PolynomialA3D();
|
|
|
|
} else {
|
|
|
|
std::cerr << "something wrong";
|
|
|
|
exit(37);
|
|
|
|
}
|
|
|
|
} else if ( attractorType == "logistic" ) {
|
|
|
|
myAttractor = new Logistic(dimension);
|
|
|
|
|
|
|
|
} else if ( attractorType == "unravel" ) {
|
|
|
|
if ( dimension == 3 ) {
|
|
|
|
myAttractor = new Unravel3D();
|
|
|
|
} else {
|
|
|
|
std::cerr << "somtheing wrong";
|
|
|
|
exit(37);
|
|
|
|
}
|
|
|
|
} else {
|
|
|
|
std::cout << "'" << attractorType << "' not recognized" << std::endl;
|
|
|
|
exit(3);
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
// read parameters
|
|
|
|
const unsigned int numberOfParameters = myAttractor->getNumberOfParameters();
|
|
|
|
|
|
|
|
for ( unsigned int i = 0; i < numberOfParameters; i++ ) {
|
|
|
|
stfu::node attractorParameters = attractor.getChild("parameters");
|
|
|
|
(*myAttractor)[i] = atof(attractorParameters.getValue(i).c_str());
|
|
|
|
std::cout << " Parameter " << i << " set to " << (*myAttractor)[i] << ", ";
|
|
|
|
}
|
|
|
|
|
|
|
|
std::cout << std::endl << " Reading file complete" << std::endl;
|
|
|
|
}
|
|
|
|
|
|
|
|
Attractor::~Attractor(){
|
|
|
|
delete myAttractor;
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
// this should probably done in the projector section
|
|
|
|
void Attractor::init_range() {
|
|
|
|
// stabilize attractor
|
|
|
|
for ( unsigned int i = 0; i < 100000; i++ ) {
|
|
|
|
iterate();
|
|
|
|
}
|
|
|
|
|
|
|
|
// initialize projectors with dimension and first point
|
|
|
|
const unsigned int dimension = myAttractor->getDimension();
|
|
|
|
const double * point = myAttractor->vector();
|
|
|
|
for ( std::vector<Projector*>::iterator it = projectors.begin(); it != projectors.end(); it++ ) {
|
|
|
|
(*it)->extern_dim = dimension;
|
|
|
|
(*it)->intern_dim = 2;
|
|
|
|
(*it)->init(point);
|
|
|
|
}
|
|
|
|
|
|
|
|
// update ranges
|
|
|
|
for ( unsigned int i = 0; i < 100000; i++ ) {
|
|
|
|
iterate();
|
|
|
|
|
|
|
|
for ( std::vector<Projector*>::iterator it = projectors.begin(); it != projectors.end(); it++ ) {
|
|
|
|
(*it)->update_range(point);
|
|
|
|
}
|
|
|
|
}
|
|
|
|
for ( std::vector<Projector*>::iterator it = projectors.begin(); it != projectors.end(); it++ ) {
|
|
|
|
(*it)->finish_range();
|
|
|
|
}
|
|
|
|
|
|
|
|
}
|
|
|
|
|
|
|
|
bool Attractor::is_chaos() {
|
|
|
|
/*
|
|
|
|
check existence of attractor:
|
|
|
|
Escaping
|
|
|
|
Single point attractor
|
|
|
|
Lyapunov exponent
|
|
|
|
*/
|
|
|
|
/*
|
|
|
|
double sum = 0;
|
|
|
|
for ( unsigned int i = 0; i < dim; i++ ) {
|
|
|
|
const double dist = 0; //new_point[i] - point[i];
|
|
|
|
sum += dist*dist;
|
|
|
|
}
|
|
|
|
if ( sum >= 1.0e7 ) {
|
|
|
|
// big change => Escaping
|
|
|
|
return false;
|
|
|
|
}
|
|
|
|
if ( sum <= 1.0e-7 ) {
|
|
|
|
// small change => singularity
|
|
|
|
return false;
|
|
|
|
}
|
|
|
|
return true;
|
|
|
|
*/
|
|
|
|
return true;
|
|
|
|
}
|
|
|
|
|
|
|
|
void Attractor::iterate() {
|
|
|
|
(*myAttractor)();
|
|
|
|
}
|
|
|
|
|
|
|
|
void Attractor::plot() {
|
|
|
|
for ( std::vector<Projector *>::iterator it = projectors.begin(); it != projectors.end(); it++ ) {
|
|
|
|
const double * point = myAttractor->vector();
|
|
|
|
(*it)->plot(point);
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
IO & control
|
|
|
|
*/
|
|
|
|
void Attractor::output() {
|
|
|
|
const unsigned int dimension = myAttractor->getDimension();
|
|
|
|
const double * point = myAttractor->vector();
|
|
|
|
|
|
|
|
for ( unsigned int i = 0; i < dimension; i++ ) {
|
|
|
|
std::cout << point[i] << " ";
|
|
|
|
}
|
|
|
|
std::cout << std::endl;
|
|
|
|
}
|