diff options
Diffstat (limited to 'chaos/src/chaos_search.hpp')
| -rw-r--r-- | chaos/src/chaos_search.hpp | 190 |
1 files changed, 190 insertions, 0 deletions
diff --git a/chaos/src/chaos_search.hpp b/chaos/src/chaos_search.hpp new file mode 100644 index 0000000..fded6bc --- /dev/null +++ b/chaos/src/chaos_search.hpp @@ -0,0 +1,190 @@ +// +// +// chaos~ +// Copyright (C) 2005 Tim Blechmann +// +// This program 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. +// +// This program 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 this program; see the file COPYING. If not, write to +// the Free Software Foundation, Inc., 59 Temple Place - Suite 330, +// Boston, MA 02111-1307, USA. + +#include "chaos_base.hpp" + +#define MAXDIMENSION 5 + +template <class system> class chaos_search + : public flext_base +{ + FLEXT_HEADER(chaos_search, flext_base); + +public: + + /* local data for system, output and interpolation */ + system m_system; /* the system */ + + + data_t min[MAXDIMENSION]; /* minimal coordinates */ + data_t max[MAXDIMENSION]; /* maximal coordinates */ + data_t final[MAXDIMENSION]; /* initial coordinates for outlet */ + + data_t ly; /* lyapunov exponent */ + + int m_transient_steps; /* steps before starting the analysis */ + int m_asymptotic_steps; /* steps for the analysis */ + + + void get_tsteps(int &i) + { + i = m_transient_steps; + } + + void set_tsteps(int i) + { + if (i > 0) + m_transient_steps = i; + else + m_transient_steps = 0; + } + + void get_asteps(int &i) + { + i = m_asymptotic_steps; + } + + void set_asteps(int &i) + { + if (i > 0) + m_asymptotic_steps = i; + else + m_asymptotic_steps = 0; + } + + void print_results(void) + { + /* - send parameters to 1 + - send initial coordinates to 2 + - send minimal coordinates to 3 + - send lyapunov exponent to 4 + */ + + for (std::map<const t_symbol*,int>::iterator it = m_system.attr_ind.begin(); + it != m_system.attr_ind.end(); ++it) + { + post("key %s", it->first->s_name); + post("value %f", m_system.get_data(it->second)); + } + } + + void m_search(); + + FLEXT_CALLBACK(m_bang); + FLEXT_CALLVAR_I(get_tsteps, set_tsteps); + FLEXT_CALLVAR_I(get_asteps, set_asteps); + FLEXT_THREAD(m_search); +}; + + +/* create constructor / destructor */ + +#define CHAOS_SEARCH_INIT(SYSTEM, ATTRIBUTES) \ +FLEXT_HEADER(SYSTEM##_search, chaos_search<SYSTEM>) \ + \ +SYSTEM##_search(int argc, t_atom* argv ) \ +{ \ + int size = m_system.get_num_eq(); \ + \ + \ + m_asymptotic_steps = 10000; \ + m_transient_steps = 100; \ + \ + AddOutList("parameters"); \ + AddOutList("initial coordinates"); \ + AddOutList("minimal coordinates"); \ + AddOutList("maximal coordinates"); \ + AddOutFloat("lyapunov exponent"); \ + \ + FLEXT_ADDATTR_VAR("transient_steps", get_tsteps, set_tsteps); \ + FLEXT_ADDATTR_VAR("steps", get_asteps, set_asteps); \ + ATTRIBUTES; \ + FLEXT_ADDMETHOD_(0,"search", m_search); \ +} \ + \ + \ +FLEXT_ATTRVAR_I(m_transient_steps); \ +FLEXT_ATTRVAR_I(m_asymptotic_steps); + + + +template <class system> +void chaos_search<system>::m_search() +{ + int dimensions = m_system.get_num_eq(); + + ly = 0; + data_t diff_old = 0.1; + data_t last[MAXDIMENSION]; + + /* transient dynamics */ + for (int i = 0; i != m_transient_steps; ++i) + { + m_system.m_step(); + m_system.m_bash_denormals(); + m_system.m_verify(); + } + + for (int i = 0; i != dimensions; ++i) + { + last[i] = min[i] = max[i] = m_system.m_data[i]; + } + + /* now we start the analysis */ + + for (int i = 0; i != m_asymptotic_steps; ++i) + { + m_system.m_step(); + m_system.m_bash_denormals(); + m_system.m_verify(); + + data_t diff = 0; + for (int j = 0; j != dimensions; ++j) + { + /* update min/max */ + data_t datum = m_system.m_data[j]; + if (datum > max[j]) + max[j] = datum; + else if (datum < min[j]) + min[j] = datum; + + /* sum up diff */ + diff += (last[j] - datum) * (last[j] - datum); + + last[j] = datum; + } + diff = sqrt(diff); + + if (diff < 0) + diff = -diff; + + ly += log(diff / diff_old); + diff_old = diff; + + /* todo: maybe some overflow checking */ + if (diff == 0) + break; + } + + + ly /= m_asymptotic_steps; + + print_results(); +} |