From 013f13943f44a903ff30bda3c247cd5397f9b41e5c68a524ea1329ef3b09585f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Stefan=20Br=C3=BCns?= Date: Sun, 5 Jan 2025 17:22:33 +0000 Subject: [PATCH] - Fix some build issues with Octave 9.x, add * 0001-Fix-const-correctness-invalid-used-of-non-const-fort.patch * 0001-Fix-element-wise-plus-operator.patch OBS-URL: https://build.opensuse.org/package/show/science/octave-forge-tisean?expand=0&rev=5 --- .gitattributes | 23 + .gitignore | 1 + ...tness-invalid-used-of-non-const-fort.patch | 124 + 0001-Fix-element-wise-plus-operator.patch | 115 + octave-forge-tisean.changes | 19 + octave-forge-tisean.spec | 71 + tisean-0.2.3.tar.gz | 3 + tisean-drop-error_state-use.patch | 3725 +++++++++++++++++ 8 files changed, 4081 insertions(+) create mode 100644 .gitattributes create mode 100644 .gitignore create mode 100644 0001-Fix-const-correctness-invalid-used-of-non-const-fort.patch create mode 100644 0001-Fix-element-wise-plus-operator.patch create mode 100644 octave-forge-tisean.changes create mode 100644 octave-forge-tisean.spec create mode 100644 tisean-0.2.3.tar.gz create mode 100644 tisean-drop-error_state-use.patch diff --git a/.gitattributes b/.gitattributes new file mode 100644 index 0000000..9b03811 --- /dev/null +++ b/.gitattributes @@ -0,0 +1,23 @@ +## Default LFS +*.7z filter=lfs diff=lfs merge=lfs -text +*.bsp filter=lfs diff=lfs merge=lfs -text +*.bz2 filter=lfs diff=lfs merge=lfs -text +*.gem filter=lfs diff=lfs merge=lfs -text +*.gz filter=lfs diff=lfs merge=lfs -text +*.jar filter=lfs diff=lfs merge=lfs -text +*.lz filter=lfs diff=lfs merge=lfs -text +*.lzma filter=lfs diff=lfs merge=lfs -text +*.obscpio filter=lfs diff=lfs merge=lfs -text +*.oxt filter=lfs diff=lfs merge=lfs -text +*.pdf filter=lfs diff=lfs merge=lfs -text +*.png filter=lfs diff=lfs merge=lfs -text +*.rpm filter=lfs diff=lfs merge=lfs -text +*.tbz filter=lfs diff=lfs merge=lfs -text +*.tbz2 filter=lfs diff=lfs merge=lfs -text +*.tgz filter=lfs diff=lfs merge=lfs -text +*.ttf filter=lfs diff=lfs merge=lfs -text +*.txz filter=lfs diff=lfs merge=lfs -text +*.whl filter=lfs diff=lfs merge=lfs -text +*.xz filter=lfs diff=lfs merge=lfs -text +*.zip filter=lfs diff=lfs merge=lfs -text +*.zst filter=lfs diff=lfs merge=lfs -text diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..57affb6 --- /dev/null +++ b/.gitignore @@ -0,0 +1 @@ +.osc diff --git a/0001-Fix-const-correctness-invalid-used-of-non-const-fort.patch b/0001-Fix-const-correctness-invalid-used-of-non-const-fort.patch new file mode 100644 index 0000000..7ea3220 --- /dev/null +++ b/0001-Fix-const-correctness-invalid-used-of-non-const-fort.patch @@ -0,0 +1,124 @@ +From 1c672392dffcc5a9a124241ac966d9a81102dccf Mon Sep 17 00:00:00 2001 +From: =?UTF-8?q?Stefan=20Br=C3=BCns?= +Date: Wed, 1 Jan 2025 20:24:11 +0100 +Subject: [PATCH] Fix const-correctness, invalid used of non-const fortran_vec + +fortran_vec() is a non-const method, replace with the const data() method. +--- + src/__lzo_gm__.cc | 4 ++-- + src/__lzo_run__.cc | 2 +- + src/__lzo_test__.cc | 8 ++++---- + src/routines_c/find_multi_neighbors.cc | 2 +- + src/routines_c/tsa.h | 2 +- + 5 files changed, 9 insertions(+), 9 deletions(-) + mode change 100755 => 100644 src/__lzo_gm__.cc + mode change 100755 => 100644 src/__lzo_run__.cc + +diff --git a/src/__lzo_gm__.cc b/src/__lzo_gm__.cc +old mode 100755 +new mode 100644 +index 40b5065..c8b6476 +--- a/src/__lzo_gm__.cc ++++ b/src/__lzo_gm__.cc +@@ -42,7 +42,7 @@ void make_fit(const Matrix& series, octave_idx_type dim, + octave_idx_type LENGTH = series.rows (); + for (octave_idx_type i=0;i box (dim_vector(NMAX,NMAX)); + +diff --git a/src/__lzo_run__.cc b/src/__lzo_run__.cc +old mode 100755 +new mode 100644 +index 03ff253..14bca5c +--- a/src/__lzo_run__.cc ++++ b/src/__lzo_run__.cc +@@ -136,7 +136,7 @@ void make_zeroth(const Matrix &series, TISEAN_rand &generator, + for (octave_idx_type d=0;d &box, +- long *list,double **x, ++ long *list,const double **x, + octave_idx_type bs,octave_idx_type dim, + octave_idx_type emb,octave_idx_type del, + double eps, unsigned long *flist) +diff --git a/src/routines_c/tsa.h b/src/routines_c/tsa.h +index b044728..33d51b5 100644 +--- a/src/routines_c/tsa.h ++++ b/src/routines_c/tsa.h +@@ -66,7 +66,7 @@ extern octave_idx_type exclude_interval(octave_idx_type,long,long, + + extern octave_idx_type find_multi_neighbors(const Matrix &, + const MArray &, +- long *,double **, ++ long *,const double **, + octave_idx_type, octave_idx_type, + octave_idx_type, octave_idx_type, + double,unsigned long *); +-- +2.47.1 + diff --git a/0001-Fix-element-wise-plus-operator.patch b/0001-Fix-element-wise-plus-operator.patch new file mode 100644 index 0000000..49d610a --- /dev/null +++ b/0001-Fix-element-wise-plus-operator.patch @@ -0,0 +1,115 @@ +From 91eb92651920d306d7281b718d6c903e03e62691 Mon Sep 17 00:00:00 2001 +From: =?UTF-8?q?Stefan=20Br=C3=BCns?= +Date: Wed, 1 Jan 2025 22:14:59 +0100 +Subject: [PATCH] Fix element-wise plus operator + +The ".+" notation was deprecated in 7.x and removed with 9.x, it was +equivalent with the regular "+" operator. + +The operation broadcasts both index ranges, and then sums them +elementwise to create indices representing a sliding window. +--- + inst/av_d2.m | 2 +- + inst/c2g.m | 6 +++--- + inst/c2t.m | 4 ++-- + inst/endtoend.m | 2 +- + inst/upo.m | 2 +- + inst/upoembed.m | 4 ++-- + 6 files changed, 10 insertions(+), 10 deletions(-) + +diff --git a/inst/av_d2.m b/inst/av_d2.m +index 3c897cd..e44f021 100644 +--- a/inst/av_d2.m ++++ b/inst/av_d2.m +@@ -149,7 +149,7 @@ function output = av_d2 (d2_c2d_c1_out, varargin) + # Create smoothing functions + # Smoothes single column + smooth_column = @(col) sum (col((-aver:aver) ... +- .+(aver+1:(length(col)-aver)).'),2) ... ++ +(aver+1:(length(col)-aver)).'),2) ... + ./(2*aver + 1); + # Smoothes one cell + smooth_cell = @(cell_mat) [smooth_column(cell_mat{1}(:,1)), ... +diff --git a/inst/c2g.m b/inst/c2g.m +index c0dc02e..964a34d 100644 +--- a/inst/c2g.m ++++ b/inst/c2g.m +@@ -123,9 +123,9 @@ function output = c2g (d2_out) + + # Create column vectors instead of using loop + k_id = 1:length(tmp.c2)-1; +- f = exp((emat(k_id+1).*cmat(k_id).-emat(k_id).*cmat(k_id+1)) +- ./(emat(k_id+1).-emat(k_id))); +- d = (cmat(k_id+1).-cmat(k_id))./(emat(k_id+1).-emat(k_id)); ++ f = exp((emat(k_id+1).*cmat(k_id)-emat(k_id).*cmat(k_id+1)) ++ ./(emat(k_id+1)-emat(k_id))); ++ d = (cmat(k_id+1)-cmat(k_id))./(emat(k_id+1)-emat(k_id)); + a = emat(k_id); + b = emat(k_id+1); + +diff --git a/inst/c2t.m b/inst/c2t.m +index 4cb0d43..6f0b55d 100644 +--- a/inst/c2t.m ++++ b/inst/c2t.m +@@ -109,8 +109,8 @@ function output = c2t (d2_c1_out) + emat = log (corr_sums(idx,1)); + cmat = log (corr_sums(idx,2)); + +- b = (emat(2:end) .* cmat(1:end-1) .- emat(1:end-1) .* cmat(2:end)) ... +- ./ (emat(2:end) .- emat(1:end-1)); ++ b = (emat(2:end) .* cmat(1:end-1) - emat(1:end-1) .* cmat(2:end)) ... ++ ./ (emat(2:end) - emat(1:end-1)); + a = (cmat(2:end) - cmat(1:end-1)) ./ (emat(2:end) - emat(1:end-1)); + + cint = (exp (b) ./ a) ... +diff --git a/inst/endtoend.m b/inst/endtoend.m +index 2a9154d..e05a2ef 100644 +--- a/inst/endtoend.m ++++ b/inst/endtoend.m +@@ -209,7 +209,7 @@ S= reshape (S, [rows(S),1,columns(S)]); + endwhile + etot = original_col_S; + for nj = 0:(rows (S) - nmaxp) +- [x,s] = mismatch (S((1:nmaxp).+nj,1,:)); ++ [x,s] = mismatch (S((1:nmaxp)+nj,1,:)); + xj(1+nj) = x; + sj(1+nj) = s; + endfor +diff --git a/inst/upo.m b/inst/upo.m +index 4c17bf0..898fc46 100644 +--- a/inst/upo.m ++++ b/inst/upo.m +@@ -80,7 +80,7 @@ + ## n-the orbit you need to: + ## @example + ## +-## nth_orbit_data = orbit_data(sum(olens(1:n-1)).+(1:olens(n))); ++## nth_orbit_data = orbit_data(sum(olens(1:n-1))+(1:olens(n))); + ## + ## @end example + ## @item acc +diff --git a/inst/upoembed.m b/inst/upoembed.m +index c14b160..3143629 100644 +--- a/inst/upoembed.m ++++ b/inst/upoembed.m +@@ -126,7 +126,7 @@ values of olens"); + ## Create delay vectors + j = (m:-1:1); + delay_vec = @(x) orbit_data(sum(olens(1:idx(x)-1))+ ... +- mod(((1:olens(idx(x))+1).').-... ++ mod(((1:olens(idx(x))+1).')-... + (j-1).*delay -1 +m.*olens(idx(x)),... + olens(idx(x)))+1); + output = arrayfun (delay_vec,(1:length(idx)).','UniformOutput', false); +@@ -138,7 +138,7 @@ values of olens"); + # for k=1:olens(i)+1; + # for j=m:-1:1; + # output(sum(olens(1:i-1)+1)+k,m+1-j) = ... +- # orbit_data(sum(olens(1:i-1))+mod(k.-(j-1).*delay-1+m*olens(i),olens(i))+1); ++ # orbit_data(sum(olens(1:i-1))+mod(k-(j-1).*delay-1+m*olens(i),olens(i))+1); + # endfor + # endfor + # endfor +-- +2.47.1 + diff --git a/octave-forge-tisean.changes b/octave-forge-tisean.changes new file mode 100644 index 0000000..ca44fd2 --- /dev/null +++ b/octave-forge-tisean.changes @@ -0,0 +1,19 @@ +------------------------------------------------------------------- +Wed Jan 1 19:25:41 UTC 2025 - Stefan Brüns + +- Fix some build issues with Octave 9.x, add + * 0001-Fix-const-correctness-invalid-used-of-non-const-fort.patch + * 0001-Fix-element-wise-plus-operator.patch + +------------------------------------------------------------------- +Wed Apr 19 08:51:06 UTC 2023 - Atri Bhattacharya + +- Add tisean-drop-error_state-use.patch -- Drop the use of + error_state to support octave >= 8 + (https://savannah.gnu.org/bugs/index.php?61583). + +------------------------------------------------------------------- +Mon Aug 24 13:30:56 UTC 2015 - dmitry_r@opensuse.org + +- Initial package, version 0.2.3 + diff --git a/octave-forge-tisean.spec b/octave-forge-tisean.spec new file mode 100644 index 0000000..e2a6247 --- /dev/null +++ b/octave-forge-tisean.spec @@ -0,0 +1,71 @@ +# +# spec file for package octave-forge-tisean +# +# Copyright (c) 2025 SUSE LLC +# +# All modifications and additions to the file contributed by third parties +# remain the property of their copyright owners, unless otherwise agreed +# upon. The license for this file, and modifications and additions to the +# file, is the same license as for the pristine package itself (unless the +# license for the pristine package is not an Open Source License, in which +# case the license is the MIT License). An "Open Source License" is a +# license that conforms to the Open Source Definition (Version 1.9) +# published by the Open Source Initiative. + +# Please submit bugfixes or comments via https://bugs.opensuse.org/ +# + + +%define octpkg tisean +Name: octave-forge-%{octpkg} +Version: 0.2.3 +Release: 0 +Summary: Nonlinear Time Series Analysis +License: GPL-3.0-or-later +Group: Productivity/Scientific/Math +URL: https://gnu-octave.github.io/packages/tisean/ +Source0: https://downloads.sourceforge.net/project/octave/Octave%20Forge%20Packages/Individual%20Package%20Releases/%{octpkg}-%{version}.tar.gz +# PATCH-FIX-UPSTREAM tisean-drop-error_state-use.patch badshah400@gmail.com -- Drop the use of error_state to support octave >= 8 (https://savannah.gnu.org/bugs/index.php?61583) +Patch0: tisean-drop-error_state-use.patch +Patch1: 0001-Fix-const-correctness-invalid-used-of-non-const-fort.patch +Patch2: 0001-Fix-element-wise-plus-operator.patch +BuildRequires: gcc-c++ +BuildRequires: gcc-fortran +BuildRequires: octave-devel +Requires: octave-cli >= 4.0.0 +Requires: octave-forge-signal >= 1.3.0 + +%description +TISEAN stands for TIme SEries ANalysis. +This is part of Octave-Forge project. + +%prep +%setup -q -c %{name}-%{version} +pushd %{octpkg}-%{version} +%autopatch -p1 +# Fix missing namespace +find src/ -iname \*.cc -exec sed -i -e 's@set_warning_state\s*(@octave::\0@g' '{}' \; +popd +%octave_pkg_src + +%build +%octave_pkg_build + +%install +%octave_pkg_install + +%check +%octave_pkg_test + +%post +%octave --eval "pkg rebuild" + +%postun +%octave --eval "pkg rebuild" + +%files +%defattr(-,root,root) +%{octpackages_dir}/%{octpkg}-%{version} +%{octlib_dir}/%{octpkg}-%{version} + +%changelog diff --git a/tisean-0.2.3.tar.gz b/tisean-0.2.3.tar.gz new file mode 100644 index 0000000..03761a9 --- /dev/null +++ b/tisean-0.2.3.tar.gz @@ -0,0 +1,3 @@ +version https://git-lfs.github.com/spec/v1 +oid sha256:bbc94ea6af663092863a437f88215e84d3a59a60b9f51f6ffd9dbe90606a8259 +size 715494 diff --git a/tisean-drop-error_state-use.patch b/tisean-drop-error_state-use.patch new file mode 100644 index 0000000..1d16f28 --- /dev/null +++ b/tisean-drop-error_state-use.patch @@ -0,0 +1,3725 @@ +# HG changeset patch +# User Markus Mützel +# Date 1638189809 -3600 +# Mon Nov 29 13:43:29 2021 +0100 +# Node ID e40a599d68cf3061d04a9dac30e36751ed20acb2 +# Parent 71f2c8fde0c59f4942dd22da89bba162e8ae6c97 + +# Updated by Marius Schamschula +# Date 2023-03-10 + +Remove usage of `error_state` (bug #61583). + +* src/pcregexp.cc (many files): Remove usage of `error_state`. It was +unconditionally set to 0 since about 6 years ago and will finally be removed in +Octave 8. + +diff -r 71f2c8fde0c5 -r e40a599d68cf src/__boxcount__.cc +--- a/src/__boxcount__.cc Mon Aug 26 12:51:20 2019 -0400 ++++ b/src/__boxcount__.cc Mon Nov 29 13:43:29 2021 +0100 +@@ -194,89 +194,85 @@ + + octave_idx_type length=LENGTH-(maxembed-1)*DELAY; + +- if ( ! error_state) ++ // Calculate output ++ double heps = EPSMAX*EPSFAKTOR; ++ octave_idx_type epsi_old = 0; ++ for (octave_idx_type k=0;k histo (maxembed * dimension, 0.0); ++ start_box(series, which_dims, histo, maxembed, dimension, DELAY, ++ length, epsi, Q); + +- octave_idx_type epsi_test; +- do { +- heps /= EPSFAKTOR; +- epsi_test=(octave_idx_type)(1./heps); +- } while (epsi_test <= epsi_old); ++ if (Q != 1.0) ++ for (std::vector::iterator it = histo.begin (); ++ it != histo.end (); it++) ++ *it = log(*it)/(1.0-Q); ++ ++ histo_list.push_back (histo); ++ } ++ ++ // Create and assign output ++ dim_vector dv (maxembed, dimension); ++ string_vector keys; ++ keys.append (std::string("dim")); ++ keys.append (std::string("entropy")); ++ octave_map output (dv, keys); + +- octave_idx_type epsi = epsi_test; +- epsi_old = epsi; +- deps[k] = heps; ++ for (octave_idx_type i=0;i histo (maxembed * dimension, 0.0); +- start_box(series, which_dims, histo, maxembed, dimension, DELAY, +- length, epsi, Q); +- +- if (Q != 1.0) +- for (std::vector::iterator it = histo.begin (); +- it != histo.end (); it++) +- *it = log(*it)/(1.0-Q); +- +- histo_list.push_back (histo); ++ std::list>::const_iterator it_hist; ++ it_hist = histo_list.cbegin (); ++ for (octave_idx_type j=0;jhist[i],histo_el->hist[i]); ++ entropy_out(j,0) = deps[j]*maxinterval; ++ entropy_out(j,1) = (*it_hist)[i]; ++ entropy_out(j,2) = (*it_hist)[i]; ++ } ++ else ++ { ++ // old fprintf(fHq,"%e %e %e\n",deps[j]*maxinterval, ++ // histo_el->hist[i], ++ // histo_el->hist[i]-histo_el->hist[i-1]); ++ entropy_out(j,0) = deps[j]*maxinterval; ++ entropy_out(j,1) = (*it_hist)[i]; ++ entropy_out(j,2) = (*it_hist)[i] ++ - (*it_hist)[i-1]; ++ } ++ it_hist++; + } + +- // Create and assign output +- dim_vector dv (maxembed, dimension); +- string_vector keys; +- keys.append (std::string("dim")); +- keys.append (std::string("entropy")); +- octave_map output (dv, keys); +- +- for (octave_idx_type i=0;i>::const_iterator it_hist; +- it_hist = histo_list.cbegin (); +- for (octave_idx_type j=0;jhist[i],histo_el->hist[i]); +- entropy_out(j,0) = deps[j]*maxinterval; +- entropy_out(j,1) = (*it_hist)[i]; +- entropy_out(j,2) = (*it_hist)[i]; +- } +- else +- { +- // old fprintf(fHq,"%e %e %e\n",deps[j]*maxinterval, +- // histo_el->hist[i], +- // histo_el->hist[i]-histo_el->hist[i-1]); +- entropy_out(j,0) = deps[j]*maxinterval; +- entropy_out(j,1) = (*it_hist)[i]; +- entropy_out(j,2) = (*it_hist)[i] +- - (*it_hist)[i-1]; +- } +- it_hist++; +- } ++ output.assign (idx_vector(which_dims[i][1]), ++ idx_vector(which_dims[i][0]), ++ tmp); ++ } + +- tmp.setfield ("entropy",entropy_out); +- +- output.assign (idx_vector(which_dims[i][1]), +- idx_vector(which_dims[i][0]), +- tmp); +- } +- +- +- retval(0) = output; +- } ++ retval(0) = output; + } + return retval; + } +diff -r 71f2c8fde0c5 -r e40a599d68cf src/__c1__.cc +--- a/src//__c1__.cc.orig 2015-08-14 17:25:52.000000000 -0500 ++++ b/src//__c1__.cc 2023-03-09 18:58:35.000000000 -0600 +@@ -78,65 +78,61 @@ + octave_idx_type iverb = verbose; + + +- if (! error_state) +- { ++ octave_idx_type lines_read = input.rows (); //nmax in d1() ++ octave_idx_type columns_read = input.columns (); + +- octave_idx_type lines_read = input.rows (); //nmax in d1() +- octave_idx_type columns_read = input.columns (); + +- +- dim_vector dv (maxdim - mindim + 1, 1); +- string_vector keys; +- keys.append (std::string("dim")); +- keys.append (std::string("c1")); +- octave_map output (dv, keys); +- +- // Seed the rand() function for d1() +- F77_XFCN (rand, RAND, (sqrt(seed))); +- +- for (octave_idx_type m = mindim; m <= maxdim; m++) +- { +- octave_scalar_map tmp (keys); +- tmp.setfield ("dim", m); +- +- // Creat c1 output +- Matrix c1_out ((octave_idx_type) ((0 - log (1./(lines_read +- -(m-1) * delay)) + +- (log (2.) /resolution)) +- / (log (2.) /resolution)) +- , 2); +- +- double pr = 0.0; +- octave_idx_type current_row = 0; +- for (double pl = log (1./(lines_read - (m-1)*delay)); +- pl <= 0.0; pl += log (2.) / resolution) +- { +- double pln = pl; +- double rln; +- +- F77_XFCN (d1, D1, +- (lines_read, columns_read, lines_read, +- input.fortran_vec (), delay, m, cmin, +- pr, pln, rln, tmin, kmax, iverb)); +- +- if (pln != pr) +- { +- pr = pln; +- c1_out(current_row,0) = exp (rln); +- c1_out(current_row,1) = exp (pln); +- current_row += 1; +- } +- +- } +- // Resize output +- c1_out.resize (current_row, 2); +- tmp.setfield ("c1", c1_out); +- +- output.assign (idx_vector(m-mindim), tmp); +- } +- +- retval(0) = output; +- } +- } +- return retval; ++ dim_vector dv (maxdim - mindim + 1, 1); ++ string_vector keys; ++ keys.append (std::string("dim")); ++ keys.append (std::string("c1")); ++ octave_map output (dv, keys); ++ ++ // Seed the rand() function for d1() ++ F77_XFCN (rand, RAND, (sqrt(seed))); ++ ++ for (octave_idx_type m = mindim; m <= maxdim; m++) ++ { ++ octave_scalar_map tmp (keys); ++ tmp.setfield ("dim", m); ++ ++ // Creat c1 output ++ Matrix c1_out ((octave_idx_type) ((0 - log (1./(lines_read ++ -(m-1) * delay)) + ++ (log (2.) /resolution)) ++ / (log (2.) /resolution)) ++ , 2); ++ ++ double pr = 0.0; ++ octave_idx_type current_row = 0; ++ for (double pl = log (1./(lines_read - (m-1)*delay)); ++ pl <= 0.0; pl += log (2.) / resolution) ++ { ++ double pln = pl; ++ double rln; ++ ++ F77_XFCN (d1, D1, ++ (lines_read, columns_read, lines_read, ++ input.fortran_vec (), delay, m, cmin, ++ pr, pln, rln, tmin, kmax, iverb)); ++ ++ if (pln != pr) ++ { ++ pr = pln; ++ c1_out(current_row,0) = exp (rln); ++ c1_out(current_row,1) = exp (pln); ++ current_row += 1; ++ } ++ ++ } ++ // Resize output ++ c1_out.resize (current_row, 2); ++ tmp.setfield ("c1", c1_out); ++ ++ output.assign (idx_vector(m-mindim), tmp); ++ } ++ ++ retval(0) = output; ++ } ++ return retval; + } +diff -r 71f2c8fde0c5 -r e40a599d68cf src/__d2__.cc +--- a/src/__d2__.cc Mon Aug 26 12:51:20 2019 -0400 ++++ b/src/__d2__.cc Mon Nov 29 13:43:29 2021 +0100 +@@ -365,251 +365,247 @@ + time_t lasttime; + time(&lasttime); + +- if (! error_state) ++ bool imin_too_large = false; ++ bool pause_calc = false; ++ // Calculate the outputs ++ for (octave_idx_type n = 1 + counter; n < nmax && !imin_too_large ++ && !pause_calc; n++) + { +- +- bool imin_too_large = false; +- bool pause_calc = false; +- // Calculate the outputs +- for (octave_idx_type n = 1 + counter; n < nmax && !imin_too_large +- && !pause_calc; n++) ++ counter += 1; ++ bool smaller = 0; ++ octave_idx_type sn = scr[n-1]; ++ double epsinv = 1.0 / EPSMAX; ++ octave_idx_type x,y; ++ if (DIM > 1) + { +- counter += 1; +- bool smaller = 0; +- octave_idx_type sn = scr[n-1]; +- double epsinv = 1.0 / EPSMAX; +- octave_idx_type x,y; +- if (DIM > 1) +- { +- x=(octave_idx_type)(series[0][sn]*epsinv)&(NMAX-1); +- y=(octave_idx_type)(series[1][sn]*epsinv)&(NMAX-1); +- } +- else +- { +- x=(octave_idx_type)(series[0][sn]*epsinv)&(NMAX-1); +- y=(octave_idx_type)(series[0][sn+DELAY]*epsinv)&(NMAX-1); +- } +- list[sn]=box[x][y]; +- box[x][y]=sn; +- listc1[sn]=boxc1[x]; +- boxc1[x]=sn; ++ x=(octave_idx_type)(series[0][sn]*epsinv)&(NMAX-1); ++ y=(octave_idx_type)(series[1][sn]*epsinv)&(NMAX-1); ++ } ++ else ++ { ++ x=(octave_idx_type)(series[0][sn]*epsinv)&(NMAX-1); ++ y=(octave_idx_type)(series[0][sn+DELAY]*epsinv)&(NMAX-1); ++ } ++ list[sn]=box[x][y]; ++ box[x][y]=sn; ++ listc1[sn]=boxc1[x]; ++ boxc1[x]=sn; + +- octave_idx_type i=imin; +- while (found[maxembed][i] >= MAXFOUND) +- { +- smaller=1; +- if (++i > (HOWOFTEN-1)) +- break; +- } +- if (smaller) +- { +- imin=i; +- if (imin <= (HOWOFTEN-1)) +- { +- EPSMAX = epsm[imin]; +- double epsinv = 1.0/EPSMAX; +- for (octave_idx_type i1=0;i1 1) +- { +- x=(octave_idx_type)(series[0][sn]*epsinv) +- &(NMAX-1); +- y=(octave_idx_type)(series[1][sn]*epsinv) +- &(NMAX-1); +- } +- else +- { +- x=(octave_idx_type)(series[0][sn]*epsinv) +- &(NMAX-1); +- y=(octave_idx_type)(series[0][sn+DELAY]*epsinv) +- &(NMAX-1); +- } +- list[sn]=box[x][y]; +- box[x][y]=sn; +- listc1[sn]=boxc1[x]; +- boxc1[x]=sn; +- } +- } +- } +- ++ octave_idx_type i=imin; ++ while (found[maxembed][i] >= MAXFOUND) ++ { ++ smaller=1; ++ if (++i > (HOWOFTEN-1)) ++ break; ++ } ++ if (smaller) ++ { ++ imin=i; + if (imin <= (HOWOFTEN-1)) + { +- octave_idx_type lnorm=n; +- if (MINDIST > 0) ++ EPSMAX = epsm[imin]; ++ double epsinv = 1.0/EPSMAX; ++ for (octave_idx_type i1=0;i1=0)?sn-MINDIST:0; +- octave_idx_type n2=(sn+MINDIST 1) +- make_c2_dim(series, found, list, box, DIM, imin, EPSMAX1, +- EPSMAX, EPSMIN, EMBED, DELAY, MINDIST, +- HOWOFTEN, scr[n]); +- make_c2_1(series, found, listc1, boxc1, imin, EPSMAX1, +- EPSMAX, EPSMIN, MINDIST, HOWOFTEN, scr[n]); +- for (octave_idx_type i=imin;i WHEN) || (n == (nmax-1)) || +- (imin > (HOWOFTEN-1)) || (counter % it_pause == 0)) +- { +- time(&lasttime); +- +- if (imin > (HOWOFTEN-1)) ++ for (octave_idx_type i1=0;i1 1) ++ { ++ x=(octave_idx_type)(series[0][sn]*epsinv) ++ &(NMAX-1); ++ y=(octave_idx_type)(series[1][sn]*epsinv) ++ &(NMAX-1); ++ } ++ else ++ { ++ x=(octave_idx_type)(series[0][sn]*epsinv) ++ &(NMAX-1); ++ y=(octave_idx_type)(series[0][sn+DELAY]*epsinv) ++ &(NMAX-1); ++ } ++ list[sn]=box[x][y]; ++ box[x][y]=sn; ++ listc1[sn]=boxc1[x]; ++ boxc1[x]=sn; + } +- pause_calc = true; + } + } + +- // Create vars output +- octave_scalar_map vars; +- +- +- // Create vars output +- // old fprintf(fstat,"Center points treated so far= %ld\n",n); +- vars.setfield ("treated", counter); +- // old fprintf(fstat,"Maximum epsilon in the moment= %e\n", +- // epsm[imin]); +- vars.setfield ("eps", epsm[imin]); +- +- if (counter < nmax - 1 && imin_too_large == false) ++ if (imin <= (HOWOFTEN-1)) + { +- vars.setfield ("counter", counter); +- vars.setfield ("found", found_matrix); +- vars.setfield ("norm", norm_matrix); +- vars.setfield ("boxc1", boxc1_matrix); +- vars.setfield ("box", box_matrix); +- vars.setfield ("list", list_matrix); +- vars.setfield ("listc1", listc1_matrix); +- vars.setfield ("imin", imin); +- vars.setfield ("EPSMAX1",EPSMAX1); +- vars.setfield ("EPSMAX", EPSMAX); +- vars.setfield ("EPSMIN", EPSMIN); +- } +- +- // Create values output +- dim_vector dv (DIM * EMBED, 1); +- string_vector keys; +- keys.append (std::string("dim")); +- keys.append (std::string("c2")); +- keys.append (std::string("d2")); +- keys.append (std::string("h2")); +- octave_map values (dv, keys); +- +- for (octave_idx_type i=0;i 0) + { +- eps /= epsfactor; ++ octave_idx_type sn=scr[n]; ++ octave_idx_type n1=(sn-MINDIST>=0)?sn-MINDIST:0; ++ octave_idx_type n2=(sn+MINDIST 0) && (found[i][j] > 0.0) +- && (found[i][j-1] > 0.0)) +- { +- // old fprintf(fout,"%e %e\n",eps, +- // log(found[i][j-1]/found[i][j]/norm[j-1] +- // *norm[j]) +- // /log(epsfactor)); +- d2_out(d2_row,0) = eps; +- d2_out(d2_row,1) = log(found[i][j-1]/found[i][j] +- /norm[j-1]*norm[j]) +- /log(epsfactor); +- d2_row += 1; +- } +- +- // Calculate h2 output +- if (i < 1) +- { +- if (found[0][j] > 0.0) +- { +- // old fprintf(fout,"%e %e\n",eps, +- // -log(found[0][j]/norm[j])); +- h2_out(h2_row,0) = eps; +- h2_out(h2_row,1) = -log(found[0][j]/norm[j]); +- h2_row += 1; +- } +- } +- else +- { +- if ((found[i-1][j] > 0.0) && (found[i][j] > 0.0)) +- { +- // old fprintf(fout,"%e %e\n",eps, +- // log(found[i-1][j]/found[i][j])); +- h2_out(h2_row,0) = eps; +- h2_out(h2_row,1) = log(found[i-1][j] +- /found[i][j]); +- h2_row += 1; +- } +- } +- +- // Calculate c2 output +- if (norm[j] > 0.0) +- { +- // old fprintf(fout,"%e %e\n",eps, +- // found[i][j]/norm[j]); +- c2_out(c2_row,0) = eps; +- c2_out(c2_row,1) = found[i][j]/norm[j]; +- c2_row += 1; +- } +- } +- // Prepare d2 output +- d2_out.resize (d2_row, 2); +- tmp.setfield ("d2", d2_out); +- // Prepare h2 output +- h2_out.resize (h2_row, 2); +- tmp.setfield ("h2", h2_out); +- // Prepare c2 output +- c2_out.resize (c2_row, 2); +- tmp.setfield ("c2", c2_out); +- +- values.assign (idx_vector(i), tmp); ++ if (EMBED*DIM > 1) ++ make_c2_dim(series, found, list, box, DIM, imin, EPSMAX1, ++ EPSMAX, EPSMIN, EMBED, DELAY, MINDIST, ++ HOWOFTEN, scr[n]); ++ make_c2_1(series, found, listc1, boxc1, imin, EPSMAX1, ++ EPSMAX, EPSMIN, MINDIST, HOWOFTEN, scr[n]); ++ for (octave_idx_type i=imin;i WHEN) || (n == (nmax-1)) || ++ (imin > (HOWOFTEN-1)) || (counter % it_pause == 0)) ++ { ++ time(&lasttime); ++ ++ if (imin > (HOWOFTEN-1)) ++ { ++ // old exit(0); ++ imin_too_large = true; ++ } ++ pause_calc = true; ++ } ++ } ++ ++ // Create vars output ++ octave_scalar_map vars; ++ ++ ++ // Create vars output ++ // old fprintf(fstat,"Center points treated so far= %ld\n",n); ++ vars.setfield ("treated", counter); ++ // old fprintf(fstat,"Maximum epsilon in the moment= %e\n", ++ // epsm[imin]); ++ vars.setfield ("eps", epsm[imin]); ++ ++ if (counter < nmax - 1 && imin_too_large == false) ++ { ++ vars.setfield ("counter", counter); ++ vars.setfield ("found", found_matrix); ++ vars.setfield ("norm", norm_matrix); ++ vars.setfield ("boxc1", boxc1_matrix); ++ vars.setfield ("box", box_matrix); ++ vars.setfield ("list", list_matrix); ++ vars.setfield ("listc1", listc1_matrix); ++ vars.setfield ("imin", imin); ++ vars.setfield ("EPSMAX1",EPSMAX1); ++ vars.setfield ("EPSMAX", EPSMAX); ++ vars.setfield ("EPSMIN", EPSMIN); + } + ++ // Create values output ++ dim_vector dv (DIM * EMBED, 1); ++ string_vector keys; ++ keys.append (std::string("dim")); ++ keys.append (std::string("c2")); ++ keys.append (std::string("d2")); ++ keys.append (std::string("h2")); ++ octave_map values (dv, keys); ++ ++ for (octave_idx_type i=0;i 0) && (found[i][j] > 0.0) ++ && (found[i][j-1] > 0.0)) ++ { ++ // old fprintf(fout,"%e %e\n",eps, ++ // log(found[i][j-1]/found[i][j]/norm[j-1] ++ // *norm[j]) ++ // /log(epsfactor)); ++ d2_out(d2_row,0) = eps; ++ d2_out(d2_row,1) = log(found[i][j-1]/found[i][j] ++ /norm[j-1]*norm[j]) ++ /log(epsfactor); ++ d2_row += 1; ++ } ++ ++ // Calculate h2 output ++ if (i < 1) ++ { ++ if (found[0][j] > 0.0) ++ { ++ // old fprintf(fout,"%e %e\n",eps, ++ // -log(found[0][j]/norm[j])); ++ h2_out(h2_row,0) = eps; ++ h2_out(h2_row,1) = -log(found[0][j]/norm[j]); ++ h2_row += 1; ++ } ++ } ++ else ++ { ++ if ((found[i-1][j] > 0.0) && (found[i][j] > 0.0)) ++ { ++ // old fprintf(fout,"%e %e\n",eps, ++ // log(found[i-1][j]/found[i][j])); ++ h2_out(h2_row,0) = eps; ++ h2_out(h2_row,1) = log(found[i-1][j] ++ /found[i][j]); ++ h2_row += 1; ++ } ++ } ++ ++ // Calculate c2 output ++ if (norm[j] > 0.0) ++ { ++ // old fprintf(fout,"%e %e\n",eps, ++ // found[i][j]/norm[j]); ++ c2_out(c2_row,0) = eps; ++ c2_out(c2_row,1) = found[i][j]/norm[j]; ++ c2_row += 1; ++ } ++ } ++ // Prepare d2 output ++ d2_out.resize (d2_row, 2); ++ tmp.setfield ("d2", d2_out); ++ // Prepare h2 output ++ h2_out.resize (h2_row, 2); ++ tmp.setfield ("h2", h2_out); ++ // Prepare c2 output ++ c2_out.resize (c2_row, 2); ++ tmp.setfield ("c2", c2_out); ++ ++ values.assign (idx_vector(i), tmp); ++ } ++ ++ ++ // Assign outputs ++ retval(0) = values; ++ retval(1) = vars; + } ++ + return retval; + } +diff -r 71f2c8fde0c5 -r e40a599d68cf src/__delay__.cc +--- a/src/__delay__.cc Mon Aug 26 12:51:20 2019 -0400 ++++ b/src/__delay__.cc Mon Nov 29 13:43:29 2021 +0100 +@@ -60,56 +60,52 @@ + + OCTAVE_LOCAL_BUFFER (octave_idx_type, inddelay, alldim); + +- if (! error_state) +- { +- octave_idx_type rundel=0; +- octave_idx_type runmdel=0; +- +- unsigned int delsum; +- for (octave_idx_type i = 0; i < indim; i++) +- { +- delsum = 0; +- inddelay[rundel++] = delsum; +- +- for (octave_idx_type j = 1; j < formatdelay(i); j++) +- { +- delsum += delaylist(runmdel++); +- inddelay[rundel++] = delsum; +- } +- } +- +- octave_idx_type maxemb = 0; +- for (octave_idx_type i = 0; i < alldim; i++) +- maxemb = (maxemb < inddelay[i])? inddelay[i] : maxemb; ++ octave_idx_type rundel=0; ++ octave_idx_type runmdel=0; + +- octave_idx_type outdim = 0; +- for (octave_idx_type i = 0; i < indim; i++) +- outdim += formatdelay(i); +- +- octave_idx_type out_rows = (length > maxemb) ? length - maxemb : 0; +- +- Matrix series (out_rows, outdim); +- unsigned int embsum; +- for (octave_idx_type i = maxemb; i < length; i++) +- { +- rundel = 0; +- embsum = 0; ++ unsigned int delsum; ++ for (octave_idx_type i = 0; i < indim; i++) ++ { ++ delsum = 0; ++ inddelay[rundel++] = delsum; + +- for (octave_idx_type j = 0; j < indim; j++) +- { +- octave_idx_type emb = formatdelay(j); +- +- for (octave_idx_type k = 0; k < emb; k++) +- series(i-maxemb, embsum+k) = data(i-inddelay[rundel++], j); +- +- // previously fprintf(stdout,"%e ",series[j][i-inddelay[rundel++]]); +- embsum += emb; +- } +- // previously fprintf(stdout,"\n"); ++ for (octave_idx_type j = 1; j < formatdelay(i); j++) ++ { ++ delsum += delaylist(runmdel++); ++ inddelay[rundel++] = delsum; + } +- retval(0) = series; + } + ++ octave_idx_type maxemb = 0; ++ for (octave_idx_type i = 0; i < alldim; i++) ++ maxemb = (maxemb < inddelay[i])? inddelay[i] : maxemb; ++ ++ octave_idx_type outdim = 0; ++ for (octave_idx_type i = 0; i < indim; i++) ++ outdim += formatdelay(i); ++ ++ octave_idx_type out_rows = (length > maxemb) ? length - maxemb : 0; ++ ++ Matrix series (out_rows, outdim); ++ unsigned int embsum; ++ for (octave_idx_type i = maxemb; i < length; i++) ++ { ++ rundel = 0; ++ embsum = 0; ++ ++ for (octave_idx_type j = 0; j < indim; j++) ++ { ++ octave_idx_type emb = formatdelay(j); ++ ++ for (octave_idx_type k = 0; k < emb; k++) ++ series(i-maxemb, embsum+k) = data(i-inddelay[rundel++], j); ++ ++ // previously fprintf(stdout,"%e ",series[j][i-inddelay[rundel++]]); ++ embsum += emb; ++ } ++ // previously fprintf(stdout,"\n"); ++ } ++ retval(0) = series; + } + return retval; + } +diff -r 71f2c8fde0c5 -r e40a599d68cf src/__false_nearest__.cc +--- a/src/__false_nearest__.cc Mon Aug 26 12:51:20 2019 -0400 ++++ b/src/__false_nearest__.cc Mon Nov 29 13:43:29 2021 +0100 +@@ -179,78 +179,75 @@ + check_alloc(vcomp=(unsigned int*)malloc(sizeof(int)*(maxdim))); + check_alloc(vemb=(unsigned int*)malloc(sizeof(int)*(maxdim))); + +- if ( ! error_state) +- { +- for (i=0;i= minn) { ++ make_correction(input,n,nfound); ++ ok[n]=epscount; ++ if (epscount == 1) ++ resize_eps=1; ++ allfound++; + } +- epsilon=mineps; +- all_done=0; +- epscount=1; +- allfound=0; +- if (verbosity) +- octave_stdout << "Starting iteration " << iter << "\n"; +- while(!all_done) { +- mmb(input, epsilon); +- all_done=1; +- for (n=emb_offset;n= minn) { +- make_correction(input,n,nfound); +- ok[n]=epscount; +- if (epscount == 1) +- resize_eps=1; +- allfound++; ++ else ++ all_done=0; ++ } ++ if (verbosity) ++ octave_stdout << "Corrected " << allfound << " points with epsilon= " ++ << epsilon*d_max_max << "\n"; ++ if (std::isinf(epsilon*d_max_max)) ++ { ++ error_with_id ("Octave:tisean", "cannot reduce noise on input data"); ++ return retval; ++ } ++ epsilon *= epsfac; ++ epscount++; ++ } ++ if (verbosity) ++ octave_stdout << "Start evaluating the trend\n"; ++ ++ epsilon=mineps; ++ allfound=0; ++ for (i=1;i 2*(dim*embed+1)) +- { +- make_fit (series, found, error_array, +- i,dim, embed, delay, STEP, actfound); +- // Checking if the fit was correct +- // If any errors were raised: end function +- if (error_state) +- { +- return retval; +- } +- pfound++; +- avfound += (double)(actfound-1); +- for (octave_idx_type j=0;j 1) +- { +- double sumerror=0.0; +- for (octave_idx_type j=0;j 2*(dim*embed+1)) ++ { ++ make_fit (series, found, error_array, ++ i,dim, embed, delay, STEP, actfound); ++ pfound++; ++ avfound += (double)(actfound-1); ++ for (octave_idx_type j=0;j 1) ++ { ++ double sumerror=0.0; ++ for (octave_idx_type j=0;j input_max(dim-1)) +- ? input_max(0) : input_max(dim-1); +- maximum_epsilon *= EPSF; ++ // Calculate the maximum epsilon that makes sense ++ // On the basis of 'i' and 'j' from put_in_boxes () ++ NDArray input_max = input.max (); ++ double maximum_epsilon = (input_max(0) > input_max(dim-1)) ++ ? input_max(0) : input_max(dim-1); ++ maximum_epsilon *= EPSF; + +- // Create output +- Matrix output (FLENGTH, dim); +- for (octave_idx_type i=0;i maximum_epsilon) + { +- // If epsilon became too large +- // there is no sense in continuing +- if (epsilon > maximum_epsilon) ++ error_with_id ("Octave:tisean", "The neighbourhood size" ++ " became too large during search," ++ " no sense continuing"); ++ return retval; ++ } ++ ++ epsilon*=EPSF; ++ put_in_boxes(series, LENGTH, list, box, epsilon, dim, embed, ++ DELAY); ++ octave_idx_type actfound; ++ actfound=hfind_neighbors (series, cast, found, list, box, ++ epsilon, dim, embed, DELAY); ++ if (actfound >= MINN) ++ { ++ if (!do_zeroth) ++ make_fit(series, cast, found, dim, embed, DELAY, ++ actfound, newcast); ++ else ++ make_zeroth(series, found, dim, actfound,newcast); ++ ++ for (octave_idx_type j=0;j= MINN) ++ done=1; ++ for (octave_idx_type j=0;j 2.0) || (newcast[j] < -1.0)) + { ++ error_with_id("Octave:tisean","forecast failed, " ++ "escaping data region"); + return retval; + } +- +- for (octave_idx_type j=0;j 2.0) || (newcast[j] < -1.0)) +- { +- error_with_id("Octave:tisean","forecast failed, " +- "escaping data region"); +- return retval; +- } +- } +- double *swap=cast[0]; +- for (octave_idx_type j=0;j input_max(dim-1)) ++ ? input_max(0) : input_max(dim-1); ++ maximum_epsilon *= 1.1; ++ ++ // Calculate lyapunov exponents ++ for (double eps=eps0;!alldone;eps*=1.1) + { + +- // Calculate the maximum epsilon that makes sense +- // On the basis of 'i' and 'j' from put_in_boxes () +- NDArray input_max = input.max (); +- double maximum_epsilon = (input_max(0) > input_max(dim-1)) +- ? input_max(0) : input_max(dim-1); +- maximum_epsilon *= 1.1; +- +- // Calculate lyapunov exponents +- for (double eps=eps0;!alldone;eps*=1.1) ++ // If epsilon became too large ++ // there is no sense in continuing ++ if (eps > maximum_epsilon) + { +- +- // If epsilon became too large +- // there is no sense in continuing +- if (eps > maximum_epsilon) +- { +- error_with_id ("Octave:tisean", "The neighbourhood size" +- " became too large during search," +- " no sense continuing"); +- return retval; +- } +- +- put_in_boxes(series, box, list, eps, length, dim, delay, steps); +- alldone=1; +- for (octave_idx_type n=0;n<=maxlength;n++) +- { +- if (!done[n]) +- done[n]=make_iterate(series, box, list, found, lyap, eps, +- length, dim, delay, steps, mindist, +- n); +- alldone &= done[n]; +- } +- if (verbose) +- printf("epsilon: %e already found: %ld\n",eps*max_val, +- found[0]); ++ error_with_id ("Octave:tisean", "The neighbourhood size" ++ " became too large during search," ++ " no sense continuing"); ++ return retval; + } + +- // Create output +- Matrix output (steps+1, 2); +- octave_idx_type count = 0; +- for (octave_idx_type i=0;i<=steps;i++) +- if (found[i]) +- { +- // old fprintf(file,"%d %e\n",i,lyap[i]/found[i]/2.0); +- output(count,0) = i; +- output(count,1) = lyap[i]/found[i]/2.0; +- count += 1; +- } ++ put_in_boxes(series, box, list, eps, length, dim, delay, steps); ++ alldone=1; ++ for (octave_idx_type n=0;n<=maxlength;n++) ++ { ++ if (!done[n]) ++ done[n]=make_iterate(series, box, list, found, lyap, eps, ++ length, dim, delay, steps, mindist, ++ n); ++ alldone &= done[n]; ++ } ++ if (verbose) ++ printf("epsilon: %e already found: %ld\n",eps*max_val, ++ found[0]); ++ } + +- // Resize output to match number of found points +- if (count < output.numel ()) +- output.resize (count, 2); ++ // Create output ++ Matrix output (steps+1, 2); ++ octave_idx_type count = 0; ++ for (octave_idx_type i=0;i<=steps;i++) ++ if (found[i]) ++ { ++ // old fprintf(file,"%d %e\n",i,lyap[i]/found[i]/2.0); ++ output(count,0) = i; ++ output(count,1) = lyap[i]/found[i]/2.0; ++ count += 1; ++ } + +- // Assign output +- retval(0) = output; +- } ++ // Resize output to match number of found points ++ if (count < output.numel ()) ++ output.resize (count, 2); ++ ++ // Assign output ++ retval(0) = output; + } + return retval; + } +diff -r 71f2c8fde0c5 -r e40a599d68cf src/__lyap_spec__.cc +--- a/src/__lyap_spec__.cc Mon Aug 26 12:51:20 2019 -0400 ++++ b/src/__lyap_spec__.cc Mon Nov 29 13:43:29 2021 +0100 +@@ -220,13 +220,6 @@ + Matrix solved = mat.solve (vec); + double *solved_arr = solved.fortran_vec (); + +- // If errors were raised (a singular matrix was encountered), +- // there is no sense in countinuing +- if (error_state) +- { +- return ; +- } +- + double new_vec = solved_arr[0]; + for (octave_idx_type i=1;i<=alldim;i++) + dynamics[d][i-1] = solved_arr[i]; +@@ -419,149 +412,138 @@ + } + // end old indexes = make_multi_index(); + +- if (!error_state) +- { +- +- // Promote warnings connected with singular matrixes to errors +- set_warning_state ("Octave:nearly-singular-matrix","error"); +- set_warning_state ("Octave:singular-matrix","error"); ++ // Promote warnings connected with singular matrixes to errors ++ set_warning_state ("Octave:nearly-singular-matrix","error"); ++ set_warning_state ("Octave:singular-matrix","error"); + +- // Prepare data for running first time +- if (count == 0) +- { +- averr_matrix.fill (0.0); ++ // Prepare data for running first time ++ if (count == 0) ++ { ++ averr_matrix.fill (0.0); + +- if (epsset) +- epsmin /= maxinterval; ++ if (epsset) ++ epsmin /= maxinterval; + +- // old rnd_init(0x098342L); +- TISEAN_rand generator (0x098342L); +- for (octave_idx_type i=0;i<10000;i++) +- generator.rnd_long(); +- for (octave_idx_type i=0;i +- ::max (); +- } +- gram_schmidt(alldim, delta,lfactor); ++ // old rnd_init(0x098342L); ++ TISEAN_rand generator (0x098342L); ++ for (octave_idx_type i=0;i<10000;i++) ++ generator.rnd_long(); ++ for (octave_idx_type i=0;i ++ ::max (); ++ } ++ gram_schmidt(alldim, delta,lfactor); + +- avneig = 0.0; +- aveps = 0.0; +- } ++ avneig = 0.0; ++ aveps = 0.0; ++ } + +- // Create output +- Matrix lyap_exp (1, 1 + alldim); +- time_t lasttime; +- time(&lasttime); +- bool pause_calc = false; +- for (octave_idx_type i = count + (EMBED-1) * DELAY; +- i < ITERATIONS && !pause_calc; i++) ++ // Create output ++ Matrix lyap_exp (1, 1 + alldim); ++ time_t lasttime; ++ time(&lasttime); ++ bool pause_calc = false; ++ for (octave_idx_type i = count + (EMBED-1) * DELAY; ++ i < ITERATIONS && !pause_calc; i++) ++ { ++ count++; ++ make_dynamics(series, box, indexes, epsmin, epsset, EPSSTEP, ++ EMBED, MINNEIGHBORS, LENGTH, DIMENSION, count, ++ avneig, aveps, dynamics, averr, i); ++ make_iteration(DIMENSION, alldim, dynamics, delta); ++ gram_schmidt(alldim, delta,lfactor); ++ for (octave_idx_type j=0;j OUT) || (i == (ITERATIONS-1)) ++ || (count % it_pause == 0) ) + { +- count++; +- make_dynamics(series, box, indexes, epsmin, epsset, EPSSTEP, +- EMBED, MINNEIGHBORS, LENGTH, DIMENSION, count, +- avneig, aveps, dynamics, averr, i); +- // If there was an error +- // (matrix singularity or not enough neighbors) +- // No sense continuing +- if (error_state) ++ time(&lasttime); ++ ++ // Create spectrum output ++ // old fprintf(stdout,"%ld ",count); ++ lyap_exp(0,0) = count; ++ for (octave_idx_type j=0;j OUT) || (i == (ITERATIONS-1)) +- || (count % it_pause == 0) ) +- { +- time(&lasttime); +- +- // Create spectrum output +- // old fprintf(stdout,"%ld ",count); +- lyap_exp(0,0) = count; +- for (octave_idx_type j=0;j box (dim_vector(NMAX,NMAX)); + +- if ( ! error_state) +- { +- // Estimate maximum possible output size +- octave_idx_type output_rows = (octave_idx_type) +- ((log(EPS1) - log(EPS0)) / log (EPSF)); +- output_rows += 2; ++ // Estimate maximum possible output size ++ octave_idx_type output_rows = (octave_idx_type) ++ ((log(EPS1) - log(EPS0)) / log (EPSF)); ++ output_rows += 2; + +- // Create output +- Matrix output (output_rows,dim+4); +- octave_idx_type count = 0; ++ // Create output ++ Matrix output (output_rows,dim+4); ++ octave_idx_type count = 0; + +- for (double epsilon=EPS0;epsilon 2*(dim*embed+1)) +- { +- make_fit (input, dim, i, actfound, STEP, found, +- error_array); +- pfound++; +- avfound += (double)(actfound-1); +- for (octave_idx_type j=0;j 1) ++ octave_idx_type actfound; ++ actfound = find_multi_neighbors(input,box,list,hser, ++ NMAX,dim,embed,delay, ++ epsilon,hfound); ++ actfound = exclude_interval(actfound,i-causal+1, ++ i+causal+(embed-1)*delay-1, ++ hfound,found); ++ if (actfound > 2*(dim*embed+1)) + { +- double sumerror=0.0; +- for (octave_idx_type j=0;j 1) ++ { ++ double sumerror=0.0; ++ for (octave_idx_type j=0;j box (dim_vector(NMAX,NMAX)); + +- if ( ! error_state) +- { +- for (octave_idx_type j=0;j indexes (dim_vector (alldim, 2)); +- for (octave_idx_type i=0;i= MINN) { +- if (setsort) { +- epsilon0 += epsilon; +- count++; +- sort(input, found, abstand, cast, actfound, (embed-1)*DELAY); +- actfound=MINN; +- } +- make_zeroth(input, generator, setnoise, var, Q, found, +- actfound, newcast); +- +- for (octave_idx_type j=0;j indexes (dim_vector (alldim, 2)); ++ for (octave_idx_type i=0;i= MINN) { ++ if (setsort) { ++ epsilon0 += epsilon; ++ count++; ++ sort(input, found, abstand, cast, actfound, (embed-1)*DELAY); ++ actfound=MINN; ++ } ++ make_zeroth(input, generator, setnoise, var, Q, found, ++ actfound, newcast); ++ ++ for (octave_idx_type j=0;j box (dim_vector(NMAX, NMAX)); + + // Compute forecast error +- if ( ! error_state) +- { +- for (octave_idx_type i=0;i= MINN) ++ { ++ if (setsort) ++ { ++ sort(input, found, abstand, embed, DELAY, MINN, ++ actfound, hser); ++ actfound=MINN; ++ } ++ for (octave_idx_type j=1;j<=STEP;j++) { ++ make_fit(input,dim,hi,actfound,j,found,error_array,diffs); ++ } ++ done[i]=1; ++ } ++ alldone &= done[i]; ++ } ++ } + +- // Compute estimates +- octave_idx_type actfound; +- octave_idx_type hi; +- while (!alldone) { +- alldone=1; +- epsilon*=EPSF; +- make_multi_box(input,box,list,LENGTH-(long)STEP,NMAX,dim, +- embed,DELAY,epsilon); +- for (octave_idx_type i=(embed-1)*DELAY;i 1) ++ { ++ ind_forecast_err.resize(clength - (embed-1)*DELAY, dim); ++ for (octave_idx_type i=(embed-1)*DELAY;i= MINN) +- { +- if (setsort) +- { +- sort(input, found, abstand, embed, DELAY, MINN, +- actfound, hser); +- actfound=MINN; +- } +- for (octave_idx_type j=1;j<=STEP;j++) { +- make_fit(input,dim,hi,actfound,j,found,error_array,diffs); +- } +- done[i]=1; +- } +- alldone &= done[i]; +- } +- } +- +- // Create relative forecast error output +- Matrix rel_forecast_err (STEP, dim + 1); +- for (octave_idx_type i=0;i 1) +- { +- ind_forecast_err.resize(clength - (embed-1)*DELAY, dim); +- for (octave_idx_type i=(embed-1)*DELAY;i 0) +- make_cast(forecast, series, coding, results, LENGTH, CLENGTH, N, +- DIM, DELAY, maxencode, pars, std_dev); ++ // Create forecast ++ NDArray forecast (dim_vector (0,0)); ++ if (CLENGTH > 0) ++ make_cast(forecast, series, coding, results, LENGTH, CLENGTH, N, ++ DIM, DELAY, maxencode, pars, std_dev); + +- // Assign outputs +- retval(0) = free_par; +- retval(1) = fit_norm; +- retval(2) = coeffs; +- retval(3) = sample_err; +- retval(4) = forecast; +- } ++ // Assign outputs ++ retval(0) = free_par; ++ retval(1) = fit_norm; ++ retval(2) = coeffs; ++ retval(3) = sample_err; ++ retval(4) = forecast; + } + return retval; + } +diff -r 71f2c8fde0c5 -r e40a599d68cf src/__rbf__.cc +--- a/src/__rbf__.cc Mon Aug 26 12:51:20 2019 -0400 ++++ b/src/__rbf__.cc Mon Nov 29 13:43:29 2021 +0100 +@@ -152,169 +152,161 @@ + for (octave_idx_type j=0;j= MINN) { ++ for (octave_idx_type j=1;j<=STEP;j++) ++ error_array[j-1] += make_fit (series1, series2, found, ++ i,actfound,j); ++ done[i]=1; ++ } ++ alldone &= done[i]; ++ } ++ } ++ ++ // Create output ++ Matrix output (STEP,1); ++ for (octave_idx_type i=0;i= MINN) { +- for (octave_idx_type j=1;j<=STEP;j++) +- error_array[j-1] += make_fit (series1, series2, found, +- i,actfound,j); +- done[i]=1; +- } +- alldone &= done[i]; +- } +- } +- +- // Create output +- Matrix output (STEP,1); +- for (octave_idx_type i=0;i 1)) +- { +- transposed = 1; +- in_out1 = in_out1.transpose(); +- } +- +- int lines_read = in_out1.numel(); +- NDArray in_out2 (Matrix (lines_read, 1)); +- +- F77_XFCN (ts_lazy, TS_LAZY, +- (m, rv, imax, lines_read, +- in_out1.fortran_vec(), in_out2.fortran_vec())); +- +- // Transpose the output to resemble the input +- if (transposed) +- { +- in_out1 = in_out1.transpose(); +- in_out2 = in_out2.transpose(); +- } +- +- retval(0) = in_out1; +- retval(1) = in_out2; +- } +- } +- return retval; ++ // If vector is in 1 row: transpose (we will transpose the output to fit) ++ transposed = 0; ++ ++ if ((rows == 1) && (cols > 1)) ++ { ++ transposed = 1; ++ in_out1 = in_out1.transpose(); ++ } ++ ++ int lines_read = in_out1.numel(); ++ NDArray in_out2 (Matrix (lines_read, 1)); ++ ++ F77_XFCN (ts_lazy, TS_LAZY, ++ (m, rv, imax, lines_read, ++ in_out1.fortran_vec(), in_out2.fortran_vec())); ++ ++ // Transpose the output to resemble the input ++ if (transposed) ++ { ++ in_out1 = in_out1.transpose(); ++ in_out2 = in_out2.transpose(); ++ } ++ ++ retval(0) = in_out1; ++ retval(1) = in_out2; ++ } ++ return retval; + } + + /* +diff -r 71f2c8fde0c5 -r e40a599d68cf src/mutual.cc +--- a/src/mutual.cc Mon Aug 26 12:51:20 2019 -0400 ++++ b/src/mutual.cc Mon Nov 29 13:43:29 2021 +0100 +@@ -193,60 +193,57 @@ + int32NDArray h2 (dim_vector(partitions, partitions)); + OCTAVE_LOCAL_BUFFER (long, array, length); + +- if (! error_state) +- { +- // Load array ++ // Load array + +- // Rescale data and load array +- // NOTE: currently supports only vectors so (dim == 1) always +- if (dim == 1){ +- double mint, interval; +- rescale_data(input,0,length,&mint,&interval); ++ // Rescale data and load array ++ // NOTE: currently supports only vectors so (dim == 1) always ++ if (dim == 1){ ++ double mint, interval; ++ rescale_data(input,0,length,&mint,&interval); + +- for (octave_idx_type i=0;i= (long)length) +- corrlength=length-1; ++ double shannon = make_cond_entropy (0, h1, h11, h2, array, length, ++ partitions); ++ if (corrlength >= (long)length) ++ corrlength=length-1; + +- // Construct the output +- Matrix delay (corrlength+1,1); +- // To save memory +- int minf_len = 1; ++ // Construct the output ++ Matrix delay (corrlength+1,1); ++ // To save memory ++ int minf_len = 1; + +- if (nargout > 1) +- minf_len = corrlength+1; +- Matrix mutual_inf (minf_len,1); ++ if (nargout > 1) ++ minf_len = corrlength+1; ++ Matrix mutual_inf (minf_len,1); + +- // Assign output +- delay(0,0) = 0; +- mutual_inf(0,0) = shannon; +- for (octave_idx_type tau=1;tau<=corrlength;tau++) { ++ // Assign output ++ delay(0,0) = 0; ++ mutual_inf(0,0) = shannon; ++ for (octave_idx_type tau=1;tau<=corrlength;tau++) { + +- // fprintf(stdout,"%ld %e %e\n",tau,condent,condent/log((double)partitions)); +- delay(tau,0) = tau; +- if (nargout > 1) +- { +- mutual_inf(tau,0) = make_cond_entropy(tau, h1, h11, h2, array, +- length, partitions); +- } ++ // fprintf(stdout,"%ld %e %e\n",tau,condent,condent/log((double)partitions)); ++ delay(tau,0) = tau; ++ if (nargout > 1) ++ { ++ mutual_inf(tau,0) = make_cond_entropy(tau, h1, h11, h2, array, ++ length, partitions); + } ++ } + +- if (transposed) +- { +- delay = delay.transpose(); +- if (nargout > 1) +- mutual_inf = mutual_inf.transpose(); +- } +- retval(0) = delay; +- retval(1) = mutual_inf; ++ if (transposed) ++ { ++ delay = delay.transpose(); ++ if (nargout > 1) ++ mutual_inf = mutual_inf.transpose(); + } ++ retval(0) = delay; ++ retval(1) = mutual_inf; + } + return retval; + } +--- a/src/__lfo_test__.cc.orig 2015-08-14 17:25:52.000000000 -0500 ++++ b/src/__lfo_test__.cc 2023-03-09 19:39:33.000000000 -0600 +@@ -212,13 +212,6 @@ + Matrix solved_vec = mat.solve (vec); + double *solved_vec_arr = solved_vec.fortran_vec (); + +- // If errors were raised (a singular matrix was encountered), +- // there is no sense in countinueing +- if (error_state) +- { +- return ; +- } +- + newcast[i]=foreav[i]; + for (octave_idx_type j=0;j input_max(COMP-1)) +- ? input_max(0) : input_max(COMP-1); +- maximum_epsilon *= EPSF; +- +- // Calculate output +- bool alldone=0; +- while (!alldone) { +- alldone=1; +- +- // If epsilon became too large there is no sense in continuing +- if (epsilon > maximum_epsilon) +- { +- error_with_id ("Octave:tisean", "The neighbourhood size became" +- " too large during search, no sense" +- " continuing"); +- return retval; +- } +- epsilon*=EPSF; +- put_in_boxes(series, LENGTH, STEP, COMP, hdim, epsilon, list, box); +- for (octave_idx_type i=(EMBED-1)*DELAY;i MINN) +- { +- make_fit(series, indexes, found, STEP, DIM, COMP, +- actfound,i,newcast); +- // Checking if the fit was correct +- // If any errors were raised: end function +- if (error_state) +- { +- return retval; +- } +- for (octave_idx_type j=0;j input_max(COMP-1)) ++ ? input_max(0) : input_max(COMP-1); ++ maximum_epsilon *= EPSF; ++ ++ // Calculate output ++ bool alldone=0; ++ while (!alldone) { ++ alldone=1; ++ ++ // If epsilon became too large there is no sense in continuing ++ if (epsilon > maximum_epsilon) ++ { ++ error_with_id ("Octave:tisean", "The neighbourhood size became" ++ " too large during search, no sense" ++ " continuing"); ++ return retval; ++ } ++ epsilon*=EPSF; ++ put_in_boxes(series, LENGTH, STEP, COMP, hdim, epsilon, list, box); ++ for (octave_idx_type i=(EMBED-1)*DELAY;i MINN) ++ { ++ make_fit(series, indexes, found, STEP, DIM, COMP, ++ actfound,i,newcast); ++ for (octave_idx_type j=0;j